angelfp / wake-t Goto Github PK

The README.md currently shows only the Wake-T logo.

It would be nice to add, at least:

• A description of the code and its working principle.
• Installation instructions.
• Any useful badges such as PyPI, License and Travis status.

Currently, the equations of motion for the particle beams assume that the transverse wakefield experienced by the particles is Wx = Ex - cBy and Wy = E_y + cBx. This is based on the assumption that the beam particles propagate with v_z = c. This assumption should be removed for making the code more general, and the openpmd output should include the transverse components of E and B separately, and not combined in a single W.

In the Baxevanis quasistatic solver, the plasma is modeled as a single column of particles which moves from the front to the back of the simulation grid. Currently, the past positions of the particles as the column is propagated are not stored, meaning that at any moment only a single column of particles exists. This makes it a bit less trivial to deposit the charge in the grid (and in the future also the susceptibility, which is needed for the laser envelope solver).

This is because the current implementation of the charge deposition method charge_distribution_cyl expects a full 3D particle distribution as input, and deposits all its charge at once in a new array with the same size as the grid (plus guard cells). Therefore, if we were to deposit the charge of the plasma particles with this method, only a vertical stripe of charge would be deposited. Since charge_distribution_cyl generates a new charge array every time it's called, it is not currently possible to obtain a full 2D charge density map of the plasma.

To fix this, we could modify charge_distribution_cyl to accept as an input an already existing charge array, so that at every step of the plasma solver we simply add more charge to this existing array. Conceptually this would look something like this in the main loop of the plasma solver (the calculate_wakefields method):

# From the beginning we already have these variables:
xi_max  # longitudinal position front of the simulation box (xi is the longitudinal coordinate in the comoving frame)
r, pr  # numpy arrays containing the transverse position (radius) and momentum of the plasma particles in the column.
n_xi, n_r  # number of grid points along each direction
dxi  # The longitudinal step size

# Initialize charge array (where the charge of plasma particle will be deposited)
rho = np.zeros((nz+4, nr+4))  # Includes guard cells outside of simulation boundaries (2 at each boundary)

# Deposit charge of initial column
rho = charge_distribution_cyl(xi_max, r, rho, other_params)

# Main loop
for step in np.arange(n_xi):
    xi = xi_max - dxi * step

    # Evolve plasma to next xi step.
    r, pr = evolve_plasma(r, pr, q, xi, dxi, laser_params, beam_source, s_d)

    # Deposit charge of updated plasma column
    rho = charge_distribution_cyl(xi, r, rho, other_params)

By doing this, when the main loop is done (i.e. when the column of plasma particles has been propagated until the leftmost edge of the simulation domain) we will have a 2D rho array with the whole plasma density distribution.

We can divide this task in a few steps:

• Modify charge_distribution_cyl to accept as an input the charge array, so that it does not have to generate one every time.
• If an array is given as an input, do not set its values to zero, but simply deposit more charge into it.
• Implement the charge depostion as described in the loop above.
• Once this is working, implement method for the susceptibility deposition, and add to main loop in the same way.
• Finally, both of these 2D fields can be added to the openPMD output files generated by the simulation so that we can easily visualize them.

Ideally, we would have a test where we create a lasy laser, save it to file, and read it into Wake-T. However, this is not yet possible with lasy (mtr data not yet supported in the openpmd-viewer).

When saving the bunch to file, the propagation distance info is lost. This is because wake_t.utilities.bunch_saving.save_bunch_to_file utilises ParticleBunch.get_bunch_matrix, where only the position in the xi plane is saved. So, you do not have any info on position in the beam line when reading bunch from file. I saw that there is get_openpmd_diagnostics_data in ParticleBunch class. Could this be used to save the data when saving to openpmd?

As I understand it, Wake-T includes space-charge effects. In principle it can be used to model beamlines that include plasma lenses, quadrupoles, etc. We were considering to use it for modelling electron bunch refocusing in a solenoid, would this work? Is there a particular example to start from? Particularly for defining the solenoid.

Inspecting the list of BLAST codes, I see some have various degrees of overlap with Wake-T functionality. More specifically,

• HiPACE++ can model both beam driven and laser-driven plasma acceleration stages quasi-statically
• ImpactX can do particle tracking, including space charge and various beam elements

What does Wake-T fit into the BLAST code universe? For instance, is its relation to ImpactX similar to the relation between fbpic and WarpX?

Most of the wakefield models are implemented in normalized units, and use the plasma density at the current z position to do the normalization. This obviously fails when the density is 0, which is something that can happen especially when custom density
profiles with ramps are used.

Dear Wake-T developers,

Thanks for the excellent code. Is it possible to make a plasma density profile as given in the attached figure.

Best,
F. Koester

Right now, when using the get_from_file method to read a particle species, the user has to specify the species_name argument even if there is only one species. Also, this name is not being preserved by Wake-T unless it is manually specified in the name argument. It would be much more convenient if the default behavior was to automatically use the name of the stored species for both arguments.

The method to track a particle bunch through an element currently takes several arguments. This is particularly true for the plasma stages, where a large number of optional arguments have been added over time to adapt to different tracking modes.

Many of these options could be moved to the initial definition of the elements, such that the tracking method does not take any argument except for the bunch.

Tasks to perform for plasma elements:

• Have a single tracking method and not one for the numerical solver and another for the analytical. Simply add an option tracking_mode to the PlasmaStage.__init__ which can take the 'analytical' or 'numerical' values.
• Change the name of current option mode for wakefield_model, which takes a string with the name of the wakefield model to use. Move this option to PlasmaStage.__init__.

In general:

• The tracking method of all elements must have the same name and arguments.
• After this, a general track method could be defined which takes as input a bunch and a beamline (a list with beamline elements) which calls the tracking method of each element and returns the new bunch distribution.

Dear All,

I am looking for a clarification for "ActivePlasmaLens" module.
If I want the plasma lens to include "Wakefields = True" then how "ActivePlasmaLens" parameters should be described.

For example, if I use

drift_1 = Drift(length=0.1) # bunch enter in the APL after drifting for 10 cm in vacuum
drift_2 = Drift(length=0.3) # bunch focusing is detected at 30 cm after coming out of APL

apl = ActivePlasmaLens(length=0.2, for_strength=400, wakefields=True, density=1.0e23, wakefield_model='quasi static_2d')

beamline = Beamline([drift_1, apl, drift_2])

now my question is 5 more parameters are required as positional arguments in "apl". How I should define them in "apl" ?

should I do:

xi_min=10e-2, # entrance of apl
xi_max=30.0e-2, # exit of apl
r_max=500.0e-6, # radial extent
r_min=-500e-6, # radial extent

and could you suggest the criteria of choosing values for n_r and n_xi.

Thank you.

While the charge deposition can be performed with both linear and cubic particle shapes, the field interpolation is currently only linear. An option for field gathering with cubic shapes would be convenient, and would allow for smoother beam phase spaces at lower resolutions.

Currently, any particle distribution read with get_from_file is automatically centered in z around the average of the distribution. This is confusing and should be avoided by default.

Proper support for multi-beam simulations should be added. This is key specially for beam-driven plasma acceleration.

Dear All,

I have a question about calculation of the "focusing strength" of an Active Plasma Lens (APL).

We give directly the value of the focusing strength of the APL in the code. In general, for APL, in my understanding, it should depend on at least two parameters:

(1) Plasma density and its length in the APL, and
(2) Applied voltage or current.

Both these should contribute to the net azimuthal magnetic field. Flow of current depends on the movement of the plasma electrons towards the electrodes (if we consider capillary discharge).

Could you please, help me to understand if we fix the length of APL and plasma density inside the APL then how to determine the focusing strength.

When I follow the reference PRL, 115, 184802 (2015), the focusing strength only depends on the current and radius, but no plasma density dependence.

Thanks a lot.

Similarly to how it's done in APtools, a new workflow should be implemented for automatically publishing new releases to PyPI.

Does it make sense to transition the Wake-T input to the PICMI standard?

Dear @AngelFP,

Is it possible to simulate an electron beam crossing the wakefield propagating perpendicular to its direction, i.e., wakefield is propagating along "z" direction and electron beam is propagating along "r" direction ?

I want to use an electron beam to probe the wakefield and track this electron beam interacting with the wakefield. Do you think it can be done with wake-T ?

Thank you.

Best regards,
Zhang

Hi @AngelFP ,

can we use Wake-T to simulate electron beam transportation through quadrupole, chicane etc. ?

Thank you.

The JIT compilation of numba methods currently takes about ~10s. It would be convenient to enable the numba caching option in order to avoid this step every time a new simulation is started.

There is an issue with numerical precision when calculating the a_i, b_i coefficients, possibly arising from the recursive nature of the problem, which leads to large numbers. Precision problems arise when the number of plasma particles and the radial extent of the plasma increase beyond a certain (not yet well-known) limit. This results in a noisy azimuthal magnetic field b_theta and wrong plasma evolution.

Related PR: #23

The API reference in the documentation is currently mostly empty and not well formatted.

Although it works correctly, the sqrt(2) factor difference in a_0 between the linear and circular polarization should not be there. This difference should be in the ponderomotive force, not on the amplitude.

Hi, I wanted to report an artifact I found when I analyzed an LWFA simulation.
It first became visible in the gradient of the transverse focusing force inside the plasma bubble at the position of the accelerated particle beam.

Image of the normalized focusing force $$\frac{1}{E_0}\frac{\mathrm{d}(E_x - c B_y)}{\mathrm{d}r}$$

Zoom-in:

Here is a Python file (as TXT) that reproduces the run.
run_wake-t.py.txt

It changes with the numerical resolution in $r$. See here an average over all slices at all beam positions.

Here with 10x radial resolution

These spikes flatten out when the number of plasma particles is increased. (The original number of particles was 2)
With ppc=4 it is still visible, with ppc=8 and higher it is basically gone.

This artifact can cause numerical emittance growth in the beam.

Already discussed with @AngelFP and @RemiLehe. This issue is for tracking purposes and can be closed with the implementation of either documentation for users or maybe a correction.

I have 2 issues that arise for my use case:

When initialising a ParticleBunch instance from parameters, (eg get_gaussian_bunch_from_twiss), the naming of the input parameters do not match the naming of the output obtained from analyze_bunch.
Consequently, one cannot simply run

bunch_params = analyze_bunch(ExistingBunch)
NewBunch = get_gaussian_bunch_from_twiss(**bunch_params)

I wrote a function which rewrites the names to match the required input to fix the issue or me, but I guess it could just be renamed in wake-t..?

In analyze_bunch, basically all data (eg energy spread) has to be calculated (even when the bunch was initialised from parameters, not a file) because that initial info is not stored... Saving input parameters in ParticleBunch and/or adding routines to calc them would also resolve my first issue.

There are currently no unit tests checking the correctness and possible bugs in the code. A comprehensive set of tests should be added, as well as the capability of automatically executing them for new PRs.

See: openPMD/openPMD-api#1092

Either revert to old methods for setting dt, time, etc. or wait for fix in the openPMD-api.

In some cases a progress bar is not desirable, for example when running many Wake-T simulations in batches and capturing the output to file. It would be convenient to have a simple argument to turn off progress bars

There is currently a bug where the code will crash when dz_fields = 0

Traceback (most recent call last):
  File "C:\Users\ferran\Repos\Wake-T\examples\track_plasma_qs2d.py", line 45, in <module>
    bunch_list = plasma.track(bunch, opmd_diag=opmd_diag)
  File "C:\Users\ferran\Repos\Wake-T\wake_t\beamline_elements\field_element.py", line 105, in track
    tracker = Tracker(
  File "C:\Users\ferran\Repos\Wake-T\wake_t\tracking\tracker.py", line 82, in __init__
    field.adjust_dt(self.t_final)
  File "C:\Users\ferran\Repos\Wake-T\wake_t\fields\numerical_field.py", line 74, in adjust_dt
    n_updates = np.ceil(t_final / self.dt_update)
ZeroDivisionError: float division by zero

Changes in the upcoming 0.5.0 release break the Wake-T adaptors implemented in Ocelot. They should be updated.

With the current implementation, the calculation of plasma wakefields (and laser envelope evolution) can only be triggered when the track method of the plasma element is called. This therefore requires a ParticleBunch. However, for many applications (e.g. laser guiding optimization), a beam is not needed. It would be interesting to add a different method than track for these kind of simulations.

Since its origin, the only output in Wake-T after tracking a particle distribution through a beamline (or beamline element) has been a simple list of n_out ParticleBunches. To make Wake-T more useful and resemble better what a user typically expects from a plasma simulation code, a proper output to disk should be implemented. This will be done following the openPMD standard by taking advantage of the openPMD-api.

This issue keeps track and compiles all the required and wished features of the diagnostics (as well as bugs/issues which appear during implementation).

Required features:

• Output of particle distributions at each step, as defined by n_out, in:
  • PlasmaStage.
  • PlasmaRamp.
  • PlasmaLens.
  • TMElement.
• Output of all fields used internally (Ez, Wr, rho, CSR wakes...) at each step, as defined by n_out.
  • PlasmaStage.
  • PlasmaRamp.
  • PlasmaLens.
  • TMElement.
• Possibility of specifying name and path of diagnostics directory.
• Write diagnostics of whole Beamline to the same folder (this might include plasma elements, drifts, magnets, etc. which require different attributes and datasets in the output).
• In the plasma elements, output individual components of E and B and not just the combined wakefield W.

Nice to have:

• Option to choose field or particle diagnostics only.
• Or even better, option to specify individual fields and particle quantities to write.
• Option to output plasma particles in quasistatic_2d model.
• Option of writing diagnostics of each element to a separate directory.

Bugs/issues to solve:

• A name parameter should be added to ParticleBunch. Right now all bunches are called electron_bunch in the diagnostics.
• Check that the position field attribute it correctly set (cell centered vs node centered).
• A cleaner way of setting the field dimensions and other parameter could be implemented.
• ED_PIC metadata is currently missing in particle species. Check if/what to add because Wake-T is not a PIC code.
• Same for fields, it is not clear what to write in some ED_PIC metadata.
• No output is currently generated at t=0.
• The attribute iterationFormat is not correct.
• Some fields are transposed in the output and do not agree with the axis order specified in attributes such as axisLabels.
• gridGlobalOffset is not correct. It should add (and not substract) xi_min. Also, it is incorrectly named 'local_offset' in some parts of the code.
• Set iteration encoding to file based.

This would allow using the openpmd-viewer ParticleTracker.

The time step for the beam particle pusher is currently automatically chosen to be 1/10 of the betatron period in the plasma stage. An option should be added to allow the user to change this.

Current code is not compliant with the PEP8 standards:

• Format current code according to PEP8.
• Add automatic checks with Travis and flake8 for new commits.

Many methods can benefit from parallelization, in particular the beam tracking, but are currently only serial. Given that most are compiled with numba, adding parallelism should be relatively straightforward.

Hello, you might consider adding Wake-T to the community-driven list of Accelerator Simulation Codes available at https://github.com/campa-consortium/accelerator_simulation_codes <3

Sometimes they appear as cm^{-3} and others are used as m^{-3}.

When ParticleBunch is imported, it will fail if matplotlib is not installed.
I think the culprit is this line in particle_bunch.py :
from aptools.plotting.quick_diagnostics import full_phase_space
Since there is a import matplotlib instance in aptools.plotting.quick_diagnostics.