Electronic Boltzmann physics [WIP]

bolt/src/electronic_boltzmann/advection_terms.py

@@ -0,0 +1,32 @@
+"""
+Here we define the advection terms for the 
+nonrelativistic Boltzmann equation.
+"""
+
+
+def A_q(p1, p2, p3, params):
+    """Return the terms A_q1, A_q2."""
+    return (p1, p2)
+
+# If necessary, additional terms as a function of the arguments
+# passed to A_p may be used:s
+# For instance:
+def T1(q1, q2, p1, p2, p3):
+    return(q1*q2)
+
+# This can then be called inside A_p if needed:
+# F1 = (params.char....)(E1 + ....) + T1(q1, q2, p1, p2, p3)
+
+def A_p(q1, q2, p1, p2, p3,
+        E1, E2, E3, B1, B2, B3,
+        params
+       ):
+    """Return the terms A_p1, A_p2 and A_p3."""
+    F1 =   (params.charge_electron / params.mass_particle) \
+         * (E1 + p2 * B3 - p3 * B2)
+    F2 =   (params.charge_electron / params.mass_particle) \
+         * (E2 - p1 * B3 + p3 * B1)
+    F3 =   (params.charge_electron / params.mass_particle) \
+         * (E3 - p2 * B1 + p1 * B2)
+
+    return (F1, F2, F3)

bolt/src/electronic_boltzmann/collision_operator.py

@@ -0,0 +1,56 @@
+"""Contains the function which returns the Source/Sink term."""
+
+
+import numpy as np
+import arrayfire as af
+
+# Using af.broadcast, since p1, p2, p3 are of size (1, 1, Np1*Np2*Np3)
+# All moment quantities are of shape (Nq1, Nq2)
+# By wrapping with af.broadcast, we can perform batched operations
+# on arrays of different sizes.
+@af.broadcast
+def f0(p1, p2, p3, n, T, p1_bulk, p2_bulk, p3_bulk, params):
+    """Return the Local MB distribution."""
+    m = params.mass_particle
+    k = params.boltzmann_constant
+
+    if (params.p_dim == 3):
+        f0 = n * (m / (2 * np.pi * k * T))**(3 / 2)  \
+               * af.exp(-m * (p1 - p1_bulk)**2 / (2 * k * T)) \
+               * af.exp(-m * (p2 - p2_bulk)**2 / (2 * k * T)) \
+               * af.exp(-m * (p3 - p3_bulk)**2 / (2 * k * T))
+
+    elif (params.p_dim == 2):
+        f0 = n * (m / (2 * np.pi * k * T)) \
+               * af.exp(-m * (p1 - p1_bulk)**2 / (2 * k * T)) \
+               * af.exp(-m * (p2 - p2_bulk)**2 / (2 * k * T))
+
+    else:
+        f0 = n * af.sqrt(m / (2 * np.pi * k * T)) \
+               * af.exp(-m * (p1 - p1_bulk)**2 / (2 * k * T))
+
+    af.eval(f0)
+    return (f0)
+
+
+def BGK(f, q1, q2, p1, p2, p3, moments, params):
+    """Return BGK operator -(f-f0)/tau."""
+    n = moments('density')
+
+    p1_bulk = moments('mom_p1_bulk') / n
+    p2_bulk = moments('mom_p2_bulk') / n
+    p3_bulk = moments('mom_p3_bulk') / n
+
+    T =   (1 / params.p_dim) \
+        * (  moments('energy') 
+           - n * p1_bulk**2
+           - n * p2_bulk**2
+           - n * p3_bulk**2
+          ) / n
+
+    C_f = -(  f
+            - f0(p1, p2, p3, n, T, p1_bulk, p2_bulk, p3_bulk, params)
+           ) / params.tau(q1, q2, p1, p2, p3)
+
+    af.eval(C_f)
+    return(C_f)

bolt/src/electronic_boltzmann/moment_defs.py

@@ -0,0 +1,19 @@
+moment_exponents = dict(density      = [0, 0, 0],
+                        mom_p1_bulk  = [1, 0, 0],
+                        mom_p2_bulk  = [0, 1, 0],
+                        mom_p3_bulk  = [0, 0, 1],
+                        energy       = [2, 2, 2],
+                        q_q1         = [[1, 0, 0], [2, 2, 2]],
+                        q_q2         = [[0, 1, 0], [2, 2, 2]],
+                        q_q3         = [[0, 0, 1], [2, 2, 2]]
+                        )
+
+moment_coeffs = dict(density     = [1, 0, 0],
+                     mom_p1_bulk = [1, 0, 0],
+                     mom_p2_bulk = [0, 1, 0],
+                     mom_p3_bulk = [0, 0, 1],
+                     energy      = [1, 1, 1],
+                     q_q1        = [[1, 0, 0], [1, 1, 1]],
+                     q_q2        = [[0, 1, 0], [1, 1, 1]],
+                     q_q3        = [[0, 0, 1], [1, 1, 1]]
+                     )