Use params dictionary

Airy_2d/Airy2d.py

@@ -137,6 +137,8 @@ def main(args):
     shift = source_shift + r_w
     chi_rad = math.radians(chi_deg)
 
+    params = dict(W_y=w_0, k=k1, M=M, W=W)
+
     # --------------------------------------------------------------------------
     # placement of the dielectric interface within the computational cell
     # --------------------------------------------------------------------------
@@ -172,12 +174,16 @@ def Delta_x(alpha):
     # --------------------------------------------------------------------------
     # beam profile distribution (field amplitude) at the waist of the beam
     # --------------------------------------------------------------------------
-    def Gauss(r, W_y=w_0):
+    def Gauss(r, params):
         """Gauss profile."""
+        W_y = params['W_y']
+
         return math.exp(-(r.y / W_y)**2)
 
-    def Ai_inc(r, W_y=w_0, M=M, W=W):
+    def Ai_inc(r, params):
         """Incomplete Airy function."""
+        W_y, M, W = params['W_y'], params['M'], params['W']
+
         (result,
          real_tol,
          imag_tol) = complex_quad(lambda xi:
@@ -186,8 +192,8 @@ def Ai_inc(r, W_y=w_0, M=M, W=W):
         return result
 
     if test_output:
-        print("w_0:", w_0)
-        print("Airy function 1:", Ai_inc(mp.Vector3(1, -0.3, 1), w_0, 0, 4))
+        print("w_0:", params['W_y'])
+        print("Airy function 1:", Ai_inc(mp.Vector3(1, -0.3, 1), params))
 
     # --------------------------------------------------------------------------
     # spectrum amplitude distribution
@@ -196,24 +202,28 @@ def Heaviside(x):
         """Theta (Heaviside step) function."""
         return 0 if x < 0 else 1
 
-    def f_Gauss(k_y, W_y=w_0):
+    def f_Gauss(k_y, params):
         """Gaussian spectrum amplitude."""
+        W_y = params['W_y']
+
         return math.exp(-(k_y*W_y/2)**2)
 
-    def f_Airy(k_y, W_y=w_0, M=M, W=W):
+    def f_Airy(k_y, params):
         """Airy spectrum amplitude."""
+        W_y, M, W = params['W_y'], params['M'], params['W']
+
         return W_y*cmath.exp(1.0j*(-1/3)*(k_y*W_y)**3) \
             * Heaviside(W_y*k_y - (M-W)) * Heaviside((M+W) - (W_y*k_y))
 
     if test_output:
-        print("Airy spectrum:", f_Airy(0.2, w_0, 0, 4))
+        print("Airy spectrum:", f_Airy(0.2, params))
 
     # --------------------------------------------------------------------------
     # plane wave decomposition
     # (purpose: calculate field amplitude at light source position if not
     #           coinciding with beam waist)
     # --------------------------------------------------------------------------
-    def psi(r, f, x):
+    def psi(r, f, x, params):
         """Field amplitude function."""
         try:
             getattr(psi, "called")
@@ -230,7 +240,8 @@ def phase(k_y, x, y):
             (result,
              real_tol,
              imag_tol) = complex_quad(lambda k_y:
-                                      f(k_y) * cmath.exp(1j*phase(k_y, x, r.y)),
+                                      f(k_y, params) *
+                                      cmath.exp(1j*phase(k_y, x, r.y)),
                                       -k1, k1, limit=100)
         except Exception as e:
             print(type(e).__name__ + ":", e)
@@ -246,7 +257,7 @@ def phase(k_y, x, y):
         r = mp.Vector3(0, y, z)
 
         print()
-        print("psi :", psi(r, f_Airy, x))
+        print("psi :", psi(r, f_Airy, x, params))
         sys.exit()
 
     # --------------------------------------------------------------------------
@@ -275,10 +286,10 @@ def phase(k_y, x, y):
                          component=mp.Ez if s_pol else mp.Ey,
                          size=mp.Vector3(0, 9, 0),
                          center=mp.Vector3(source_shift, 0, 0),
-                         #amp_func=lambda r: Gauss(r, w_0)
-                         #amp_func=lambda r: Ai_inc(r, w_0, M, W)
-                         #amp_func=lambda r: psi(r, lambda k_y: f_Gauss(k_y, w_0), shift)
-                         amp_func=lambda r: psi(r, lambda k_y: f_Airy(k_y, w_0, M, W), shift)
+                         #amp_func=lambda r: Gauss(r, params)
+                         #amp_func=lambda r: Ai_inc(r, params)
+                         #amp_func=lambda r: psi(r, f_Gauss, shift, params)
+                         amp_func=lambda r: psi(r, f_Airy, shift, params)
                         )
                ]