Add new tests for particle absorption on EB with EM solver

Examples/Tests/CMakeLists.txt

@@ -15,6 +15,7 @@ add_subdirectory(electrostatic_sphere)
 add_subdirectory(electrostatic_sphere_eb)
 add_subdirectory(embedded_boundary_cube)
 add_subdirectory(embedded_boundary_diffraction)
+add_subdirectory(embedded_boundary_em_particle_absorption)
 add_subdirectory(embedded_boundary_python_api)
 add_subdirectory(embedded_boundary_rotated_cube)
 add_subdirectory(embedded_circle)

Examples/Tests/embedded_boundary_em_particle_absorption/CMakeLists.txt

@@ -0,0 +1,39 @@
+# Add tests (alphabetical order) ##############################################
+#
+
+if(WarpX_EB)
+    add_warpx_test(
+        test_3d_embedded_boundary_em_particle_absorption_sh_factor_1  # name
+        3  # dims
+        1  # nprocs
+        inputs_test_3d_sh_factor_1  # inputs
+        "analysis.py"  # analysis
+        "analysis_default_regression.py --path diags/diag1"  # checksum
+        OFF  # dependency
+    )
+endif()
+
+
+if(WarpX_EB)
+    add_warpx_test(
+        test_2d_embedded_boundary_em_particle_absorption_sh_factor_1  # name
+        2  # dims
+        1  # nprocs
+        inputs_test_2d_sh_factor_1  # inputs
+        "analysis.py"  # analysis
+        "analysis_default_regression.py --path diags/diag1"  # checksum
+        OFF  # dependency
+    )
+endif()
+
+if(WarpX_EB)
+    add_warpx_test(
+        test_rz_embedded_boundary_em_particle_absorption_sh_factor_1  # name
+        RZ  # dims
+        1  # nprocs
+        inputs_test_rz_sh_factor_1  # inputs
+        "analysis.py"  # analysis
+        "analysis_default_regression.py --path diags/diag1"  # checksum
+        OFF  # dependency
+    )
+endif()

Examples/Tests/embedded_boundary_em_particle_absorption/analysis.py

@@ -0,0 +1,40 @@
+#!/usr/bin/env python
+
+# TODO: update description: explain two particles separate, absorbed by the boundary
+# with EM solver, there can be a static spurious divE (spurious charge build-up) that remains at the position where particle was absorbed
+# this checks this is not the case. Since divE also has a (physical) component that propagate along the EB (due to EM waves reflecting),
+# we average in time to remove this component and only check the static part.
+# TODO: Adjust tolerances
+"""
+This analysis script checks for any spurious charge build-up at the embedded boundary, when particles are removed in 3D.
+It averages the divergence of the electric field (divE) over the last 20 time steps and compares the results with a specified tolerance.
+"""
+
+from openpmd_viewer import OpenPMDTimeSeries
+
+ts = OpenPMDTimeSeries("./diags/diag1/")
+
+divE_stacked = ts.iterate(
+    lambda iteration: ts.get_field("divE", iteration=iteration)[0]
+)
+start_avg_iter = 32
+end_avg_iter = 100
+divE_avg = divE_stacked[start_avg_iter:end_avg_iter].mean(axis=0)
+
+dim = ts.fields_metadata["divE"]["geometry"]
+if dim == "3dcartesian":
+    tolerance = 5e-11
+elif dim == "2dcartesian":
+    tolerance = 3e-10
+elif dim == "thetaMode":
+    tolerance = 3e-11
+
+
+def check_tolerance(array, tolerance):
+    assert abs(array).max() <= tolerance, (
+        f"Test did not pass: the max error {abs(array).max()} exceeded the tolerance of {tolerance}."
+    )
+    print("All elements of are within the tolerance.")
+
+
+check_tolerance(divE_avg, tolerance)

Examples/Tests/embedded_boundary_em_particle_absorption/inputs_base_2d

@@ -0,0 +1,45 @@
+max_step = 100
+
+geometry.dims = 2
+
+amr.n_cell = 32 32
+amr.max_level = 0
+amr.blocking_factor = 8
+amr.max_grid_size = 256
+geometry.prob_lo     = -10 -10
+geometry.prob_hi     =  10  10
+
+
+# Boundary condition
+boundary.field_lo = pec pec
+boundary.field_hi = pec pec
+
+algo.charge_deposition = standard
+algo.field_gathering = energy-conserving
+warpx.use_filter = 0
+warpx.const_dt = 1.203645751e-09
+warpx.eb_implicit_function = "(x**2 + y**2 + z**2 - 56.25)"
+
+particles.species_names = electron positron
+
+electron.charge = -q_e
+electron.mass = m_e
+electron.injection_style = "SingleParticle"
+electron.single_particle_pos = 0.0   0.0 0.0
+electron.single_particle_u = 1.e20   0.0 0.3e20  # gamma*beta
+electron.single_particle_weight = 1.0
+
+positron.charge = q_e
+positron.mass = m_e
+positron.injection_style = "SingleParticle"
+positron.single_particle_pos = 0.0   0.0 0.0
+positron.single_particle_u = -1.e20   0.0 -0.3e20  # gamma*beta
+positron.single_particle_weight = 1.0
+
+
+# Diagnostics
+diagnostics.diags_names = diag1
+diag1.intervals = 1
+diag1.diag_type = Full
+diag1.fields_to_plot = divE rho
+diag1.format = openpmd

Examples/Tests/embedded_boundary_em_particle_absorption/inputs_base_3d

@@ -0,0 +1,43 @@
+max_step = 100
+amr.n_cell = 32 32 32
+amr.max_level = 0
+amr.blocking_factor = 8
+amr.max_grid_size = 256
+geometry.dims = 3
+geometry.prob_lo     = -10 -10 -10
+geometry.prob_hi     =  10  10  10
+
+# Boundary condition
+boundary.field_lo = pec pec pec
+boundary.field_hi = pec pec pec
+
+algo.charge_deposition = standard
+algo.field_gathering = energy-conserving
+warpx.cfl = 1.0
+warpx.use_filter = 0
+
+warpx.eb_implicit_function = "(x**2 + y**2 + z**2 - 56.25)"
+
+particles.species_names = electron positron
+
+electron.charge = -q_e
+electron.mass = m_e
+electron.injection_style = "SingleParticle"
+electron.single_particle_pos = 0.0   0.0 0.0
+electron.single_particle_u = 1.e20   0.0 0.3e20  # gamma*beta
+electron.single_particle_weight = 1.0
+
+positron.charge = q_e
+positron.mass = m_e
+positron.injection_style = "SingleParticle"
+positron.single_particle_pos = 0.0   0.0 0.0
+positron.single_particle_u = -1.e20   0.0 -0.3e20  # gamma*beta
+positron.single_particle_weight = 1.0
+
+
+# Diagnostics
+diagnostics.diags_names = diag1
+diag1.intervals = 1
+diag1.diag_type = Full
+diag1.fields_to_plot = divE rho
+diag1.format = openpmd

Examples/Tests/embedded_boundary_em_particle_absorption/inputs_base_rz

@@ -0,0 +1,43 @@
+max_step = 100
+
+geometry.dims = RZ
+
+amr.n_cell = 32 32
+amr.max_level = 0
+amr.blocking_factor = 8
+amr.max_grid_size = 256
+geometry.prob_lo     = 0 -10
+geometry.prob_hi     =  10  10
+
+# Boundary condition
+boundary.field_lo = none pec
+boundary.field_hi = pec pec
+
+algo.charge_deposition = standard
+algo.field_gathering = energy-conserving
+warpx.use_filter = 0
+
+warpx.eb_implicit_function = "(x**2 + y**2 + z**2 - 56.25)"
+
+particles.species_names = electron positron
+
+electron.charge = -q_e
+electron.mass = m_e
+electron.injection_style = "SingleParticle"
+electron.single_particle_pos = 1.0 0.0 0.0
+electron.single_particle_u = 0.0  0.0 0.3e20  # gamma*beta
+electron.single_particle_weight = 1.0
+
+positron.charge = q_e
+positron.mass = m_e
+positron.injection_style = "SingleParticle"
+positron.single_particle_pos = 1.0 0.0 0.0
+positron.single_particle_u = 0.0  0.0 -0.3e20  # gamma*beta
+positron.single_particle_weight = 1.0
+
+# Diagnostics
+diagnostics.diags_names = diag1
+diag1.intervals = 1
+diag1.diag_type = Full
+diag1.fields_to_plot = divE rho
+diag1.format = openpmd

Examples/Tests/embedded_boundary_em_particle_absorption/inputs_test_2d_sh_factor_1

@@ -0,0 +1,6 @@
+# base input parameters
+FILE = inputs_base_2d
+
+# test input parameters
+# Order of particle shape factors
+algo.particle_shape = 1

Examples/Tests/embedded_boundary_em_particle_absorption/inputs_test_3d_sh_factor_1

@@ -0,0 +1,6 @@
+# base input parameters
+FILE = inputs_base_3d
+
+# test input parameters
+# Order of particle shape factors
+algo.particle_shape = 1

Examples/Tests/embedded_boundary_em_particle_absorption/inputs_test_rz_sh_factor_1

@@ -0,0 +1,6 @@
+# base input parameters
+FILE = inputs_base_rz
+
+# test input parameters
+# Order of particle shape factors
+algo.particle_shape = 1

Regression/Checksum/benchmarks_json/test_2d_embedded_boundary_em_particle_absorption_sh_factor_1.json

@@ -0,0 +1,24 @@
+{
+  "electron": {
+    "particle_momentum_x": 0.0,
+    "particle_momentum_y": 0.0,
+    "particle_momentum_z": 0.0,
+    "particle_position_x": 0.0,
+    "particle_position_y": 0.0,
+    "particle_position_z": 0.0,
+    "particle_weight": 0.0
+  },
+  "lev=0": {
+    "divE": 1.6305300553737757e-07,
+    "rho": 0.0
+  },
+  "positron": {
+    "particle_momentum_x": 0.0,
+    "particle_momentum_y": 0.0,
+    "particle_momentum_z": 0.0,
+    "particle_position_x": 0.0,
+    "particle_position_y": 0.0,
+    "particle_position_z": 0.0,
+    "particle_weight": 0.0
+  }
+}

Regression/Checksum/benchmarks_json/test_3d_embedded_boundary_em_particle_absorption_sh_factor_1.json

@@ -0,0 +1,24 @@
+{
+  "electron": {
+    "particle_momentum_x": 0.0,
+    "particle_momentum_y": 0.0,
+    "particle_momentum_z": 0.0,
+    "particle_position_x": 0.0,
+    "particle_position_y": 0.0,
+    "particle_position_z": 0.0,
+    "particle_weight": 0.0
+  },
+  "lev=0": {
+    "divE": 8.149809003960999e-07,
+    "rho": 0.0
+  },
+  "positron": {
+    "particle_momentum_x": 0.0,
+    "particle_momentum_y": 0.0,
+    "particle_momentum_z": 0.0,
+    "particle_position_x": 0.0,
+    "particle_position_y": 0.0,
+    "particle_position_z": 0.0,
+    "particle_weight": 0.0
+  }
+}

Regression/Checksum/benchmarks_json/test_rz_embedded_boundary_em_particle_absorption_sh_factor_1.json

@@ -0,0 +1,24 @@
+{
+  "electron": {
+    "particle_momentum_x": 0.0,
+    "particle_momentum_y": 0.0,
+    "particle_momentum_z": 0.0,
+    "particle_position_x": 0.0,
+    "particle_position_y": 0.0,
+    "particle_position_z": 0.0,
+    "particle_weight": 0.0
+  },
+  "lev=0": {
+    "divE": 5.999732410984964e-08,
+    "rho": 0.0
+  },
+  "positron": {
+    "particle_momentum_x": 0.0,
+    "particle_momentum_y": 0.0,
+    "particle_momentum_z": 0.0,
+    "particle_position_x": 0.0,
+    "particle_position_y": 0.0,
+    "particle_position_z": 0.0,
+    "particle_weight": 0.0
+  }
+}