How to specify the propagation direction in Mie scattering of a gold nanosphere?

[rolly-ng]

Hi, I would like to repeat the following example as given in the discussion and the tutorial.
https://github.com/NanoComp/meep/discussions/1911#discussioncomment-1999005
https://github.com/NanoComp/meep/issues/1530#issuecomment-1002192047
https://meep.readthedocs.io/en/latest/Scheme_Tutorials/Basics/#mie-scattering-of-a-lossless-dielectric-sphere

Here is my ctl file,

(define-param r 0.025) ;; radius of sphere

(define wvl-min (* 2 pi r 2))
(define wvl-max (* 2 pi r 20))

(define frq-min (/ wvl-max))
(define frq-max (/ wvl-min))
(define frq-cen (* 0.5 (+ frq-min frq-max)))
(define dfrq (- frq-max frq-min))
(define nfrq 1024)

;; at least 8 pixels per smallest wavelength, i.e. (floor (/ 8 wvl-min))
(set-param! resolution 250)

(define dpml (* 0.25 wvl-min))
(define dair (* 0.25 wvl-min))

(define boundary-layers (list (make pml (thickness dpml))))
(set! pml-layers boundary-layers)

;;(define symm (list (make mirror-sym (direction Y))
;; (make mirror-sym (direction Z) (phase -1))))
;;(set! symmetries symm)

(define s (* 2 (+ dpml dair r)))
(define cell (make lattice (size s s s)))
(set! geometry-lattice cell)

;; (is-integrated? true) necessary for any planewave source extending into PML
(define pw-src (make source
(src (make gaussian-src (frequency frq-cen) (fwidth dfrq) (is-integrated? true)))
(center (+ (* -0.5 s) dpml) 0 0)
(size 0 s s)
(component Ez)))
(set! sources (list pw-src))

(set! k-point (vector3 0 0 0))

(define dpad (* 0.5 dair))

(define box-x1 (add-flux frq-cen dfrq nfrq
(make flux-region (center (- r dpad) 0 0) (size 0 (* 2 (+ r dpad)) (* 2 (+ r dpad))))))
(define box-x2 (add-flux frq-cen dfrq nfrq
(make flux-region (center (+ r dpad) 0 0) (size 0 (* 2 (+ r dpad)) (* 2 (+ r dpad))))))
(define box-y1 (add-flux frq-cen dfrq nfrq
(make flux-region (center 0 (- r dpad) 0) (size (* 2 (+ r dpad)) 0 (* 2 (+ r dpad))))))
(define box-y2 (add-flux frq-cen dfrq nfrq
(make flux-region (center 0 (+ r dpad) 0) (size (* 2 (+ r dpad)) 0 (* 2 (+ r dpad))))))
(define box-z1 (add-flux frq-cen dfrq nfrq
(make flux-region (center 0 0 (- r dpad)) (size (* 2 (+ r dpad)) (* 2 (+ r dpad)) 0))))
(define box-z2 (add-flux frq-cen dfrq nfrq
(make flux-region (center 0 0 (+ r dpad)) (size (* 2 (+ r dpad)) (* 2 (+ r dpad)) 0))))

(run-sources+ 100)

(display-fluxes box-x1)

(save-flux "box-x1-flux" box-x1)
(save-flux "box-x2-flux" box-x2)
(save-flux "box-y1-flux" box-y1)
(save-flux "box-y2-flux" box-y2)
(save-flux "box-z1-flux" box-z1)
(save-flux "box-z2-flux" box-z2)

(reset-meep)

(set! geometry (list
(make sphere
(material Au)
(radius r)
(center 0))))

(set! geometry-lattice cell)

(set! pml-layers boundary-layers)

;;(set! symmetries symm)

(set! sources (list pw-src))

(set! k-point (vector3 0))

(define box-x1 (add-flux frq-cen dfrq nfrq
(make flux-region (center (- r dpad) 0 0) (size 0 (* 2 (+ r dpad)) (* 2 (+ r dpad))))))
(define box-x2 (add-flux frq-cen dfrq nfrq
(make flux-region (center (+ r dpad) 0 0) (size 0 (* 2 (+ r dpad)) (* 2 (+ r dpad))))))
(define box-y1 (add-flux frq-cen dfrq nfrq
(make flux-region (center 0 (- r dpad) 0) (size (* 2 (+ r dpad)) 0 (* 2 (+ r dpad))))))
(define box-y2 (add-flux frq-cen dfrq nfrq
(make flux-region (center 0 (+ r dpad) 0) (size (* 2 (+ r dpad)) 0 (* 2 (+ r dpad))))))
(define box-z1 (add-flux frq-cen dfrq nfrq
(make flux-region (center 0 0 (- r dpad)) (size (* 2 (+ r dpad)) (* 2 (+ r dpad)) 0))))
(define box-z2 (add-flux frq-cen dfrq nfrq
(make flux-region (center 0 0 (+ r dpad)) (size (* 2 (+ r dpad)) (* 2 (+ r dpad)) 0))))

(load-minus-flux "box-x1-flux" box-x1)
(load-minus-flux "box-x2-flux" box-x2)
(load-minus-flux "box-y1-flux" box-y1)
(load-minus-flux "box-y2-flux" box-y2)
(load-minus-flux "box-z1-flux" box-z1)
(load-minus-flux "box-z2-flux" box-z2)

(run-sources+ 1000)

(display-fluxes box-x1 box-x2 box-y1 box-y2 box-z1 box-z2)

The result is reasonable with a peak at 0.3 which is 530nm.

However, if I change the propagation direction to the Z direction, that is
(define pw-src (make source
(src (make gaussian-src (frequency frq-cen) (fwidth dfrq) (is-integrated? true)))
(center 0 0 (+ (* -0.5 s) dpml))
(size s s 0)
(component Ey)))

I have incorrect flux1, which is the empty geometry without the Au nanosphere, and these are the few lines of my flux1,

Correct reading with component Ez:
...
run 0 finished at t = 103.492 (51746 timesteps)
flux1:, 0.3183098861837905, 2.1315036348733428e-10
flux1:, 0.3211102664141465, 2.2554234889234997e-10
flux1:, 0.32391064664450236, 2.3859559185299824e-10
...

Incorrect reading with component Ey:
...
run 0 finished at t = 103.492 (51746 timesteps)
flux1:, 0.3183098861837905, 1.6209939379924598e-25
flux1:, 0.3211102664141465, -7.144624787788959e-25
flux1:, 0.32391064664450236, -7.108705197512597e-25
...

But, my flux2 seems correct in both cases?
With component Ez,
...
run 1 finished at t = 1003.4920000000001 (501746 timesteps)
flux2:, 0.3183098861837905, 8.795193507679312e-13, -1.865852489222941e-15, -1.838384679228736e-15, 1.0378996540777725e-16, -2.2768462855503596e-15, 3.388196953727363e-17
flux2:, 0.3211102664141465, 9.319988233661141e-13, -1.9745950992522867e-15, -1.966660001284201e-15, 1.137226468387698e-16, -2.432314903349079e-15, 3.715067430317086e-17
...
With component Ey,
...
run 1 finished at t = 1003.4920000000001 (501746 timesteps)
flux2:, 0.3183098861837905, -1.8383846695536913e-15, 1.0378996537395662e-16, -2.276846285258645e-15, 3.3881969540043714e-17, 8.795193507739996e-13, -1.8658524892331504e-15
flux2:, 0.3211102664141465, -1.9666599911590894e-15, 1.1372264684522782e-16, -2.432314903313549e-15, 3.715067427083348e-17, 9.31998823373213e-13, -1.9745950992002957e-15
...

So It seems changing the direction will have impact on the empty case reading only? What did I miss in my code?

Thanks,
rolly

[rolly-ng]

Just realized that this line has to be adjusted to save the Z direction flux.

(display-fluxes box-z1)

Problem solved and thanks for your attention!