Dear Steven and Meep users,
I am writing a ctl code to model a cone optical fiber tip with center
aperture (diameter 100nm) and ring gratings on the side wall of the Au
coating layer. So this is a perfect cylindrical symmetrical geometry.
The source wavelength is 405 nm with circular polarization. The output
shows there is light comes out from the side gratings. The grating gap
is 50 nm width. But nothing from the center aperture around (r=0), even
the center aperture size is bigger (100nm). THERE SHOULD BE SOMETHING
COMING OUT FROM THE CENTER. For example if the lattice size is 50x50
(50 nm characteristic length, resolution is 10), the output hdf5 matrix
size will be 501x500, instead of 500x500. So is this extra 1 unit coming
from r=0? Angular momentum "m" is 3 or 4. But if I put m=1, the
calculation will be unstable (number will be inf).
How do I ignore the r=0?
Quote from meep reference:
|"m| [|number|]
For |CYLINDRICAL| simulations, specifies that the angular ?
dependence of the fields is of the form /e/^/i//m/? (default is
|m=0|). If the simulation cell includes the origin /r/ = 0, then |m|
must be an integer."
If I want to get more accurate data from the center axis (r=0), how do I
setup it?
Attached is my ctl file and the outputs plot of eps structure and
denergy output. Side one angle is 15 degree. The light source locates
at the right side of the page. Tip center locates at the top center of
the page. So don't be confused.
Due to the size limitation, I posted these plots at link below.
http://picasaweb.google.com/btld123/MeepOutput#
Thank you!
Codes:
;Date Jan 09, 2009
;3D simulation of the smaller tip.
;System memory has been upgraded to 6GB, and running 64 bit linux.
;The memory will be upgraded to 8GB.
;Make a smaller model of the laser Tip to reduce calculation time and
memory requirement.
;Top to bottom size is 1.5um.
;Change coordinate to be Cylindrical instead of Cartisian
;Collect data at certain point in time domain to do a spectral analysis
;FFT will be running by Matlab.
;Geomery detail see note page 88.
;The setup in program TipGratingsV1.ctl is a little bit off for the cone
angle.
;The cone angle in this case is 15 degree, which is smaller than the
actual angle
;setup by the cone dimensions.
;The TipGratingsV2.ctl fixed this problem. It also fixed the wrong
caluclation that is based
;on outer cone. It should use the inner cone to do the calculation.
(reset-meep)
(define-param UnitL 50); characteristic unit length equals to 50 nm
(define-param WaveL 405);free space wavelength of the incident lightwave
in nm.
(define PI 3.1415928);define constant pi
(define fcen (/ UnitL WaveL));free space radian frequency as defined in Meep
(define-param S_r 50);r lattice dimension
;(define-param S_p (* 2 pi));phi lattice dimension
(define-param S_z 50);z lattice dimension
(define C_L 8.5026);cone length in z direction
(define C1_R1 32.7322);inner cone upper radius in r direction
(define C1_R2 1);inner cone lower radius in r direction. B/c in meep,
the r=0 axis
;is the axis of symmetry. r=0 should not be included in the calculation.
Otherwise
;there will be a blank stripe in the center, which doesnot reflect the
real situation.
;So for a 100 nm aperture (diameter), the setup in script should be 2
for radius given
;the 50nm characteristic length. KEY POINTS: Why there is not field
output at the center
;aperture?????
(define C_t_up 2); metal coating thickness at upper opening
(define C_t_down 2); metal coating thickness at lower bottom
(define-param C2_R1 (+ C1_R1 C_t_up));outter cone upper radius in r
direction
(define-param C2_R2 (+ C1_R2 C_t_down));outter cone lower radius in r
direction
(define-param P_t 10);define the PML thickness
(define C_C (- (/ (- S_z C_L) 2) P_t));cone center point in z axis
(define-param Src_Cr (/ C1_R1 2));center of the source in r axis
(define-param Src_Cz (- (/ S_z 2) P_t));center of the source in z axis
(define-param Src_Size (- C1_R1 2));size of the source in r direction
(set! default-material air)
(set! dimensions CYLINDRICAL)
(set-param! resolution 10); 10 pixels per unit
(set-param! m 4); m factor for angular momentum
(set-param! Courant 0.25)
(set! geometry-lattice (make lattice (size S_r no-size S_z)))
; The geometry of this structure is a bit tricky to specify. First the
outer cone needs to be
; drawn as a solid of Au. Then the outer most ring needs to be done
first. The outer most ring
; is drawn first by an air cylinder with the same radius as the ring
outer edged. Then another
; cylinder with the same location and height but with the radius of the
ring inner edge is drawn
; as a solid of gold. The same goes with the other rings until it
reaches the bottom aperture openning.
; The last step is drawing the inner cone filled with the Silicon dioxide.
; First the paramters that specify the cylinder locations are defined here.
; See note book drawing on page 91. The matlab program for calculating
it is stored at the HP
; notebook's /Desktop/Meep/work/Lsp.m. There was a difference in the
matlab code and the drawing on
; lab notebook that the ring starts at d2 in the outer cone.
; Detail table is listed on page 92 of the notebook.
(define dy1 0.6224)
(define dy2 1.4171)
(define dx1 4.3228)
(define dx2 5.2887)
(define-param d2 (+ C1_R2 dx1))
(define-param d3 (+ C1_R2 dx2))
(define-param d4 (+ d2 dx2))
(define-param d5 (+ d3 dx2))
(define-param d6 (+ d4 dx2))
(define-param d7 (+ d5 dx2))
(define-param d8 (+ d6 dx2))
(define-param d9 (+ d7 dx2))
(define-param d10 (+ d8 dx2))
(define-param d11 (+ d9 dx2))
(define-param H2 dy1)
(define-param H3 (+ H2 dy2))
(define-param H4 (+ H3 dy2))
(define-param H5 (+ H4 dy2))
(define-param H6 (+ H5 dy2))
; Ring 2
(define-param C_h2 (- C_L H2))
(define-param Cz_H2 (- (/ (- S_z (- C_L H2)) 2) P_t))
(define-param Rin_2 d2)
(define-param Rout_2 d3)
; Ring 3
(define-param C_h3 (- C_L H3))
(define-param Cz_H3 (- (/ (- S_z (- C_L H3)) 2) P_t))
(define-param Rin_3 d4)
(define-param Rout_3 d5)
; Ring 4
(define-param C_h4 (- C_L H4))
(define-param Cz_H4 (- (/ (- S_z (- C_L H4)) 2) P_t))
(define-param Rin_4 d6)
(define-param Rout_4 d7)
; Ring 5
(define-param C_h5 (- C_L H5))
(define-param Cz_H5 (- (/ (- S_z (- C_L H5)) 2) P_t))
(define-param Rin_5 d8)
(define-param Rout_5 d9)
; Ring 6
(define-param C_h6 (- C_L H6))
(define-param Cz_H6 (- (/ (- S_z (- C_L H6)) 2) P_t))
(define-param Rin_6 d10)
(define-param Rout_6 d11)
(set! eps-averaging? false)
;If true (the default), then subpixel averaging is used when
initializing the dielectric function (see the Farjadpour et al.
reference in Citing Meep). The input variables subpixel-maxeval (default
100000) and subpixel-tol (default 1.0e-4) specify the maximum number of
function evaluations and the integration tolerance for subpixel
averaging. Increasing/decreasing these, respectively, will cause a more
accurate (but slower) computation of the average ? (with diminishing
returns for the actual FDTD error).
(define omega_d 4.43)
(define Au ; definition of material dispersion
(make dielectric (epsilon 5.967)
(polarizations
(make polarizability
(omega 1e-8) (gamma 3.33e-02)
(delta-epsilon (* (* omega_d omega_d) 1e+16)))
(make polarizability
(omega 1.36) (gamma 0.22)
(delta-epsilon 1.09)))))
(define silica
(make dielectric (epsilon (* 1.55 1.55))))
(set! geometry
(list
(make cone (center 0 0 C_C) (material Au) (radius C2_R1)
(radius2 C2_R2) (axis 0 0 -1) (height C_L))
(make cylinder (center 0 0 Cz_H6) (material air) (radius Rout_6)
(axis 0 0 1) (height C_h6))
(make cylinder (center 0 0 Cz_H6) (material Au) (radius Rin_6)
(axis 0 0 1) (height C_h6))
(make cylinder (center 0 0 Cz_H5) (material air) (radius Rout_5)
(axis 0 0 1) (height C_h5))
(make cylinder (center 0 0 Cz_H5) (material Au) (radius Rin_5)
(axis 0 0 1) (height C_h5))
(make cylinder (center 0 0 Cz_H4) (material air) (radius Rout_4)
(axis 0 0 1) (height C_h4))
(make cylinder (center 0 0 Cz_H4) (material Au) (radius Rin_4)
(axis 0 0 1) (height C_h4))
(make cylinder (center 0 0 Cz_H3) (material air) (radius Rout_3)
(axis 0 0 1) (height C_h3))
(make cylinder (center 0 0 Cz_H3) (material Au) (radius Rin_3)
(axis 0 0 1) (height C_h3))
(make cylinder (center 0 0 Cz_H2) (material air) (radius Rout_2)
(axis 0 0 1) (height C_h2))
(make cylinder (center 0 0 Cz_H2) (material Au) (radius Rin_2)
(axis 0 0 1) (height C_h2))
;The last cylinder is needed to make sure the aperture is opened and not
blocked by numerical error.
(make cylinder (center 0 0 0) (material air) (radius C1_R2)
(axis 0 0 1) (height S_z))
(make cone (center 0 0 C_C) (material silica) (radius C1_R1)
(radius2 C1_R2) (axis 0 0 -1) (height C_L)
)))
(set! pml-layers (list
(make pml (thickness P_t))))
(set! sources (list
(make source (src (make continuous-src (frequency fcen)))
(component Ep) (amplitude 0+1i) (center Src_Cr 0 Src_Cz) (size C1_R1 0 0))
(make source (src (make continuous-src (frequency fcen)))
(component Er) (amplitude 1) (center Src_Cr 0 Src_Cz) (size C1_R1 0 0))))
(use-output-directory)
(define-param D_20nm (- (/ S_z 2) P_t C_L 0.4))
(define-param D_50nm (- (/ S_z 2) P_t C_L 1))
(define-param D_100nm (- (/ S_z 2) P_t C_L 2))
(define-param D_200nm (- (/ S_z 2) P_t C_L 4))
(define-param D_500nm (- (/ S_z 2) P_t C_L 10))
(run-until 500
(at-beginning output-epsilon)
(after-time 500
(to-appended "Er_20nm"
(at-every 0.5
(in-volume (volume (center (/ S_r 2) 0 D_20nm) (size S_r 0 0))
output-dpwr)))
(to-appended "Er_50nm"
(at-every 0.5
(in-volume (volume (center (/ S_r 2) 0 D_50nm) (size S_r 0 0))
output-dpwr)))
(to-appended "Er_100nm"
(at-every 0.5
(in-volume (volume (center (/ S_r 2) 0 D_100nm) (size S_r 0 0))
output-dpwr)))
(to-appended "Er_200nm"
(at-every 0.5
(in-volume (volume (center (/ S_r 2) 0 D_200nm) (size S_r 0 0))
output-dpwr)))
(to-appended "Er_500nm"
(at-every 0.5
(in-volume (volume (center (/ S_r 2) 0 D_500nm) (size S_r 0 0))
output-dpwr))))
(at-every 0.5 output-dpwr))
_______________________________________________
meep-discuss mailing list
[email protected]
http://ab-initio.mit.edu/cgi-bin/mailman/listinfo/meep-discuss