1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
|
---
# Near to Far Field Spectra
---
We demonstrate Meep's [near-to-far-field transformation](../Scheme_User_Interface.md#near-to-far-field-spectra) feature using two examples. There are three steps to using the near-to-far-field feature. First, we need to define the "near" surface(s) as a set of surfaces capturing *all* outgoing radiation in the desired direction(s). Second, we run the simulation using a pulsed source (or possibly, the frequency-domain solver) to allow Meep to accumulate the Fourier transforms on the near surface(s). Third, we have Meep compute the far fields at any desired points with the option to save the far fields to an HDF5 file.
[TOC]
### Radiation Pattern of an Antenna
In this example, we compute the [radiation pattern](https://en.wikipedia.org/wiki/Radiation_pattern) of an antenna. This involves an electric-current point dipole source as the emitter in vacuum. We will compute the radiation pattern for three different polarizations of the input source. The source is placed in the middle of the 2d cell which is surrounded by PMLs. The near fields are obtained on a bounding box positioned just outside of the PML. The far fields are computed at equally-spaced points along the circumference of a circle having a radius many (i.e., 1000) times larger than the source wavelength and lying outside of the cell. The simulation geometry is shown in the following schematic.
<center>
</center>
The simulation script is in [examples/antenna-radiation.ctl](https://github.com/stevengj/meep/blob/master/scheme/examples/antenna-radiation.ctl).
```scm
(set-param! resolution 50)
(define-param sxy 4)
(define-param dpml 1)
(set! geometry-lattice (make lattice (size (+ sxy (* 2 dpml)) (+ sxy (* 2 dpml)) no-size)))
(set! pml-layers (list (make pml (thickness dpml))))
(define-param fcen 1.0)
(define-param df 0.4)
(define-param src-cmpt Ez)
(set! sources (list (make source (src (make gaussian-src (frequency fcen) (fwidth df))) (center 0) (component src-cmpt))))
(if (= src-cmpt Ex)
(set! symmetries (list (make mirror-sym (direction Y)))))
(if (= src-cmpt Ey)
(set! symmetries (list (make mirror-sym (direction X)))))
(if (= src-cmpt Ez)
(set! symmetries (list (make mirror-sym (direction X)) (make mirror-sym (direction Y)))))
(define nearfield-box
(add-near2far fcen 0 1
(make near2far-region (center 0 (* 0.5 sxy)) (size sxy 0))
(make near2far-region (center 0 (* -0.5 sxy)) (size sxy 0) (weight -1))
(make near2far-region (center (* 0.5 sxy) 0) (size 0 sxy))
(make near2far-region (center (* -0.5 sxy) 0) (size 0 sxy) (weight -1))))
(define flux-box
(add-flux fcen 0 1
(make flux-region (center 0 (* 0.5 sxy)) (size sxy 0))
(make flux-region (center 0 (* -0.5 sxy)) (size sxy 0) (weight -1))
(make flux-region (center (* 0.5 sxy) 0) (size 0 sxy))
(make flux-region (center (* -0.5 sxy) 0) (size 0 sxy) (weight -1))))
(run-sources+ (stop-when-fields-decayed 50 src-cmpt (vector3 0 0) 1e-8))
(print "flux:, " (list-ref (get-fluxes flux-box) 0) "\n")
(define-param r (/ 1000 fcen)) ; 1000 wavelengths out from the source
(define-param npts 100) ; number of points in [0,2*pi) range of angles
(map (lambda (n)
(let ((ff (get-farfield nearfield-box (vector3 (* r (cos (* 2 pi (/ n npts)))) (* r (sin (* 2 pi (/ n npts)))) 0))))
(print "farfield:, " (number->string n) ", " (number->string (* 2 pi (/ n npts))))
(map (lambda (m)
(print ", " (number->string (list-ref ff m))))
(arith-sequence 0 1 6))
(print "\n")))
(arith-sequence 0 1 npts))
```
We use the `get-farfield` routine to compute the far fields by looping over a set of 100 points along the circumference of the circle with radius 1 mm. We compute the far fields at a wavelength of 1 μm for three different polarizations of the current source by setting the `src-cmpt` parameter to E$_x$, E$_y$, and E$_z$ in separate runs. The output consists of eight columns containing for each point: index identifier (integer), angle (radians), and six field components (E$_x$, E$_y$, E$_z$, H$_x$, H$_y$, H$_z$). Note that the far fields are always complex even though the near fields are real (as in this example). We also compute the flux from the source using the same bounding box.
The script is run and the output piped to a file using the following shell commands. The far field results are extracted from the output and placed in a separate file.
```sh
meep src-cmpt=Ez antenna-radiation.ctl |tee source_Jz_farfields.out
grep farfield: source_Jz_farfields.out |cut -d , -f2- > source_Jz_farfields.dat
```
From the far fields at each point $\mathbf{r}$, we can compute the radial flux: $\sqrt{P_x^2+P_y^2}$, where P$_x$ and P$_y$ are the components of the Poynting vector $\mathbf{P}(\mathbf{r})=(P_x,P_y,P_z)=\mathrm{Re}\, \mathbf{E}(\mathbf{r})^*\times\mathbf{H}(\mathbf{r})$. Since this is a 2d simulation, $P_z$ is always 0. We plot the radial flux normalized by its maximum value over the entire interval to obtain a range of values between 0 and 1. These are shown below in the linearly-scaled, polar-coordinate plots. As expected, the J$_x$ and J$_y$ sources produce [dipole](https://en.wikipedia.org/wiki/Electric_dipole_moment) radiation patterns while J$_z$ has a monopole pattern. These plots were generated using the following Octave/Matlab script.
```matlab
d = dlmread("source_Jz_farfields.dat",",");
Ex = conj(d(:,3));
Ey = conj(d(:,4));
Ez = conj(d(:,5));
Hx = d(:,6);
Hy = d(:,7);
Hz = d(:,8);
Px = real(Ey.*Hz-Ez.*Hy);
Py = real(Ez.*Hx-Ex.*Hz);
Pz = real(Ex.*Hy-Ey.*Hx);
Pr = sqrt(Px.^2+Py.^2);
angles = d(:,2)
polar(angles,Pr/max(Pr),'b-');
set(gca, 'xtick', [0 0.5 1.0]);
```
By [Poynting's theorem](https://en.wikipedia.org/wiki/Poynting%27s_theorem), the flux spectrum which is obtained by integrating around a closed surface should be the same whether it is calculated from the near or far fields (unless there are sources or absorbers in between), with slight differences due to discretization errors. The integral of the radial flux along the circumference of a circle with radius 1000 μm is obtained via:
```matlab
r = 1000 % circle radius (same units as Meep simulation)
disp(sprintf("flux:, %0.16f",sum(Pr)*2*pi*r/length(Pr)))
```
The far-field flux for the J$_z$ source is `2.457249`. The near-field flux, shown in the simulation output in the line prefixed by `flux:,`, is `2.456196`. This is a ratio of `0.999571`. Similarly, for the J$_x$ source, the far- and near-field flux values are `1.227260` and `1.227786` which is a ratio of `0.999571`. This ratio will converge to one as the resolution is increased.
<center>
</center>
### Far-Field Intensity of a Cavity
For this demonstration, we will compute the far-field spectra of a resonant cavity mode in a holey waveguide; a structure we had explored in [Tutorial/Resonant Modes and Transmission in a Waveguide Cavity](Resonant_Modes_and_Transmission_in_a_Waveguide_Cavity.md). The script is in [examples/cavity-farfield.ctl](https://github.com/stevengj/meep/blob/master/scheme/examples/cavity-farfield.ctl). The structure is shown at the bottom of the left image below.
To set this up, we simply remove the last portion of [examples/holey-wvg-cavity.ctl](https://github.com/stevengj/meep/blob/master/scheme/examples/holey-wvg-cavity.ctl), beginning right after the line:
```scm
(set! symmetries
(list (make mirror-sym (direction Y) (phase -1))
(make mirror-sym (direction X) (phase -1))))
```
and insert the following lines:
```scm
(define-param d1 0.2)
(define nearfield
(add-near2far fcen 0 1
(make near2far-region (center 0 (+ (* 0.5 w) d1))
(size (- sx (* 2 dpml)) 0))
(make near2far-region (center (+ (* -0.5 sx) dpml) (+ (* 0.5 w) (* 0.5 d1)))
(size 0 d1)
(weight -1.0))
(make near2far-region (center (- (* 0.5 sx) dpml) (+ (* 0.5 w) (* 0.5 d1)))
(size 0 d1))))
```
We are creating a "near" bounding surface, consisting of three separate regions surrounding the cavity, that captures <i>all</i> outgoing waves in the top-half of the cell. Note that the *x*-normal surface on the left has a `weight` of -1 corresponding to the direction of the *outward normal* vector relative to the *x* direction so that the far-field spectra is correctly computed from the outgoing fields, similar to the flux and force features. The parameter `d1` is the distance between the edge of the waveguide and the bounding surface, as shown in the schematic above, and we will demonstrate that changing this parameter does not change the far-field spectra which we compute at a single frequency corresponding to the cavity mode.
We then time step the fields until, at a random point, they have sufficiently decayed away as the cell is surrounded by PMLs, and output the far-field spectra over a rectangular area that lies <i>outside</i> of the cell:
```scm
(run-sources+ (stop-when-fields-decayed 50 Hz (vector3 0.12 -0.37) 1e-8))
(define-param d2 20)
(define-param h 4)
(output-farfields nearfield
(string-append "spectra-" (number->string d1) "-" (number->string d2) "-" (number->string h))
(volume (center 0 (+ (* 0.5 w) d2 (* 0.5 h))) (size (- sx (* 2 dpml)) h)) resolution)
```
The first item to note is that the far-field region is located <i>outside</i> of the cell, although in principle it can be located anywhere. The second is that the far-field spectra can be interpolated onto a spatial grid that has any given resolution but in this example we used the same resolution as the simulation. Note that the simulation itself used purely real fields but the output, given its analytical nature, contains complex fields. Finally, given that the far-field spectra is derived from the Fourier-transformed fields which includes an arbitrary constant factor, we should expect an overall scale and phase difference in the results obtained using the near-to-far-field feature with those from a corresponding simulation involving the full computational volume. The key point is that the results will be qualitatively but not quantitatively identical. The data will be written out to an HDF5 file having a filename prefix with the values of the three main parameters. This file will includes the far-field spectra for all six field components, including real and imaginary parts.
We run the above modified control file and in post-processing create an image of the real and imaginary parts of H$_z$ over the far-field region which is shown in insets (a) above. For comparison, we compute the steady-state fields using a larger cell that contains within it the far-field region. This involves a continuous source and complex fields. Results are shown in figure (b) above. The difference in the relative phases among any two points within each of the two field spectra is zero, which can be confirmed numerically. Also, as would be expected, it can be shown that increasing `d1` does not change the far-field spectra as long as the results are sufficiently converged. This indicates that discretization effects are irrelevant.
In general, it is tricky to interpret the overall scale and phase of the far fields, because it is related to the scaling of the Fourier transforms of the near fields. It is simplest to use the `near2far` feature in situations where the overall scaling is irrelevant, e.g. when you are computing a ratio of fields in two simulations, or a fraction of the far field in some region, etcetera.
|
