Running a Simulation#

Now that we have tuned the simulation’s parameter it is time to run it. Because this is another notebook we copy the first cell of the Pre Simulation Testing notebook and then run the simulation with only three lines of code.

The YAML input file can be found at input_file and this notebook at notebook

[1]:
# Import the usual libraries
import os

# Import sarkas
from sarkas.processes import Simulation

# Create the file path to the YAML input file
input_file_name = os.path.join('input_files', 'yukawa_mks_p3m.yaml')

Then we run our simulation with only three lines of code

[2]:
sim = Simulation(input_file_name)
sim.setup(read_yaml=True)
sim.run()






 __            _
/ _\ __ _ _ __| | ____ _ ___
\ \ / _` | '__| |/ / _` / __|
_\ \ (_| | |  |   < (_| \__ \
\__/\__,_|_|  |_|\_\__,_|___/



An open-source pure-python molecular dynamics suite for non-ideal plasmas.





* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
                                   Simulation
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *

Job ID:  yocp
Job directory:  Simulations/yocp_pppm

Equilibration dumps directory:
 Simulations/yocp_pppm/Simulation/Equilibration/dumps
Production dumps directory:
 Simulations/yocp_pppm/Simulation/Production/dumps

Equilibration Thermodynamics file:
 Simulations/yocp_pppm/Simulation/Equilibration/EquilibrationEnergy_yocp.csv
Production Thermodynamics file:
 Simulations/yocp_pppm/Simulation/Production/ProductionEnergy_yocp.csv

PARTICLES:
Total No. of particles =  10000
No. of species =  1
Species ID: 0
        Name: H
        No. of particles = 10000
        Number density = 1.620000e+32 [N/m^3]
        Atomic weight = 1.0002 [a.u.]
        Mass = 1.673000e-27 [kg]
        Mass density = 2.710260e+05 [kg/m^3]
        Charge number/ionization degree = 1.0000
        Charge = 1.602177e-19 [C]
        Temperature = 1.450565e+04 [K] = 1.250000e+00 [eV]
        Debye Length = 6.530052e-13 [1/m]
        Plasma Frequency = 1.675504e+16 [rad/s]

SIMULATION AND INITIAL PARTICLE BOX:
Units:  mks
Wigner-Seitz radius = 1.137973e-11 [m]
No. of non-zero box dimensions =  3
Box side along x axis = 3.472931e+01 a_ws = 3.952104e-10 [m]
Box side along y axis = 3.472931e+01 a_ws = 3.952104e-10 [m]
Box side along z axis = 3.472931e+01 a_ws = 3.952104e-10 [m]
Box Volume = 6.172840e-29 [m^3]
Initial particle box side along x axis = 3.472931e+01 a_ws = 3.952104e-10 [m]
Initial particle box side along y axis = 3.472931e+01 a_ws = 3.952104e-10 [m]
Initial particle box side along z axis = 3.472931e+01 a_ws = 3.952104e-10 [m]
Initial particle box Volume = 6.172840e-29 [m^3]
Boundary conditions: periodic

ELECTRON PROPERTIES:
Number density: n_e = 1.620000e+32 [N/m^3]
Wigner-Seitz radius: a_e = 1.137973e-11 [m]
Temperature: T_e = 1.450565e+07 [K] = 1.250000e+03 [eV]
de Broglie wavelength: lambda_deB = 1.957093e-11 [m]
Thomas-Fermi length: lambda_TF = 2.272532e-11 [m]
Fermi wave number: k_F = 1.686470e+11 [1/m]
Fermi Energy: E_F = 1.083628e+03 [eV]
Relativistic parameter: x_F = 6.512461e-02 --> E_F = 1.082482e+03 [eV]
Degeneracy parameter: Theta = 1.153532e+00
Coupling: r_s = 0.215046,  Gamma_e = 0.101230
Warm Dense Matter Parameter: W = 1.9838e-01
Chemical potential: mu = -2.8605e-01 k_B T_e = -3.2996e-01 E_F

POTENTIAL:  yukawa
electron temperature = 1.4506e+07 [K] = 1.2500e+03 [eV]
kappa = 0.5008
Gamma_eff = 101.23

ALGORITHM:  pppm
Charge assignment order: 6
FFT aliases: [3, 3, 3]
Mesh: 64 x 64 x 64
Ewald parameter alpha = 0.6220 / a_ws = 5.465900e+10 [1/m]
Mesh width = 0.5426, 0.5426, 0.5426 a_ws
           = 6.1752e-12, 6.1752e-12, 6.1752e-12 [m]
Mesh size * Ewald_parameter (h * alpha) = 0.3375, 0.3375, 0.3375
                                        ~ 1/2, 1/2, 1/2
rcut = 5.5100 a_ws = 6.270200e-11 [m]
No. of PP cells per dimension =  6,  6,  6
No. of particles in PP loop =   1078
No. of PP neighbors per particle =    167
PM Force Error = 5.800215e-06
PP Force Error = 2.804557e-06
Tot Force Error = 6.442673e-06

THERMOSTAT:
Type: berendsen
First thermostating timestep, i.e. relaxation_timestep = 50
Berendsen parameter tau: 1.000 [timesteps]
Berendsen relaxation rate: 1.000 [1/timesteps]
Thermostating temperatures:
Species ID 0: T_eq = 1.450565e+04 [K] = 1.250000e+00 [eV]

INTEGRATOR:
Type: verlet
Time step = 2.000000e-18 [s]
Total plasma frequency = 1.675504e+16 [rad/s]
w_p dt = 0.0335 ~ 1/29

Equilibration:
No. of equilibration steps = 5000
Total equilibration time = 1.0000e-14 [s] ~ 167 w_p T_eq
snapshot interval step = 10
snapshot interval time = 2.0000e-17 [s] = 0.3351 w_p T_snap
Total number of snapshots = 500

Production:
No. of production steps = 5000
Total production time = 1.0000e-14 [s] ~ 167 w_p T_prod
snapshot interval step = 10
snapshot interval time = 2.0000e-17 [s] = 0.3351 w_p T_snap
Total number of snapshots = 500


------------------------ Initialization Times ------------------------


Potential Initialization Time: 0 hrs 0 min 6 sec

Particles Initialization Time: 0 sec 1 msec 271 usec 628 nsec

Total Simulation Initialization Time: 0 hrs 0 min 6 sec


--------------------------- Equilibration ----------------------------


Equilibration Time: 0 hrs 32 min 45 sec


----------------------------- Production -----------------------------


Production Time: 0 hrs 32 min 16 sec

Total Time: 1 hrs 5 min 4 sec

The verbose output of the simulation gives the same information as PreProcess. The first few lines tell us where the simulations’ dumps and the thermodynamics files are be stored. Notice that these are diffrent than what was printed by PreProcessing.

NOTE: This simulation was run on a 2019 Dell XPS 8930 with Intel Core i7-8700K @ 3.70Ghz and 48GB of RAM running Ubuntu 18.04. Actual run times will be different on your computer.

Dumps#

Each dump is a .npz file containing the following information for each of phase:

Equilibration:
 id = particles id's numbers,
 names = species name of each particle,
 pos = Particles positions N x 3 array,
 vel = Particles velocities N x 3 array,
 acc = Particles acceleration N x 3 array,
 virial = Virial matrix N x 3 x 3 array,
 time = Time in seconds

Production:
 id = particles id's numbers,
 names = species name of each particle,
 pos = Particles positions N x 3 array,
 vel = Particles velocities N x 3 array,
 acc = Particles acceleration N x 3 array,
 cntr = Number of times each particle has been folded back into the simulation box, N x 3 array,
 rdf_hist = Instantaneous Histogram of the Radial Distribution Function,
 virial = Virial matrix N x 3 x 3 array,
 time = Time in seconds

This information is saved by the .dump() method of the InputOutput.

The same method saves Energy and Temperature information in the Thermodynamics files also. This can be accessed by sim.io.eq_energy_filename and sim.io.prod_energy_filename. The information saved is the following:

Time, Total Energy, Total Kinetic Energy, Potential Energy, Total Temperature

In addition, in the case of multispecies plasmas, the kinetic energy and the temperature of each species will be saved.

Once the simulation is complete it is time to analyze the data. Follow this link.

[ ]: