From b8e49e1a6c99c9aae0cfbdf43cf920bdc6d1e00d Mon Sep 17 00:00:00 2001 From: Olga Shapoval Date: Fri, 17 Jan 2025 19:05:06 -0800 Subject: [PATCH 01/13] Added CI to test secondary ion emission in RZ. --- Examples/Tests/CMakeLists.txt | 1 + .../secondary_ion_emission/CMakeLists.txt | 14 + .../Tests/secondary_ion_emission/analysis.py | 62 +++++ .../analysis_default_regression.py | 1 + ...ts_test_rz_secondary_ion_emission_picmi.py | 255 ++++++++++++++++++ .../test_rz_secondary_ion_emission_picmi.json | 26 ++ 6 files changed, 359 insertions(+) create mode 100644 Examples/Tests/secondary_ion_emission/CMakeLists.txt create mode 100644 Examples/Tests/secondary_ion_emission/analysis.py create mode 120000 Examples/Tests/secondary_ion_emission/analysis_default_regression.py create mode 100644 Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py create mode 100644 Regression/Checksum/benchmarks_json/test_rz_secondary_ion_emission_picmi.json diff --git a/Examples/Tests/CMakeLists.txt b/Examples/Tests/CMakeLists.txt index c4713123ae6..7d67666bb0b 100644 --- a/Examples/Tests/CMakeLists.txt +++ b/Examples/Tests/CMakeLists.txt @@ -70,6 +70,7 @@ add_subdirectory(resampling) add_subdirectory(restart) add_subdirectory(restart_eb) add_subdirectory(rigid_injection) +add_subdirectory(secondary_ion_emission) add_subdirectory(scraping) add_subdirectory(effective_potential_electrostatic) add_subdirectory(silver_mueller) diff --git a/Examples/Tests/secondary_ion_emission/CMakeLists.txt b/Examples/Tests/secondary_ion_emission/CMakeLists.txt new file mode 100644 index 00000000000..e6e38138a08 --- /dev/null +++ b/Examples/Tests/secondary_ion_emission/CMakeLists.txt @@ -0,0 +1,14 @@ +# Add tests (alphabetical order) ############################################## +# + +if(WarpX_EB) + add_warpx_test( + test_rz_secondary_ion_emission_picmi # name + RZ # dims + 1 # nprocs + inputs_test_rz_secondary_ion_emission_picmi.py # inputs + "analysis.py diags/diag1/" # analysis + "analysis_default_regression.py --path diags/diag1/" # checksum + OFF # dependency + ) +endif() diff --git a/Examples/Tests/secondary_ion_emission/analysis.py b/Examples/Tests/secondary_ion_emission/analysis.py new file mode 100644 index 00000000000..7f43b8031b8 --- /dev/null +++ b/Examples/Tests/secondary_ion_emission/analysis.py @@ -0,0 +1,62 @@ +#!/usr/bin/env python +""" +This script tests the last coordinates of the emitted secondary electrons. +The EB sphere is centered on O and has a radius of 0.2. +The proton is initially at: (0,0,-0.25) and moves with a velocity: +(0.1e6,0,1.5e6) with a time step of dt = 0.000000015. +The simulation uses a fixed random seed (np.random.seed(10025015)) +to ensure the emission of secodnary electrons. +An input file inputs_test_rz_secondary_ion_emission_picmi.py is used. +""" + +import sys + +import numpy as np +import yt +from openpmd_viewer import OpenPMDTimeSeries + +yt.funcs.mylog.setLevel(0) + +# Open plotfile specified in command line +filename = sys.argv[1] +ts = OpenPMDTimeSeries(filename) + +it = ts.iterations +x, y, z = ts.get_particle(["x", "y", "z"], species="electrons", iteration=it[-1]) +print("x", x) +print("y", y) +print("z", z) +# Analytical results calculated +x_analytic = [0.004028, 0.003193] +y_analytic = [-0.0001518, -0.0011041] +z_analytic = [-0.19967, -0.19926] + +N_sec_e = np.size(z) # number of the secondary electrons + +assert N_sec_e == 2, ( + "Test did not pass: for this set up we expect 2 secondary electrons emitted" +) + +tolerance = 1e-3 + +for i in range(0, N_sec_e): + print("\n") + print(f"Electron # {i}:") + print("NUMERICAL coordinates of the emitted electrons:") + print("x=%5.5f, y=%5.5f, z=%5.5f" % (x[i], y[i], z[i])) + print("\n") + print("ANALYTICAL coordinates of the point of contact:") + print("x=%5.5f, y=%5.5f, z=%5.5f" % (x_analytic[i], y_analytic[i], z_analytic[i])) + + diff_x = np.abs((x[i] - x_analytic[i]) / x_analytic[i]) + diff_y = np.abs((y[i] - y_analytic[i]) / y_analytic[i]) + diff_z = np.abs((z[i] - z_analytic[i]) / z_analytic[i]) + + print("\n") + print("percentage error for x = %5.4f %%" % (diff_x * 100)) + print("percentage error for y = %5.4f %%" % (diff_y * 100)) + print("percentage error for z = %5.4f %%" % (diff_z * 100)) + + assert (diff_x < tolerance) and (diff_y < tolerance) and (diff_z < tolerance), ( + "Test particle_boundary_interaction did not pass" + ) diff --git a/Examples/Tests/secondary_ion_emission/analysis_default_regression.py b/Examples/Tests/secondary_ion_emission/analysis_default_regression.py new file mode 120000 index 00000000000..d8ce3fca419 --- /dev/null +++ b/Examples/Tests/secondary_ion_emission/analysis_default_regression.py @@ -0,0 +1 @@ +../../analysis_default_regression.py \ No newline at end of file diff --git a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py new file mode 100644 index 00000000000..5849ba4c59e --- /dev/null +++ b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py @@ -0,0 +1,255 @@ +#!/usr/bin/env python3 +# --- Input file for secondary-ion emission testing in RZ via callback function. +import numpy as np +from scipy.constants import e, elementary_charge, m_e, proton_mass + +from pywarpx import callbacks, particle_containers, picmi + +########################## +# numerics parameters +########################## + +dt = 0.000000015 + +# --- Nb time steps +Te = 0.0259 # in eV +dist_th = np.sqrt(Te * elementary_charge / m_e) + +max_steps = 3 +diagnostic_interval = 1 + +# --- grid +nr = 64 +nz = 64 + +rmin = 0.0 +rmax = 2 +zmin = -2 +zmax = 2 +delta_H = 0.4 +E_HMax = 250 + +np.random.seed(10025015) +########################## +# numerics components +########################## + +grid = picmi.CylindricalGrid( + number_of_cells=[nr, nz], + n_azimuthal_modes=1, + lower_bound=[rmin, zmin], + upper_bound=[rmax, zmax], + lower_boundary_conditions=["none", "dirichlet"], + upper_boundary_conditions=["dirichlet", "dirichlet"], + lower_boundary_conditions_particles=["none", "reflecting"], + upper_boundary_conditions_particles=["absorbing", "reflecting"], +) + +solver = picmi.ElectrostaticSolver( + grid=grid, method="Multigrid", warpx_absolute_tolerance=1e-7 +) + +embedded_boundary = picmi.EmbeddedBoundary( + implicit_function="-(x**2+y**2+z**2-radius**2)", radius=0.2 +) + +########################## +# physics components +########################## + +i_dist = picmi.ParticleListDistribution( + x=0.0, y=0.0, z=-0.25, ux=0.1e6, uy=0.0, uz=1.50e6, weight=1 +) +electrons = picmi.Species( + particle_type="electron", # Specify the particle type + name="electrons", # Name of the species +) + +ions = picmi.Species( + name="ions", + mass=proton_mass, + charge=e, + initial_distribution=i_dist, + warpx_save_particles_at_eb=1, +) + +########################## +# diagnostics +########################## + +field_diag = picmi.FieldDiagnostic( + name="diag1", + grid=grid, + period=diagnostic_interval, + data_list=["Er", "Ez", "phi", "rho"], # , "rho_electrons"], + warpx_format="openpmd", +) + +part_diag = picmi.ParticleDiagnostic( + name="diag1", + period=diagnostic_interval, + species=[ions, electrons], + warpx_format="openpmd", +) + +########################## +# simulation setup +########################## +# num_procs = [1, 1] +# warpx_numprocs=num_procs, +sim = picmi.Simulation( + solver=solver, + time_step_size=dt, + max_steps=max_steps, + warpx_embedded_boundary=embedded_boundary, + warpx_amrex_the_arena_is_managed=1, +) + +sim.add_species( + electrons, + layout=picmi.GriddedLayout(n_macroparticle_per_cell=[0, 0, 0], grid=grid), +) + +sim.add_species( + ions, + layout=picmi.GriddedLayout(n_macroparticle_per_cell=[10, 1, 1], grid=grid), +) + +sim.add_diagnostic(part_diag) +sim.add_diagnostic(field_diag) + +sim.initialize_inputs() +sim.initialize_warpx() + +########################## +# python particle data access +########################## + + +def concat(list_of_arrays): + if len(list_of_arrays) == 0: + # Return a 1d array of size 0 + return np.empty(0) + else: + return np.concatenate(list_of_arrays) + + +def sigma_nascap(energy_kEv, delta_H, E_HMax): + """ + Compute sigma_nascap for each element in the energy array using a loop. + + Parameters: + - energy: ndarray or list, energy values in KeV + - delta_H: float, parameter for the formula + - E_HMax: float, parameter for the formula in KeV + + Returns: + - numpy array, computed probability sigma_nascap + """ + sigma_nascap = np.array([]) + # Loop through each energy value + for energy in energy_kEv: + if energy > 0.0: + sigma = ( + delta_H + * (E_HMax + 1.0) + / (E_HMax * 1.0 + energy) + * np.sqrt(energy / 1.0) + ) + else: + sigma = 0.0 + sigma_nascap = np.append(sigma_nascap, sigma) + return sigma_nascap + + +def secondary_emission(): + buffer = particle_containers.ParticleBoundaryBufferWrapper() # boundary buffer + # STEP 1: extract the different parameters of the boundary buffer (normal, time, position) + lev = 0 # level 0 (no mesh refinement here) + n = buffer.get_particle_boundary_buffer_size("ions", "eb") + elect_pc = particle_containers.ParticleContainerWrapper("electrons") + + if n != 0: + r = concat(buffer.get_particle_boundary_buffer("ions", "eb", "x", lev)) + theta = concat(buffer.get_particle_boundary_buffer("ions", "eb", "theta", lev)) + z = concat(buffer.get_particle_boundary_buffer("ions", "eb", "z", lev)) + x = r * np.cos(theta) # from RZ coordinates to 3D coordinates + y = r * np.sin(theta) + ux = concat(buffer.get_particle_boundary_buffer("ions", "eb", "ux", lev)) + uy = concat(buffer.get_particle_boundary_buffer("ions", "eb", "uy", lev)) + uz = concat(buffer.get_particle_boundary_buffer("ions", "eb", "uz", lev)) + w = concat(buffer.get_particle_boundary_buffer("ions", "eb", "w", lev)) + nx = concat(buffer.get_particle_boundary_buffer("ions", "eb", "nx", lev)) + ny = concat(buffer.get_particle_boundary_buffer("ions", "eb", "ny", lev)) + nz = concat(buffer.get_particle_boundary_buffer("ions", "eb", "nz", lev)) + delta_t = concat( + buffer.get_particle_boundary_buffer("ions", "eb", "deltaTimeScraped", lev) + ) + energy_ions = 0.5 * proton_mass * w * (ux**2 + uy**2 + uz**2) + energy_ions_in_kEv = energy_ions / (1.602176634e-19 * 1000) + sigma_nascap_ions = sigma_nascap(energy_ions_in_kEv, delta_H, E_HMax) + # Loop over all ions in the EB buffer + for i in range(0, n): + sigma = sigma_nascap_ions[i] + sigma_int = int(sigma) + rn = np.random.uniform(sigma_int, sigma_int + 1 + np.finfo(float).eps) + if rn < sigma: + Ne_sec = ( + sigma_int + 1 + ) # number of the secondary electrons to be emitted + for j in [0, Ne_sec - 1]: + xe = np.array([]) + ye = np.array([]) + ze = np.array([]) + we = np.array([]) + delta_te = np.array([]) + uxe = np.array([]) + uye = np.array([]) + uze = np.array([]) + + # Random thermal momenta distribution + ux_th = np.random.normal(0, dist_th) + uy_th = np.random.normal(0, dist_th) + uz_th = np.random.normal(0, dist_th) + + un_th = nx[i] * ux_th + ny[i] * uy_th + nz[i] * uz_th + + if un_th < 0: + ux_th_reflect = ( + -2 * un_th * nx[i] + ux_th + ) # for a "mirror reflection" u(sym)=-2(u.n)n+u + uy_th_reflect = -2 * un_th * ny[i] + uy_th + uz_th_reflect = -2 * un_th * nz[i] + uz_th + + uxe = np.append(uxe, ux_th_reflect) + uye = np.append(uye, uy_th_reflect) + uze = np.append(uze, uz_th_reflect) + else: + uxe = np.append(uxe, ux_th) + uye = np.append(uye, uy_th) + uze = np.append(uze, uz_th) + + xe = np.append(xe, x[i]) + ye = np.append(ye, y[i]) + ze = np.append(ze, z[i]) + we = np.append(we, w[i]) + delta_te = np.append(delta_te, delta_t[i]) + + elect_pc.add_particles( + x=xe + (dt - delta_te) * uxe, + y=ye + (dt - delta_te) * uye, + z=ze + (dt - delta_te) * uze, + ux=uxe, + uy=uye, + uz=uze, + w=we, + ) + buffer.clear_buffer() # reinitialise the boundary buffer + + +# using the new particle container modified at the last step +callbacks.installafterstep(secondary_emission) +########################## +# simulation run +########################## +sim.step(max_steps) # the whole process is done "max_steps" times diff --git a/Regression/Checksum/benchmarks_json/test_rz_secondary_ion_emission_picmi.json b/Regression/Checksum/benchmarks_json/test_rz_secondary_ion_emission_picmi.json new file mode 100644 index 00000000000..92e3b6ec0cc --- /dev/null +++ b/Regression/Checksum/benchmarks_json/test_rz_secondary_ion_emission_picmi.json @@ -0,0 +1,26 @@ +{ + "electrons": { + "particle_momentum_x": 6.863524723051401e-26, + "particle_momentum_y": 1.0324224192517369e-25, + "particle_momentum_z": 5.621885683115495e-26, + "particle_position_x": 0.0072217326490580675, + "particle_position_y": 0.001255871169011761, + "particle_position_z": 0.39892796754417836, + "particle_weight": 2.0 + }, + "ions": { + "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": { + "Er": 3.996898574068735e-08, + "Ez": 3.77948124311196e-08, + "phi": 2.359076749421693e-08, + "rho": 4.6171309216540875e-15 + } +} \ No newline at end of file From 02de61d6fdc0e0011919fe0fed81a34c8f848998 Mon Sep 17 00:00:00 2001 From: Olga Shapoval Date: Tue, 21 Jan 2025 09:19:06 -0800 Subject: [PATCH 02/13] Clean-up --- .../inputs_test_rz_secondary_ion_emission_picmi.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py index 5849ba4c59e..ed299483b0e 100644 --- a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py +++ b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py @@ -186,7 +186,7 @@ def secondary_emission(): buffer.get_particle_boundary_buffer("ions", "eb", "deltaTimeScraped", lev) ) energy_ions = 0.5 * proton_mass * w * (ux**2 + uy**2 + uz**2) - energy_ions_in_kEv = energy_ions / (1.602176634e-19 * 1000) + energy_ions_in_kEv = energy_ions / (e * 1000) sigma_nascap_ions = sigma_nascap(energy_ions_in_kEv, delta_H, E_HMax) # Loop over all ions in the EB buffer for i in range(0, n): @@ -197,7 +197,7 @@ def secondary_emission(): Ne_sec = ( sigma_int + 1 ) # number of the secondary electrons to be emitted - for j in [0, Ne_sec - 1]: + for i_elec in [0, Ne_sec - 1]: xe = np.array([]) ye = np.array([]) ze = np.array([]) From 25412977f129724ce353ec30158d43a18434cd2e Mon Sep 17 00:00:00 2001 From: Olga Shapoval Date: Tue, 21 Jan 2025 10:07:42 -0800 Subject: [PATCH 03/13] More clean-up --- .../inputs_test_rz_secondary_ion_emission_picmi.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py index ed299483b0e..a086441b92c 100644 --- a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py +++ b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py @@ -197,7 +197,7 @@ def secondary_emission(): Ne_sec = ( sigma_int + 1 ) # number of the secondary electrons to be emitted - for i_elec in [0, Ne_sec - 1]: + for _ in range(Ne_sec): xe = np.array([]) ye = np.array([]) ze = np.array([]) From a8a9b68c6b280faab5558e666798943a877c41c5 Mon Sep 17 00:00:00 2001 From: Olga Shapoval Date: Tue, 21 Jan 2025 16:54:48 -0800 Subject: [PATCH 04/13] Further clean-up --- .../Tests/secondary_ion_emission/analysis.py | 22 +++++++++---------- ...ts_test_rz_secondary_ion_emission_picmi.py | 4 +--- 2 files changed, 12 insertions(+), 14 deletions(-) diff --git a/Examples/Tests/secondary_ion_emission/analysis.py b/Examples/Tests/secondary_ion_emission/analysis.py index 7f43b8031b8..83a4a3bbd2f 100644 --- a/Examples/Tests/secondary_ion_emission/analysis.py +++ b/Examples/Tests/secondary_ion_emission/analysis.py @@ -43,20 +43,20 @@ print("\n") print(f"Electron # {i}:") print("NUMERICAL coordinates of the emitted electrons:") - print("x=%5.5f, y=%5.5f, z=%5.5f" % (x[i], y[i], z[i])) + print(f"x={x[i]:5.5f}, y={y[i]:5.5f}, z={z[i]:5.5f}") print("\n") print("ANALYTICAL coordinates of the point of contact:") - print("x=%5.5f, y=%5.5f, z=%5.5f" % (x_analytic[i], y_analytic[i], z_analytic[i])) + print(f"x={x_analytic[i]:5.5f}, y={y_analytic[i]:5.5f}, z={z_analytic[i]:5.5f}") - diff_x = np.abs((x[i] - x_analytic[i]) / x_analytic[i]) - diff_y = np.abs((y[i] - y_analytic[i]) / y_analytic[i]) - diff_z = np.abs((z[i] - z_analytic[i]) / z_analytic[i]) + rel_err_x = np.abs((x[i] - x_analytic[i]) / x_analytic[i]) + rel_err_y = np.abs((y[i] - y_analytic[i]) / y_analytic[i]) + rel_err_z = np.abs((z[i] - z_analytic[i]) / z_analytic[i]) print("\n") - print("percentage error for x = %5.4f %%" % (diff_x * 100)) - print("percentage error for y = %5.4f %%" % (diff_y * 100)) - print("percentage error for z = %5.4f %%" % (diff_z * 100)) + print(f"Relative percentage error for x = {rel_err_x * 100:5.4f} %") + print(f"Relative percentage error for y = {rel_err_y * 100:5.4f} %") + print(f"Relative percentage error for z = {rel_err_z * 100:5.4f} %") - assert (diff_x < tolerance) and (diff_y < tolerance) and (diff_z < tolerance), ( - "Test particle_boundary_interaction did not pass" - ) + assert ( + (rel_err_x < tolerance) and (rel_err_y < tolerance) and (rel_err_x < tolerance) + ), "Test particle_boundary_interaction did not pass" diff --git a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py index a086441b92c..7c320a40069 100644 --- a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py +++ b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py @@ -95,8 +95,7 @@ ########################## # simulation setup ########################## -# num_procs = [1, 1] -# warpx_numprocs=num_procs, + sim = picmi.Simulation( solver=solver, time_step_size=dt, @@ -147,7 +146,6 @@ def sigma_nascap(energy_kEv, delta_H, E_HMax): - numpy array, computed probability sigma_nascap """ sigma_nascap = np.array([]) - # Loop through each energy value for energy in energy_kEv: if energy > 0.0: sigma = ( From 0a21435ece68562c1dd15494929acdb08f8259e5 Mon Sep 17 00:00:00 2001 From: Olga Shapoval Date: Wed, 22 Jan 2025 10:27:30 -0800 Subject: [PATCH 05/13] Skipped pre-commit changes. --- .../inputs_test_rz_secondary_ion_emission_picmi.py | 13 +------------ 1 file changed, 1 insertion(+), 12 deletions(-) diff --git a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py index 7c320a40069..ac3df540a53 100644 --- a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py +++ b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py @@ -145,18 +145,7 @@ def sigma_nascap(energy_kEv, delta_H, E_HMax): Returns: - numpy array, computed probability sigma_nascap """ - sigma_nascap = np.array([]) - for energy in energy_kEv: - if energy > 0.0: - sigma = ( - delta_H - * (E_HMax + 1.0) - / (E_HMax * 1.0 + energy) - * np.sqrt(energy / 1.0) - ) - else: - sigma = 0.0 - sigma_nascap = np.append(sigma_nascap, sigma) + sigma_nascap = np.array([(delta_H * (E_HMax + 1.0) / (E_HMax + energy) * np.sqrt(energy) if energy > 0.0 else 0.0) for energy in energy_kEv]) return sigma_nascap From b3730bf0ef5214afff6c902e812305c3c0ee7688 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 22 Jan 2025 18:29:27 +0000 Subject: [PATCH 06/13] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- .../inputs_test_rz_secondary_ion_emission_picmi.py | 11 ++++++++++- 1 file changed, 10 insertions(+), 1 deletion(-) diff --git a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py index ac3df540a53..c64f72ef50c 100644 --- a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py +++ b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py @@ -145,7 +145,16 @@ def sigma_nascap(energy_kEv, delta_H, E_HMax): Returns: - numpy array, computed probability sigma_nascap """ - sigma_nascap = np.array([(delta_H * (E_HMax + 1.0) / (E_HMax + energy) * np.sqrt(energy) if energy > 0.0 else 0.0) for energy in energy_kEv]) + sigma_nascap = np.array( + [ + ( + delta_H * (E_HMax + 1.0) / (E_HMax + energy) * np.sqrt(energy) + if energy > 0.0 + else 0.0 + ) + for energy in energy_kEv + ] + ) return sigma_nascap From 8801226792541d7e8617c2a036452cc5550025fc Mon Sep 17 00:00:00 2001 From: Olga Shapoval Date: Wed, 22 Jan 2025 10:36:20 -0800 Subject: [PATCH 07/13] Skipped pre-commit changes --- .../inputs_test_rz_secondary_ion_emission_picmi.py | 11 +---------- 1 file changed, 1 insertion(+), 10 deletions(-) diff --git a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py index c64f72ef50c..ac3df540a53 100644 --- a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py +++ b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py @@ -145,16 +145,7 @@ def sigma_nascap(energy_kEv, delta_H, E_HMax): Returns: - numpy array, computed probability sigma_nascap """ - sigma_nascap = np.array( - [ - ( - delta_H * (E_HMax + 1.0) / (E_HMax + energy) * np.sqrt(energy) - if energy > 0.0 - else 0.0 - ) - for energy in energy_kEv - ] - ) + sigma_nascap = np.array([(delta_H * (E_HMax + 1.0) / (E_HMax + energy) * np.sqrt(energy) if energy > 0.0 else 0.0) for energy in energy_kEv]) return sigma_nascap From 0d5330ec189d840b61dcab00d1e6da3488896cef Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 22 Jan 2025 18:37:01 +0000 Subject: [PATCH 08/13] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- .../inputs_test_rz_secondary_ion_emission_picmi.py | 11 ++++++++++- 1 file changed, 10 insertions(+), 1 deletion(-) diff --git a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py index ac3df540a53..c64f72ef50c 100644 --- a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py +++ b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py @@ -145,7 +145,16 @@ def sigma_nascap(energy_kEv, delta_H, E_HMax): Returns: - numpy array, computed probability sigma_nascap """ - sigma_nascap = np.array([(delta_H * (E_HMax + 1.0) / (E_HMax + energy) * np.sqrt(energy) if energy > 0.0 else 0.0) for energy in energy_kEv]) + sigma_nascap = np.array( + [ + ( + delta_H * (E_HMax + 1.0) / (E_HMax + energy) * np.sqrt(energy) + if energy > 0.0 + else 0.0 + ) + for energy in energy_kEv + ] + ) return sigma_nascap From 79393c1c4be3a1c69a3608cac2f8381e5ab8bf17 Mon Sep 17 00:00:00 2001 From: Olga Shapoval Date: Wed, 22 Jan 2025 13:17:02 -0800 Subject: [PATCH 09/13] Reverted changes in sigma_nascap function --- .../inputs_test_rz_secondary_ion_emission_picmi.py | 14 +++++++++++++- 1 file changed, 13 insertions(+), 1 deletion(-) diff --git a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py index ac3df540a53..1671f74fc22 100644 --- a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py +++ b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py @@ -145,7 +145,19 @@ def sigma_nascap(energy_kEv, delta_H, E_HMax): Returns: - numpy array, computed probability sigma_nascap """ - sigma_nascap = np.array([(delta_H * (E_HMax + 1.0) / (E_HMax + energy) * np.sqrt(energy) if energy > 0.0 else 0.0) for energy in energy_kEv]) + sigma_nascap = np.array([]) + # Loop through each energy value + for energy in energy_kEv: + if energy > 0.0: + sigma = ( + delta_H + * (E_HMax + 1.0) + / (E_HMax * 1.0 + energy) + * np.sqrt(energy / 1.0) + ) + else: + sigma = 0.0 + sigma_nascap = np.append(sigma_nascap, sigma) return sigma_nascap From 14e9ca625162424d7b02cc685840e576111cdbaa Mon Sep 17 00:00:00 2001 From: Olga Shapoval Date: Thu, 23 Jan 2025 15:00:33 -0800 Subject: [PATCH 10/13] Fixed typo introduced during the clean-up --- Examples/Tests/secondary_ion_emission/analysis.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Examples/Tests/secondary_ion_emission/analysis.py b/Examples/Tests/secondary_ion_emission/analysis.py index 83a4a3bbd2f..16dda3d70e5 100644 --- a/Examples/Tests/secondary_ion_emission/analysis.py +++ b/Examples/Tests/secondary_ion_emission/analysis.py @@ -58,5 +58,5 @@ print(f"Relative percentage error for z = {rel_err_z * 100:5.4f} %") assert ( - (rel_err_x < tolerance) and (rel_err_y < tolerance) and (rel_err_x < tolerance) + (rel_err_x < tolerance) and (rel_err_y < tolerance) and (rel_err_z < tolerance) ), "Test particle_boundary_interaction did not pass" From f74bb9e30502c5f5a97288dabd851eb1b103a6c5 Mon Sep 17 00:00:00 2001 From: Olga Shapoval Date: Tue, 28 Jan 2025 20:32:07 -0800 Subject: [PATCH 11/13] Clean-up --- .../Tests/secondary_ion_emission/analysis.py | 14 +- ...ts_test_rz_secondary_ion_emission_picmi.py | 135 ++++++++++-------- .../test_rz_secondary_ion_emission_picmi.json | 20 +-- 3 files changed, 96 insertions(+), 73 deletions(-) diff --git a/Examples/Tests/secondary_ion_emission/analysis.py b/Examples/Tests/secondary_ion_emission/analysis.py index 16dda3d70e5..2479a4f18a3 100644 --- a/Examples/Tests/secondary_ion_emission/analysis.py +++ b/Examples/Tests/secondary_ion_emission/analysis.py @@ -27,9 +27,17 @@ print("y", y) print("z", z) # Analytical results calculated -x_analytic = [0.004028, 0.003193] -y_analytic = [-0.0001518, -0.0011041] -z_analytic = [-0.19967, -0.19926] +# x_analytic = [0.004028, 0.003193] +# y_analytic = [-0.0001518, -0.0011041] +# z_analytic = [-0.19967, -0.19926] + +x_analytic = [-0.091696, 0.011599] +y_analytic = [-0.002282, -0.0111624] +z_analytic = [-0.200242, -0.201728] + +# x_analytic = [-0.09169647, 0.01159922] +# y_analytic = [-0.00228188, -0.01116237] +# z_analytic = [-0.20024175, -0.20172786] N_sec_e = np.size(z) # number of the secondary electrons diff --git a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py index 1671f74fc22..5b6248da33c 100644 --- a/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py +++ b/Examples/Tests/secondary_ion_emission/inputs_test_rz_secondary_ion_emission_picmi.py @@ -1,5 +1,13 @@ #!/usr/bin/env python3 -# --- Input file for secondary-ion emission testing in RZ via callback function. +# This is the script that tests secondary ion emission when ions hit an embedded boundary +# with a specified secondary emission yield of delta_H = 0.4. Specifically, a callback +# function at each time step ensures that the correct number of secondary electrons is +# emitted when ions impact the embedded boundary, following the given secondary emission +# model defined in sigma_nescap function. This distribution depends on the ion's energy and +# suggests that for an ion incident with 1 keV energy, an average of 0.4 secondary +# electrons will be emitted. +# Simulation is initialized with four ions with i_dist distribution and spherical +# embedded boundary given by implicit function. import numpy as np from scipy.constants import e, elementary_charge, m_e, proton_mass @@ -9,7 +17,7 @@ # numerics parameters ########################## -dt = 0.000000015 +dt = 0.000000075 # --- Nb time steps Te = 0.0259 # in eV @@ -56,10 +64,21 @@ ########################## # physics components ########################## - i_dist = picmi.ParticleListDistribution( - x=0.0, y=0.0, z=-0.25, ux=0.1e6, uy=0.0, uz=1.50e6, weight=1 + x=[ + 0.025, + 0.0, + -0.1, + -0.14, + ], + y=[0.0, 0.0, 0.0, 0], + z=[-0.26, -0.29, -0.25, -0.23], + ux=[0.18e6, 0.1e6, 0.15e6, 0.21e6], + uy=[0.0, 0.0, 0.0, 0.0], + uz=[8.00e5, 7.20e5, 6.40e5, 5.60e5], + weight=[1, 1, 1, 1], ) + electrons = picmi.Species( particle_type="electron", # Specify the particle type name="electrons", # Name of the species @@ -67,7 +86,7 @@ ions = picmi.Species( name="ions", - mass=proton_mass, + particle_type="proton", charge=e, initial_distribution=i_dist, warpx_save_particles_at_eb=1, @@ -81,7 +100,7 @@ name="diag1", grid=grid, period=diagnostic_interval, - data_list=["Er", "Ez", "phi", "rho"], # , "rho_electrons"], + data_list=["Er", "Ez", "phi", "rho"], warpx_format="openpmd", ) @@ -190,60 +209,56 @@ def secondary_emission(): # Loop over all ions in the EB buffer for i in range(0, n): sigma = sigma_nascap_ions[i] - sigma_int = int(sigma) - rn = np.random.uniform(sigma_int, sigma_int + 1 + np.finfo(float).eps) - if rn < sigma: - Ne_sec = ( - sigma_int + 1 - ) # number of the secondary electrons to be emitted - for _ in range(Ne_sec): - xe = np.array([]) - ye = np.array([]) - ze = np.array([]) - we = np.array([]) - delta_te = np.array([]) - uxe = np.array([]) - uye = np.array([]) - uze = np.array([]) - - # Random thermal momenta distribution - ux_th = np.random.normal(0, dist_th) - uy_th = np.random.normal(0, dist_th) - uz_th = np.random.normal(0, dist_th) - - un_th = nx[i] * ux_th + ny[i] * uy_th + nz[i] * uz_th - - if un_th < 0: - ux_th_reflect = ( - -2 * un_th * nx[i] + ux_th - ) # for a "mirror reflection" u(sym)=-2(u.n)n+u - uy_th_reflect = -2 * un_th * ny[i] + uy_th - uz_th_reflect = -2 * un_th * nz[i] + uz_th - - uxe = np.append(uxe, ux_th_reflect) - uye = np.append(uye, uy_th_reflect) - uze = np.append(uze, uz_th_reflect) - else: - uxe = np.append(uxe, ux_th) - uye = np.append(uye, uy_th) - uze = np.append(uze, uz_th) - - xe = np.append(xe, x[i]) - ye = np.append(ye, y[i]) - ze = np.append(ze, z[i]) - we = np.append(we, w[i]) - delta_te = np.append(delta_te, delta_t[i]) - - elect_pc.add_particles( - x=xe + (dt - delta_te) * uxe, - y=ye + (dt - delta_te) * uye, - z=ze + (dt - delta_te) * uze, - ux=uxe, - uy=uye, - uz=uze, - w=we, - ) - buffer.clear_buffer() # reinitialise the boundary buffer + # Ne_sec is number of the secondary electrons to be emitted + Ne_sec = int(sigma + np.random.uniform()) + for _ in range(Ne_sec): + xe = np.array([]) + ye = np.array([]) + ze = np.array([]) + we = np.array([]) + delta_te = np.array([]) + uxe = np.array([]) + uye = np.array([]) + uze = np.array([]) + + # Random thermal momenta distribution + ux_th = np.random.normal(0, dist_th) + uy_th = np.random.normal(0, dist_th) + uz_th = np.random.normal(0, dist_th) + + un_th = nx[i] * ux_th + ny[i] * uy_th + nz[i] * uz_th + + if un_th < 0: + ux_th_reflect = ( + -2 * un_th * nx[i] + ux_th + ) # for a "mirror reflection" u(sym)=-2(u.n)n+u + uy_th_reflect = -2 * un_th * ny[i] + uy_th + uz_th_reflect = -2 * un_th * nz[i] + uz_th + + uxe = np.append(uxe, ux_th_reflect) + uye = np.append(uye, uy_th_reflect) + uze = np.append(uze, uz_th_reflect) + else: + uxe = np.append(uxe, ux_th) + uye = np.append(uye, uy_th) + uze = np.append(uze, uz_th) + + xe = np.append(xe, x[i]) + ye = np.append(ye, y[i]) + ze = np.append(ze, z[i]) + we = np.append(we, w[i]) + delta_te = np.append(delta_te, delta_t[i]) + + elect_pc.add_particles( + x=xe + (dt - delta_te) * uxe, + y=ye + (dt - delta_te) * uye, + z=ze + (dt - delta_te) * uze, + ux=uxe, + uy=uye, + uz=uze, + w=we, + ) + buffer.clear_buffer() # reinitialise the boundary buffer # using the new particle container modified at the last step diff --git a/Regression/Checksum/benchmarks_json/test_rz_secondary_ion_emission_picmi.json b/Regression/Checksum/benchmarks_json/test_rz_secondary_ion_emission_picmi.json index 92e3b6ec0cc..cfc84819e97 100644 --- a/Regression/Checksum/benchmarks_json/test_rz_secondary_ion_emission_picmi.json +++ b/Regression/Checksum/benchmarks_json/test_rz_secondary_ion_emission_picmi.json @@ -1,11 +1,11 @@ { "electrons": { - "particle_momentum_x": 6.863524723051401e-26, - "particle_momentum_y": 1.0324224192517369e-25, - "particle_momentum_z": 5.621885683115495e-26, - "particle_position_x": 0.0072217326490580675, - "particle_position_y": 0.001255871169011761, - "particle_position_z": 0.39892796754417836, + "particle_momentum_x": 5.621885683102775e-26, + "particle_momentum_y": 1.2079178196118306e-25, + "particle_momentum_z": 1.2496342823828099e-25, + "particle_position_x": 0.10329568998704057, + "particle_position_y": 0.013444257249267193, + "particle_position_z": 0.4019696082583948, "particle_weight": 2.0 }, "ions": { @@ -18,9 +18,9 @@ "particle_weight": 0.0 }, "lev=0": { - "Er": 3.996898574068735e-08, - "Ez": 3.77948124311196e-08, - "phi": 2.359076749421693e-08, - "rho": 4.6171309216540875e-15 + "Er": 1.772547702166409e-06, + "Ez": 2.2824957684716966e-06, + "phi": 4.338168233265556e-07, + "rho": 1.933391680367631e-15 } } \ No newline at end of file From 196ad45775d3b427bbf56320c8da42918b85935b Mon Sep 17 00:00:00 2001 From: Remi Lehe Date: Thu, 30 Jan 2025 10:30:48 -0800 Subject: [PATCH 12/13] Apply suggestions from code review --- .../Tests/secondary_ion_emission/analysis.py | 24 +++++-------------- 1 file changed, 6 insertions(+), 18 deletions(-) diff --git a/Examples/Tests/secondary_ion_emission/analysis.py b/Examples/Tests/secondary_ion_emission/analysis.py index 2479a4f18a3..384a1ce41ec 100644 --- a/Examples/Tests/secondary_ion_emission/analysis.py +++ b/Examples/Tests/secondary_ion_emission/analysis.py @@ -1,12 +1,11 @@ #!/usr/bin/env python """ -This script tests the last coordinates of the emitted secondary electrons. -The EB sphere is centered on O and has a radius of 0.2. -The proton is initially at: (0,0,-0.25) and moves with a velocity: -(0.1e6,0,1.5e6) with a time step of dt = 0.000000015. -The simulation uses a fixed random seed (np.random.seed(10025015)) -to ensure the emission of secodnary electrons. -An input file inputs_test_rz_secondary_ion_emission_picmi.py is used. +This script checks that electron secondary emission (implemented by a callback function) works as intended. + +In this test, four ions hit a spherical embedded boundary, and produce secondary +electrons with a probability of `0.4`. We thus expect ~2 electrons to be produced. +This script tests the number of electrons emitted and checks that their position is +close to the embedded boundary. """ import sys @@ -23,22 +22,11 @@ it = ts.iterations x, y, z = ts.get_particle(["x", "y", "z"], species="electrons", iteration=it[-1]) -print("x", x) -print("y", y) -print("z", z) -# Analytical results calculated -# x_analytic = [0.004028, 0.003193] -# y_analytic = [-0.0001518, -0.0011041] -# z_analytic = [-0.19967, -0.19926] x_analytic = [-0.091696, 0.011599] y_analytic = [-0.002282, -0.0111624] z_analytic = [-0.200242, -0.201728] -# x_analytic = [-0.09169647, 0.01159922] -# y_analytic = [-0.00228188, -0.01116237] -# z_analytic = [-0.20024175, -0.20172786] - N_sec_e = np.size(z) # number of the secondary electrons assert N_sec_e == 2, ( From 4b142248a697173feded0f1a8f9f1b77c5652947 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 30 Jan 2025 18:30:57 +0000 Subject: [PATCH 13/13] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- Examples/Tests/secondary_ion_emission/analysis.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/Examples/Tests/secondary_ion_emission/analysis.py b/Examples/Tests/secondary_ion_emission/analysis.py index 384a1ce41ec..8c2ed5b4af6 100644 --- a/Examples/Tests/secondary_ion_emission/analysis.py +++ b/Examples/Tests/secondary_ion_emission/analysis.py @@ -2,9 +2,9 @@ """ This script checks that electron secondary emission (implemented by a callback function) works as intended. -In this test, four ions hit a spherical embedded boundary, and produce secondary -electrons with a probability of `0.4`. We thus expect ~2 electrons to be produced. -This script tests the number of electrons emitted and checks that their position is +In this test, four ions hit a spherical embedded boundary, and produce secondary +electrons with a probability of `0.4`. We thus expect ~2 electrons to be produced. +This script tests the number of electrons emitted and checks that their position is close to the embedded boundary. """