Following the introduction on the first day, the second day is focused on explaining how to up-scale existing simulation workflows. The tutorial is based on the molecular dynamics simulation code LAMMPS. Still the same protocols can also be executed with DFT codes like VASP and S/PHI/nX.
Following the tutorial from the first day, the energy volume curve is calculated using a for-loop
over 11 different strains. Start by importing numpy
as well as the pyiron Project
class and create a new project named scaleup
.
import numpy as np
from pyiron import Project
pr = Project("scaleup")
As a first step the lattice constant is optimized for the specific lattice constant of the selected interatomic potential using LAMMPS. By using auto-completion (tab-completion) only the job_name
and the element in this case iron Fe
have to be selected manually. A cubic cell is created by adding the option cubic=True
when creating the bulk structure:
job_mini = pr.create.job.Lammps(job_name="job_mini")
job_mini.structure = pr.create.structure.ase.bulk("Fe", cubic=True)
potential = job_mini.list_potentials()[0]
potential
'1997--Ackland-G-J--Fe--LAMMPS--ipr1'
Optimise the lattice constant of the cubic simulation cell using the calc_minimize()
function with the additional parameter pressue=0.0
. The potential is assinged to the job object and afterwards the calculation is executed using run()
:
job_mini.calc_minimize(pressure=0)
job_mini.potential = potential
job_mini.run()
The job job_mini was saved and received the ID: 2
The lattice constant is extracted from the final structure of the minimization calculation using the get_structure()
function:
alat_guess = job_mini.get_structure().cell[0,0]
np.round(alat_guess, 3)
2.866
Based on the initial guess of the lattice constant 11 strains ranging from $90$% to $110$% are applied to calculate the energy volume curve. The calculation are named based on the applied strain and the required structures are again created using the pr.create.structure.ase.bulk()
function. All calculation are executed directly inline.
for strain in np.linspace(0.9, 1.1, 11):
job_name = "lmp_" + str(np.round(strain,2)).replace(".", "_")
job = pr.create.job.Lammps(job_name=job_name)
print(alat_guess * strain**(1/3))
job.structure = pr.create.structure.ase.bulk("Fe", a=alat_guess * strain**(1/3), )
job.potential = potential
job.run()
2.7675753208676444 The job lmp_0_9 was saved and received the ID: 3 2.7879258702511898 The job lmp_0_92 was saved and received the ID: 4 2.8079835805272397 The job lmp_0_94 was saved and received the ID: 5 2.8277587643206537 The job lmp_0_96 was saved and received the ID: 6 2.8472611651817745 The job lmp_0_98 was saved and received the ID: 7 2.866499999891875 The job lmp_1_0 was saved and received the ID: 8 2.885483996845223 The job lmp_1_02 was saved and received the ID: 9 2.90422143094099 The job lmp_1_04 was saved and received the ID: 10 2.9227201553630375 The job lmp_1_06 was saved and received the ID: 11 2.9409876305783165 The job lmp_1_08 was saved and received the ID: 12 2.9590309508440615 The job lmp_1_1 was saved and received the ID: 13
Overview of the existing calculation in the current project using job_table()
:
pr.job_table()
id | status | chemicalformula | job | subjob | projectpath | project | timestart | timestop | totalcputime | computer | hamilton | hamversion | parentid | masterid | |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
0 | 2 | finished | Fe2 | job_mini | /job_mini | /home/jovyan/ | scaleup/ | 2021-03-22 23:51:52.778780 | 2021-03-22 23:51:53.409059 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | None |
1 | 3 | finished | Fe | lmp_0_9 | /lmp_0_9 | /home/jovyan/ | scaleup/ | 2021-03-22 23:51:54.545651 | 2021-03-22 23:51:55.084443 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | None |
2 | 4 | finished | Fe | lmp_0_92 | /lmp_0_92 | /home/jovyan/ | scaleup/ | 2021-03-22 23:51:55.922362 | 2021-03-22 23:51:56.454064 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | None |
3 | 5 | finished | Fe | lmp_0_94 | /lmp_0_94 | /home/jovyan/ | scaleup/ | 2021-03-22 23:51:57.354493 | 2021-03-22 23:51:57.921399 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | None |
4 | 6 | finished | Fe | lmp_0_96 | /lmp_0_96 | /home/jovyan/ | scaleup/ | 2021-03-22 23:51:58.821766 | 2021-03-22 23:51:59.324438 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | None |
5 | 7 | finished | Fe | lmp_0_98 | /lmp_0_98 | /home/jovyan/ | scaleup/ | 2021-03-22 23:52:00.450136 | 2021-03-22 23:52:01.021797 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | None |
6 | 8 | finished | Fe | lmp_1_0 | /lmp_1_0 | /home/jovyan/ | scaleup/ | 2021-03-22 23:52:02.045621 | 2021-03-22 23:52:02.621343 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | None |
7 | 9 | finished | Fe | lmp_1_02 | /lmp_1_02 | /home/jovyan/ | scaleup/ | 2021-03-22 23:52:03.866085 | 2021-03-22 23:52:04.411873 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | None |
8 | 10 | finished | Fe | lmp_1_04 | /lmp_1_04 | /home/jovyan/ | scaleup/ | 2021-03-22 23:52:05.616177 | 2021-03-22 23:52:06.127510 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | None |
9 | 11 | finished | Fe | lmp_1_06 | /lmp_1_06 | /home/jovyan/ | scaleup/ | 2021-03-22 23:52:07.527880 | 2021-03-22 23:52:08.085600 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | None |
10 | 12 | finished | Fe | lmp_1_08 | /lmp_1_08 | /home/jovyan/ | scaleup/ | 2021-03-22 23:52:09.550049 | 2021-03-22 23:52:10.132919 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | None |
11 | 13 | finished | Fe | lmp_1_1 | /lmp_1_1 | /home/jovyan/ | scaleup/ | 2021-03-22 23:52:11.452091 | 2021-03-22 23:52:12.000244 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | None |
The simulation results are collected with a simple for loop by iterating over the jobs in the project using iterjobs()
and appending the final volume of the calculation and the total energy to a set of lists:
vol_lst, eng_lst = [], []
for job in pr.iter_jobs():
if "lmp_" in job.job_name:
vol_lst.append(job["output/generic/volume"][-1])
eng_lst.append(job["output/generic/energy_tot"][-1])
Finally the simulation results of the results are visualised using matplotlib
to plot the energy volume curve:
import matplotlib.pyplot as plt
plt.plot(vol_lst, eng_lst)
plt.xlabel("Volume")
plt.ylabel("Energy")
Text(0, 0.5, 'Energy')
Using for-loops
to create large datasets is limited. A more scalable approach is using already existing classes and combine those like building blocks. For the example of calculating an energy volume curve the Murnaghan
class implements the creation of different strained structures as well as the management of the related calculation. In particular these job classes or GenericMasters
support various modes of communication with the HPC queuing system. For example submitting individual calculation or executing all calculation in a single job either in serial or in parallel. Still for this presentation all jobs are executed inline.
As a start a template job is created which uses the same interatomic potential as the previous calculation and the same already optimised crystal structure:
job_fe = pr.create.job.Lammps(job_name="job_fe")
job_fe.structure = job_mini.get_structure()
job_fe.potential = potential
The same template job is used for multiple calculation in the following, therefore it is copied for the different applications using the copy_template()
function:
job_murn = job_fe.copy_template(
new_job_name="job_murn"
)
From the newly created copy a new job is created of type Murnaghan
. By creating the job from the previous job template the relation of the job and it is location on the filesystem is defined. The input of the Murnaghan
class can again be edited using the edge notation:
murn = job_murn.create_job(
job_type=pr.job_type.Murnaghan,
job_name="murn"
)
murn.input
Parameter | Value | Comment | |
---|---|---|---|
0 | num_points | 11 | number of sample points |
1 | fit_type | polynomial | ['polynomial', 'birch', 'birchmurnaghan', 'murnaghan', 'pouriertarantola', 'vinet'] |
2 | fit_order | 3 | order of the fit polynom |
3 | vol_range | 0.1 | relative volume variation around volume defined by ref_ham |
When the run method of the Murnaghan
class is called, the individual calculation are started and after the successful completion the aggregated information is stored in the Murnaghan
object:
murn.run()
The job murn was saved and received the ID: 14 The job strain_0_9 was saved and received the ID: 15 The job strain_0_92 was saved and received the ID: 16 The job strain_0_94 was saved and received the ID: 17 The job strain_0_96 was saved and received the ID: 18 The job strain_0_98 was saved and received the ID: 19 The job strain_1_0 was saved and received the ID: 20 The job strain_1_02 was saved and received the ID: 21 The job strain_1_04 was saved and received the ID: 22 The job strain_1_06 was saved and received the ID: 23 The job strain_1_08 was saved and received the ID: 24 The job strain_1_1 was saved and received the ID: 25 job_id: 15 finished job_id: 16 finished job_id: 17 finished job_id: 18 finished job_id: 19 finished job_id: 20 finished job_id: 21 finished job_id: 22 finished job_id: 23 finished job_id: 24 finished job_id: 25 finished
In addition to the data aggregation methods the Murnaghan
class also provides plotting methods for commonly used plots:
murn.plot()
Following the same principles also the Phonopy
interface is implemented in pyiron. This enables using Phonopy
in combination with all simulation codes implemented in pyiron. The required steps are more or less identical. Again starting by copying the job template using the copy_template()
function. Starting with a bulk calculation as a reference:
job_bulk = job_fe.copy_template(
new_job_name="job_bulk"
)
job_bulk.run()
The job job_bulk was saved and received the ID: 26
Following the bluk calculation the Phonopy
job is created by again copying the job template and creating a Phonopy
job from the job object:
job_phono = job_fe.copy_template(
new_job_name="job_phono"
)
Again the input object of the Phonopy
job provides access to the most commonly used input parameters:
phono = job_phono.create_job(
job_type=pr.job_type.PhonopyJob,
job_name="phono"
)
phono.input
Parameter | Value | Comment | |
---|---|---|---|
0 | interaction_range | 10.000000 | Minimal size of supercell, Ang |
1 | factor | 15.633302 | Frequency unit conversion factor (default for VASP) |
2 | displacement | 0.010000 | atoms displacement, Ang |
3 | dos_mesh | 20.000000 | mesh size for DOS calculation |
4 | primitive_matrix | NaN |
Once more the job is executed using the run()
function:
phono.run()
The job phono was saved and received the ID: 27 The job job_phono_0 was saved and received the ID: 28
And the density of states can be plotted directly from the Phonopy
job object:
phono.plot_dos()
<AxesSubplot:title={'center':'Phonon DOS vs Energy'}, xlabel='Frequency [THz]', ylabel='DOS'>
In addition to the calculation of the density of states Phonopy
also enables the calculation of the free energy with the harmonic and quasi-harmonic approximation. Starting with the harmonic approximation, the free energy is calculated for temperatures ranging from $0$K to $800$K:
temperatures = np.linspace(0, 800, 41)
bulk_thermal_properties = phono.get_thermal_properties(temperatures=temperatures)
The free energy over temperature is plotted by combining the vibrational and the bulk energy:
plt.plot(temperatures, job_bulk.output.energy_pot[-1] + bulk_thermal_properties.free_energies)
plt.xlabel("Temperature [K]")
plt.ylabel("Free energy ($U+F_{vib}$) [eV]");
The major advantage of the GenericMaster
classes is that they can not only be used to combine regular calculation, but also GenericMaster
classes themselves. Like building blocks the different pyiron classes can be combined and nested. For the quasi-harmonic approximation the volume dependent free energy is calculated and the optimal volume for a given temperature is identifed by fitting. Starting once more by copying the template job using the copy_template()
function:
job_strain = job_fe.copy_template(
new_job_name="job_strain"
)
Aftewards the same template job is used to create both a Murnaghan
job and a Phonopy
job:
murn_strain = job_strain.create_job(
job_type=pr.job_type.Murnaghan,
job_name="murn_strain"
)
phono_strain = job_strain.create_job(
job_type=pr.job_type.PhonopyJob,
job_name="phono_strain"
)
The Phonopy
job is then used to create a QuasiHarmonicJob
which calculates the free energy for each of the different volumes:
quasi_strain = phono_strain.create_job(
job_type=pr.job_type.QuasiHarmonicJob,
job_name="quasi_strain"
)
The input of the QuasiHarmonicJob
is adjusted to range the temperature up to $800$K using $41$ steps:
quasi_strain.input["temperature_end"] = 800
quasi_strain.input["temperature_steps"] = 41
After all calculation are setup they can be executed by again calling the run()
commands on both the Murnaghan
job and the QuasiHarmonicJob
job:
murn_strain.run()
quasi_strain.run()
The job murn_strain was saved and received the ID: 29 The job strain_0_9 was saved and received the ID: 30 The job strain_0_92 was saved and received the ID: 31 The job strain_0_94 was saved and received the ID: 32 The job strain_0_96 was saved and received the ID: 33 The job strain_0_98 was saved and received the ID: 34 The job strain_1_0 was saved and received the ID: 35 The job strain_1_02 was saved and received the ID: 36 The job strain_1_04 was saved and received the ID: 37 The job strain_1_06 was saved and received the ID: 38 The job strain_1_08 was saved and received the ID: 39 The job strain_1_1 was saved and received the ID: 40 job_id: 30 finished job_id: 31 finished job_id: 32 finished job_id: 33 finished job_id: 34 finished job_id: 35 finished job_id: 36 finished job_id: 37 finished job_id: 38 finished job_id: 39 finished job_id: 40 finished The job quasi_strain was saved and received the ID: 41 The job strain_0_9 was saved and received the ID: 42 The job job_strain_0 was saved and received the ID: 43 The job strain_0_92 was saved and received the ID: 44 The job job_strain_0 was saved and received the ID: 45 The job strain_0_94 was saved and received the ID: 46 The job job_strain_0 was saved and received the ID: 47 The job strain_0_96 was saved and received the ID: 48 The job job_strain_0 was saved and received the ID: 49 The job strain_0_98 was saved and received the ID: 50 The job job_strain_0 was saved and received the ID: 51 The job strain_1_0 was saved and received the ID: 52 The job job_strain_0 was saved and received the ID: 53 The job strain_1_02 was saved and received the ID: 54 The job job_strain_0 was saved and received the ID: 55 The job strain_1_04 was saved and received the ID: 56 The job job_strain_0 was saved and received the ID: 57 The job strain_1_06 was saved and received the ID: 58 The job job_strain_0 was saved and received the ID: 59 The job strain_1_08 was saved and received the ID: 60 The job job_strain_0 was saved and received the ID: 61 The job strain_1_1 was saved and received the ID: 62 The job job_strain_0 was saved and received the ID: 63
After the successful completion of all calculation the results can be visualised using matplotlib
starting by creating a color map:
import matplotlib
cmap = matplotlib.cm.get_cmap('coolwarm')
The free energy is plotted over volume and temperature, starting at a low temperature of $0$K up to the $800$K:
# Iterate over the the output of the QuasiHarmonicJob to plot the temperature dependent free energy over volume
for i, [t, free_energy, v] in enumerate(
zip(quasi_strain["output/temperatures"].T,
quasi_strain["output/free_energy"].T,
quasi_strain["output/volumes"].T)):
color = cmap(i/len(quasi_strain["output/temperatures"].T))
# Add the energy of the Murnaghan Job to the vibrational free energy
plt.plot(v, free_energy + murn_strain["output/energy"], color=color)
# Add labels to the plot
plt.xlabel("Volume")
plt.ylabel("Free Energy")
# Add a color bar to visualise the temperature dependence
temperatures = quasi_strain["output/temperatures"]
normalize = matplotlib.colors.Normalize(vmin=temperatures.min(), vmax=temperatures.max())
scalarmappaple = matplotlib.cm.ScalarMappable(norm=normalize, cmap=cmap)
scalarmappaple.set_array(temperatures)
cbar = plt.colorbar(scalarmappaple)
cbar.set_label("Temperature")
In addition to the free energy over temperature and volume, the optimised volume over temperature allows to quantify the volume expansion:
# Use the optimise_volume() function of the QuasiHarmonicJob
v0_lst, free_eng_lst, entropy_lst, cv_lst = quasi_strain.optimise_volume(
# It requires the output energy of the energy volume curve as additional input
bulk_eng=murn_strain["output/energy"]
)
temperature_output_lst = quasi_strain["output/temperatures"][0]
plt.plot(temperature_output_lst, v0_lst)
plt.xlabel("Temperature")
plt.ylabel("Volume")
Text(0, 0.5, 'Volume')
Finally based on the volume expansion over temperature the quasi-harmonic free energy is calculated:
plt.plot(temperature_output_lst, free_eng_lst)
plt.xlabel("Temperature")
plt.ylabel("Free Energy")
Text(0, 0.5, 'Free Energy')
This example demonstrates how pyiron can be used to calculate free energies.
From a more technical perspective the GenericMaster
classes are only one kind of tool integrated in pyiron to up-scale existing simulation protocols. While they provide a lot of flexibility like nesting multiple levels of GenericMaster
objects, the development of GenericMaster
classes requires a general understanding of pyiron. A more flexible way to up-scale an existing workflow is using Scriptjob
objects which can execute existing python scripts and also existing jupyter notebooks. So once a new simulation workflow is implemented in a jupyter notebook this jupyter notebook can be used as an input for a Scriptjob
to loop over a control parameter like looping over existing interatomic potentials. In this example we create a simple Scriptjob
which just contains a python script:
python_script = """\
from pyiron import Project
pr = Project("script")
job = pr.create.job.Lammps("lmp")
job.structure = pr.create.structure.ase.bulk("Fe")
job.run()
"""
The python script which creates a single LAMMPS calculation is written to a python file named script.py
.
with open("script.py", "w") as f:
f.writelines(python_script)
Afterwards a ScriptJob
object is created from the project object and the python file script.py
is attached as script_path
before the job is executed using the run()
function. Again the ScriptJob
object behaves like any other job object in pyiron:
job_script = pr.create.job.ScriptJob("script")
job_script.script_path = "script.py"
job_script.run()
The job script was saved and received the ID: 64
By looking at the job table it is visible that not only the ScriptJob
itself but also the job executed from the python script named lmp
is listed in the project:
pr.job_table()
id | status | chemicalformula | job | subjob | projectpath | project | timestart | timestop | totalcputime | computer | hamilton | hamversion | parentid | masterid | |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
0 | 2 | finished | Fe2 | job_mini | /job_mini | /home/jovyan/ | scaleup/ | 2021-03-22 23:51:52.778780 | 2021-03-22 23:51:53.409059 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | NaN |
1 | 3 | finished | Fe | lmp_0_9 | /lmp_0_9 | /home/jovyan/ | scaleup/ | 2021-03-22 23:51:54.545651 | 2021-03-22 23:51:55.084443 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | NaN |
2 | 4 | finished | Fe | lmp_0_92 | /lmp_0_92 | /home/jovyan/ | scaleup/ | 2021-03-22 23:51:55.922362 | 2021-03-22 23:51:56.454064 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | NaN |
3 | 5 | finished | Fe | lmp_0_94 | /lmp_0_94 | /home/jovyan/ | scaleup/ | 2021-03-22 23:51:57.354493 | 2021-03-22 23:51:57.921399 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | NaN |
4 | 6 | finished | Fe | lmp_0_96 | /lmp_0_96 | /home/jovyan/ | scaleup/ | 2021-03-22 23:51:58.821766 | 2021-03-22 23:51:59.324438 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | NaN |
... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
59 | 61 | finished | Fe128 | job_strain_0 | /job_strain_0 | /home/jovyan/ | scaleup/quasi_strain_hdf5/strain_1_08_hdf5/ | 2021-03-22 23:54:24.027821 | 2021-03-22 23:54:24.729845 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | 60.0 |
60 | 62 | finished | Fe2 | strain_1_1 | /strain_1_1 | /home/jovyan/ | scaleup/quasi_strain_hdf5/ | 2021-03-22 23:54:28.338953 | 2021-03-22 23:54:34.997655 | 6.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1#1/1 | PhonopyJob | 0.0.1 | None | 41.0 |
61 | 63 | finished | Fe128 | job_strain_0 | /job_strain_0 | /home/jovyan/ | scaleup/quasi_strain_hdf5/strain_1_1_hdf5/ | 2021-03-22 23:54:31.821842 | 2021-03-22 23:54:32.536800 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | 62.0 |
62 | 64 | finished | None | script | /script | /home/jovyan/ | scaleup/ | 2021-03-22 23:54:50.762741 | 2021-03-22 23:54:55.087347 | 4.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Script | 0.1 | None | NaN |
63 | 65 | finished | Fe | lmp | /lmp | /home/jovyan/ | scaleup/script_hdf5/script/script/ | 2021-03-22 23:54:53.820840 | 2021-03-22 23:54:54.309320 | 0.0 | pyiron@jupyter-pyiron-2dswimm-2dworkshop-2d2021-2dw6x4leku#1 | Lammps | 0.1 | None | NaN |
64 rows × 15 columns
Therefore while offering less flexibility the ScriptJob
class is a simple way to create reusable snippets if certain parts of a simulation protocol are commonly repeated in different projects.
The third way to accelerate existing simulation protocols in pyiron is by using the so-called interactive mode. This interactive mode is enabled by simply setting the flag job.server.run_mode.interactive = True
and it enables communication with serial and parallel simulation codes during run time by leveraging various communication interfaces. For example:
INTERACTIVE
tag set in the INCAR
file which is unfortunatley not officially documented.pyiron implements a generic interface for these different communication standards, compareable to ipi-code and again enables the user to interact with the codes interactively by keeping the simulation code running in the background. In combination with the Murnaghan
job a single LAMMPS
calculation can be used as quantum engine to evaluate 11 different strains.
job_murn = job_fe.copy_template(
new_job_name="job_murn_int"
)
job_murn.server.run_mode.interactive = True
murn = job_murn.create_job(
job_type=pr.job_type.Murnaghan,
job_name="murn_int"
)
murn.run()
murn.plot()
The job murn_int was saved and received the ID: 66 The job murn_int_job_murn_int was saved and received the ID: 67
After the general introduction to pyiron on the first day the second session introduced multiple methods to up-scale existing simulation protocols:
GenericMaster
classes which can be combined like building blocks.ScriptJob
class which enables the easy reuse of existing simulation protocols.Additional exercises:
ScriptJob
for a Murnaghan calculation, by using the script_job.input
to set the input for the ScriptJob
and then inside the ScriptJob
use pr.get_external_input()
to retrieve the input.job.structure
object and calling job.run()
after each change to evaluate the structure in interactive mode. By setting job.interactive_enforce_structure_reset = True
the structure is always reset for each step, while the flag is typically set to false for applications like calculating molecular dynamics trajectories or thermodynamic integration.