Skip to content

Commit 48d7ebb

Browse files
mesonepigrecomesonepigreco
committed
Added a documentation of the gold function
1 parent 6b3c9a3 commit 48d7ebb

1 file changed

Lines changed: 184 additions & 1 deletion

File tree

‎UserGuide/start.rst‎

Lines changed: 184 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,12 +1,195 @@
11
Quick start
22
===========
33

4+
In this chapter we provide ready to use examples to setup your first SSCHA calculation.
5+
6+
The free energy of gold: a simulation in the NVT ensemble
7+
---------------------------------------------------------
8+
9+
This simple tutorial explains how to setup a SSCHA calculation starting just from the structure, in this case a cif file we downloaded from the [Materials Project](https://materialsproject.org/materials/mp-81/) database.
10+
11+
You can find there a lot of structures ready to use for your SSCHA runs.
12+
13+
For the purpouse of this tutorial, we are going to use the EMT force field, so that the calculation can be run in a laptop without the need of a supercomputer.
14+
We explain in a later section how to couple the SSCHA with a cluster to submit the same calculation fully ab-initio.
15+
16+
Starting from the Gold structure in the primitive cell, to run the SSCHA we need:
17+
- Compute the harmonic phonons (dynamical matrix)
18+
- Remove imaginary frequencies (if any)
19+
- Run the SSCHA
20+
21+
We prepared an input file in the form of a python script (tested with python-sscha version 1.2) which makes all these passages automatically.
22+
23+
.. code-block:: python
24+
25+
# Import the sscha code
26+
import sscha, sscha.Ensemble, sscha.SchaMinimizer, sscha.Relax, sscha.Utilities
27+
28+
# Import the cellconstructor library to manage phonons
29+
import cellconstructor as CC, cellconstructor.Phonons
30+
import cellconstructor.Structure, cellconstructor.calculators
31+
32+
# Import the force field of Gold
33+
import ase, ase.calculators
34+
from ase.calculators.emt import EMT
35+
36+
# Import numerical and general pourpouse libraries
37+
import numpy as np, matplotlib.pyplot as plt
38+
import sys, os
39+
40+
41+
"""
42+
Here we load the primitive cell of Gold from a cif file.
43+
And we use CellConstructor to compute phonons from finite differences.
44+
The phonons are computed on a q-mesh 4x4x4
45+
"""
46+
47+
gold_structure = CC.Structure.Structure()
48+
gold_structure.read_generic_file("Au.cif")
49+
50+
# Get the force field for gold
51+
calculator = EMT()
52+
53+
# Relax the gold structure (useless since for symmetries it is already relaxed)
54+
relax = CC.calculators.Relax(gold_structure, calculator)
55+
gold_structure_relaxed = relax.static_relax()
56+
57+
# Compute the harmonic phonons
58+
# NOTE: if the code is run with mpirun, the calculation goes in parallel
59+
gold_harmonic_dyn = CC.Phonons.compute_phonons_finite_displacements(gold_structure_relaxed, calculator, supercell = (4,4,4))
60+
61+
# Impose the symmetries and
62+
# save the dynamical matrix in the quantum espresso format
63+
gold_harmonic_dyn.Symmetrize()
64+
gold_harmonic_dyn.save_qe("harmonic_dyn")
65+
66+
67+
# If the dynamical matrix has imaginary frequencies, remove them
68+
gold_harmonic_dyn.ForcePositiveDefinite()
69+
70+
"""
71+
gold_harmonic_dyn is ready to start the SSCHA calculation.
72+
73+
Now let us initialize the ensemble, and the calculation at 300 K.
74+
We will run a NVT calculation, using 100 configurations at each step
75+
"""
76+
77+
TEMPERATURE = 300
78+
N_CONFIGS = 50
79+
MAX_ITERATIONS = 20
80+
81+
# Initialize the random ionic ensemble
82+
ensemble = sscha.Ensemble.Ensemble(gold_harmonic_dyn, TEMPERATURE)
83+
84+
# Initialize the free energy minimizer
85+
minim = sscha.SchaMinimizer.SSCHA_Minimizer(ensemble)
86+
minim.set_minimization_step(0.01)
87+
88+
# Initialize the NVT simulation
89+
relax = sscha.Relax.SSCHA(minim, calculator, N_configs = N_CONFIGS,
90+
max_pop = MAX_ITERATIONS)
91+
92+
# Define the I/O operations
93+
# To save info about the free energy minimization after each step
94+
ioinfo = sscha.Utilities.IOInfo()
95+
ioinfo.SetupSaving("minim_info")
96+
relax.setup_custom_functions(custom_function_post = ioinfo.CFP_SaveAll)
97+
98+
99+
# Run the NVT simulation (save the stress to compute the pressure)
100+
relax.relax(get_stress = True)
101+
102+
# If instead you want to run a NPT simulation, use
103+
# The target pressure is given in GPa.
104+
#relax.vc_relax(target_press = 0)
105+
106+
# You can also run a mixed simulation (NVT) but with variable lattice parameters
107+
#relax.vc_relax(fix_volume = True)
108+
109+
# Now we can save the final dynamical matrix
110+
# And print in stdout the info about the minimization
111+
relax.minim.finalize()
112+
relax.minim.dyn.save_qe("sscha_T{}_dyn".format(TEMPERATURE))
113+
114+
115+
116+
While the input may seem long, it is heavily commented, but lets go through it step by step.
117+
At the very beginning, we simply import the sscha libraries, cellconstructor, the math libraries and the force field. This is done in python with the `import` statemets.
118+
119+
The first real part of the code is:
120+
121+
.. code-block:: python
122+
123+
gold_structure = CC.Structure.Structure()
124+
gold_structure.read_generic_file("Au.cif")
125+
126+
# Get the force field for gold
127+
calculator = EMT()
128+
129+
# Relax the gold structure (useless since for symmetries it is already relaxed)
130+
relax = CC.calculators.Relax(gold_structure, calculator)
131+
gold_structure_relaxed = relax.static_relax()
132+
133+
Here we initialize a cellconstructor structure from the cif file downloaded from the material database (*Au.cif*). We initialize the EMT calculator from ASE, and relax the structure.
134+
135+
In the case of Gold the relaxation is useless, as it is a FCC structure with Fm-3m symmetry group and 1 atom per primitive cell. This means the atomic positions have no degrees of freedom, thus the relaxation will end before even start.
136+
137+
In the next part of the code, we perform the harmonic phonon calculation using cellconstructor and a finite displacement approach:
138+
139+
.. code-block:: python
140+
141+
gold_harmonic_dyn = CC.Phonons.compute_phonons_finite_displacements(gold_structure_relaxed, calculator, supercell = (4,4,4))
142+
143+
# Impose the symmetries and
144+
# save the dynamical matrix in the quantum espresso format
145+
gold_harmonic_dyn.Symmetrize()
146+
gold_harmonic_dyn.save_qe("harmonic_dyn")
147+
148+
149+
# If the dynamical matrix has imaginary frequencies, remove them
150+
gold_harmonic_dyn.ForcePositiveDefinite()
151+
152+
153+
The method `compute_phonons_finite_displacements` is documented in the CellConstructor guide. It requires the structure (in this case `gold_structure_relaxed`), the force-field (`calculator`) and the supercell for the calculation. In this case we use a 4x4x4 (equivalent to 64 atoms). This may not be sufficient to converge all the properties, especially at very high temperature, but it is just a start.
154+
155+
Note that `compute_phonons_finite_displacements` works in parallel with MPI, therefore, if the script is executed with `mpirun -np 16 python myscript.py` it will split the calculations of the finite displacements across 16 processors. You need to have mpi4py installed.
156+
157+
After computing the harmonic phonons in gold_harmonic_dyn, we impose the correct symmetrization and the acousitic sum rule with the `Symmetrize` method, and save the result in the quantum ESPRESSO format with `save_qe`.
158+
This should not be the case for Gold, however, if we have a structure which has imaginary phonon frequencies, we need to get rid of them before starting the SSCHA. This is achieved with `ForcePositiveDefinite` (see CellConstructor documentation for more details on how these methods work).
159+
160+
161+
**Now we are ready to submit the SSCHA calculation in the NVT ensemble!**.
162+
The important parameters are:
163+
- The temperature
164+
- The number of random configurations in the ensemble
165+
- The maximum number of iterations
166+
167+
168+
These parameters are almost self-explaining. However, we give a brief overview of how the SSCHA works to help you understand which are the best one for your case.
169+
While MD or MC calculation represent the equilibrium probability distribution over time of the system by updating a single structure, the SSCHA encodes the whole probability distribution as an analytical function. Therefore, to compute properties, we can generate on the fly the ionic configurations that represent the equilibrium distributions.
170+
The number of random configuration is exactly how many ionic configuration we generate to compute the properties (Free energy and Stress tensors)
171+
172+
173+
174+
175+
176+
177+
178+
179+
180+
181+
182+
183+
184+
185+
186+
4187

5188
To quickly start using the code, we recommend using the jupyter notebooks with examples we provide in the Tutorials directory of the source code.
6189

7190
Tutorials are organized as follows:
8191

9-
1. Setup the first calculation: PbTe tutorial. Here you learn how to set up a SSCHA calculation starting just with the structure (we provide a .cif file of the PbTe at high temperature). The tutorial will guide you step by step. You will learn how to: prepare the starting data needed for the SSCHA calculation, generate a random ensemble, save the ensemble and prepare input files for your favorite ab-initio code, read back the energies and the forces inside SSCHA, run a SSCHA minimization. You will also learn how to use ASE and the Cluster module to automatize the calculation of the ensemble and submit it to a HPC system.
192+
1. Setup from the structure and manual submission: PbTe tutorial. Here you learn how to set up a SSCHA calculation starting just with the structure (we provide a .cif file of the PbTe at high temperature). The tutorial will guide you step by step. You will learn how to: prepare the starting data needed for the SSCHA calculation, generate a random ensemble, save the ensemble and prepare input files for your favorite ab-initio code, read back the energies and the forces inside SSCHA, run a SSCHA minimization. You will also learn how to use ASE and the Cluster module to automatize the calculation of the ensemble and submit it to a HPC system.
10193
2. Automatic relaxation with a force field: SnTe_ToyModel. Here, we show how to use a force-field for a SSCHA calculation, running everything on your computer. We also will explain how to calculate the free energy hessian for second-order phase transitions, and study a phase transition as a function of temperature.
11194
3. Variable cell relaxation: LaH10 tutorial. Here you learn how to perform an automatic calculation with a variable cell. You will exploit the quantum effect to search the high-temperature superconductive phase (Fm-3m) of LaH10 below 200 GPa, starting from a distorted structure.
12195
4. Hessian matrix calculation for second-order phase transitions: H3S tutorial. Here you reproduce the calculation of the Hessian of the free energy to assert the stability of the H3S phase.

0 commit comments

Comments
 (0)