HomerReid / scuff-em

A comprehensive and full-featured computational physics suite for boundary-element analysis of electromagnetic scattering, fluctuation-induced phenomena (Casimir forces and radiative heat transfer), nanophotonics, RF device engineering, electrostatics, and more. Includes a core library with C++ and python APIs as well as many command-line applications.
http://www.homerreid.com/scuff-em
GNU General Public License v2.0
126 stars 50 forks source link

Numerical instability in "quasistatic" calculations #82

Open gevero opened 8 years ago

gevero commented 8 years ago

Hi Homer

I am performing scattering and local field calculations with for the structure depicted below, using the scuff-scatter command line utility. All the files and inputs for the simulation are attacched below.

bowtie3

As you can se in the plot reported below, I monitored the local field at the (0,0,z) coordinates, i.e. in the middle of the gap sweeping along the triangle thickness, for different wavelengths:

Best

spectral_cross_section

LambdaFile.in.txt input_bowtie3.args.txt bowtie3.geo.txt bowtie3.EPfile.txt bowtie3.msh.txt e_ag_rashke_b_sc.edb.txt bowtie3.scuffgeo.txt

gevero commented 8 years ago

Hi homer

I tried to tackle the problem with a couple of different approaches:

The second approach seems to work brilliantly. Structure that gave highly unstable and unreliable results now converge robustly to physically meaningful results. Therefore I guess that the problem in the integration depends on the absolute scale of the modeled structure, and not on the relative scale with respect to the incident wavelength. If you can give me a nod of approval, I am closing the issue.

Thanks a lot

Giovanni

HomerReid commented 8 years ago

Yes, your solution was the right one. There is no intrinsic length scale in SCUFF-EM, so a mesh file describing a triangle of length L=1 that is simulated at a wavelength of Lambda=1 could equally well describe a one-micron structure at 1-micron wavelength, or a 1-mm structure a 1-mm wavelength. The panel-integration algorithms work best when dimensionless numbers like L/Lambda have numerical values within a few orders of magnitude of 1, i.e. in the range [0.001:1000]. The warning messages you observed are indicating that this condition is violated by your structure at your frequency.

The only place where an intrinsic length scale is assumed is in evaluating frequency-dependent material properties, and the length scale assumed in this case is microns, i.e. the assumption is that L=1 means that the object in question is 1 micron long. In API programs or python scripts this assumption can be changed by calling e.g. MatProp::SetLengthUnit(1.0e-3) (to set the default length unit to millimeters). For the command-line codes there is no direct functionality exposed to enable this (maybe there should be...) but the workaround you adopted, of modifying your description of the frequency dependence of the material, is fully equivalent and should work fine.

Regarding the failure of SCUFF-STATIC: Could it be simply that this involved the use of the wrong value of \epsilon? SCUFF-STATIC determines the DC permittivity of your structure by evaluating your permittivity function at Omega=0. Your structure exhibits a sharp frequency dependence despite being in the quasistatic regime, which means what you are really observing is a sharp permittivity dependence: each calculation you do at a given wavelength is essentially a DC calculation with a different value of the DC permittivity. You should get equivalent results by doing SCUFF_STATIC calculations with the permittivity of your structure set to epsilon(omega) for each omega value you are using. If you find otherwise, let me know. (Or, if you don't want to bother since you already got things working your way, that's fine too.)

SCUFF_STATIC will be faster and require much less memory for a given structure than doing a full finite-frequency calculation, so in particular it would enable you to go up to finer meshing resolution, but at the expense of some hassle involving separately setting CONST_EPS_XX+XXi values in the .scuffgeo file for each frequency you want to simulate, although it would be easy to write a little script to automate that.