Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Implement proper integration weights in cylindrical coordinates #2467

Open
smartalecH opened this issue Apr 10, 2023 · 3 comments · May be fixed by #2470
Open

Implement proper integration weights in cylindrical coordinates #2467

smartalecH opened this issue Apr 10, 2023 · 3 comments · May be fixed by #2470
Assignees
Labels

Comments

@smartalecH
Copy link
Collaborator

As discussed in #2452, the integration weights in cylindrical coordinates are implemented identically as those in cartesian coordinates (which is incorrect). The procedure to fix this is outlined in the code itself:

meep/src/loop_in_chunks.cpp

Lines 112 to 113 in 1fe3899

Integration Weights in Cylindrical Coordinates
FIXME: implement this below?

but it's not trivial to implement, particularly with cases at or crossing r=0 (which happens a lot...)

Anecdotally, this doesn't seem to be an issue if R is rather large. However, for smaller cells, the error compounds substantially (and doesn't "easily" disappear with resolution).

@smartalecH smartalecH self-assigned this Apr 10, 2023
@smartalecH smartalecH added the bug label Apr 11, 2023
@smartalecH
Copy link
Collaborator Author

@oskooi, i've got a draft for a fix to this. You mentioned you have a test problem. Can you link me to it?

@oskooi
Copy link
Collaborator

oskooi commented Apr 12, 2023

The test that I had in mind is to demonstrate that the normalized radiated flux $P(r) / r^2$ is independent of the location of the point source for $r \approx 0$. Currently, $P(r) / r^2 \rightarrow C$ for a constant $C$ only for $r \rightarrow \infty$ as shown in the plots of #2108 (comment).

To set this up, we can simply modify the tutorial script python/examples/point_dipole_cyl.py. The function led_flux returns the radiated flux (flux_air) and the total emitted flux (flux_src) separately as a 2-tuple:

for m in ms:
flux_air, flux_src = led_flux(
layer_thickness,
dipole_height,
rp,
m,
)
flux_air_tot += flux_air if m == 0 else 2 * flux_air
flux_src_tot += flux_src if m == 0 else 2 * flux_src

The key metric to check would be flux_air_tot / rpos**2 vs. rpos for several values of rpos near 0.

Note that this test does not use the LDOS feature.

@oskooi
Copy link
Collaborator

oskooi commented Apr 12, 2023

Also note that in Tutorial/Nonaxisymmetric Dipole Sources, I intentionally used "large" values of the source position (rpos > 3.5) to ensure the extraction-efficiency results were reasonably accurate:

rpos = [3.5, 6.7, 9.5]
for rp in rpos:

Hopefully once this bug is fixed, we can modify this tutorial to include point-source locations closer to $r = 0$.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

Successfully merging a pull request may close this issue.

2 participants