WarpX
WarpX copied to clipboard
Periodic boundaries eat particles in single precision
We noticed that particle counts are being depleted in a single precision (fields and particles) simulation with only a reflecting and periodic boundaries. We have seen similar issues of particle depletion in our other simulation runs/plots, posted/referenced below.
We've replicated this issue using the cpp input file attached.
Originally posted by @roelof-groenewald in https://github.com/ECP-WarpX/WarpX/issues/3231#issuecomment-1203174466
# This is a 2D electrostatic simulation containing a circle embedded boundary at
# the center of the grid domain, pre-seeded with equal density of argon ions
# and electrons.
max_step = 25000
warpx.const_dt = 3.99e-13
warpx.do_electrostatic = labframe
warpx.self_fields_required_precision = 1e-3
warpx.self_fields_verbosity = 0
warpx.eb_potential(x,y,z,t) = -10
warpx.eb_implicit_function = -((x-0.00005)**2+(z-0.00005)**2-1e-05**2)
warpx.self_fields_absolute_tolerance = 0.02
algo.load_balance_intervals = 5
algo.load_balance_costs_update = timers
algo.load_balance_efficiency_ratio_threshold = 1.001
algo.load_balance_with_sfc = 0
algo.load_balance_knapsack_factor = 2
amr.n_cell = 128 128
amr.max_grid_size = 16
amr.max_level = 0
geometry.dims = 2
geometry.prob_lo = 0.0 0.0
geometry.prob_hi = 0.0001 0.0001
boundary.field_lo = pec periodic
boundary.field_hi = pec periodic
boundary.particle_lo = reflecting periodic
boundary.particle_hi = reflecting periodic
algo.particle_shape = 1
particles.species_names = electrons ar_ions
electrons.species_type = electron
electrons.injection_style = NRandomPerCell
electrons.initialize_self_fields = 1
electrons.num_particles_per_cell = 5
electrons.density = 10.2e16
electrons.profile = constant
electrons.momentum_distribution_type = maxwell_boltzmann
electrons.theta = (kb*30000/(m_e*clight^2))
electrons.save_particles_at_zlo = 1
electrons.save_particles_at_zhi = 1
electrons.save_particles_at_eb = 1
ar_ions.species_type = argon
ar_ions.charge = q_e
ar_ions.injection_style = NRandomPerCell
ar_ions.initialize_self_fields = 1
ar_ions.num_particles_per_cell = 5
ar_ions.profile = constant
ar_ions.density = 10.2e16
ar_ions.momentum_distribution_type = gaussian
ar_ions.ux_m = 0.0
ar_ions.uy_m = 0.0
ar_ions.uz_m = 0.0
ar_ions.ux_th = 2.6285641070031447e-06
ar_ions.uy_th = 2.6285641070031447e-06
ar_ions.uz_th = 2.6285641070031447e-06
ar_ions.save_particles_at_zlo = 1
ar_ions.save_particles_at_zhi = 1
ar_ions.save_particles_at_eb = 1
diagnostics.diags_names = diag1
diag1.diag_type = Full
diag1.format = plotfile
diag1.intervals = 1
diag1.fields_to_plot = phi rho_electrons rho_ar_ions```
One thought I had in looking at the code is that the arithmetic in AMReX's enforcePeriodic
feels a bit sketchy. It does iterated addition/subtraction of the domain length from the particle's position, until it's in bounds. That can potentially accumulate error, but I don't know where particles would then get lost
In investigating this, there are a number of distinct issues in play. I'm digging into the one that's of most immediate concern for Modern Electron's usage, which appears with a switch from DP fields/DP particles to DP fields/SP particles. Going to SP/SP exposes much more dramatic physical differences that we're not prepared to isolate or tackle yet.
Bottom line: Work on this issue can probably be put on hold at least temporarily, if there are other things for you to focus on