Skip to content

Instantly share code, notes, and snippets.

@oskooi
Last active January 20, 2023 06:23
Show Gist options
  • Save oskooi/9be8dafb58d60562d17ae8cfb0771344 to your computer and use it in GitHub Desktop.
Save oskooi/9be8dafb58d60562d17ae8cfb0771344 to your computer and use it in GitHub Desktop.
Determines whether the fields in a cylindrical-coordinate simulation with m=±1 are purely real
import argparse
import numpy as np
import meep as mp
def dipole_emission(res: float, dpml_r: float, dpml_z: float, s: float,
m: int, fcen: float) -> float:
"""Computes the complex fields from an Er point source at r=z=0 in vacuum."""
cell_size = mp.Vector3(s+dpml_r,0,s+2*dpml_z)
pml_layers = [
mp.PML(dpml_r, direction=mp.R),
mp.PML(dpml_z, direction=mp.Z),
]
sources = [
mp.Source(
src=mp.CustomSource(src_func=lambda t: np.sin(2*np.pi*fcen*t)),
center=mp.Vector3(),
component=mp.Er,
),
]
sim = mp.Simulation(
resolution=res,
cell_size=cell_size,
dimensions=mp.CYLINDRICAL,
m=m,
sources=sources,
boundary_layers=pml_layers,
)
sim.run(until=40.6)
er = sim.get_array(
component=mp.Er,
center=mp.Vector3(0.5*s,0,0),
size=mp.Vector3(s,0,s),
cmplx=True,
)
return er
if __name__ == '__main__':
parser = argparse.ArgumentParser()
parser.add_argument(
'-res',
type=float,
default=50.,
help='resolution (default: 50.0)',
)
parser.add_argument(
'-dpml_r',
type=float,
default=1.,
help='PML thickness in r direction (default: 1.0)',
)
parser.add_argument(
'-dpml_z',
type=float,
default=1.,
help='PML thickness in z direction (default: 1.0)',
)
parser.add_argument(
'-s',
type=float,
default=5.,
help='thickess of non-PML region (default: 5.0)',
)
args = parser.parse_args()
er0 = dipole_emission(args.res, args.dpml_r, args.dpml_z, args.s, -1, +1.)
er1 = dipole_emission(args.res, args.dpml_r, args.dpml_z, args.s, -1, -1.)
er2 = dipole_emission(args.res, args.dpml_r, args.dpml_z, args.s, +1, +1.)
er3 = dipole_emission(args.res, args.dpml_r, args.dpml_z, args.s, +1, -1.)
def is_real(x):
return np.all(np.conj(x) == x)
# print out the first 20 elements of the Er arrays
# for inspection of their imaginary parts
print(f"er0:, {er0[:20]}")
print(f"er1:, {er1[:20]}")
print(f"er2:, {er2[:20]}")
print(f"er3:, {er3[:20]}")
print(f"real? er0: {is_real(er0)}, er1: {is_real(er1)}, "
f"er2: {is_real(er2)}, er3: {is_real(er3)}")
@oskooi
Copy link
Author

oskooi commented Jan 20, 2023

output from running this script:

run 0 finished at t = 40.61 (4061 timesteps)
er0:, [[ 4.97362786e-04+0.j -2.87743027e-03+0.j -2.64652356e-03+0.j ...
  -2.61486309e-04+0.j -2.71377079e-04+0.j -2.77701179e-04+0.j]
 [ 3.31736140e-04+0.j -2.63761179e-03+0.j -2.45164624e-03+0.j ...
  -2.50864191e-04+0.j -2.62602309e-04+0.j -2.70692758e-04+0.j]
 [-2.06848158e-04+0.j -2.35198918e-03+0.j -2.20251263e-03+0.j ...
  -2.39333991e-04+0.j -2.52613174e-04+0.j -2.62506871e-04+0.j]
 ...
 [-1.08385315e-03+0.j  3.44963506e-03+0.j  3.16224061e-03+0.j ...
  -3.06395092e-05+0.j -5.70990524e-05+0.j -8.21864684e-05+0.j]
 [-8.38156376e-04+0.j  3.65826090e-03+0.j  3.36381841e-03+0.j ...
  -1.70343626e-05+0.j -4.36051665e-05+0.j -6.90490534e-05+0.j]
 [-4.65454692e-04+0.j  3.79258515e-03+0.j  3.52540904e-03+0.j ...
  -3.76449054e-06+0.j -3.02802036e-05+0.j -5.57899726e-05+0.j]]
er1:, [[-4.97362786e-04+0.j  2.87743027e-03+0.j  2.64652356e-03+0.j ...
   2.61486309e-04+0.j  2.71377079e-04+0.j  2.77701179e-04+0.j]
 [-3.31736140e-04+0.j  2.63761179e-03+0.j  2.45164624e-03+0.j ...
   2.50864191e-04+0.j  2.62602309e-04+0.j  2.70692758e-04+0.j]
 [ 2.06848158e-04+0.j  2.35198918e-03+0.j  2.20251263e-03+0.j ...
   2.39333991e-04+0.j  2.52613174e-04+0.j  2.62506871e-04+0.j]
 ...
 [ 1.08385315e-03+0.j -3.44963506e-03+0.j -3.16224061e-03+0.j ...
   3.06395092e-05+0.j  5.70990524e-05+0.j  8.21864684e-05+0.j]
 [ 8.38156376e-04+0.j -3.65826090e-03+0.j -3.36381841e-03+0.j ...
   1.70343626e-05+0.j  4.36051665e-05+0.j  6.90490534e-05+0.j]
 [ 4.65454692e-04+0.j -3.79258515e-03+0.j -3.52540904e-03+0.j ...
   3.76449054e-06+0.j  3.02802036e-05+0.j  5.57899726e-05+0.j]]
er2:, [[ 4.97362786e-04+0.j -2.87743027e-03+0.j -2.64652356e-03+0.j ...
  -2.61486309e-04+0.j -2.71377079e-04+0.j -2.77701179e-04+0.j]
 [ 3.31736140e-04+0.j -2.63761179e-03+0.j -2.45164624e-03+0.j ...
  -2.50864191e-04+0.j -2.62602309e-04+0.j -2.70692758e-04+0.j]
 [-2.06848158e-04+0.j -2.35198918e-03+0.j -2.20251263e-03+0.j ...
  -2.39333991e-04+0.j -2.52613174e-04+0.j -2.62506871e-04+0.j]
 ...
 [-1.08385315e-03+0.j  3.44963506e-03+0.j  3.16224061e-03+0.j ...
  -3.06395092e-05+0.j -5.70990524e-05+0.j -8.21864684e-05+0.j]
 [-8.38156376e-04+0.j  3.65826090e-03+0.j  3.36381841e-03+0.j ...
  -1.70343626e-05+0.j -4.36051665e-05+0.j -6.90490534e-05+0.j]
 [-4.65454692e-04+0.j  3.79258515e-03+0.j  3.52540904e-03+0.j ...
  -3.76449054e-06+0.j -3.02802036e-05+0.j -5.57899726e-05+0.j]]
er3:, [[-4.97362786e-04+0.j  2.87743027e-03+0.j  2.64652356e-03+0.j ...
   2.61486309e-04+0.j  2.71377079e-04+0.j  2.77701179e-04+0.j]
 [-3.31736140e-04+0.j  2.63761179e-03+0.j  2.45164624e-03+0.j ...
   2.50864191e-04+0.j  2.62602309e-04+0.j  2.70692758e-04+0.j]
 [ 2.06848158e-04+0.j  2.35198918e-03+0.j  2.20251263e-03+0.j ...
   2.39333991e-04+0.j  2.52613174e-04+0.j  2.62506871e-04+0.j]
 ...
 [ 1.08385315e-03+0.j -3.44963506e-03+0.j -3.16224061e-03+0.j ...
   3.06395092e-05+0.j  5.70990524e-05+0.j  8.21864684e-05+0.j]
 [ 8.38156376e-04+0.j -3.65826090e-03+0.j -3.36381841e-03+0.j ...
   1.70343626e-05+0.j  4.36051665e-05+0.j  6.90490534e-05+0.j]
 [ 4.65454692e-04+0.j -3.79258515e-03+0.j -3.52540904e-03+0.j ...
   3.76449054e-06+0.j  3.02802036e-05+0.j  5.57899726e-05+0.j]]
real? er0: True, er1: True, er2: True, er3: True

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment