NanoComp / meep

free finite-difference time-domain (FDTD) software for electromagnetic simulations
GNU General Public License v2.0
1.22k stars 622 forks source link

LDOS in 2D #23

Closed amrit-poudel closed 6 years ago

amrit-poudel commented 8 years ago

I tried comparing LDOS computed in meep with analytical value in vacuum, but could not get the two results agree.

Analytically, LDOS in 2D = 2\omega/(\pi c^2) \ Tr[ Im [G] ], G = dyadic Greens function of the vacuum. In 2D, Im[G] = 1/4, for all three diagonal components.

I use a continuous-wave point source at a fixed frequency to compute the LDOS at that frequency and at that location(I could use other temporal sources to get LDOS at several frequencies in a single FDTD run, but for now, I just want to compare the numerical with analytical results and for simplicity I choose a continuous-wave source).

I am using meep as a c++ library and use dft_ldos object with update and ldos functions to compute the LDOS.

The result is off by an order of magnitude.

What is the unit of LDOS in 2D in meep?

stevengj commented 8 years ago

See http://ab-initio.mit.edu/wiki/index.php/Meep_Tutorial/Local_density_of_states

stevengj commented 8 years ago

(Don't use a CW source, use a gaussian pulse. You want something that decays in time so that you can stop the simulation after a finite time and have the (discrete-time) Fourier transform. For a CW source, the accumulated Fourier transform grows linearly with time.)

amrit-poudel commented 8 years ago

Thank you for your suggestion. I will use a Gaussian source.

Just so I understand, does meep perform discrete Fourier transform at each FDTD time step and then sums up the result to obtain the final value?

stevengj commented 8 years ago

Yes, it accumulates a DTFT of the fields at each timestep, as described here: http://ab-initio.mit.edu/wiki/index.php/Meep_Introduction#Transmission.2Freflection_spectra

amrit-poudel commented 8 years ago

Thanks again.

I tried computing LDOS in 2D in vacuum using a point dipole source with unit amplitude whose temporal part is a Gaussian. As far as I understand, the Gaussian source in meep = exp(i * 2_pi_f * t) * Gaussian envelope. I chose linear frequency, f = 0.00997/2/pi (in meep units) and fwidth =f.

I ran the simulation for 20 * period, where period = 1/f. By this time, source obviously dies out, so the DTFT should just work fine.

Then I compute LDOS (in a frequency range: f-fwidth/2 to f+fwidth/2) numerically, but this does not agree with the analytical value given by: LDOS = 2 * w/(\pi) * Tr[ Im [G] ] = 0.00476, where I used Tr[Im[G]] = 3/4. G is dimensionless in 2D.

I am confident about my computational setup (width of pml, choice of frequency, computational domain, etc) since I can reproduce the analytical value of the Greens function using a cw-source. Yet, for some reason, I cannot get correct result for LDOS. I am not quite sure what I am missing.

The analytical formula of LDOS (and the corresponding factor of 2, due to electric+magnetic, and \pi, etc) is consistent with the formula given in chapter 4 of your online text book.

UPDATE: Please ignore my comment. It works perfectly.

amrit-poudel commented 8 years ago

Sorry to reopen this issue, but I noticed that LDOS computed by meep seems to depend on the location of the point source, which should not be the case in vacuum. I find this behavior very strange. The numerical LDOS agrees with analytical result at some locations but deviates significantly at others (nearby locations).

stevengj commented 8 years ago

Does it converge with resolution?

amrit-poudel commented 8 years ago

I noticed that dft_ldos.ldos() function returns different results on a single processor vs multiple processors.

The ratio of LDOS on multiple processors is correct, but individual LDOS is totally incorrect on multiple processors. On a single core processor, however, both individual and the ratio of LDOS are correct.

Could this be a bug or is there a different way to run ldos simulation on multiple nodes/processors that I am unaware of?

I compute ldos on multiple nodes/processors using meep-mpi:

initialize mpi(argc, argv);

Call dft_ldos.ldos() function to compute ldos.

I am using meep v 1.3.

Is there a quicker way to fix this bug other than running simulation on a single core?

erikshipton commented 7 years ago

I ran into problems with LDOS using multiple processors and a random phase of the point source. When using multiple processors I get a different LDOS every time I run the same exact problem. On a single processor the same LDOS is calculated every run.

I'm also getting a different near2far result when using the random phased point source on multiple processors each time I run the same exact problem. On a single processor I get the same near2far each run.

oskooi commented 7 years ago

This is a bug in dft_ldos which occurs whenever a point source is placed between two chunks (e.g., placing a source at the origin and using more than one processor). The serial version, as you have noticed, is unaffected since the non-PML regions are not divided into chunks.

The issue is related to an over-counting of Jsum in dfot_ldos::update as a single source is duplicated within two chunks during the source initialization (i.e., src_vol_chunkloop in sources.cpp). A possible fix would involve determining that the source has been duplicated and then to properly scale Jsum. This check would need to happen at each time step which is inefficient.

An alternative approach may be to add a new flag member to dft_ldos which is set during the source initialization to indicate that a source duplication has occurred.