Basic Tutorial
Make sure to first go through the Configuration
steps.
Since DFControl is aimed at improving the day to day quality of life of a material's scientist/anyone running DFT codes, we will look at a simple demonstration of how by creating and submitting some Quantum-Espresso calculations on Si starting from a cif file specifying the structure.
using DFControl
First we download the cif file, extract the Structure
and assign the right pseudos to it. In this case Si (F d -3 m :1) from http://www.crystallography.net/cod/9011998.cif
using Downloads
cif_file = Downloads.download("http://www.crystallography.net/cod/9011998.cif", "Si.cif")
structure = Structure(cif_file)
set_pseudos!(structure, load(local_server(), PseudoSet("pbesol")))
***Warning : Site occupancies not found, assuming all occupancies = 1.
┌ Info: Structure extracted from Si.cif
│ cell parameters:
│ a = (5.43094 Å, 0.0 Å, 0.0 Å)
│ b = (0.0 Å, 5.43094 Å, 0.0 Å)
│ c = (0.0 Å, 0.0 Å, 5.43094 Å)
│ nat = 8
└ elements = [:Si]
This assumes that the "pbesol"
pseudopotential set was installed during the configuration step.
Next we specify the executables with which to run the QE calculations. We assume here that mpirun
is installed and in the user's PATH, and that QE is installed, and to be found in the /opt/qe/bin
, change this according to your own setup. The first argument to the constructor can be used as a label to later retrieve the executable after it was saved.
pw_exec = Exec(name="pw", path="/opt/qe/bin/pw.x",flags=Dict("-nk" => 4))
Exec("pw", "/opt/qe/bin/pw.x", Dict("-nk" => 4), String[], true)
Then we create the first calculation for our job, we name it scf, which will be used to reference it later. We also pass the executables to be used and additional flags to be set to the constructor. Afterwards we set the kpoints to be used in the scf calculation.
scf_calculation = Calculation("scf", :calculation => "scf"; exec = pw_exec)
set_kpoints!(scf_calculation, (6, 6, 6, 1, 1, 1))
Calculation{QE}:
name = scf
exec = pw.x
run = true
data = [:k_points]
flags:
&control
verbosity => high
calculation => scf
The code recognizes internally that this 6-Tuple corresponds to a K_POINTS (automatic)
block in QE. Alternatively (leading to an identical final result):
scf_calculation = Calculation("scf", :calculation => "scf"; exec = pw_exec,
data = [InputData(:k_points, :automatic,
(6, 6, 6, 1, 1, 1))])
Calculation{QE}:
name = scf
exec = pw.x
run = true
data = [:k_points]
flags:
&control
calculation => scf
We can now define our job: job = Job("Si", structure, [scfcalculation], :ecutwfc => 20, :convthr => 1e-6; dir="job") Additional calculations would be be added to the list [scf_calculation]
. The flag => value pairs will set the specified flags to that value for all calculations in the job that allow recognize that flag, so it's ideal for things like cutoffs and smearing etc.
job = Job("Si", structure, [scf_calculation], :ecutwfc => 40.0, :occupations => "smearing", :degauss=>0.01, :conv_thr => 1e-6, :nbnd => 18;
dir = dir, server=gethostname(), environment="default")
We are now ready to submit the job, which will run in the current working directory
submit(job)
Calculation{QE}:
name = bands
exec = pw.x
run = true
data = [:k_points]
flags:
&control
verbosity => high
calculation => bands
&system
ecutwfc => 20.0
&electrons
conv_thr => 1.0e-6
This will generate and save all calculation files, and the corresponding job script (job.sh
), to the server specified in job.server
, and then the job will be submitted on the server.
After the job finishes the outputs can be parsed through
outputdata(job)
Dict{String, Dict{Symbol, Any}} with 1 entry:
"scf" => Dict(:converged=>true, :scf_iteration=>[1, 2, 3, 4, 5, 6], :highest_…
this returns a dictionary with the input names as keys and a Dict with parsed outputs as values. i.e.
outputdata(job)["scf"]
Dict{Symbol, Any} with 14 entries:
:converged => true
:scf_iteration => [1, 2, 3, 4, 5, 6]
:highest_occupied => 6.0183
:timing => DFControl.TimingData[TimingData("init_run", 75 millisec…
:total_energy => [-90.6902, -90.6986, -90.7021, -90.7024, -90.7024, -90.…
:n_KS_states => 16
:accuracy => [0.210495, 0.00940295, 0.00048403, 6.398e-5, 1.21e-6, 2…
:bands => DFControl.Band[Band(10 k_points, eigvals: -5.7383 eV -…
:scf_steps => 1
:fermi => 6.0183
:initial_structure => Structure(volume: 160.1866745508759 Å^3, atoms: 8 Si)
:n_scf => 0
:finished => true
:n_electrons => 32
takes the outputdata of the input that we previously named "scf" when creating it.
Now that the scf calculation finished succesfully, the next step is usually to have a look at the bandstructure. For this we generate a bands calculation, using the scf calculation as a template and a generated high symmetry k-point path with 20 kpoints per segment.
bands_calc = Calculations.gencalc_bands(job["scf"], Structures.high_symmetry_kpath(job.structure, 20))
Calculation{QE}:
name = bands
exec = pw.x
run = true
data = [:k_points]
flags:
&control
verbosity => high
calculation => bands
&system
ecutwfc => 20.0
&electrons
conv_thr => 1.0e-6
Observe the :calculation => "bands"
, and automatic setting of the :verbosity => "high"
flags. We now push!
this calculation to the job queue.
push!(job, bands_calc)
2-element Vector{Calculation}:
Calculation{QE}(name: scf)
Calculation{QE}(name: bands)
However, since we know the scf succeeded there is no need to rerun it. To un-schedule it we do
job["scf"].run = false
false
Printing the job will now highlight the scheduled calculations differently from the non-scheduled ones
job
+-------------------------------------------JOB-------------------------------------------+
| name: Si |
| dir: /home/runner/work/DFControl.jl/DFControl.jl/src/../docs/src/assets/job |
| version: 0 |
| server: fv-az396-476, alive |
| versions: |
| last submission: 2023-09-14T14:10:18 |
| state: Unknown |
+-----------------------------------------------------------------------------------------+
(scheduled, not scheduled)
scf
bands
Seeing that all is right we submit the job again
submit(job)
+-------------------------------------------JOB-------------------------------------------+
| name: Si |
| dir: /home/runner/work/DFControl.jl/DFControl.jl/src/../docs/src/assets/job |
| version: 0 |
| server: fv-az396-476, alive |
| versions: |
| last submission: 2023-09-14T14:10:18 |
| state: Unknown |
+-----------------------------------------------------------------------------------------+
(scheduled, not scheduled)
scf
bands
We can access the bands through
bands = readbands(job);
or
bands = outputdata(job)["bands"][:bands];
They can be plotted too
fermi = readfermi(job) # = outputdata(job)["scf"][:fermi]
using Plots
plot(bands; fermi = fermi)
Since more info (such as the structure) is available in the job, plotting the job leads to a richer plot
plot(job, -10, 1)
As can be seen in the Advanced Usage, additional information generated by additional calculations will be picked up by DFControl in order to create richer plots.