Radiation pattern for dipole emitter in cylindrical coordinates

[kalosu]

Hi there again friends from the MEEP community,

I wanted to ask for your help for extracting the radiation pattern of a dipole emitter while running my simulations in cylindrical coordinates.

To start, I am trying to evaluate the radiation pattern for a dipole in vacuum sitting at the origin of my simulation domain. I use a near to far-field transformation via get_farfield and Near2FarRegion to compute the E and H components as a function of a given angle of observation (I check only in between 0 and 90 degrees with respect to the r axis)

For calculating the Poynting flux (and the radiation pattern $P(\theta)$) I have followed the tutorial given in https://meep.readthedocs.io/en/latest/Python_Tutorials/Cylindrical_Coordinates/

Nevertheless, the far-field radiation pattern that I get is not what I would expect. Or am I not understanding it correctly and for a dipole point sitting at the origin in vacuum we would not get the same far-field pattern as for when the problem is solved in cartesian coordinates?

Notice that I am just interested for the moment on the radiation pattern on the rz plane (essentially would be the same as for cartesian coordinates right??)

I have compared the cylindrical coordinates radiation pattern to the one obtained for a 2D cartesian simulation (with the source polarized along x)

As you can see from below plot, the results from the cylindrical case do not match at all the ones obtained for the 2D cartesian case. (for the cylindrical case I am just evaluating from 0 to 90 degrees)

ificant difference (I have used the same simulation parameters, the only thing that I have changed is the implementation in between cylindrical and cartesian coordinates but the other global parameters were the same)

You can also take a look to the domains that I used in where the blue lines indicate the position of my Near2FarRegion regions.

Could it be that I am missing some scaling factor or something similar? Any thoughts on what could be the problem in this case?

[oskooi]

There is an unresolved bug involving an Er point source at $r = 0$ with $|m| = 1$ as described in #2538 (comment). A simple workaround that you can try is to move the Er source a couple of voxels away from $r = 0$. See also #2644.

[kalosu]

Hi @oskooi , thanks for the feedback and details.

I have just found that there is a new tutorial for computing the far field radiation pattern for a dipole emitter embedded within a dielectric disk (https://meep.readthedocs.io/en/latest/Python_Tutorials/Near_to_Far_Field_Spectra/)

In that tutorial I see that only the m=-1 mode is considered even though the dipole is not located at the origin. Is this just for the purpose of the example or is the single m=-1 mode sufficient in this case? (as far as I understand it, higher m modes are needed to correctly represent the total field components right??)

I have tried changing the position of my dipole (I am now just trying for a dipole in vacuum) along r and I include more m modes in the expansion. Nevertheless I cannot get the far-field pattern to converge to the pattern of a dipole emitter in vacuum. Any idea about what could I test in this case?

[oskooi]

In that tutorial I see that only the m=-1 mode is considered even though the dipole is not located at the origin. Is this just for the purpose of the example or is the single m=-1 mode sufficient in this case? (as far as I understand it, higher m modes are needed to correctly represent the total field components right??)

In Tutorial/Radiation Pattern of a Disc in Cylindrical Coordinates, the dipole was intentionally placed at $r \neq 0$ to avoid the bug that I mentioned. Also, the example involved a single $m = -1$ although as you noted the total flux requires summing the power over $m$.

I have tried changing the position of my dipole (I am now just trying for a dipole in vacuum) along r and I include more m modes in the expansion. Nevertheless I cannot get the far-field pattern to converge to the pattern of a dipole emitter in vacuum. Any idea about what could I test in this case?

As a reference, to obtain the radiation pattern $P(\theta)$ for a dipole at $r \neq 0$, you need to follow the second of two approaches described in Tutorial/Cylindrical Coordinates/Nonaxisymmetric Dipole Sources. This is discussed in the paragraph which starts with the sentence "However, if you want to compute $P(\theta) = \int P(\theta, \phi) d\phi$ requires three steps." Note that the summation in step 2 involves terms $P_{m=0}(\theta)$ and $2 P_{m}(\theta)$. Did you remember to use the correct weights?

[kalosu]

Hi @oskooi,

Thanks for the feedback.

I have tried using the second approach from the Nonaxisymmetric dipole sources tutorial for computing $P(\theta)$ for a dipole at $r\neq0$ using 6 m mode terms (from -5 to 0). I have set the dipole at $r=1.5\Delta r$ as suggested in the comments from the references you shared with me before.

However, I still cannot obtain the far-field pattern for a dipole emitter in vacuum. You can clearly see the discrepancy between the result from a 3D cartesian simulation and the 2D cylindrical one in the following plot (As you know, for the cylindrical coordinates case, the zero $\theta$ angle points along z, whereas in the cartesian simulations the zero $\theta$ angle points along x ). Disregarding the direction of the zero $\theta$, you can see that for the cylindrical coordinates case, there is some high residue flux towards large angles which should not be the case for a dipole embedded in vacuum right?

As for the way in which I compute the $P(\theta)$ pattern, I first calculate the Poynting flux and take the dot product with a unit vector defined over a circle (following the same definition for having $\theta$ zero pointing towards z):

for n in range(npts):
    ff = sim.get_farfield(n2f_obj,mp.Vector3(r*np.sin(angles[n]),0,r*np.cos(angles[n])))
    if int(m_val)==0:
        E[:,n]= [ff[j] for j in range(3)]
        H[:,n]= [ff[j+3] for j in range(3)]
    # if int(m_val)==-1:
    else: 
        E[:,n]= [ff[j] for j in range(3)]
        H[:,n]= [ff[j+3] for j in range(3)]
        # E[:,n] = [np.conj(ff[j]) for j in range(3)]
        # H[:,n] = [ff[j+3] for j in range(3)]

n_vec = np.zeros((3,npts))
n_vec[0,:] = np.sin(angles)
n_vec[2,:] = np.cos(angles)
cross_eval = np.cross(E,np.conjugate(H),axis=0)
if int(m_val)==0:
    P_eval += (np.sum(cross_eval*n_vec,axis=0))
else: 
    P_eval += 2*np.sum(cross_eval*n_vec,axis=0)

Then once the loops over the m modes are finished, I calculate the total $P(\theta)$ patterns by taking the real part of P_eval and multiplying it by $2 \pi$

if mp.am_master():
    P_eval = 2*np.pi*np.real(P_eval)

Any idea about what could be happening and what could I do?

[kalosu]

Hi again,

I have managed to get a somehow closer shape also considering step number 3 from the second approach mentioned in the Nonaxisymmetric dipole sources tutorial.

In specific, I calculate the total emitted power by accounting for the $r^2 sin(\theta) d\theta$ factors and then compute the ratio of the power radiated towards a given solid angle $\theta'$ and the total power. You can observe how the radiation pattern now looks closer to what is expected

However, the crossing points are not really correct (the position of the curve as a function of $r$ and $\theta$). I have tried increasing the resolution and the number of m modes being evaluated, without observing any change in the pattern.

Any idea about what could be happening here?

Also, why would I need to follow up to step 3 in order to extract the radiation pattern? Would it not be sufficient to compute $P(\theta)$ up to what is indicated in step number 2 for obtaining the radiation pattern??

[kalosu]

Hi again,

I have given it a try by using the DFT fields instead of a near2far-field transformation. For this, I extended the cylindrical domain along r for a sufficiently large distance and set a z value of 5 with the dipole emitter positioned at z=0.

Similarly as before, I set the dipole at $r=1.5\Delta r$ and considered 6 m modes (from m=0 to -5)

I extracted then all field components (er,ep,ez,hr,hp,hz) and estimated the Poynting flux for each m mode mode separately.

I then projected the points along r to angles by calculating the angle made between each point along r and the distance to the DFT monitor I used (2.5 for z and 30 for the maximum r value)

Following you can take a look to the results again for a dipole in vacuum. These results were obtained by using all m modes Poynting flux, excepting m=0. For comparison, I have included the radiation pattern from the analytical expression for a dipole source $(\propto \frac{e^{kr}}{r} sin(\theta) )$ and the results obtained from a 3D cartesian simulation. You can clearly see that there is a mismatch in the obtained far-field radiation pattern similarly to the N2F case which I presented in my previous comment.

If I include the contribution from the m=0 term, then things get worse:

The problem from above case (by including the m=0 contribution) comes from the form of the Poynting flux for m=0. There seems to be a large overshoot towards r=0 for m=0, differently from the other m mode cases.

Below I show the estimated Poynting flux terms from all used modes including m=0

You can clearly observe that there is a large disparity in the flux values for pz and pr in the case of m=0, with respect to the rest of flux terms.
I can also show the same plot but without considering the m=0 contributions:

In this case, the major contributions come from m=-1 and m=-2. For higher order modes, the contributions are really small (in the order of 1e^-13).

Finally, I have compared the field components extracted from different m values (from m=0 to 3). I took the absolute value of the extracted complex amplitudes and present the comparisons in the following plot

Clearly, the magnitude of the m=0 field components is one or two orders of magnitude larger than for the other modes. Therefore, this would explain why there is a large overshoot towards r=0 for when the m=0 terms are included in the Poynting flux estimation.

What do you think? Are these results for m=0 consistent? Should these not be somehow similar in magnitude to the m=1 and m=2 mode results? Any idea about what could be going on and about possible solutions that I could use, in order to correctly use the cylindrical coordinates solver? I can see that in this case, adding the position r offset to the dipole did not really help in solving the issues related to the case in where r=0.

Any comments and help will be appreciated.

[kalosu]

Hi there,

I have given it another try, but instead of computing the far-field flux over a circle, I just defined a flat plane extending along r (perpendicular to z in cylindrical coordinates).

Similarly to above, I can show you the results from evaluating the flux crossing through the imaginary 1D line in where I just consider the z components of the Poynting flux vector. First I show the results obtained by neglecting the contributions from the m=0 mode.

The other two curves shown in the plot correspond to the reference analytical expression and to the results obtained from a 2D cartesian coordinates simulation. You can clearly see that the results obtained from the 2D simulation and the analytical expression are quite similar and for the cylindrical coordinates result there is a small mismatch at the center.

Now, if I include the m=0 contributions:

Clearly, the difference is even more pronounced.

We can compare the individual Poynting components. For a reference, I display the Poynting components obtained from the cartesian 2D simulation. (In this case, the y component is the one that I have used for above plots). The equivalent contributions for the m=0 and m=-1 modes are shown as well (in cylindrical coordinates with $\phi=0$ we have $r=x$, $phi=y$ and $z=z$ based on the documentation). Finally, on the left top plot, we can see the individual Poynting z contributions from the 3 considered modes (0,-1,-2).

As you can see, the increasing contribution from the m=0 Poynting z term makes the difference with respect to the reference shape to increase.

What could I be missing here? Or is there any problem at the moment for computing the contributions from m=0 even if I position the emitter at a $r\neq 0$ position?

If you have any insights or comments, it will be highly appreciated.