aboutsummaryrefslogtreecommitdiff
path: root/Examples/Modules/RigidInjection
diff options
context:
space:
mode:
Diffstat (limited to 'Examples/Modules/RigidInjection')
-rwxr-xr-xExamples/Modules/RigidInjection/analysis_rigid_injection_BoostedFrame.py53
-rwxr-xr-xExamples/Modules/RigidInjection/analysis_rigid_injection_LabFrame.py70
-rw-r--r--Examples/Modules/RigidInjection/inputs88
-rw-r--r--Examples/Modules/RigidInjection/inputs.BoostedFrame49
-rw-r--r--Examples/Modules/RigidInjection/inputs.LabFrame41
5 files changed, 213 insertions, 88 deletions
diff --git a/Examples/Modules/RigidInjection/analysis_rigid_injection_BoostedFrame.py b/Examples/Modules/RigidInjection/analysis_rigid_injection_BoostedFrame.py
new file mode 100755
index 000000000..497a30097
--- /dev/null
+++ b/Examples/Modules/RigidInjection/analysis_rigid_injection_BoostedFrame.py
@@ -0,0 +1,53 @@
+#! /usr/bin/env python
+
+'''
+Analysis script of a WarpX simulation of rigid injection in a boosted frame.
+
+A Gaussian electron beam starts from -5 microns, propagates rigidly up to
+20 microns after which it expands due to emittance only (the focal position is
+20 microns). The beam width is measured after ~50 microns, and compared with
+the theory (with a 5% error allowed).
+
+The simulation runs in a boosted frame, and the analysis is done in the lab
+frame, i.e., on the back-transformed diagnostics.
+'''
+
+import sys, os, yt, glob
+import numpy as np
+import scipy.constants as scc
+import read_raw_data
+yt.funcs.mylog.setLevel(0)
+
+# filename = sys.argv[1]
+
+def get_particle_field(snapshot, species, field):
+ fn = snapshot + '/' + species
+ files = glob.glob(os.path.join(fn, field + '_*'))
+ files.sort()
+ all_data = np.array([])
+ for f in files:
+ data = np.fromfile(f)
+ all_data = np.concatenate((all_data, data))
+ return all_data
+
+# Read data from back-transformed diagnostics
+snapshot = './lab_frame_data/snapshot00001'
+header = './lab_frame_data/Header'
+allrd, info = read_raw_data.read_lab_snapshot(snapshot, header)
+z = np.mean( get_particle_field(snapshot, 'beam', 'z') )
+w = np.std ( get_particle_field(snapshot, 'beam', 'x') )
+
+# initial parameters
+z0 = 20.e-6
+w0 = 1.e-6
+theta0 = np.arcsin(0.1)
+
+# Theoretical beam width after propagation if rigid ON
+wth = np.sqrt( w0**2 + (z-z0)**2*theta0**2 )
+error = np.abs((w-wth)/wth)
+
+# Print error and assert small error
+print("Beam position: " + str(z))
+print("Beam width : " + str(w))
+print("error: " + str(error))
+assert( error < 0.03 )
diff --git a/Examples/Modules/RigidInjection/analysis_rigid_injection_LabFrame.py b/Examples/Modules/RigidInjection/analysis_rigid_injection_LabFrame.py
new file mode 100755
index 000000000..86214ad72
--- /dev/null
+++ b/Examples/Modules/RigidInjection/analysis_rigid_injection_LabFrame.py
@@ -0,0 +1,70 @@
+#! /usr/bin/env python
+
+'''
+Analysis script of a WarpX simulation of rigid injection.
+
+A Gaussian electron beam starts from -5 microns, propagates rigidly up to
+20 microns after which it expands due to emittance only (the focal position is
+20 microns). The beam width is measured after ~50 microns, and compared with
+the theory (with a 5% error allowed).
+
+As a help to the user, the script also compares beam width to the theory in
+case rigid injection is OFF (i.e., the beam starts expanding from -5 microns),
+in which case a warning is raised.
+'''
+
+import sys
+import yt
+import numpy as np
+import scipy.constants as scc
+yt.funcs.mylog.setLevel(0)
+
+filename = sys.argv[1]
+
+# WarpX headers include more data when rigid injection is used,
+# which gives an error with the last yt release.
+# To avoid this issue, the three last lines of WarpXHeader are removed if
+# needed.
+def remove_rigid_lines(plotfile, nlines_if_rigid):
+ header_name = plotfile + '/WarpXHeader'
+ f = open(header_name, 'r')
+ file_lines = f.readlines()
+ nlines = len(file_lines)
+ f.close()
+ if nlines == nlines_if_rigid:
+ f = open(header_name, 'w')
+ f.writelines(file_lines[:-3])
+ f.close()
+
+# Remove rigid injection header lines
+remove_rigid_lines(filename, 18)
+# Read beam parameters
+ds = yt.load( filename )
+ad = ds.all_data()
+# Beam longitudinal position
+z = np.mean(ad['beam', 'particle_position_y'].v)
+# Beam width
+w = np.std(ad['beam', 'particle_position_x'].v)
+
+# initial parameters
+z0 = 20.e-6
+z0_no_rigid = -5.e-6
+w0 = 1.e-6
+theta0 = np.arcsin(0.1)
+
+# Theoretical beam width after propagation if rigid OFF
+# Inform the user if rigid injection simply off (just to be kind)
+wth_no_rigid = np.sqrt( w0**2 + (z-z0_no_rigid)**2*theta0**2 )
+error_no_rigid = np.abs((w-wth_no_rigid)/wth_no_rigid)
+if ( error_no_rigid < 0.05):
+ print("error no rigid: " + str(error_no_rigid))
+ print("Looks like the beam defocuses as if rigid injection were OFF")
+
+# Theoretical beam width after propagation if rigid ON
+wth = np.sqrt( w0**2 + (z-z0)**2*theta0**2 )
+error = np.abs((w-wth)/wth)
+# Print error and assert small error
+print("Beam position: " + str(z))
+print("Beam width : " + str(w))
+print("error: " + str(error))
+assert( error < 0.05 )
diff --git a/Examples/Modules/RigidInjection/inputs b/Examples/Modules/RigidInjection/inputs
deleted file mode 100644
index 3c0d52ed0..000000000
--- a/Examples/Modules/RigidInjection/inputs
+++ /dev/null
@@ -1,88 +0,0 @@
-# Maximum number of time steps
-max_step = 60
-
-# number of grid points
-amr.n_cell = 32 32 32
-
-# Maximum allowable size of each subdomain in the problem domain;
-# this is used to decompose the domain for parallel calculations.
-amr.max_grid_size = 16
-
-# Maximum level in hierarchy (for now must be 0, i.e., one level in total)
-amr.max_level = 0
-
-amr.plot_int = 1 # How often to write plotfiles. "<= 0" means no plotfiles.
-
-# Geometry
-geometry.coord_sys = 0 # 0: Cartesian
-geometry.is_periodic = 1 1 0 # Is periodic?
-geometry.prob_lo = -2. -2. -4. # physical domain
-geometry.prob_hi = 2. 2. 4.
-
-# Verbosity
-warpx.verbose = 1
-
-# Algorithms
-
-# interpolation
-interpolation.nox = 3
-interpolation.noy = 3
-interpolation.noz = 3
-
-# CFL
-warpx.cfl = 1.0
-
-# Information about the particle species
-particles.nspecies = 1
-particles.species_names = electrons
-particles.rigid_injected_species = electrons
-
-
-
-#
-# The electron species information
-#
-
-electrons.charge = -q_e
-electrons.mass = m_e
-electrons.injection_style = "gaussian_beam"
-electrons.x_rms = 0.1
-electrons.y_rms = 0.1
-electrons.z_rms = 0.1
-electrons.x_m = 0.
-electrons.y_m = 0.
-electrons.z_m = -1.5
-electrons.npart = 1000
-electrons.q_tot = -8.010883097437485e-07
-
-electrons.profile = "constant"
-electrons.density = 1
-electrons.momentum_distribution_type = "gaussian"
-electrons.ux_m = 0.
-electrons.uy_m = 0.
-electrons.uz_m = 3.
-electrons.ux_th = 0.01
-electrons.uy_th = 0.01
-electrons.uz_th = 0.01
-
-electrons.xmin = -2
-electrons.xmax = 2
-electrons.ymin = -2
-electrons.ymax = 2
-electrons.zmin = -2
-electrons.zmax = 2
-
-electrons.zinject_plane = 0.
-electrons.projected = true
-electrons.focused = false
-
-warpx.do_pml = 0
-
-# Moving window
-warpx.do_moving_window = 0
-warpx.moving_window_dir = z
-warpx.moving_window_v = 1.0 # in units of the speed of light
-
-# Boosted frame
-warpx.gamma_boost = 1.5
-warpx.boost_direction = z
diff --git a/Examples/Modules/RigidInjection/inputs.BoostedFrame b/Examples/Modules/RigidInjection/inputs.BoostedFrame
new file mode 100644
index 000000000..c7a60f14f
--- /dev/null
+++ b/Examples/Modules/RigidInjection/inputs.BoostedFrame
@@ -0,0 +1,49 @@
+# stop_time = 1.5e-13
+
+warpx.zmax_plasma_to_compute_max_step = 50.e-6
+warpx.gamma_boost = 5.
+warpx.boost_direction = z
+warpx.do_boosted_frame_diagnostic = 1
+warpx.num_snapshots_lab = 2
+warpx.dt_snapshots_lab = 1.8679589331096515e-13
+
+amr.n_cell = 256 512 # 32 64
+amr.max_grid_size = 512
+amr.blocking_factor = 16
+amr.max_level = 0
+amr.plot_int = 100000
+geometry.coord_sys = 0 # 0: Cartesian
+geometry.is_periodic = 1 0 # Is periodic?
+geometry.prob_lo = -50.e-6 -10.e-6
+geometry.prob_hi = 50.e-6 0.e-6
+warpx.cfl = .999
+warpx.do_pml = 1
+warpx.do_moving_window = 1
+warpx.moving_window_dir = z
+warpx.moving_window_v = 1.0 # in units of the speed of light
+warpx.serialize_ics = 1
+particles.nspecies = 1
+particles.species_names = beam
+particles.rigid_injected_species = beam
+beam.charge = -q_e
+beam.mass = m_e
+beam.injection_style = "gaussian_beam"
+beam.x_rms = 1.e-6
+beam.y_rms = 1.e-6
+beam.z_rms = .5e-6
+beam.x_m = 0.
+beam.y_m = 0.
+beam.z_m = -5.e-6
+beam.npart = 10000
+beam.q_tot = -1.e-20
+beam.momentum_distribution_type = "gaussian"
+beam.ux_m = 0.0
+beam.uy_m = 0.0
+beam.uz_m = 1000.
+beam.ux_th = 100.
+beam.uy_th = 100.
+beam.uz_th = 0.
+beam.zinject_plane = 20.e-6
+beam.rigid_advance = true
+beam.projected = true
+beam.focused = false
diff --git a/Examples/Modules/RigidInjection/inputs.LabFrame b/Examples/Modules/RigidInjection/inputs.LabFrame
new file mode 100644
index 000000000..730b46e01
--- /dev/null
+++ b/Examples/Modules/RigidInjection/inputs.LabFrame
@@ -0,0 +1,41 @@
+stop_time = 1.5e-13
+amr.n_cell = 32 64
+amr.max_grid_size = 256
+amr.blocking_factor = 16
+amr.max_level = 0
+amr.plot_int = 10000
+geometry.coord_sys = 0 # 0: Cartesian
+geometry.is_periodic = 1 0 # Is periodic?
+geometry.prob_lo = -50.e-6 -10.e-6
+geometry.prob_hi = 50.e-6 0.e-6
+warpx.cfl = .999
+warpx.do_pml = 1
+warpx.do_moving_window = 1
+warpx.moving_window_dir = z
+warpx.moving_window_v = 1.0 # in units of the speed of light
+warpx.serialize_ics = 1
+particles.nspecies = 1
+particles.species_names = beam
+particles.rigid_injected_species = beam
+beam.charge = -q_e
+beam.mass = m_e
+beam.injection_style = "gaussian_beam"
+beam.x_rms = 1.e-6
+beam.y_rms = 1.e-6
+beam.z_rms = .5e-6
+beam.x_m = 0.
+beam.y_m = 0.
+beam.z_m = -5.e-6
+beam.npart = 2000
+beam.q_tot = -1.e-20
+beam.momentum_distribution_type = "gaussian"
+beam.ux_m = 0.0
+beam.uy_m = 0.0
+beam.uz_m = 1000.
+beam.ux_th = 100.
+beam.uy_th = 100.
+beam.uz_th = 0.
+beam.zinject_plane = 20.e-6
+beam.rigid_advance = true
+beam.projected = true
+beam.focused = false