SmileiPIC / Smilei

Particle-in-cell code for plasma simulation
https://smileipic.github.io/Smilei
344 stars 120 forks source link

Accept input from 'LASY' laser python module #760

Open Nkehoe-QUB opened 2 weeks ago

Nkehoe-QUB commented 2 weeks ago

The problem and the context

Accept laser initialisation from the LASY Python module.

https://lasydoc.readthedocs.io/en/stable/index.html

The solution

Failed attempts

Trying to implement into existing laser blocks isn't straight forward or, possibly, doable.

beck-llr commented 2 weeks ago

Hi. It is actually doable with python calls in the namelist (I have already done it) but it is not automated and have to be adjusted to each case. A full automated support would be great but unfortunately none of us has the time to work on that in the short term at least.

Here are some snippet that I used that might help. Recover the laser envelope and mesh properties:

#Read values from propgation and decomposition
f = h5py.File('lasy_laser_He_05_000002.h5', 'r')
dset = f['data']['0']['meshes']['laserEnvelope']
data = dset[0,:,:] 
#Time and space steps
datadt,datadr = f['data']['0']['meshes']['laserEnvelope'].attrs.get("gridSpacing")
#Number of points in time and r
datant,datanr = data.shape

#According to the polarization, the data stores Ex in the openPMD standard. That is Ey in the Smilei standard.
#The corresponding B component is Bz = Ey/c
dataB   = data/cst.c # B envelope

Reconstruction of the field:

### Data interpolator ###
#Time and radius axis of the data
datat = np.linspace(0, datant*datadt, datant)
datar = np.linspace(0, datanr*datadr, datanr)
#Complex envelope interpolator
B_interp = RegularGridInterpolator((datat, datar), dataB)
#Bz component reconstruction from the complex envelope: Bz = Re( B exp(-iw_0t) )
nprect = np.vectorize(cmath.rect)
def Bz(t,r):
    return (B_interp((t,r))*nprect(1, -t)).real.astype(complex)

def br_mode0(r,t):
    return 0.j
def br_mode1(r,t):
    return 1j * Bz(t,r)
def bt_mode0(r,t):
    return 0.j
def bt_mode1(r,t):
    return Bz(t,r)

And the laser block becomes:

Laser(
    box_side         = "xmin",
    space_time_profile_AM = [br_mode0, bt_mode0, br_mode1, bt_mode1],
)

Hope this helps. It is not directly using the Lasy python module but it is a lazy workaround :-) Just be very careful about the normalizations.