Geometry optimisation#
In this example we will perform a geometry optimisation of silicon using VASP. We will use the VASP code and the AiiDA plugin for VASP.
Setting up the environment#
The code block below configures the environment for this example. Please see the single point example for more details.
from aiida_vasp.utils.temp_profile import *
print(load_temp_profile())
# Uncomment the below line to create a localhost Computer if you have not done so
comp = orm.Computer('localhost', 'localhost', transport_type='core.local', scheduler_type='core.direct')
comp.store()
# Some configuration may be needed for first-time user
from pathlib import Path
import os
comp.set_workdir('/tmp/aiida_run/')
comp.configure()
vasp_path = !which mock-vasp
vasp_code = orm.InstalledCode(comp, vasp_path[0], default_calc_job_plugin='vasp.vasp')
vasp_code.label ='mock-vasp'
vasp_code.store()
os.environ['MOCK_VASP_REG_BASE'] = str((Path() / 'mock_registry').absolute())
os.environ['MOCK_VASP_UPLOAD_PREFIX'] = 'relax'
# Upload the POTCAR files
from aiida_vasp.data.potcar import PotcarData, PotcarFileData
from pathlib import Path
print(PotcarData.upload_potcar_family(str(Path('potcars').absolute()), "PBE.EXAMPLE", "PBE.EXAMPLE"))
# Setting up the silicon structure
from ase.build import bulk
si = bulk('Si', 'diamond', 5.4)
si_node = orm.StructureData(ase=si)
Profile<uuid='df4506db86b148339ed6f6790c12ac07' name='myprofile'>
(1, 1, 1)
Running the relaxation#
Similar to the single point calculation tutorial, we will use a BuilderUpdater to setup the inputs
for the VaspRelaxWorkChain.
from aiida import orm
from aiida_vasp.workchains.v2 import VaspRelaxUpdater
upd = VaspRelaxUpdater().apply_preset(si_node, code='mock-vasp@localhost')
upd.builder.vasp.potential_family = 'PBE.EXAMPLE'
upd.builder
Process class: VaspRelaxWorkChain
Inputs:
metadata:
label: ''
relax_settings:
algo: cg
clean_reuse: true
convergence_absolute: false
convergence_max_iterations: 5
convergence_mode: last
convergence_on: true
convergence_positions: 0.1
convergence_shape_angles: 0.1
convergence_shape_lengths: 0.1
convergence_volume: 0.01
energy_cutoff: null
force_cutoff: 0.03
hybrid_calc_bootstrap: false
hybrid_calc_bootstrap_wallclock: 3600
keep_magnetization: false
keep_sp_workdir: false
perform: true
positions: true
reuse: false
shape: true
steps: 60
volume: true
structure: Si
vasp:
code: mock-vasp@localhost
dynamics: {}
kpoints_spacing: 0.05
metadata: {}
options:
max_wallclock_seconds: 3600
resources:
num_machines: 1
tot_num_mpiprocs: 1
withmpi: false
parameters:
incar:
algo: fast
ediff: 2.0e-06
encut: 250
gga: ps
ismear: 0
ispin: 2
lasph: true
lcharg: false
lreal: false
lvhar: true
lwave: false
nedos: 100
nelm: 200
nelmin: 4
nwrite: 1
prec: normal
sigma: 0.05
potential_family: PBE.EXAMPLE
potential_mapping:
Si: Si
settings: {}
Here we can see that the input nodes are different from that of the single point calculation (workflow).
Nodes such as the settings and parameters go into a vasp input names rather than at the root level.
This is because the VaspRelaxWorkChain exposes the inputs of the VaspWorkChain in the vasp namespace.
Also note that the structure input is still at the root level.
This is because the atomic structure is an essential input for the relaxation (and another calculations).
Note that there is a also a relax_settings input, which contains control parameters for the relaxation.
You can see the available options by inspecting the field of the builder:
upd.get_input_help('relax_settings')
algo: str
Default: The algorithm to use for relaxation
energy_cutoff: Optional
Default: The cut off energy difference when the relaxation is stopped (e.g. EDIFF)
force_cutoff: float
Default: The maximum force when the relaxation is stopped (e.g. EDIFFG)
steps: int
Default: Number of relaxation steps to perform (eg. NSW)
positions: bool
Default: If True, perform relaxation of the atomic positions
shape: bool
Default: If True, perform relaxation of the cell shape
volume: bool
Default: If True, perform relaxation of the cell volume
convergence_on: bool
Default: If True, perform convergence checks within the workchain
convergence_absolute: bool
Default: If True, use absolute values where possible when performing convergence checks
convergence_max_iterations: int
Default: Maximum iterations for convergence checking
convergence_positions: float
Default: The cutoff value for the convergence check on positions in Angstram. A negative value by pass the check.
convergence_volume: float
Default: The cutoff value for the convergence check on volume between the two structures. A negative value by pass the check.
convergence_shape_lengths: float
Default: The cutoff value for the convergence check on the lengths of the unit cell vectors, between input and the outputs structure. A negative value by pass the check.
convergence_shape_angles: float
Default: The cutoff value for the convergence check on the angles of the unit cell vectors, between input and the outputs structure. A negative value by pass the check.
convergence_mode: str
Default: Mode of the convergence check for positions. 'inout' for checking input/output structure, or 'last' to check only the change of the last step.
reuse: bool
Default: Whether reuse the previous calculation by copying over the remote folder
clean_reuse: bool
Default: Whether to perform a final cleaning of the reused calculations
keep_sp_workdir: bool
Default: Whether to keep the workdir of the final singlepoint calculation
perform: bool
Default: Do not perform any relaxation if set to 'False'
hybrid_calc_bootstrap: bool
Default: Whether to bootstrap hybrid calculation by performing standard DFT first
hybrid_calc_bootstrap_wallclock: int
Default: Wall time limit in second for the bootstrap calculation
keep_magnetization: bool
Default: Whether to keep magnetization from the previous calculation if possible
or simply hit <Shift>+<Tab> after upd.builder.relax_settings in the Notebook code cell.
We can now run the relaxation using the same run_get_node method as in the single point example.
results = upd.run_get_node()
09/13/2025 06:23:30 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|run_relax]: launching VaspWorkChain<19> iterations #1
09/13/2025 06:23:30 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [19|VaspWorkChain|run_process]: launching VaspCalculation<22> iteration #1
09/13/2025 06:23:33 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [19|VaspWorkChain|results]: work chain completed after 1 iterations
09/13/2025 06:23:33 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [19|VaspWorkChain|on_terminated]: cleaned remote folders of calculations: 22
09/13/2025 06:23:33 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|analyze_convergence]: Checking the convergence of the relaxation.
/home/docs/checkouts/readthedocs.org/user_builds/aiida-vasp-plugin/envs/v4.1.0/lib/python3.11/site-packages/aiida_vasp/workchains/v2/relax.py:845: RuntimeWarning: invalid value encountered in divide
delta.relative.pos_lengths = np.array([np.linalg.norm(vector) for vector in site_vectors]) / a_lengths
09/13/2025 06:23:33 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|check_positions_convergence]: max site position change is 0.01057, tolerance is 0.1 - OK
09/13/2025 06:23:33 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|check_volume_convergence]: cell volume changed by 0.01336, tolerance is 0.01 - NOT OK
09/13/2025 06:23:33 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|check_shape_convergence]: cell lengths changed by max 0.004474, tolerance is 0.1 - OK
09/13/2025 06:23:33 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|check_shape_convergence]: cell angles changed by max 1.184e-16, tolerance is 0.1 - OK
09/13/2025 06:23:33 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|analyze_convergence]: Maximum force in the structure 0 - OK
09/13/2025 06:23:33 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|analyze_convergence]: Relaxation did not converge, restarting the relaxation.
09/13/2025 06:23:34 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|run_relax]: launching VaspWorkChain<28> iterations #2
09/13/2025 06:23:34 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [28|VaspWorkChain|run_process]: launching VaspCalculation<31> iteration #1
09/13/2025 06:23:36 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [28|VaspWorkChain|results]: work chain completed after 1 iterations
09/13/2025 06:23:36 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [28|VaspWorkChain|on_terminated]: cleaned remote folders of calculations: 31
09/13/2025 06:23:37 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|analyze_convergence]: Checking the convergence of the relaxation.
/home/docs/checkouts/readthedocs.org/user_builds/aiida-vasp-plugin/envs/v4.1.0/lib/python3.11/site-packages/aiida_vasp/workchains/v2/relax.py:845: RuntimeWarning: invalid value encountered in divide
delta.relative.pos_lengths = np.array([np.linalg.norm(vector) for vector in site_vectors]) / a_lengths
09/13/2025 06:23:37 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|check_positions_convergence]: max site position change is 0, tolerance is 0.1 - OK
09/13/2025 06:23:37 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|check_volume_convergence]: cell volume changed by 0, tolerance is 0.01 - OK
09/13/2025 06:23:37 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|check_shape_convergence]: cell lengths changed by max 0, tolerance is 0.1 - OK
09/13/2025 06:23:37 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|check_shape_convergence]: cell angles changed by max 0, tolerance is 0.1 - OK
09/13/2025 06:23:37 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|analyze_convergence]: Maximum force in the structure 0 - OK
09/13/2025 06:23:37 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|analyze_convergence]: Relaxation is converged, finishing with a final static calculation.
09/13/2025 06:23:37 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|store_relaxed]: attaching the node StructureData<35> as 'relax.structure'
09/13/2025 06:23:37 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|run_static_calculation]: performing a final calculation using the relaxed structure.
09/13/2025 06:23:38 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [16|VaspRelaxWorkChain|run_static_calculation]: launching VaspWorkChain<38> iterations #3
09/13/2025 06:23:38 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [38|VaspWorkChain|run_process]: launching VaspCalculation<41> iteration #1
09/13/2025 06:23:40 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [38|VaspWorkChain|results]: work chain completed after 1 iterations
09/13/2025 06:23:40 AM <842> aiida.orm.nodes.process.workflow.workchain.WorkChainNode: [REPORT] [38|VaspWorkChain|on_terminated]: cleaned remote folders of calculations: 41
The VaspRelaxWorkChain also has a misc output node, which is in fact the output of the last VaspWorkChain it launched.
results.node.outputs.misc.get_dict()
{'total_energies': {'energy_extrapolated': -10.62045741,
'energy_extrapolated_electronic': -10.62045741},
'stress': [[0.07859658, 0.0, -0.0],
[-0.0, 0.07859658, 0.0],
[0.0, 0.0, 0.07859658]],
'forces': [[0.0, 0.0, -0.0], [-0.0, -0.0, 0.0]],
'fermi_level': 6.03669271,
'band_properties': {'cbm': 6.2549,
'vbm': 5.7938,
'is_direct_gap': False,
'band_gap': 0.4611},
'version': '6.2.0',
'magnetization': [-0.0004373],
'site_magnetization': {'sphere': {'x': {'site_moment': {},
'total_magnetization': {}},
'y': {'site_moment': {}, 'total_magnetization': {}},
'z': {'site_moment': {}, 'total_magnetization': {}}},
'full_cell': [-0.0004373]},
'run_stats': {'mem_usage_base': 30000.0,
'mem_usage_nonl-proj': 478.0,
'mem_usage_fftplans': 1184.0,
'mem_usage_grid': 2809.0,
'mem_usage_one-center': 12.0,
'mem_usage_wavefun': 2194.0,
'total_cpu_time_used': 3.241,
'user_time': 3.221,
'system_time': 0.019,
'elapsed_time': 3.566,
'maximum_memory_used': 100124.0,
'average_memory_used': None},
'run_status': {'nelm': 200,
'nsw': 0,
'last_iteration_index': [1, 13],
'finished': True,
'ionic_converged': None,
'electronic_converged': True,
'consistent_nelm_breach': False,
'contains_nelm_breach': False,
'nbands': 9},
'notifications': []}
The final relaxed structure is stored in the relax.structure output node.
relaxed_node = results.node.outputs.relax.structure
print(f"Volume before relaxations: {si_node.get_cell_volume():3f} A^3")
print(f"Volume after relaxations: {relaxed_node.get_cell_volume():3f} A^3")
Volume before relaxations: 39.366000 A^3
Volume after relaxations: 40.087704 A^3
Why use VaspRelaxWorkChain?#
Surely VASP has the option to run relaxation itself (ISIF=3, IBRION=2, NSW=60), so why do we need a VaspRelaxWorkChain in the first place?
The reason is tha VaspRelaxWorkChain run more strict check and verifies if the structure is really fully relaxed.
Also, since VASP run variable cell relaxation in the fixed-basis mode, the effective plane wave cut off energies changes with the cell volume inside a single calculation.
This means multiple calculations are needed to get fully (and self-consistently) relaxed structures.
In addition, the VaspRelaxWorkChain will run a single point calculation in the end to get the correct energy that is consistent with the plane wave cut-off energy (ENCUT) originally set.
The called and called_descendants property can be used to inspect the sub-processes launched by the VaspRelaxWorkChain
print("Directly called workchains:", results.node.called)
print("All called processes:", results.node.called_descendants)
Directly called workchains: [<WorkChainNode: uuid: e0babb15-94f8-48e5-b72d-f0604da406ad (pk: 19) (aiida.workflows:vasp.v2.vasp)>, <WorkChainNode: uuid: d016ce82-fcb2-4160-9af3-2ffc5208bda2 (pk: 28) (aiida.workflows:vasp.v2.vasp)>, <WorkChainNode: uuid: 9ea3af72-06f9-42a8-8edf-6cd8f8680cb3 (pk: 38) (aiida.workflows:vasp.v2.vasp)>]
All called processes: [<WorkChainNode: uuid: e0babb15-94f8-48e5-b72d-f0604da406ad (pk: 19) (aiida.workflows:vasp.v2.vasp)>, <CalcJobNode: uuid: 24ab8269-daf8-4ddf-a02f-c3d1baf4a695 (pk: 22) (aiida.calculations:vasp.vasp)>, <WorkChainNode: uuid: d016ce82-fcb2-4160-9af3-2ffc5208bda2 (pk: 28) (aiida.workflows:vasp.v2.vasp)>, <CalcJobNode: uuid: 5bd5c6ca-67d4-431d-9b2c-89a3783a2ba7 (pk: 31) (aiida.calculations:vasp.vasp)>, <WorkChainNode: uuid: 9ea3af72-06f9-42a8-8edf-6cd8f8680cb3 (pk: 38) (aiida.workflows:vasp.v2.vasp)>, <CalcJobNode: uuid: 81bae8f6-8e11-4d06-8ef6-12112679c089 (pk: 41) (aiida.calculations:vasp.vasp)>]