aboutsummaryrefslogtreecommitdiff
path: root/Docs/source/usage/python.rst
blob: 4052dc0896f09258749e841b3eb94baa58843370 (plain) (blame)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
.. _usage-picmi:

Python (PICMI)
==============

WarpX uses the `PICMI standard <https://github.com/picmi-standard/picmi>`__ for its Python input files.
Python version 3.6 or newer is required.

Example input files can be found in :ref:`the examples section <usage-examples>`.
The examples support running in both modes by commenting and uncommenting the appropriate lines.

.. _usage-picmi-parameters:

Parameters
----------

Simulation and grid setup
^^^^^^^^^^^^^^^^^^^^^^^^^

The `Simulation` object is the central object in a PICMI script.
It defines the simulation time, field solver, registered species, etc.

.. autoclass:: picmistandard.PICMI_Simulation
    :members: step, add_species, add_laser, write_input_file

Field solvers define the updates of electric and magnetic fields.

.. autoclass:: picmistandard.PICMI_ElectromagneticSolver

.. autoclass:: picmistandard.PICMI_ElectrostaticSolver

Grid define the geometry and discretization.

.. autoclass:: picmistandard.PICMI_Cartesian3DGrid

.. autoclass:: picmistandard.PICMI_Cartesian2DGrid

.. autoclass:: picmistandard.PICMI_Cartesian1DGrid

.. autoclass:: picmistandard.PICMI_CylindricalGrid

For convenience, the PICMI interface defines the following constants,
which can be used directly inside any PICMI script. The values are in SI units.

- ``picmi.constants.c``: The speed of light in vacuum.
- ``picmi.constants.ep0``: The vacuum permittivity :math:`\epsilon_0`
- ``picmi.constants.mu0``: The vacuum permeability :math:`\mu_0`
- ``picmi.constants.q_e``: The elementary charge (absolute value of the charge of an electron).
- ``picmi.constants.m_e``: The electron mass
- ``picmi.constants.m_p``: The proton mass

Additionally to self-consistent fields from the field solver, external fields can be applied.

.. autoclass:: picmistandard.PICMI_ConstantAppliedField

.. autoclass:: picmistandard.PICMI_AnalyticAppliedField

.. autoclass:: picmistandard.PICMI_Mirror

Diagnostics can be used to output data.

.. autoclass:: picmistandard.PICMI_ParticleDiagnostic

.. autoclass:: picmistandard.PICMI_FieldDiagnostic

.. autoclass:: picmistandard.PICMI_ElectrostaticFieldDiagnostic

Lab-frame diagnostics diagnostics are used when running boosted-frame simulations.

.. autoclass:: picmistandard.PICMI_LabFrameParticleDiagnostic

.. autoclass:: picmistandard.PICMI_LabFrameFieldDiagnostic

Particles
^^^^^^^^^

Species objects are a collection of particles with similar properties.
For instance, background plasma electrons, background plasma ions and an externally injected beam could each be their own particle species.

.. autoclass:: picmistandard.PICMI_Species

.. autoclass:: picmistandard.PICMI_MultiSpecies

Particle distributions can be used for to initialize particles in a particle species.

.. autoclass:: picmistandard.PICMI_GaussianBunchDistribution

.. autoclass:: picmistandard.PICMI_UniformDistribution

.. autoclass:: picmistandard.PICMI_AnalyticDistribution

.. autoclass:: picmistandard.PICMI_ParticleListDistribution

Particle layouts determine how to microscopically place macro particles in a grid cell.

.. autoclass:: picmistandard.PICMI_GriddedLayout

.. autoclass:: picmistandard.PICMI_PseudoRandomLayout

Lasers
^^^^^^

Laser profiles can be used to initialize laser pulses in the simulation.

.. autoclass:: picmistandard.PICMI_GaussianLaser

.. autoclass:: picmistandard.PICMI_AnalyticLaser

Laser injectors control where to initialize laser pulses on the simulation grid.

.. autoclass:: picmistandard.PICMI_LaserAntenna


.. _usage-picmi-run:

Running
-------

WarpX can be run in one of two modes. It can run as a preprocessor, using the
Python input file to generate an input file to be used by the C++ version, or
it can be run directly from Python.

In either mode, if using a `virtual environment <https://docs.python.org/3/tutorial/venv.html>`__, be sure to activate it before compiling and running WarpX.


Running WarpX directly from Python
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

For this, a full Python installation of WarpX is required, as described in :ref:`the install documentation <install-users>` (:ref:`developers <install-developers>`).

In order to run a new simulation:

* Create a **new directory**, where the simulation will be run.

* Add a **Python script** in the directory.

The input file should have the line ``sim.step()`` which runs the simulation.

* **Run** the script with Python:

.. code-block:: bash

   mpirun -np <n_ranks> python <python_script>

where ``<n_ranks>`` is the number of MPI ranks used, and ``<python_script>``
is the name of the script.


Extending a Simulation from Python
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

When running WarpX directly from Python it is possible to interact with the simulation
by installing ``CallbackFunctions``, which will execute a given Python function at a
specific location in the WarpX simulation loop.

.. autoclass:: pywarpx.callbacks.CallbackFunctions

Places in the WarpX loop where callbacks are available include:
``afterinit``, ``beforecollisions``, ``aftercollisions``, ``beforeEsolve``, ``afterEsolve``,
``beforedeposition``, ``afterdeposition``, ``beforestep``, ``afterstep``, ``afterdiagnostics``,
``afterrestart`` and ``oncheckpointsignal``.
See the examples in *Examples/Tests/ParticleDataPython* for references on how to use
``callbacks``.

There are several "hooks" available via the ``libwarpx`` shared library to access and manipulate
simulation objects (particles, fields and memory buffers) as well as general properties
(such as processor number). These "hooks" are accessible through the `Simulation.extension` object.

.. autofunction:: pywarpx.picmi.Simulation.extension.getNProcs

.. autofunction:: pywarpx.picmi.Simulation.extension.getMyProc

.. autofunction:: pywarpx.picmi.Simulation.extension.get_nattr

.. autofunction:: pywarpx.picmi.Simulation.extension.get_nattr_species

.. autofunction:: pywarpx.picmi.Simulation.extension.getistep

.. autofunction:: pywarpx.picmi.Simulation.extension.gett_new

.. autofunction:: pywarpx.picmi.Simulation.extension.evolve

.. autofunction:: pywarpx.picmi.Simulation.extension.finalize

.. autofunction:: pywarpx.picmi.Simulation.extension.getistep

.. autofunction:: pywarpx.picmi.Simulation.extension.gett_new

.. autofunction:: pywarpx.picmi.Simulation.extension.evolve

.. autofunction:: pywarpx.picmi.Simulation.extension.getProbLo

.. autofunction:: pywarpx.picmi.Simulation.extension.getProbHi

.. autofunction:: pywarpx.picmi.Simulation.extension.getCellSize

Particles can be added to the simulation at specific positions and with specific
attribute values:

.. autofunction:: pywarpx.picmi.Simulation.extension.add_particles

Properties of the particles already in the simulation can be obtained with various
functions.

.. autofunction:: pywarpx.picmi.Simulation.extension.get_particle_count

.. autofunction:: pywarpx.picmi.Simulation.extension.get_particle_structs

.. autofunction:: pywarpx.picmi.Simulation.extension.get_particle_arrays

The ``get_particle_structs()`` and ``get_particle_arrays()`` functions are called
by several utility functions of the form ``get_particle_{comp_name}`` where
``comp_name`` is one of ``x``, ``y``, ``z``, ``r``, ``theta``, ``id``, ``cpu``,
``weight``, ``ux``, ``uy`` or ``uz``.

The index of some specific component of the particle data can be obtained.

.. autofunction:: pywarpx.picmi.Simulation.extension.get_particle_comp_index

New components can be added via Python.

.. autofunction:: pywarpx.picmi.Simulation.extension.add_real_comp

Various diagnostics are also accessible from Python.
This includes getting the deposited or total charge density from a given species
as well as accessing the scraped particle buffer. See the example in
*Examples/Modules/ParticleBoudaryScrape* for a reference on how to interact
with scraped particle data.

.. autofunction:: pywarpx.picmi.Simulation.extension.get_species_charge_sum

.. autofunction:: pywarpx.picmi.Simulation.extension.depositChargeDensity

.. autofunction:: pywarpx.picmi.Simulation.extension.get_particle_boundary_buffer_size

.. autofunction:: pywarpx.picmi.Simulation.extension.get_particle_boundary_buffer_size

.. autofunction:: pywarpx.picmi.Simulation.extension.get_particle_boundary_buffer_structs

.. autofunction:: pywarpx.picmi.Simulation.extension.get_particle_boundary_buffer

.. autofunction:: pywarpx.picmi.Simulation.extension.clearParticleBoundaryBuffer

Using Python input as a preprocessor
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

In this case, only the pure Python version needs to be installed, as described :ref:`here <developers-gnumake-python>`.

In order to run a new simulation:

* Create a **new directory**, where the simulation will be run.

* Add a **Python script** in the directory.

The input file should have the line like ``sim.write_input_file(file_name = 'inputs_from_PICMI')``
which runs the preprocessor, generating the AMReX inputs file.

* **Run** the script with Python:

.. code-block:: bash

   python <python_script>

where ``<python_script>`` is the name of the script.
This creates the WarpX input file that you can run as normal with the WarpX executable.