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

Surface conductivity treated as PEC #61

Closed odmiller closed 8 years ago

odmiller commented 8 years ago

Hi Homer,

I am trying to model scattering from open surfaces characterized by surface conductivities (e.g. graphene). I have attached sample .scuffgeo and .msh files, which I run using e.g.

scuff-scatter --geometry amorphous.scuffgeo --Omega 1 --pwDirection 0 0 1 --pwPolarization 1 0 0 --OPFTFile out.dat --PlotSurfaceCurrents --HDF5File out.hdf

I am getting NAN's for most of the power quantities. One issue I can immediately see is that scuff treats the surface as a PEC (I can verify this also by running scuff-analyze --geometry amorphous.scuffgeo). I think this happens for the following reason:

In line 277 of RWGSurface.cc (it's possible the version I'm reading is a few weeks old), the code assigns PEC to anything for which RegionLabels[0]==0. In the comments above it, this is explained as "any object without material specification or any surface without regions specification." The surface-conductivity case falls into this if statement, even though technically a material has been specified, since RegionLabels haven't been assigned. I'm not sure whether the fix should just be to change the condition in the if statement or to assign RegionLabels when the surface conductivity is assigned, but presumably one or both of these two?

Thanks in advance, Owen

P.S. I'm having troubles attaching the scuffgeo and msh files, will send separately.

HomerReid commented 8 years ago

Thanks for this report and for the analysis of what might be going wrong.

It looks like I didn't get a .msh file with your bug report, but I used your .scuffgeo and .dat files with a mesh I had lying around for a 1-micron sphere, and I got reasonable-looking results using both OPFT and DSIPFT:

scuff-scatter --geometry amorphous.scuffgeo --OPFTFile Sphere.OPFT --DSIPFTFile Sphere.DSIPFT --DSIRadius 5.0 --DSIPoints 302 --Omega 1.0 --pwDirection 0 0 1 --pwPolarization 1 0 0

This produced the following output for OPFT:

# scuff-scatter run on hikari (10/01/15::12:38:54)
# data file columns: 
# 1 omega           (rad/sec) 
# 2 surface label 
# 3 absorbed power  (watts)
# 4 scattered power (watts)
# 5 x-force         (nanonewtons)
# 6 y-force         (nanonewtons)
# 7 z-force         (nanonewtons)
# 8 x-torque        (nanonewtons * microns)
# 9 y-torque        (nanonewtons * microns)
#10 z-torque        (nanonewtons * microns)
1 Spheroid 0.000000e+00 8.084372e-03 -3.584397e-04 5.233953e-05 3.119315e-02 2.898027e-05 1.966481e-04 -1.833466e-04 

and for DSIPFT (quoting only the last line):

1 Spheroid 1.641523e-05 8.067957e-03 -2.387127e-05 9.000145e-07 3.150633e-02 2.774428e-06 1.280833e-04 -9.213342e-05 

Note that the OPFT absorbed power is 0.0, which has to do with the fact (as you observed) that RWGSurfaces with surface conductivity are treated as PEC for some purposes, and it seems the OPFT algorithm was never properly updated to handle power absorption for surface-conductivity bodies. That should be pretty easy to remedy.

On the other hand, the scattered power and Z-force reported by OPFT and DSIPFT are very similar, suggesting that the two calculations are doing the same thing and thus that the results are trustworthy. If this is the case, then I expect the absorbed power computed by DSIPFT to be accurate.

As for the NAN you reported, I wasn't able to reproduce it on my side, but on the other hand I did use my own mesh---I wonder if it could have to do with the particular mesh you used? I'll give it a try if you can attach your out.msh file.

HomerReid commented 8 years ago

In c102338dfe78a9169a2fd6ecdbfb160d6b81140d I have updated the OPFT calculation to account for finite surface conductivity when computing absorbed power. I did a couple of tests with open surfaces and got good agreement between OPFT and DSIPFT. When you get a chance, please give it a try and let me know if it works for you.

Meanwhile, I was not able to reproduce the NAN situation. When you get a chance, can you try pulling the latest code base from the repository and seeing if that fixes the NAN problem for your geometry?

In line 277 of RWGSurface.cc (it's possible the version I'm reading is a few weeks old), the code assigns PEC to anything for which RegionLabels[0]==0. In the comments above it, this is explained as "any object without material specification or any surface without regions specification." The surface-conductivity case falls into this if statement, even though technically a material has been specified, since RegionLabels haven't been assigned. I'm not sure whether the fix should just be to change the condition in the if statement or to assign RegionLabels when the surface conductivity is assigned, but presumably one or both of these two?

Yes, thanks for looking into this. Maybe the right way to think about this is that the IsPEC flag in RWGSurface should really be called IsInfinitesimallyThin or IsOpenSurface, which would encompass both PEC and IPEC ("imperfectly PEC", i.e. with finite surface conductivity) bodies. For most purposes within the code, PEC and IPEC bodies are equivalent---for example, they each contribute the same number of basis functions to the total for a given geometry, and they don't enter into the determination of whether a given evaluation point lies inside or outside a given region, etc. For both PEC and IPEC bodies we want Regions[0] to be the exterior region and Regions[1] to be a dummy nonexistent region (simply labeled "PEC"), since there is no interior for such a body (even if it is closed---the fields vanish identically inside such a body, so the interior of a closed PEC or IPEC surface does not count as a distinct region in a SCUFF calculation). So I think the right call is to leave things as is, but just reinterpret the IsPEC flag as including the case of IPEC bodies.

odmiller commented 8 years ago

Thanks, Homer! Indeed once I update the code the NaNs are gone. Should have checked this sooner.

I agree about the IsPEC flag as well. Thanks again!