2 files changed, 131 insertions, 0 deletions
diff --git a/Examples/Tests/RepellingParticles/analysis_repelling.py b/Examples/Tests/RepellingParticles/analysis_repelling.py
new file mode 100755
index 000000000..25f1fe977
--- /dev/null
+++ b/Examples/Tests/RepellingParticles/analysis_repelling.py
@@ -0,0 +1,78 @@
+#! /usr/bin/env python
+#
+# Copyright 2021 Remi Lehe
+#
+# This file is part of WarpX.
+#
+# License: BSD-3-Clause-LBNL
+
+"""
+This script tests the interaction between 2 electron macroparticles.
+
+The macroparticles are moving towards each other and slow down due to the
+repelling force between charges of the same sign. The particles should
+be moving slow enough for relativistic effects to be negligible. In
+this case the conservation of energy gives:
+
+beta(t)^2 = beta(0)^2 + 2 * w * re * log(d(0)/d(t))
+
+where:
+re is the classical electron radius
+w is the weight of the individual macroparticles
+d is the distance between them
+beta is the velocity normalized by the speed of light
+"""
+import numpy as np
+from scipy.constants import m_e, c, physical_constants
+import sys, re
+import yt
+import glob
+import matplotlib.pyplot as plt
+yt.funcs.mylog.setLevel(0)
+
+# Check plotfile name specified in command line
+last_filename = sys.argv[1]
+filename_radical = re.findall('(.*?)\d+/*$', last_filename)[0]
+
+# Loop through files, and extract the position and velocity of both particles
+x1 = []
+x2 = []
+beta1 = []
+beta2 = []
+for filename in sorted(glob.glob(filename_radical + '*')):
+    print(filename)
+    ds = yt.load(filename)
+    ad = ds.all_data()
+
+    x1.append( float(ad[('electron1','particle_position_x')]) )
+    x2.append( float(ad[('electron2','particle_position_x')]) )
+    beta1.append( float(ad[('electron1','particle_momentum_x')])/(m_e*c) )
+    beta2.append( float(ad[('electron2','particle_momentum_x')])/(m_e*c) )
+
+# Convert to numpy array
+x1 = np.array(x1)
+x2 = np.array(x2)
+beta1 = np.array(beta1)
+beta2 = np.array(beta2)
+
+# Plot velocities, compare with theory
+w = 5.e12
+re = physical_constants['classical electron radius'][0]
+beta_th = np.sqrt( beta1[0]**2 - 2*w*re*np.log( (x2[0]-x1[0])/(x2-x1) ) )
+plt.plot( beta1, '+', label='Particle 1' )
+plt.plot( -beta2, 'x', label='Particle 2' )
+plt.plot( beta_th, '*', label='Theory' )
+plt.legend(loc=0)
+plt.xlabel('Time (a.u.)')
+plt.ylabel('Normalized velocity')
+plt.savefig('Comparison.png')
+
+# Check that the results are close to the theory
+assert np.allclose( beta1[1:], beta_th[1:], atol=0.01 )
+assert np.allclose( -beta2[1:], beta_th[1:], atol=0.01  )
+
+# Run checksum regression test
+sys.path.insert(1, '../../../../warpx/Regression/Checksum/')
+import checksumAPI
+test_name = last_filename[:-9]
+checksumAPI.evaluate_checksum(test_name, last_filename)
diff --git a/Examples/Tests/RepellingParticles/inputs_2d b/Examples/Tests/RepellingParticles/inputs_2d
new file mode 100644
index 000000000..5f4aafce5
--- /dev/null
+++ b/Examples/Tests/RepellingParticles/inputs_2d
@@ -0,0 +1,53 @@
+max_step = 200
+amr.n_cell = 128 128
+
+amr.blocking_factor = 16
+amr.max_grid_size = 1024
+amr.max_level = 0
+
+# Geometry
+geometry.coord_sys   = 0                  # 0: Cartesian
+geometry.prob_lo     = -32.e-6    -32.e-6      # physical domain
+geometry.prob_hi     =  32.e-6     32.e-6
+
+# Boundary condition
+boundary.field_lo = pec pec
+boundary.field_hi = pec pec
+
+# Algorithms
+algo.current_deposition = esirkepov
+algo.field_gathering = momentum-conserving
+algo.charge_deposition = standard
+algo.particle_pusher = vay
+algo.maxwell_solver = yee
+warpx.cfl = 0.9
+warpx.use_filter = 0
+warpx.do_dive_cleaning = 1
+
+# Particle species
+particles.species_names = electron1 electron2
+
+electron1.charge = q_e
+electron1.mass = m_e
+electron1.injection_style = "singleparticle"
+electron1.single_particle_pos = 1.e-6 0. 0.
+electron1.single_particle_vel = 0. 0. 0.
+electron1.single_particle_weight = 5.e12
+electron1.initialize_self_fields = 1
+
+electron2.charge = q_e
+electron2.mass = m_e
+electron2.injection_style = "singleparticle"
+electron2.single_particle_pos = -1.e-6 0. 0.
+electron2.single_particle_vel = 0. 0. 0.
+electron2.single_particle_weight = 5.e12
+electron2.initialize_self_fields = 1
+
+# Order of particle shape factors
+algo.particle_shape = 3
+
+# Diagnostics
+diagnostics.diags_names = diag1
+diag1.intervals = 20
+diag1.diag_type = Full
+diag1.fields_to_plot = Bx By Bz Ex Ey Ez jx jy jz divE rho F