diff options
Diffstat (limited to 'Examples/Modules')
11 files changed, 537 insertions, 22 deletions
diff --git a/Examples/Modules/boosted_diags/inputs.3d.slice b/Examples/Modules/boosted_diags/inputs.3d.slice index 08f2310cd..408b18931 100644 --- a/Examples/Modules/boosted_diags/inputs.3d.slice +++ b/Examples/Modules/boosted_diags/inputs.3d.slice @@ -109,3 +109,4 @@ slice.coarsening_ratio = 1 1 1 slice.plot_int = -1 slice.num_slice_snapshots_lab = 4 slice.dt_slice_snapshots_lab = 3.3356409519815207e-12 +slice.particle_slice_width_lab = 2.e-6 diff --git a/Examples/Modules/qed/breit_wheeler/analysis_2d_tau_init.py b/Examples/Modules/qed/breit_wheeler/analysis_2d_tau_init.py index 850ecc0fe..9168c0502 100755 --- a/Examples/Modules/qed/breit_wheeler/analysis_2d_tau_init.py +++ b/Examples/Modules/qed/breit_wheeler/analysis_2d_tau_init.py @@ -1,4 +1,4 @@ -#! /usr/bin/env python3 +#! /usr/bin/env python import yt import numpy as np import scipy.stats as st diff --git a/Examples/Modules/qed/breit_wheeler/analysis_3d_optical_depth_evolution.py b/Examples/Modules/qed/breit_wheeler/analysis_3d_optical_depth_evolution.py new file mode 100755 index 000000000..617afd6f9 --- /dev/null +++ b/Examples/Modules/qed/breit_wheeler/analysis_3d_optical_depth_evolution.py @@ -0,0 +1,119 @@ +#! /usr/bin/env python +# -*- coding: utf-8 -*- + +import yt +import numpy as np +import scipy.stats as st +import sys +import math as m +import scipy.special as spe +import scipy.integrate as integ + +# This script checks if the optical depth of photons undergoing the +# Breit Wheeler process behaves as expected. Four populations of photons +# are initialized with different momenta in a background EM field. The decrease +# of the optical depth (averaged for each population) is used to estimate the +# Breit Wheeler cross section, which is then compared with the theoretical +# formula. Relative discrepancy should be smaller than 2% +# +# References: +# 1) R. Duclous et al 2011 Plasma Phys. Control. Fusion 53 015009 +# 2) A. Gonoskov et al. 2015 Phys. Rev. E 92, 023305 +# 3) M. Lobet. PhD thesis "Effets radiatifs et d'électrodynamique +# quantique dans l'interaction laser-matière ultra-relativiste" +# URL: https://tel.archives-ouvertes.fr/tel-01314224 + + +# Tolerance +tol = 2e-2 + +# EM fields +E_f = np.array([-2433321316961438, 973328526784575, 1459992790176863]) +B_f = np.array([2857142.85714286, 4285714.28571428, 8571428.57142857]) + +# Physical constants +electron_mass = 9.10938356e-31 +elementary_charge = 1.6021766208e-19 +speed_of_light = 299792458 +reduced_plank = 1.054571800e-34 +vacuum_permittivity = 8.854187817e-12 +fine_structure_constant = 0.0072973525664 +classical_elec_radius = (1./4./np.pi/vacuum_permittivity)*( elementary_charge**2 / (electron_mass * speed_of_light**2)) +lambda_ref = 1.0e-6 +field_reference = 2.0 * np.pi * electron_mass*speed_of_light * speed_of_light / (elementary_charge*lambda_ref) +schwinger_field_SI = electron_mass**2 * speed_of_light**3/(reduced_plank*elementary_charge) +schwinger_field_norm = electron_mass*speed_of_light*lambda_ref/(2.0*reduced_plank*m.pi) +#______________ + +# Initial parameters +spec_names = ["p1", "p2", "p3", "p4"] + #Assumes average = 1 at t = 0 (ok for a large number of particles) +tau_begin_avg = np.array([1.0, 1.0, 1.0, 1.0]) +mec = electron_mass*speed_of_light +p_begin = [ + np.array([2000.0,0,0])*mec, + np.array([0.0,5000.0,0.0])*mec, + np.array([0.0,0.0,10000.0])*mec, + np.array([57735.02691896, 57735.02691896, 57735.02691896])*mec +] +#______________ + +def calc_chi_gamma(p, E, B): + p = p / electron_mass / speed_of_light + E = E / field_reference + B = B * speed_of_light / field_reference + gamma_phot = np.linalg.norm(p) + c = p/gamma_phot + loc_field = gamma_phot * np.linalg.norm( E - np.dot(c,E)*c + np.cross(c,B)) + return loc_field/schwinger_field_norm + +#Auxiliary functions +def X(chi_phot, chi_ele): + if (chi_phot > chi_ele and chi_ele != 0): + return np.power(chi_phot/(chi_ele*(chi_phot-chi_ele)), 2./3.) + else: + return 1.0e30 + +def T(chi_phot): + coeff = 1./(np.pi * np.sqrt(3.) * chi_phot * chi_phot) + inner = lambda x : integ.quad(lambda s: np.sqrt(s)*spe.kv(1./3., 2./3. * s**(3./2.)), x, np.inf)[0] + return integ.quad(lambda chi_ele: + coeff*(inner(X(chi_phot, chi_ele)) - + (2.0 - chi_phot*np.power(X(chi_phot, chi_ele), 3./2.))*spe.kv(2./3., 2./3. *X(chi_phot, chi_ele)**(3./2.)) ) + , 0, chi_phot)[0] +#__________________ + +# Breit-Wheeler total cross section +def dNBW_dt(chi_phot, energy_phot): + energy_phot = energy_phot/electron_mass/speed_of_light/speed_of_light + return ((electron_mass*(speed_of_light)**2)*fine_structure_constant/reduced_plank)*(chi_phot/energy_phot)*T(chi_phot) +#__________________ + +def check(): + filename_end = sys.argv[1] + data_set_end = yt.load(filename_end) + + sim_time = data_set_end.current_time.to_value() + + all_data_end = data_set_end.all_data() + + tau_end_avg = np.array([ + np.average(all_data_end[name, 'particle_tau']) + for name in spec_names]) + + dNBW_dt_sim = (tau_begin_avg - tau_end_avg)/sim_time + + dNBW_dt_theo = [ + dNBW_dt(calc_chi_gamma(p, E_f, B_f), np.linalg.norm(p*speed_of_light)) + for p in p_begin + ] + + discrepancy = np.abs(dNBW_dt_sim-dNBW_dt_theo)/dNBW_dt_theo + + assert(np.all(discrepancy < tol)) + +def main(): + check() + +if __name__ == "__main__": + main() diff --git a/Examples/Modules/qed/breit_wheeler/inputs.2d_test_tau_init b/Examples/Modules/qed/breit_wheeler/inputs.2d_test_tau_init index 8e3dd29ca..050414185 100644 --- a/Examples/Modules/qed/breit_wheeler/inputs.2d_test_tau_init +++ b/Examples/Modules/qed/breit_wheeler/inputs.2d_test_tau_init @@ -65,16 +65,19 @@ photons.do_qed_breit_wheeler = 1 ################################# #######QED ENGINE PARAMETERS##### -qed_bw.chi_min = 0.01 -qed_bw.ignore_tables_for_test = 1 -#qed_bw.generate_table = 1 -#qed_bw.tab_dndt_chi_min = 0.01 -#qed_bw.tab_dndt_chi_max = 100 -#qed_bw.tab_dndt_how_many = 200 +qed_bw.lookup_table_mode = "dummy_builtin" + +#qed_bw.lookup_table_mode = "generate" +#qed_bw.chi_min = 0.001 +#qed_bw.tab_dndt_chi_min = 0.1 +#qed_bw.tab_dndt_chi_max = 200 +#qed_bw.tab_dndt_how_many = 64 #qed_bw.tab_pair_chi_min = 0.01 -#qed_bw.tab_pair_chi_max = 100 -#qed_bw.tab_pair_chi_how_many = 30 -#qed_bw.tab_pair_frac_how_many = 30 -#qed_bw.save_table_in = "bw_table" -#qed_bw.load_table_from = "bw_table" +#qed_bw.tab_pair_chi_max = 200 +#qed_bw.tab_pair_chi_how_many = 2 +#qed_bw.tab_pair_frac_how_many = 2 +#qed_bw.save_table_in = "bw_micro_table" + +#qed_bw.lookup_table_mode = "load" +#qed_bw.load_table_from = "bw_micro_table" ################################# diff --git a/Examples/Modules/qed/breit_wheeler/inputs.3d_test_optical_depth_evolution b/Examples/Modules/qed/breit_wheeler/inputs.3d_test_optical_depth_evolution new file mode 100644 index 000000000..0d4cad2b1 --- /dev/null +++ b/Examples/Modules/qed/breit_wheeler/inputs.3d_test_optical_depth_evolution @@ -0,0 +1,166 @@ +################################# +####### GENERAL PARAMETERS ###### +################################# +max_step = 10 +amr.n_cell = 64 64 64 +amr.max_grid_size = 32 # maximum size of each AMReX box, used to decompose the domain +amr.blocking_factor = 8 # minimum size of each AMReX box, used to decompose the domain +amr.plot_int = 10 +geometry.coord_sys = 0 # 0: Cartesian +geometry.is_periodic = 1 1 1 # Is periodic? +geometry.prob_lo = -1.e-6 -1.e-6 -1e-6 # physical domain +geometry.prob_hi = 1.e-6 1.e-6 1e-6 +amr.max_level = 0 # Maximum level in hierarchy (1 might be unstable, >1 is not supported) + +################################# +############ NUMERICS ########### +################################# +algo.current_deposition = esirkepov +algo.charge_deposition = standard +algo.field_gathering = energy-conserving +algo.particle_pusher = boris +interpolation.nox = 3 # Particle interpolation order. Must be the same in x, y, and z +interpolation.noy = 3 +interpolation.noz = 3 +warpx.verbose = 1 +warpx.do_dive_cleaning = 0 +warpx.plot_raw_fields = 0 +warpx.plot_raw_fields_guards = 0 +warpx.use_filter = 1 +warpx.cfl = 1. # if 1., the time step is set to its CFL limit +warpx.do_pml = 1 # use Perfectly Matched Layer as boundary condition +warpx.serialize_ics = 1 + +################################# +############ PLASMA ############# +################################# +particles.nspecies = 4 # number of species +particles.species_names = p1 p2 p3 p4 +particles.photon_species = p1 p2 p3 p4 +################################# + +p1.charge = -q_e +p1.mass = m_e +p1.injection_style = "NUniformPerCell" +p1.profile = "constant" +p1.xmin = -0.5e-6 +p1.ymin = -0.5e-6 +p1.zmin = -0.5e-6 +p1.xmax = 0.5e6 +p1.ymax = 0.5e6 +p1.zmax = 0.5e6 +p1.num_particles_per_cell_each_dim = 2 2 2 +p1.density = 1e19 +p1.profile = "constant" +p1.momentum_distribution_type = "gaussian" +p1.ux_m = 2000.0 +p1.uy_m = 0.0 +p1.uz_m = 0.0 +p1.ux_th = 0. +p1.uy_th = 0. +p1.uz_th = 0. +##########QED#################### +p1.do_qed = 1 +p1.do_qed_breit_wheeler = 1 +################################# + +p2.charge = -q_e +p2.mass = m_e +p2.injection_style = "NUniformPerCell" +p2.profile = "constant" +p2.xmin = -0.5e-6 +p2.ymin = -0.5e-6 +p2.zmin = -0.5e-6 +p2.xmax = 0.5e6 +p2.ymax = 0.5e6 +p2.zmax = 0.5e6 +p2.num_particles_per_cell_each_dim = 2 2 2 +p2.density = 1e19 +p2.profile = "constant" +p2.momentum_distribution_type = "gaussian" +p2.ux_m = 0.0 +p2.uy_m = 5000.0 +p2.uz_m = 0.0 +p2.ux_th = 0. +p2.uy_th = 0. +p2.uz_th = 0. +##########QED#################### +p2.do_qed = 1 +p2.do_qed_breit_wheeler = 1 +################################# + +p3.charge = -q_e +p3.mass = m_e +p3.injection_style = "NUniformPerCell" +p3.profile = "constant" +p3.xmin = -0.5e-6 +p3.ymin = -0.5e-6 +p3.zmin = -0.5e-6 +p3.xmax = 0.5e6 +p3.ymax = 0.5e6 +p3.zmax = 0.5e6 +p3.num_particles_per_cell_each_dim = 2 2 2 +p3.density = 1e19 +p3.profile = "constant" +p3.momentum_distribution_type = "gaussian" +p3.ux_m = 0.0 +p3.uy_m = 0.0 +p3.uz_m = 10000.0 +p3.ux_th = 0. +p3.uy_th = 0. +p3.uz_th = 0. +##########QED#################### +p3.do_qed = 1 +p3.do_qed_breit_wheeler = 1 +################################# + +p4.charge = -q_e +p4.mass = m_e +p4.injection_style = "NUniformPerCell" +p4.profile = "constant" +p4.xmin = -0.5e-6 +p4.ymin = -0.5e-6 +p4.zmin = -0.5e-6 +p4.xmax = 0.5e6 +p4.ymax = 0.5e6 +p4.zmax = 0.5e6 +p4.num_particles_per_cell_each_dim = 2 2 2 +p4.density = 1e19 +p4.profile = "constant" +p4.momentum_distribution_type = "gaussian" +p4.ux_m = 57735.02691896 +p4.uy_m = 57735.02691896 +p4.uz_m = 57735.02691896 +p4.ux_th = 0. +p4.uy_th = 0. +p4.uz_th = 0. +##########QED#################### +p4.do_qed = 1 +p4.do_qed_breit_wheeler = 1 +################################# + +##########QED TABLES#################### +qed_bw.lookup_table_mode = "dummy_builtin" + +#qed_bw.lookup_table_mode = "generate" +#qed_bw.chi_min = 0.001 +#qed_bw.tab_dndt_chi_min = 0.1 +#qed_bw.tab_dndt_chi_max = 200 +#qed_bw.tab_dndt_how_many = 64 +#qed_bw.tab_pair_chi_min = 0.01 +#qed_bw.tab_pair_chi_max = 200 +#qed_bw.tab_pair_chi_how_many = 2 +#qed_bw.tab_pair_frac_how_many = 2 +#qed_bw.save_table_in = "bw_micro_table" + +#qed_bw.lookup_table_mode = "load" +#qed_bw.load_table_from = "bw_micro_table" +################################# + +### EXTERNAL E FIELD ### (3e15 * [-0.811107105653813 0.324442842261525 0.486664263392288] ) +warpx.E_external_particle = -2433321316961438 973328526784575 1459992790176863 +#### + +### EXTERNAL B FIELD ### (1e7 * [0.28571429 0.42857143 0.85714286] ) +warpx.B_external_particle = 2857142.85714286 4285714.28571428 8571428.57142857 +#### diff --git a/Examples/Modules/qed/quantum_synchrotron/analysis_2d_tau_init.py b/Examples/Modules/qed/quantum_synchrotron/analysis_2d_tau_init.py index 05b313ee6..98c95d5ea 100755 --- a/Examples/Modules/qed/quantum_synchrotron/analysis_2d_tau_init.py +++ b/Examples/Modules/qed/quantum_synchrotron/analysis_2d_tau_init.py @@ -1,4 +1,4 @@ -#! /usr/bin/env python3 +#! /usr/bin/env python import yt import numpy as np import scipy.stats as st diff --git a/Examples/Modules/qed/quantum_synchrotron/inputs.2d_test_tau_init b/Examples/Modules/qed/quantum_synchrotron/inputs.2d_test_tau_init index 03a2a5798..3fa361a2f 100644 --- a/Examples/Modules/qed/quantum_synchrotron/inputs.2d_test_tau_init +++ b/Examples/Modules/qed/quantum_synchrotron/inputs.2d_test_tau_init @@ -90,16 +90,19 @@ positrons.do_qed_quantum_sync = 1 ################################# #######QED ENGINE PARAMETERS##### -qed_qs.chi_min = 0.001 -qed_qs.ignore_tables_for_test = 1 -#qed_qs.generate_table = 1 +qed_qs.lookup_table_mode = "dummy_builtin" + +#qed_qs.lookup_table_mode = "generate" +#qed_qs.chi_min = 0.001 #qed_qs.tab_dndt_chi_min = 0.001 -#qed_qs.tab_dndt_chi_max = 100 -#qed_qs.tab_dndt_how_many = 200 -#qed_qs.tab_em_chi_min = 0.01 -#qed_qs.tab_em_chi_max = 100 +#qed_qs.tab_dndt_chi_max = 200 +#qed_qs.tab_dndt_how_many = 32 +#qed_qs.tab_em_chi_min = 0.001 +#qed_qs.tab_em_chi_max = 200 #qed_qs.tab_em_chi_how_many = 2 #qed_qs.tab_em_prob_how_many = 2 -#qed_qs.save_table_in = "qs_table" -#qed_qs.load_table_from = "qs_table" +#qed_qs.save_table_in = "qs_micro_table" + +#qed_qs.lookup_table_mode = "load" +#qed_qs.load_table_from = "qs_micro_table" ################################# diff --git a/Examples/Modules/relativistic_space_charge_initialization/analysis.py b/Examples/Modules/relativistic_space_charge_initialization/analysis.py new file mode 100755 index 000000000..aad4534f5 --- /dev/null +++ b/Examples/Modules/relativistic_space_charge_initialization/analysis.py @@ -0,0 +1,75 @@ +#! /usr/bin/env python +""" +This script checks the space-charge initialization routine, by +verifying that the space-charge field of a Gaussian beam corresponds to +the expected theoretical field. +""" +import sys +import yt +import numpy as np +import matplotlib.pyplot as plt +import scipy.constants as scc +from scipy.special import gammainc +yt.funcs.mylog.setLevel(0) + +# Parameters from the Simulation +Qtot = -1.e-20 +r0 = 2.e-6 + +# Open data file +filename = sys.argv[1] +ds = yt.load( filename ) +# Extract data +ad0 = ds.covering_grid(level=0, left_edge=ds.domain_left_edge, dims=ds.domain_dimensions) +Ex_array = ad0['Ex'].to_ndarray().squeeze() +By_array = ad0['By'].to_ndarray() + +# Extract grid coordinates +Nx, Ny, Nz = ds.domain_dimensions +xmin, ymin, zmin = ds.domain_left_edge.v +Lx, Ly, Lz = ds.domain_width.v +x = xmin + Lx/Nx*(0.5+np.arange(Nx)) +y = ymin + Ly/Ny*(0.5+np.arange(Ny)) +z = zmin + Lz/Nz*(0.5+np.arange(Nz)) + +# Compute theoretical field +x_2d, y_2d, z_2d = np.meshgrid(x, y, z, indexing='ij') +r2 = x_2d**2 + y_2d**2 +factor = Qtot/scc.epsilon_0/(2*np.pi*r2) * (1-np.exp(-r2/(2*r0**2))) +factor_z = 1./(2*np.pi)**.5/r0 * np.exp(-z_2d**2/(2*r0**2)) +Ex_th = factor*factor_z*x_2d +Ey_th = factor*factor_z*y_2d + +# Plot theory and data +def make_2d(arr): + if arr.ndim == 3: + return arr[:,Ny//2,:] + else: + return arr +plt.figure(figsize=(10,10)) +plt.subplot(221) +plt.title('Ex: Theory') +plt.imshow(make_2d(Ex_th)) +plt.colorbar() +plt.subplot(222) +plt.title('Ex: Simulation') +plt.imshow(make_2d(Ex_array)) +plt.colorbar() +plt.subplot(223) +plt.title('By: Theory') +plt.imshow(make_2d(Ex_th/scc.c)) +plt.colorbar() +plt.subplot(224) +plt.title('Bz: Simulation') +plt.imshow(make_2d(By_array)) +plt.colorbar() + +plt.savefig('Comparison.png') + +# Automatically check the results +def check(E, E_th, label): + print( 'Relative error in %s: %.3f'%( + label, abs(E-E_th).max()/E_th.max())) + assert np.allclose( E, E_th, atol=0.1*E_th.max() ) + +check( Ex_array, Ex_th, 'Ex' ) diff --git a/Examples/Modules/relativistic_space_charge_initialization/inputs b/Examples/Modules/relativistic_space_charge_initialization/inputs new file mode 100644 index 000000000..eca38e074 --- /dev/null +++ b/Examples/Modules/relativistic_space_charge_initialization/inputs @@ -0,0 +1,30 @@ +max_step = 0 +amr.n_cell = 64 64 64 +amr.max_grid_size = 32 +amr.max_level = 0 + +amr.plot_int = 1 +geometry.coord_sys = 0 # 0: Cartesian +geometry.is_periodic = 0 0 0 # Is periodic? +geometry.prob_lo = -50.e-6 -50.e-6 -50.e-6 +geometry.prob_hi = 50.e-6 50.e-6 50.e-6 + +particles.nspecies = 1 +particles.species_names = beam +beam.charge = -q_e +beam.mass = m_e +beam.injection_style = "gaussian_beam" +beam.initialize_self_fields = 1 +beam.self_fields_required_precision = 1.e-4 +beam.x_rms = 2.e-6 +beam.y_rms = 2.e-6 +beam.z_rms = 2.e-6 +beam.x_m = 0. +beam.y_m = 0. +beam.z_m = 0.e-6 +beam.npart = 20000 +beam.q_tot = -1.e-20 +beam.momentum_distribution_type = "gaussian" +beam.ux_m = 0.0 +beam.uy_m = 0.0 +beam.uz_m = 100.0 diff --git a/Examples/Modules/space_charge_initialization/analysis.py b/Examples/Modules/space_charge_initialization/analysis.py new file mode 100755 index 000000000..67039bad1 --- /dev/null +++ b/Examples/Modules/space_charge_initialization/analysis.py @@ -0,0 +1,89 @@ +#! /usr/bin/env python +""" +This script checks the space-charge initialization routine, by +verifying that the space-charge field of a Gaussian beam corresponds to +the expected theoretical field. +""" +import sys +import yt +import numpy as np +import matplotlib.pyplot as plt +import scipy.constants as scc +from scipy.special import gammainc +yt.funcs.mylog.setLevel(0) + +# Parameters from the Simulation +Qtot = -1.e-20 +r0 = 2.e-6 + +# Open data file +filename = sys.argv[1] +ds = yt.load( filename ) +# Extract data +ad0 = ds.covering_grid(level=0, left_edge=ds.domain_left_edge, dims=ds.domain_dimensions) +Ex_array = ad0['Ex'].to_ndarray().squeeze() +if ds.dimensionality == 2: + # Rename the z dimension as y, so as to make this script work for 2d and 3d + Ey_array = ad0['Ez'].to_ndarray().squeeze() +elif ds.dimensionality == 3: + Ey_array = ad0['Ey'].to_ndarray() + Ez_array = ad0['Ez'].to_ndarray() + +# Extract grid coordinates +Nx, Ny, Nz = ds.domain_dimensions +xmin, ymin, zmin = ds.domain_left_edge.v +Lx, Ly, Lz = ds.domain_width.v +x = xmin + Lx/Nx*(0.5+np.arange(Nx)) +y = ymin + Ly/Ny*(0.5+np.arange(Ny)) +z = zmin + Lz/Nz*(0.5+np.arange(Nz)) + +# Compute theoretical field +if ds.dimensionality == 2: + x_2d, y_2d = np.meshgrid(x, y, indexing='ij') + r2 = x_2d**2 + y_2d**2 + factor = (Qtot/r0)/(2*np.pi*scc.epsilon_0*r2) * (1-np.exp(-r2/(2*r0**2))) + Ex_th = x_2d * factor + Ey_th = y_2d * factor +elif ds.dimensionality == 3: + x_2d, y_2d, z_2d = np.meshgrid(x, y, z, indexing='ij') + r2 = x_2d**2 + y_2d**2 + z_2d**2 + factor = Qtot/(4*np.pi*scc.epsilon_0*r2**1.5) * gammainc(3./2, r2/(2.*r0**2)) + Ex_th = factor*x_2d + Ey_th = factor*y_2d + Ez_th = factor*z_2d + +# Plot theory and data +def make_2d(arr): + if arr.ndim == 3: + return arr[:,:,Nz//2] + else: + return arr +plt.figure(figsize=(10,10)) +plt.subplot(221) +plt.title('Ex: Theory') +plt.imshow(make_2d(Ex_th)) +plt.colorbar() +plt.subplot(222) +plt.title('Ex: Simulation') +plt.imshow(make_2d(Ex_array)) +plt.colorbar() +plt.subplot(223) +plt.title('Ey: Theory') +plt.imshow(make_2d(Ey_th)) +plt.colorbar() +plt.subplot(224) +plt.title('Ey: Simulation') +plt.imshow(make_2d(Ey_array)) +plt.colorbar() +plt.savefig('Comparison.png') + +# Automatically check the results +def check(E, E_th, label): + print( 'Relative error in %s: %.3f'%( + label, abs(E-E_th).max()/E_th.max())) + assert np.allclose( E, E_th, atol=0.1*E_th.max() ) + +check( Ex_array, Ex_th, 'Ex' ) +check( Ey_array, Ey_th, 'Ey' ) +if ds.dimensionality == 3: + check( Ez_array, Ez_th, 'Ez' ) diff --git a/Examples/Modules/space_charge_initialization/inputs b/Examples/Modules/space_charge_initialization/inputs new file mode 100644 index 000000000..79309e585 --- /dev/null +++ b/Examples/Modules/space_charge_initialization/inputs @@ -0,0 +1,29 @@ +max_step = 0 +amr.n_cell = 64 64 64 +amr.max_grid_size = 32 +amr.max_level = 0 + +amr.plot_int = 1 +geometry.coord_sys = 0 # 0: Cartesian +geometry.is_periodic = 0 0 0 # Is periodic? +geometry.prob_lo = -50.e-6 -50.e-6 -50.e-6 +geometry.prob_hi = 50.e-6 50.e-6 50.e-6 + +particles.nspecies = 1 +particles.species_names = beam +beam.charge = -q_e +beam.mass = m_e +beam.injection_style = "gaussian_beam" +beam.initialize_self_fields = 1 +beam.x_rms = 2.e-6 +beam.y_rms = 2.e-6 +beam.z_rms = 2.e-6 +beam.x_m = 0. +beam.y_m = 0. +beam.z_m = 0.e-6 +beam.npart = 20000 +beam.q_tot = -1.e-20 +beam.momentum_distribution_type = "gaussian" +beam.ux_m = 0.0 +beam.uy_m = 0.0 +beam.uz_m = 0.0 |