Quickstart

Once installation is complete you can start running Sarkas. This quickstart guide will walk you through a simple example in order to check that everything is running smoothly.

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

---

Simulation

In Jupyter notebook you can run the following commands

# Import the usual libraries
%pylab
%matplotlib inline
import os

plt.style.use('MSUstyle')

# Import sarkas
from sarkas.processes import Simulation, PostProcess, PreProcess

Using matplotlib backend: <object object at 0x7fb7385b06f0>
Populating the interactive namespace from numpy and matplotlib

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

The above commands imported the required libraries and define the file path to our input file.

Let’s now run the simulation

# pre = PreProcess(input_file_name)
# pre.setup(read_yaml=True)
# pre.run()

# Initialize the Simulation class
sim = Simulation(input_file_name)
# Setup the simulation's parameters
sim.setup(read_yaml=True)
# Run the simulation
sim.run()

 _______  _______  _______  _        _______  _______
(  ____ \(  ___  )(  ____ )| \    /\(  ___  )(  ____ \
| (    \/| (   ) || (    )||  \  / /| (   ) || (    \/
| (_____ | (___) || (____)||  (_/ / | (___) || (_____
(_____  )|  ___  ||     __)|   _ (  |  ___  |(_____  )
      ) || (   ) || (\ (   |  ( \ \ | (   ) |      ) |
/\____) || )   ( || ) \ \__|  /  \ \| )   ( |/\____) |
\_______)|/     \||/   \__/|_/    \/|/     \|\_______)



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




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

Job ID: yocp_quickstart
Job directory: SarkasSimulations/yocp_quickstart
Simulation directory:
SarkasSimulations/yocp_quickstart/Simulation

Equilibration dumps directory:
SarkasSimulations/yocp_quickstart/Simulation/Equilibration/dumps
Production dumps directory:
SarkasSimulations/yocp_quickstart/Simulation/Production/dumps

Equilibration Thermodynamics file:
SarkasSimulations/yocp_quickstart/Simulation/Equilibration/EquilibrationEnergy_yocp_quickstart.csv
Production Thermodynamics file:
SarkasSimulations/yocp_quickstart/Simulation/Production/ProductionEnergy_yocp_quickstart.csv

PARTICLES:
Total No. of particles = 1000
No. of species = 1
Species ID: 0
        Name: C
        No. of particles = 1000
        Number density = 1.136669e+23 [N/cc]
        Atomic weight = 1.201048e+01 [a.u.]
        Mass = 2.008900e-23 [g]
        Mass density = 2.283455e+00 [g/cc]
        Charge number/ionization degree = 1.9766
        Charge = 9.494010e-10 [esu]
        Temperature = 5.000000e+03 [K] = 4.308667e-01 [eV]
        Debye length =  7.322426e-10 [cm]
        Plasma frequency =  2.531585e+14 [rad/s]

ELECTRON BACKGROUND PROPERTIES:
Number density: n_e = 2.246739e+23 [N/cc]
Wigner-Seitz radius: a_e = 1.020437e-08 [cm]
Temperature: T_e = 5.000000e+03 [K] = 4.308667e-01 [eV]
de Broglie wavelength: lambda_deB = 1.054132e-07 [cm]
Thomas-Fermi length: lambda_TF = 4.702909e-09[cm]
Fermi wave number: k_F = 1.880722e+08 [1/cm]
Fermi Energy: E_F = 1.347634e+01 [eV]
Relativistic parameter: x_F = 7.262582e-03 --> E_F = 1.347617e+01 [eV]
Degeneracy parameter: Theta = 3.197207e-02
Coupling: r_s = 1.928347,  Gamma_e = 1.046578
Warm Dense Matter parameter: W = 6.3813e-02
Chemical potential: mu = 3.1251e+01 k_B T_e = 9.9916e-01 E_F


SIMULATION AND INITIAL PARTICLE BOX:
Units: cgs
No. of non-zero box dimensions = 3
Wigner-Seitz radius = 1.280636e-08 [cm]
Box side along x axis = 1.611992e+01 a_ws = 2.064375e-07 [cm]
Box side along y axis = 1.611992e+01 a_ws = 2.064375e-07 [cm]
Box side along z axis = 1.611992e+01 a_ws = 2.064375e-07 [cm]
Box Volume = 8.797633e-21 [cm^3]
Initial particle box side along x axis = 1.611992e+01 a_ws = 2.064375e-07 [cm]
Initial particle box side along y axis = 1.611992e+01 a_ws = 2.064375e-07 [cm]
Initial particle box side along z axis = 1.611992e+01 a_ws = 2.064375e-07 [cm]
Initial particle box Volume = 8.797633e-21 [cm^3]
Boundary conditions: periodic
Random Seed = 314159265

SIMULATION PHASES:
Equilibration:
        No. of equilibration steps = 10000
        Total equilibration time = 5.0000e-13 [s] ~ 20 plasma periods
        snapshot interval step = 10
        snapshot interval time = 5.0000e-16 [s] = 0.0201 plasma periods
        Total number of snapshots = 1000
Production:
        No. of production steps = 10000
        Total production time = 5.0000e-13 [s] ~ 20 plasma periods
        snapshot interval step = 10
        snapshot interval time = 5.0000e-16 [s] = 0.0201 plasma periods
        Total number of snapshots = 1000

POTENTIAL:  yukawa
screening type : thomas-fermi
screening length = 4.702909e-09 [cm]
kappa = 2.7231
Gamma_eff = 101.96

ALGORITHM: pp
rcut = 4.6852 a_ws = 6.000000e-08 [cm]
No. of PP cells per dimension = [3 3 3]
No. of particles in PP loop = 662
No. of PP neighbors per particle = 102
Tot Force Error = 1.280814e-05


THERMOSTAT:
Type: berendsen
First thermostating timestep, i.e. thermalization_timestep = 300
Berendsen parameter tau: 2.000 [timesteps]
Berendsen relaxation rate: 0.500 [1/timesteps]
Thermostating temperatures:
Species ID 0: T_eq = 5.000000e+03 [K] = 4.308667e-01 [eV]

INTEGRATOR:
Equilibration Integrator Type: verlet
Production Integrator Type: verlet
Time step = 5.000000e-17 [s]
Total plasma frequency = 2.531585e+14 [rad/s]
Total plasma period = 2.481918e-14 [s]
w_p dt = 1.2658e-02 [rad] = 2.0146e-03
Timesteps per plasma cycle = 496



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


Particles Initialization Time: 0 sec 0 msec 877 usec 143 nsec

Potential Initialization Time: 1 sec 952 msec 930 usec 544 nsec

Total Simulation Initialization Time: 1 sec 970 msec 316 usec 807 nsec


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

Equilibration Time: 0 hrs 1 min 25 sec


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

Production Time: 0 hrs 1 min 23 sec

Total Time: 0 hrs 2 min 49 sec


========================= Filesize Estimates =========================

Equilibration:
        Checkpoint filesize: 0 GB 0 MB 85 KB 670 bytes
        Checkpoint folder size: 0 GB 83 MB 662 KB 304 bytes
Production:
        Checkpoint filesize: 0 GB 0 MB 85 KB 670 bytes
        Checkpoint folder size: 0 GB 83 MB 662 KB 304 bytes

Total occupied space: 0 GB 167 MB 300 KB 608 bytes

Postprocessing

Now that our simulation is complete we need to check if the simulation was physically sound. Run the following three lines will initialize the PostProcess class and calculate the observables defined in the input file.

It will also produce a plot of the Temperature and Total Energy of the Production phase.

# Initialize the Postprocessing class
postproc = PostProcess(input_file_name)
# Read the simulation's parameters and assign attributes
postproc.setup(read_yaml=True)

     _______.     ___      .______       __  ___      ___           _______.
    /       |    /   \     |   _  \     |  |/  /     /   \         /       |
   |   (----`   /  ^  \    |  |_)  |    |  '  /     /  ^  \       |   (----`
    \   \      /  /_\  \   |      /     |    <     /  /_\  \       \   \
.----)   |    /  _____  \  |  |\  \----.|  .  \   /  _____  \  .----)   |
|_______/    /__/     \__\ | _| `._____||__|\__\ /__/     \__\ |_______/



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




********************************************************************************
                                 Postprocessing
********************************************************************************

Job ID: yocp_quickstart
Job directory: SarkasSimulations/yocp_quickstart
PostProcessing directory:
SarkasSimulations/yocp_quickstart/PostProcessing

Equilibration dumps directory:
SarkasSimulations/yocp_quickstart/Simulation/Equilibration/dumps
Production dumps directory:
SarkasSimulations/yocp_quickstart/Simulation/Production/dumps

Equilibration Thermodynamics file:
SarkasSimulations/yocp_quickstart/Simulation/Equilibration/EquilibrationEnergy_yocp_quickstart.csv
Production Thermodynamics file:
SarkasSimulations/yocp_quickstart/Simulation/Production/ProductionEnergy_yocp_quickstart.csv

from sarkas.tools.observables import RadialDistributionFunction, Thermodynamics,VelocityAutoCorrelationFunction

therm = Thermodynamics()
therm.setup(postproc.parameters, phase = "production")
therm.compute()
therm.temp_energy_plot(postproc)

=========================== Thermodynamics ===========================
Data saved in:
 SarkasSimulations/yocp_quickstart/PostProcessing/Thermodynamics/Production/Thermodynamics_yocp_quickstart.h5
Data accessible via: self.dataframe_slices, self.dataframe

Time Series Data:
No. of slices = 1
No. dumps per slice = 1001
Total time: T = 5.0050e-13 [s] ~ 20 plasma periods
Time interval step: dt = 5.0000e-16 ~ 2.0146e-02 plasma period

ACF Data:
No. of acf slices = 1
No. dumps per slice = 500
Largest time lag of the autocorrelation function: tau = 2.5000e-13 [s] ~ 10 plasma periods



Thermodynamics Time: 0 sec 75 msec 179 usec 130 nsec

You will notice that both the energy and temperature oscillates wildly. This is fine as long as the percentage deviations, in the top plots, are small. You should have a temperature deviations between -2% to ~ 4-5% while energy deviations between -2% and 1%.

---

Observables

The most common observable is the radial distribution function. This can be calculated using the following code.

rdf = RadialDistributionFunction()
rdf.setup(postproc.parameters)
rdf.compute()

==================== Radial Distribution Function ====================
Data saved in:
 SarkasSimulations/yocp_quickstart/PostProcessing/RadialDistributionFunction/Production/RadialDistributionFunction_yocp_quickstart.h5
Data accessible via: self.dataframe_slices, self.dataframe
No. bins = 250
dr = 0.0187 a_ws = 2.4000e-10 [cm]
Maximum Distance (i.e. potential.rc)= 4.6852 a_ws = 6.0000e-08 [cm]
Time Series Data:
No. of slices = 1
No. dumps per slice = 1001
Total time: T = 5.0050e-13 [s] ~ 20 plasma periods
Time interval step: dt = 5.0000e-16 ~ 2.0146e-02 plasma period

Radial Distribution Function Calculation Time: 0 sec 41 msec 513 usec 995 nsec

# Initialize the Pair Distribution Function class
ax = rdf.plot(
    scaling=postproc.parameters.a_ws,
    y = [("C-C RDF", "Mean")],
    xlabel = r'$r/a_{\rm ws}$',
    ylabel = r'$g(r)$'
)
ax.legend(["C-C RDF"])

<matplotlib.legend.Legend at 0x7fb6a7767790>

Things to check in here are:

• Does \(g(r)\) go to 1 for large \(r\) values ?

• Is there a peak at \(r \sim ~1.5 a\) ?

• Is the height of this peak about ~ 1.4?

If the answer to all these question is yes than the simulation was successfull.