From 3982b0005e82c030b2476dc5efd0ff65a18f6fa1 Mon Sep 17 00:00:00 2001 From: Eya D <81635404+EyaDammak@users.noreply.github.com> Date: Wed, 6 Mar 2024 16:25:25 -0800 Subject: [PATCH] New Python test: Particle-Boundary interaction (#4729) * enable the diagnostic of ParticleScraping in Python * Update picmi.py * Update picmi.py * new test * python update * modification of the script * [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci * Update PICMI_inputs_rz.py * json * update * Update PICMI_inputs_rz.py * Update particle_boundary_interaction.json * Update PICMI_inputs_rz.py * Update PICMI_inputs_rz.py * [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci * Update PICMI_inputs_rz.py * [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci * Fix hanging script in parallel * Make the test executable * Update analysis script * Update particle_containers.py * Update PICMI_inputs_rz.py * Update analysis.py * Update analysis.py * Update particle_containers.py * [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --------- Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com> Co-authored-by: Remi Lehe Co-authored-by: Axel Huebl --- .../PICMI_inputs_rz.py | 169 ++++++++++++++++++ .../particle_boundary_interaction/analysis.py | 51 ++++++ Python/pywarpx/particle_containers.py | 3 + .../particle_boundary_interaction.json | 18 ++ Regression/WarpX-tests.ini | 21 +++ 5 files changed, 262 insertions(+) create mode 100644 Examples/Tests/particle_boundary_interaction/PICMI_inputs_rz.py create mode 100755 Examples/Tests/particle_boundary_interaction/analysis.py create mode 100644 Regression/Checksum/benchmarks_json/particle_boundary_interaction.json diff --git a/Examples/Tests/particle_boundary_interaction/PICMI_inputs_rz.py b/Examples/Tests/particle_boundary_interaction/PICMI_inputs_rz.py new file mode 100644 index 00000000000..4017a05d413 --- /dev/null +++ b/Examples/Tests/particle_boundary_interaction/PICMI_inputs_rz.py @@ -0,0 +1,169 @@ +#!/usr/bin/env python3 +# @Eya Dammak supervised by @Remi Lehe, 2024 +# --- Input file for particle-boundary interaction testing in RZ. +# --- This input is a simple case of reflection +# --- of one electron on the surface of a sphere. +import numpy as np + +from pywarpx import callbacks, particle_containers, picmi + +########################## +# numerics parameters +########################## + +dt = 1.0e-11 + +# --- Nb time steps + +max_steps = 23 +diagnostic_interval = 1 + +# --- grid + +nr = 64 +nz= 64 + +rmin = 0.0 +rmax = 2 +zmin = -2 +zmax = 2 + +########################## +# 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 = ['absorbing', '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 +########################## + +#one particle +e_dist=picmi.ParticleListDistribution(x=0.0, y=0.0, z=-0.25, ux=0.5e10, uy=0.0, uz=1.0e10, weight=1) + +electrons = picmi.Species( + particle_type='electron', name='electrons', initial_distribution=e_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', + write_dir = '.', + warpx_file_prefix = 'particle_boundary_interaction_plt' +) + +part_diag = picmi.ParticleDiagnostic(name = 'diag1', + period = diagnostic_interval, + species = [electrons], + warpx_format = 'openpmd', + write_dir = '.', + warpx_file_prefix = 'particle_boundary_interaction_plt' +) + +########################## +# simulation setup +########################## + +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=[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 mirror_reflection(): + 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) + delta_t = concat(buffer.get_particle_boundary_buffer("electrons", 'eb', 'deltaTimeScraped', lev)) + r = concat(buffer.get_particle_boundary_buffer("electrons", 'eb', 'x', lev)) + theta = concat(buffer.get_particle_boundary_buffer("electrons", 'eb', 'theta', lev)) + z = concat(buffer.get_particle_boundary_buffer("electrons", '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("electrons", 'eb', 'ux', lev)) + uy = concat(buffer.get_particle_boundary_buffer("electrons", 'eb', 'uy', lev)) + uz = concat(buffer.get_particle_boundary_buffer("electrons", 'eb', 'uz', lev)) + w = concat(buffer.get_particle_boundary_buffer("electrons", 'eb', 'w', lev)) + nx = concat(buffer.get_particle_boundary_buffer("electrons", 'eb', 'nx', lev)) + ny = concat(buffer.get_particle_boundary_buffer("electrons", 'eb', 'ny', lev)) + nz = concat(buffer.get_particle_boundary_buffer("electrons", 'eb', 'nz', lev)) + + #STEP 2: use these parameters to inject particle from the same position in the plasma + elect_pc = particle_containers.ParticleContainerWrapper('electrons') #general particle container + + ####this part is specific to the case of simple reflection. + un=ux*nx+uy*ny+uz*nz + ux_reflect=-2*un*nx+ux #for a "mirror reflection" u(sym)=-2(u.n)n+u + uy_reflect=-2*un*ny+uy + uz_reflect=-2*un*nz+uz + elect_pc.add_particles( + x=x + (dt-delta_t)*ux_reflect, y=y + (dt-delta_t)*uy_reflect, z=z + (dt-delta_t)*uz_reflect, + ux=ux_reflect, uy=uy_reflect, uz=uz_reflect, + w=w + ) #adds the particle in the general particle container at the next step + #### Can be modified depending on the model of interaction. + + buffer.clear_buffer() #reinitialise the boundary buffer + +callbacks.installafterstep(mirror_reflection) #mirror_reflection is called at the next step + # using the new particle container modified at the last step + +########################## +# simulation run +########################## + +sim.step(max_steps) #the whole process is done "max_steps" times diff --git a/Examples/Tests/particle_boundary_interaction/analysis.py b/Examples/Tests/particle_boundary_interaction/analysis.py new file mode 100755 index 00000000000..9768751ed74 --- /dev/null +++ b/Examples/Tests/particle_boundary_interaction/analysis.py @@ -0,0 +1,51 @@ +#!/usr/bin/env python +# @Eya Dammak supervised by @Remi Lehe, 2024 +""" +This script tests the last coordinate after adding an electron. +The sphere is centered on O and has a radius of 0.2 (EB) +The electron is initially at: (0,0,-0.25) and moves with a velocity: +(0.5e10,0,1.0e10) with a time step of 1e-11. +An input file PICMI_inputs_rz.py is used. +""" +import os +import sys + +import numpy as np +from openpmd_viewer import OpenPMDTimeSeries +import yt + +yt.funcs.mylog.setLevel(0) +sys.path.insert(1, '../../../../warpx/Regression/Checksum/') +import checksumAPI + +# Open plotfile specified in command line +filename = sys.argv[1] +test_name = os.path.split(os.getcwd())[1] +checksumAPI.evaluate_checksum(test_name, filename, output_format='openpmd') + +ts = OpenPMDTimeSeries('./particle_boundary_interaction_plt') + +it=ts.iterations +x,y,z=ts.get_particle(['x','y','z'], species='electrons', iteration=it[-1]) + +# Analytical results calculated +x_analytic=0.03532 +y_analytic=0.00000 +z_analytic=-0.20531 + +print('NUMERICAL coordinates of the point of contact:') +print('x=%5.5f, y=%5.5f, z=%5.5f' % (x[0], y[0], z[0])) +print('\n') +print('ANALYTICAL coordinates of the point of contact:') +print('x=%5.5f, y=%5.5f, z=%5.5f' % (x_analytic, y_analytic, z_analytic)) + +tolerance=1e-5 + +diff_x=np.abs((x[0]-x_analytic)/x_analytic) +diff_z=np.abs((z[0]-z_analytic)/z_analytic) + +print('\n') +print("percentage error for x = %5.4f %%" %(diff_x *100)) +print("percentage error for z = %5.4f %%" %(diff_z *100)) + +assert (diff_x < tolerance) and (y[0] < 1e-8) and (diff_z < tolerance), 'Test particle_boundary_interaction did not pass' diff --git a/Python/pywarpx/particle_containers.py b/Python/pywarpx/particle_containers.py index 19b0a1ab6c4..f8e687f99bd 100644 --- a/Python/pywarpx/particle_containers.py +++ b/Python/pywarpx/particle_containers.py @@ -761,6 +761,9 @@ def get_particle_boundary_buffer(self, species_name, boundary, comp_name, level) The data for the arrays are not copied, but share the underlying memory buffer with WarpX. The arrays are fully writeable. + You can find `here https://github.com/ECP-WarpX/WarpX/blob/319e55b10ad4f7c71b84a4fb21afbafe1f5b65c2/Examples/Tests/particle_boundary_interaction/PICMI_inputs_rz.py` + an example of a simple case of particle-boundary interaction (reflection). + Parameters ---------- diff --git a/Regression/Checksum/benchmarks_json/particle_boundary_interaction.json b/Regression/Checksum/benchmarks_json/particle_boundary_interaction.json new file mode 100644 index 00000000000..e35516758f2 --- /dev/null +++ b/Regression/Checksum/benchmarks_json/particle_boundary_interaction.json @@ -0,0 +1,18 @@ +{ + "electrons": { + "particle_momentum_x": 7.082238783412152e-21, + "particle_momentum_y": 0.0, + "particle_momentum_z": 7.319015172219235e-21, + "particle_position_x": 0.03532003602713794, + "particle_position_y": 0.0, + "particle_position_z": 0.20531131195378569, + "particle_weight": 1.0 + }, + "lev=0": { + "Er": 1.3287679891086577e-06, + "Ez": 1.791737740110229e-06, + "phi": 3.3983929077087634e-07, + "rho": 7.811529217958718e-16, + "rho_electrons": 7.811529217958718e-16 + } +} diff --git a/Regression/WarpX-tests.ini b/Regression/WarpX-tests.ini index 426ecb0717c..ab422d5631a 100644 --- a/Regression/WarpX-tests.ini +++ b/Regression/WarpX-tests.ini @@ -4681,3 +4681,24 @@ doVis = 0 compareParticles = 0 doComparison = 0 analysisRoutine = Examples/Tests/gaussian_beam/analysis_focusing_beam.py + +[particle_boundary_interaction] +buildDir = . +inputFile = Examples/Tests/particle_boundary_interaction/PICMI_inputs_rz.py +runtime_params = +customRunCmd = python3 PICMI_inputs_rz.py +dim = 2 +addToCompileString = USE_PYTHON_MAIN=TRUE USE_RZ=TRUE +cmakeSetupOpts = -DWarpX_DIMS="RZ" -DWarpX_EB=ON -DWarpX_PYTHON=ON +target = pip_install +restartTest = 0 +useMPI = 1 +numprocs = 2 +useOMP = 1 +numthreads = 1 +compileTest = 0 +doVis = 0 +compareParticles = 1 +particleTypes = electrons +outputFile = particle_boundary_interaction_plt +analysisRoutine = Examples/Tests/particle_boundary_interaction/analysis.py