Simple LUME-Genesis example for Genesis 1.3 Version 2 (Genesis2)¶
In [ ]:
Copied!
# Useful for debugging
%load_ext autoreload
%autoreload 2
# Useful for debugging
%load_ext autoreload
%autoreload 2
In [ ]:
Copied!
from genesis import Genesis2
# Make genesis object with some input file
G = Genesis2("data/basic/genesis.in", verbose=True)
from genesis import Genesis2
# Make genesis object with some input file
G = Genesis2("data/basic/genesis.in", verbose=True)
In [ ]:
Copied!
# Turn on field output
G["idmpfld"] = 1
# Turn on particle output
G["idmppar"] = 1
G["npart"] = 2048
# Turn on history
G["ippart"] = 10
G["ipradi"] = 0
# Turn on field output
G["idmpfld"] = 1
# Turn on particle output
G["idmppar"] = 1
G["npart"] = 2048
# Turn on history
G["ippart"] = 10
G["ipradi"] = 0
In [ ]:
Copied!
# Run genesis with default lattice
G.run()
# Run genesis with default lattice
G.run()
In [ ]:
Copied!
G.output["data"].keys()
G.output["data"].keys()
In [ ]:
Copied!
G.write_wavefront()
G.write_wavefront()
In [ ]:
Copied!
G.archive()
G.archive()
In [ ]:
Copied!
# This is the working path
G.path
# This is the working path
G.path
In [ ]:
Copied!
# Contents of this
!ls -ahl {G.path}
# Contents of this
!ls -ahl {G.path}
In [ ]:
Copied!
# Get a hash of the input
G.fingerprint()
# Get a hash of the input
G.fingerprint()
In [ ]:
Copied!
# Lattice parameters
G.lattice["param"]
# Lattice parameters
G.lattice["param"]
In [ ]:
Copied!
# Get a list of z from the output
G.output.keys()
# Get a list of z from the output
G.output.keys()
In [ ]:
Copied!
# These are the available data
G.output["data"].keys()
# These are the available data
G.output["data"].keys()
In [ ]:
Copied!
zlist = G.output["data"]["z"]
zlist.shape
zlist = G.output["data"]["z"]
zlist.shape
In [ ]:
Copied!
# Get power. This is a 2d array of: slice, z
power = G.output["data"]["power"]
power.shape
# Get power. This is a 2d array of: slice, z
power = G.output["data"]["power"]
power.shape
In [ ]:
Copied!
# Simpler plot
import matplotlib.pyplot as plt
%config InlineBackend.figure_format = 'retina'
plt.xlabel("Z (m)")
plt.ylabel("Power (GW)")
plt.plot(zlist, power[0] / 1e9)
# Simpler plot
import matplotlib.pyplot as plt
%config InlineBackend.figure_format = 'retina'
plt.xlabel("Z (m)")
plt.ylabel("Power (GW)")
plt.plot(zlist, power[0] / 1e9)
Archiving¶
In [ ]:
Copied!
afile = G.archive()
afile = G.archive()
In [ ]:
Copied!
G2 = Genesis2()
G2.load_archive(afile)
G2 = Genesis2()
G2.load_archive(afile)
In [ ]:
Copied!
# Check that fingetprints (hash on input) are the same
G.fingerprint() == G2.fingerprint()
# Check that fingetprints (hash on input) are the same
G.fingerprint() == G2.fingerprint()
In [ ]:
Copied!
import numpy as np
# Check that all output data are the same
for k in G.output["data"]:
print(k, np.all(G.output["data"][k] == G2.output["data"][k]))
import numpy as np
# Check that all output data are the same
for k in G.output["data"]:
print(k, np.all(G.output["data"][k] == G2.output["data"][k]))
Wavefront (dfl) in openPMD-wavefront format¶
This will write the loaded dfl
data to a proper openPMD-wavefront file
In [ ]:
Copied!
wfile = G.write_wavefront()
wfile = G.write_wavefront()
In [ ]:
Copied!
# Read back
import h5py
with h5py.File(wfile, "r") as h5:
component = h5["data/000000/meshes/electricField/x"]
E = component[:]
attrs = dict(component.attrs)
# This is the conversion factor to the proper SI unit for electric field
attrs
# Read back
import h5py
with h5py.File(wfile, "r") as h5:
component = h5["data/000000/meshes/electricField/x"]
E = component[:]
attrs = dict(component.attrs)
# This is the conversion factor to the proper SI unit for electric field
attrs
In [ ]:
Copied!
# Check the shapes
E.shape
# Check the shapes
E.shape
In [ ]:
Copied!
# Plot the phase
%matplotlib inline
import matplotlib.pyplot as plt
plt.rcParams["figure.figsize"] = [8, 8]
ndat = np.angle(E[:, :, 0]) # Complex angle
plt.imshow(ndat, extent=[1e3 * G["dgrid"] * i for i in [-1, 1, -1, 1]])
plt.xlabel("x (mm)")
plt.ylabel("y (mm)")
plt.show()
# Plot the phase
%matplotlib inline
import matplotlib.pyplot as plt
plt.rcParams["figure.figsize"] = [8, 8]
ndat = np.angle(E[:, :, 0]) # Complex angle
plt.imshow(ndat, extent=[1e3 * G["dgrid"] * i for i in [-1, 1, -1, 1]])
plt.xlabel("x (mm)")
plt.ylabel("y (mm)")
plt.show()
Raw particle history¶
In [ ]:
Copied!
import numpy as np
import os
import numpy as np
import os
In [ ]:
Copied!
# gamma, phase, x, y, px/mc, py/mc
b = G.output["data"]["par"]
# gamma, phase, x, y, px/mc, py/mc
b = G.output["data"]["par"]
In [ ]:
Copied!
p = b[-1]
plt.hist2d(
p[1] % (2 * np.pi), 0.511 * p[0], bins=[100, 100], cmap=plt.get_cmap("plasma")
)
plt.xlabel("phase")
plt.ylabel("energy (MeV)")
# plt.savefig('frame.png')
plt.show()
p = b[-1]
plt.hist2d(
p[1] % (2 * np.pi), 0.511 * p[0], bins=[100, 100], cmap=plt.get_cmap("plasma")
)
plt.xlabel("phase")
plt.ylabel("energy (MeV)")
# plt.savefig('frame.png')
plt.show()
In [ ]:
Copied!
def frame(
i,
):
p = b[i]
plt.hist2d(
p[1] % (2 * np.pi), 0.511 * p[0], bins=[200, 200], cmap=plt.get_cmap("plasma")
)
plt.xlabel("phase")
plt.ylabel("energy (MeV)")
# plt.savefig('frame_'+str(i)+'.png')
frame(100)
def frame(
i,
):
p = b[i]
plt.hist2d(
p[1] % (2 * np.pi), 0.511 * p[0], bins=[200, 200], cmap=plt.get_cmap("plasma")
)
plt.xlabel("phase")
plt.ylabel("energy (MeV)")
# plt.savefig('frame_'+str(i)+'.png')
frame(100)
In [ ]:
Copied!
# Make frames
# for i in range(nbunch):
# frame(i);
# Make frames
# for i in range(nbunch):
# frame(i);
In [ ]:
Copied!
# Make movie
#!ffmpeg -framerate 10 -i frame_%d.png -c:v libx264 -c:a libfdk_aac output.mp4
# Make movie
#!ffmpeg -framerate 10 -i frame_%d.png -c:v libx264 -c:a libfdk_aac output.mp4
Final particles¶
If a final particle (.dpa file) data was made, this will process the internal arrays into a high-level ParticleGroup object.
In [ ]:
Copied!
P = G.final_particles()
P
P = G.final_particles()
P
In [ ]:
Copied!
P.plot("t", "energy")
P.plot("t", "energy")
Cleanup¶
In [ ]:
Copied!
os.remove(afile)
os.remove(wfile)
os.remove(afile)
os.remove(wfile)