update PWFA example

Docs/source/usage/examples.rst

@@ -23,7 +23,13 @@ Plasma-Based Acceleration
 
    examples/lwfa/README.rst
    examples/pwfa/README.rst
-   pwfa.rst
+
+Coming soon:
+
+* LWFA: External injection in the boosted frame
+* LWFA: Ionization injection in the lab frame using a LASY data file
+* PWFA: Self-injection in the lab frame
+* MR case?
 
 
 Laser-Plasma Interaction

Docs/source/usage/faq.rst

@@ -63,6 +63,7 @@ An in-depth discussion of the boosted frame is provided in the :ref:`moving wind
 
 What about Back-transformed diagnostics (BTD)?
 ----------------------------------------------
+.. _faq-btd:
 
 .. figure:: https://user-images.githubusercontent.com/10621396/198702232-9dd595ad-479e-4170-bd25-51e2b72cd50a.png
    :alt: [fig:BTD_features] Minkowski diagram indicating several features of the back-transformed diagnostic (BTD). The diagram explains why the first BTD begins to fill at boosted time :math:`t'=0` but this doesn't necessarily correspond to lab time :math:`t=0`, how the BTD grid-spacing is determined by the boosted time step :math:`\Delta t'`, hence why the snapshot length don't correspond to the grid spacing and length in the input script, and how the BTD snapshots complete when the effective snapshot length is covered in the boosted frame.

Docs/source/usage/python.rst

@@ -95,6 +95,7 @@ Diagnostics
 .. autoclass:: pywarpx.picmi.Checkpoint
 
 .. autoclass:: pywarpx.picmi.ReducedDiagnostic
+.. _python_lab_frame_diag:
 
 Lab-frame diagnostics diagnostics are used when running boosted-frame simulations.
 

Examples/Physics_applications/plasma_acceleration/PICMI_inputs_plasma_acceleration_3d.py

@@ -0,0 +1,172 @@
+#!/usr/bin/env python3
+
+from pywarpx import picmi
+
+c = picmi.constants.c
+q_e = picmi.constants.q_e
+m_e = picmi.constants.m_e
+m_p = picmi.constants.m_p
+ep0 = picmi.constants.ep0
+
+nx = 64
+nz = 64
+
+xmin = -100.e-6
+xmax = +100.e-6
+zmin = -149.e-6
+zmax = +1.e-6
+
+moving_window_velocity = [0., 0., c]
+gamma_boost = 10.0
+
+npc_x = 1
+npc_z = 1
+
+max_steps = 10
+
+field_bc_z='open'
+field_bc_x='open'
+grid = picmi.Cartesian3DGrid(
+    number_of_cells = [nx, nx, nz],
+    lower_bound = [xmin, xmin, zmin],
+    upper_bound = [xmax, xmax, zmax],
+    lower_boundary_conditions = [field_bc_x, field_bc_x, field_bc_z],
+    upper_boundary_conditions = [field_bc_x, field_bc_x, field_bc_z],
+    lower_boundary_conditions_particles = ['absorbing']*3,
+    upper_boundary_conditions_particles = ['absorbing']*3,
+    moving_window_velocity = moving_window_velocity,
+    warpx_blocking_factor = 32,
+    warpx_max_grid_size=128
+)
+
+smoother = picmi.BinomialSmoother(n_pass=[1]*3)
+solver = picmi.ElectromagneticSolver(
+    method='CKC',
+    grid=grid,
+    cfl=0.99,
+    source_smoother=smoother,
+    warpx_pml_ncell=10)
+
+drive_charge = -10e-11
+driver_gaussian_distribution = picmi.GaussianBunchDistribution(
+    n_physical_particles= abs(drive_charge) / q_e,
+    rms_bunch_size=[2e-6, 2e-6, 4.e-6],
+    rms_velocity=[2.*c, 2.*c, 2.e4*c],
+    centroid_position=[0., 0., -20.e-6],
+    centroid_velocity=[0., 0., 2.e5*c],
+)
+beam_charge = -3.e-13
+beam_z0 = -90e-6
+beam_gaussian_distribution = picmi.GaussianBunchDistribution(
+    n_physical_particles= abs(beam_charge) / q_e,
+    rms_bunch_size=[0.5e-6, 0.5e-6, 1.e-6],
+    rms_velocity=[2.*c, 2.*c, 200.*c],
+    centroid_position=[0., 0., beam_z0],
+    centroid_velocity=[0., 0., 4000.*c],
+)
+
+dens = 1.e23
+lramp = 8.e-3
+plasma_distribution = picmi.AnalyticDistribution(
+    density_expression=
+        f'(z<lramp)*0.5*(1-cos(pi*z/lramp))*dens+(z>lramp)*dens',
+    pi=3.141592653589793,
+    dens=dens,
+    lramp = lramp,
+    lower_bound=[-70.e-6, -70.e-6, 0.],
+    upper_bound=[70.e-6, 70.e-6, 0.2],
+    fill_in=True)
+
+beam = picmi.Species(
+    particle_type='electron',
+    name='beam',
+    initial_distribution=beam_gaussian_distribution)
+driver = picmi.Species(
+    particle_type='electron',
+    name='driver',
+    initial_distribution=driver_gaussian_distribution)
+plasma_e = picmi.Species(
+    particle_type='electron',
+    name='plasma_e',
+    initial_distribution=plasma_distribution)
+plasma_i = picmi.Species(
+    particle_type='proton',
+    name='plasma_p',
+    initial_distribution=plasma_distribution)
+
+sim = picmi.Simulation(
+    solver = solver,
+    max_steps = max_steps,
+    verbose = 1,
+    warpx_use_filter = 0,
+    gamma_boost = gamma_boost,
+    warpx_use_fdtd_nci_corr = True,
+    particle_shape='cubic'
+)
+
+sim.add_species_through_plane(
+    beam,
+    layout=picmi.PseudoRandomLayout(
+        grid=grid,
+        n_macroparticles=1000
+    ),
+    injection_plane_position=0.,
+    injection_plane_normal_vector=[0.,0.,1.]
+)
+sim.add_species_through_plane(
+    driver,
+    layout=picmi.PseudoRandomLayout(
+        grid=grid,
+        n_macroparticles=1000
+    ),
+    injection_plane_position=0.,
+    injection_plane_normal_vector=[0.,0.,1.]
+)
+sim.add_species(
+    plasma_e,
+    layout=picmi.GriddedLayout(
+        grid=grid,
+        n_macroparticle_per_cell=[npc_x, npc_x, npc_z]
+    )
+)
+sim.add_species(
+    plasma_i,
+    layout=picmi.GriddedLayout(
+        grid=grid,
+        n_macroparticle_per_cell=[npc_x, npc_x, npc_z]
+    )
+)
+
+field_diag = picmi.FieldDiagnostic(
+    name = 'Python_PlasmaAcceleration_plt',
+    grid = grid,
+    period = max_steps,
+    data_list = ['Ex', 'Ey', 'Ez', 'Jx', 'Jy', 'Jz', 'part_per_cell'],
+    # write_dir='Python_PlasmaAcceleration_plt',
+    warpx_format='openpmd',
+)
+
+part_diag = picmi.ParticleDiagnostic(
+    name = 'Python_PlasmaAcceleration_plt',
+    period = max_steps,
+    species = [driver, beam, plasma_e],
+    # write_dir='Python_PlasmaAcceleration_plt',
+    warpx_format='openpmd',
+)
+
+beamrel_red_diag = picmi.ReducedDiagnostic(
+    diag_type='BeamRelevant',
+    name='beamrel',
+    species=beam,
+    period=10)
+
+sim.add_diagnostic(field_diag)
+sim.add_diagnostic(part_diag)
+sim.add_diagnostic(beamrel_red_diag)
+
+# write_inputs will create an inputs file that can be used to run
+# with the compiled version.
+sim.write_input_file(file_name = 'inputs_from_PICMI')
+
+# Alternatively, sim.step will run WarpX, controlling it from Python
+sim.step()

Examples/Physics_applications/plasma_acceleration/README.rst

@@ -11,20 +11,40 @@ On the other hand, for physical situation with strong asymmetries (e.g., non-rou
 
 Additionally, to speed up computation, this example uses the :ref:`boosted frame method <theory-boostedframe>` to effectively model long acceleration lengths.
 
-Alternatively, an other common approximation for PWFAs is quasi-static modeling, e.g., if effects such as self-injection can be ignored.
+Alternatively, another common approximation for PWFAs is quasi-static modeling, useful if effects such as self-injection can be ignored.
 In the Beam, Plasma & Accelerator Simulation Toolkit (BLAST), `HiPACE++ <https://hipace.readthedocs.io>`__ provides such methods.
 
 .. note::
 
-   TODO: The Python (PICMI) input file should use the boosted frame method, like the ``inputs_3d_boost`` file.
+    There are two FDTD Maxwell field solvers that compute the field push implemented in WarpX: the Yee and Cole-Karkkainen solver with Cowan coefficients (CKC) solvers.
+    The later includes a modification that allows the numerical dispersion of light in vacuum to be exact, and that is why we choose CKC for the example.
+
+Lorentz boosted frame
+---------------------
+
+    WarpX simulations can be done in the laboratory or a :ref:`Lorentz-boosted <theory-boostedframe>` frame.
+    In the laboratory frame, there is typically no need to model the plasma ions species, since they are mainly stationary during the short time scales associated with the motion of plasma electrons.
+    In the boosted frame, that argument is no longer valid, as ions have relativistic velocities.
+
+    To avoid spurious effects in the boosted frame, we want the longitudinal cell size to be larger than the transverse one.
+    Translating this condition to the cell transverse (:math:`d_{x}`) and longitudinal dimensions (:math:`d_{z}`) in the laboratory frame leads to: :math:`d_{x} > (d_{z} (1+\beta_{b}) \gamma_{b})`,
+    where :math:`\beta_{b}` is the boosted frame velocity in units of :math:`c`.
+
+    To determine the total number of time steps of the simulation, we could either set the ``<zmax_plasma_to_compute_max_step>`` parameter to the end of the plasma (:math:`z_{\textrm{end}}`), or compute it using:
+
+    * boosted frame edge of the simulation box, :math:`\textrm{corner} = l_{e}/ ((1-\beta_{b}) \gamma_{b})`
+    * time of interaction in the boosted frame, :math:`T = \frac{z_{\textrm{end}}/\gamma_{b}-\textrm{corner}}{c (1+\beta_{b})}`
+    * total number of iterations, :math:`i_{\textrm{max}} = T/dt`
+
+    where :math:`l_{e}` is the position of the left edge of the simulation box (in respect to propagation direction).
 
 
 Run
 ---
 
 This example can be run **either** as:
 
-* **Python** script: ``python3 PICMI_inputs_plasma_acceleration.py`` or
+* **Python** script: ``python3 PICMI_inputs_plasma_acceleration_3d.py`` or
 * WarpX **executable** using an input file: ``warpx.3d inputs_3d_boost``
 
 For `MPI-parallel <https://www.mpi-forum.org>`__ runs, prefix these lines with ``mpiexec -n 4 ...`` or ``srun -n 4 ...``, depending on the system.
@@ -33,11 +53,7 @@ For `MPI-parallel <https://www.mpi-forum.org>`__ runs, prefix these lines with `
 
    .. tab-item:: Python: Script
 
-      .. note::
-
-         TODO: This input file should use the boosted frame method, like the ``inputs_3d_boost`` file.
-
-      .. literalinclude:: PICMI_inputs_plasma_acceleration.py
+      .. literalinclude:: PICMI_inputs_plasma_acceleration_3d.py
          :language: python3
          :caption: You can copy this file from ``Examples/Physics_applications/plasma_acceleration/PICMI_inputs_plasma_acceleration.py``.
 
@@ -47,17 +63,40 @@ For `MPI-parallel <https://www.mpi-forum.org>`__ runs, prefix these lines with `
          :language: ini
          :caption: You can copy this file from ``Examples/Physics_applications/plasma_acceleration/inputs_3d_boost``.
 
-Analyze
--------
+Visualize
+---------
+
+The plasma has a density profile that consists of a long up-ramp followed by a uniform plateau.
 
-.. note::
+With this configuration the driver expels the plasma electrons, creating a bubble depleted of plasma electrons.
+This bubble is the characteristic shape of the nonlinear 3D plasma wake.
+This creates an accelerating field :math:`E_z`, as can be seen in :numref:`fig_pwfa_wake`, and transverse focusing forces (not shown).
 
-   This section is TODO.
+.. _fig_pwfa_wake:
 
+.. figure:: https://user-images.githubusercontent.com/10621396/292573839-8d6d7296-48f8-4942-94bd-548683fd35fe.png
+   :alt: Electron beam, in black, accelerated by wakefield, in red-blue colormap, driven by drive electron beam, in purple. Plasma particles are plotted in green. Set max_steps to 800 in the input to generate data for this figure.
+   :width: 100%
 
-Visualize
----------
+   Electron beam, in black, accelerated by wakefield, in red-blue colormap, driven by drive electron beam, in purple. Plasma particles are plotted in green. Set ``max_steps`` to 800 in the input to generate data for this figure. Note that the :math:`E_z` field is plotted with an interpolation option and appears smoother than the numerical resolution would suggest.
 
-.. note::
+The beam increases in energy as it propagates in the wake, as shown in :numref:`fig_pwfa_energy`.
+There is a slight inflection point in the evolution of the beam energy, showing how the beam accelerates as the wake forms around the electron beam and then gains energy linearly in time once the wake stabilizes.
+Note that this is the mean energy in the boosted frame.
+:ref:`Lab Frame <python_lab_frame_diag>`  or :ref:`back-transformed <running-cpp-parameters-diagnostics-btd>` diagnostics can be used to analyze beam energy in the lab frame. See also the :ref:`FAQ <faq-btd>`.
+
+.. _fig_pwfa_energy:
+
+.. figure:: https://user-images.githubusercontent.com/10621396/290962801-97d994f9-d48d-4f76-a37e-d14f6781d680.png
+   :alt: Energy gain of electron beam in wake through this simulation (this is the mean beam energy in the boosted frame, NOT the lab frame). Set ``max_steps`` to 800 in the input to generate data for this figure.
+   :width: 80%
+
+   Energy gain of electron beam in wake through this simulation (this is the mean beam energy in the boosted frame, NOT the lab frame). Set ``max_steps`` to 800 in the input to generate data for this figure.
+
+The plots can be generated with the following script:
+
+.. dropdown:: Script ``plot_3d.py``
 
-   This section is TODO.
+   .. literalinclude:: plot_3d.py
+      :language: python3
+      :caption: You can copy this file from ``Examples/Physics_applications/plasma_acceleration/plot_3d.py``.

Examples/Physics_applications/plasma_acceleration/inputs_2d_boost

@@ -112,5 +112,6 @@ beam.rigid_advance = true
 
 # Diagnostics
 diagnostics.diags_names = diag1
+diag1.format = openpmd
 diag1.intervals = 500
 diag1.diag_type = Full