Dear meep users:
With the help of prof. Steven Johnson above, I now have a working code to
calculate time averaged E square field. The following is an example ctl file,
for a pair of gold particles under an illumination of plane wave source. Many
parts of the code are copied from earlier posts.
Thanks for all.
Seong Kyu Kim
; -----------------------------------------------------------------------------
; Lattice and Boundaries ------------------------------------------------------
(define-param sx 0.1) ; the length of x of the lattice
(define-param sy 0.1) ; the length of y of the lattice
(define-param sz 0.1) ; the length of z of the lattice
(define-param dpml 0.1) ; thickness of PML layers
(set! geometry-lattice (make lattice (size sx sy (+ sz (* 2 dpml)) )))
(set! pml-layers (list (make pml (thickness dpml) (direction Z))))
(set! k-point (vector3 0 0 0))
(set-param! resolution 500)
; -------------------------------------------------------------------------
; Target
(define omega_d 7.05) ; plasma frequency of Au in Drude model
(define Au ; definition of material dispersion
(make dielectric (epsilon 5.967)
(polarizations
(make polarizability
(omega 1e-20) (gamma 5.310E-02)
(delta-epsilon (* (* omega_d omega_d) 1e+40)))
(make polarizability
(omega 2.168) (gamma 3.498E-01)
(delta-epsilon 1.09)))))
; -----------------------------------------------------------------------------
(set! geometry
(list
(make sphere (center -0.025 0 0) (radius 0.02) (material Au))
(make sphere (center 0.025 0 0) (radius 0.02) (material Au))))
; -----------------------------------------------------------------------------
; Source - plane wave ---------------------------------------------------------
; pw-amp is a function that returns the amplitude exp(ik(z+z0)) at a
; given point z. (We need the z0 because current amplitude functions
; in Meep are defined relative to the center of the current source,
; whereas we want a fixed origin.) Actually, it is a function of k
; and z0 that returns a function of z ...
(define ((pw-amp k z0) z) (exp (* 0+1i (vector3-dot k (vector3+ z z0)))))
(define-param fcen 1.8); pulse center frequency
(define-param df 0.4) ; frequency bandwidth
(define-param kdir (vector3 1 0 0 ) ) ; direction of k (length is irrelevant)
(define k (vector3-scale (* 2 pi fcen)
(unit-vector3 kdir))) ; k with correct length
(set! sources (list
(make source
(src (make continuous-src (frequency fcen) (fwidth df)))
(component Ex) (center 0 0 (* -0.49 sz)) (size sx sy 0)
(amp-func (pw-amp k (vector3 0 0 (* -0.49 sz) )))) ))
; -----------------------------------------------------------------------------
; Run and outputs -------------------------------------------------------------
(set! filename-prefix false)
(set! output-single-precision? true)
(define-param T 10) ; run until
(define (my-field) (output-real-field-function "e2" (list Ex Ey Ez )
(lambda (r ex ey ez ) (+ (* ex (conj ex)) (* ey (conj ey))
(* ez (conj ez))))))
(define (sum-h5 fname step-func)
(let ((first-step? true))
(lambda (to-do)
(if first-step?
(begin ; just copy the output file to fname
(set! first-step? false)
((convert-h5 true
(string-append "h5math -e \"d1\" " fname "")
step-func) to-do))
; otherwise, add the output file to fname
((convert-h5 true
(string-append "h5math -e \"d1 + d2 /71 \" " fname " " fname
"")
step-func) to-do)))))
(run-until T (at-beginning output-epsilon)
(after-time 3 (at-every 0.1 (sum-h5 "e2_avg.h5" my-field))))_______________________________________________
meep-discuss mailing list
[email protected]
http://ab-initio.mit.edu/cgi-bin/mailman/listinfo/meep-discuss
