Skip to content

Commit

Permalink
Add the verification-pbe-v1-sirius protocol for CP2K (#312)
Browse files Browse the repository at this point in the history 
Co-authored-by: hmhoseini <[email protected]>
Co-authored-by: Aliaksandr Yakutovich <[email protected]>
  • Loading branch information
@hmhoseini @HosseinMir @yakutovicha
3 people authored Aug 30, 2023
1 parent 3e4014d commit 3015631
Show file tree
Hide file tree
Showing 3 changed files with 186 additions and 47 deletions.
119 changes: 73 additions & 46 deletions aiida_common_workflows/workflows/relax/cp2k/generator.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -38,21 +38,25 @@ def dict_merge(dct, merge_dct):
dct[k] = merge_dct[k]


def get_kinds_section(structure: StructureData, basis_pseudo: str, magnetization_tags=None):
""" Write the &KIND sections given the structure and the settings_dict"""
def get_kinds_section(structure: StructureData, basis_pseudo=None, magnetization_tags=None, use_sirius=False):
"""Write the &KIND sections given the structure and the settings_dict."""
kinds = []
with open(pathlib.Path(__file__).parent / basis_pseudo, 'rb') as fhandle:
atom_data = yaml.safe_load(fhandle)
if not use_sirius:
with open(pathlib.Path(__file__).parent / basis_pseudo, 'rb') as fhandle:
atom_data = yaml.safe_load(fhandle)
ase_structure = structure.get_ase()
symbol_tag = {
(symbol, str(tag)) for symbol, tag in zip(ase_structure.get_chemical_symbols(), ase_structure.get_tags())
}
for symbol, tag in symbol_tag:
new_atom = {
'_': symbol if tag == '0' else symbol + tag,
'BASIS_SET': atom_data['basis_set'][symbol],
'POTENTIAL': atom_data['pseudopotential'][symbol],
}
if use_sirius:
new_atom['POTENTIAL'] = f'UPF {symbol}.json'
else:
new_atom['BASIS_SET'] = atom_data['basis_set'][symbol]
new_atom['POTENTIAL'] = atom_data['pseudopotential'][symbol]
if magnetization_tags:
new_atom['MAGNETIZATION'] = magnetization_tags[tag]
kinds.append(new_atom)
Expand Down Expand Up @@ -126,6 +130,12 @@ def get_file_section():
}


def get_upf_pseudos_section(structure: StructureData, pseudo_family):
"""Provide the pseudos section to the input."""
pseudos_group = orm.load_group(pseudo_family)
return pseudos_group.get_pseudos(structure=structure)


class Cp2kCommonRelaxInputGenerator(CommonRelaxInputGenerator):
"""Input generator for the `Cp2kRelaxWorkChain`."""

Expand All @@ -150,7 +160,7 @@ def define(cls, spec):
super().define(spec)
spec.input(
'protocol',
valid_type=ChoiceType(('fast', 'moderate', 'precise', 'verification-pbe-v1')),
valid_type=ChoiceType(('fast', 'moderate', 'precise', 'verification-pbe-v1', 'verification-pbe-v1-sirius')),
default='moderate',
help='The protocol to use for the automated input generation. This value indicates the level of precision '
'of the results and computational cost that the input parameters will be selected for.',
Expand Down Expand Up @@ -190,52 +200,66 @@ def _construct_builder(self, **kwargs) -> engine.ProcessBuilder:
kpoints_distance = protocol_dict.pop('kpoints_distance', None)
kpoints = self._get_kpoints(kpoints_distance, structure, reference_workchain)
mesh, _ = kpoints.get_kpoints_mesh()
if mesh != [1, 1, 1]:
if 'sirius' in protocol:
parameters['FORCE_EVAL']['PW_DFT']['PARAMETERS']['NGRIDK'] = f'{mesh[0]} {mesh[1]} {mesh[2]}'
elif mesh != [1, 1, 1]:
builder.cp2k.kpoints = kpoints

magnetization_tags = None

# Metal or insulator.
# If metal then add smearing, unoccupied orbitals, and employ diagonalization.
if electronic_type == ElectronicType.METAL:
parameters['FORCE_EVAL']['DFT']['SCF']['SMEAR'] = {
'_': 'ON',
'METHOD': 'FERMI_DIRAC',
'ELECTRONIC_TEMPERATURE': '[K] 710.5',
}
parameters['FORCE_EVAL']['DFT']['SCF']['DIAGONALIZATION'] = {
'EPS_ADAPT': '1',
}
parameters['FORCE_EVAL']['DFT']['SCF']['MIXING'] = {
'METHOD': 'BROYDEN_MIXING',
'ALPHA': '0.1',
'BETA': '1.5',
}
parameters['FORCE_EVAL']['DFT']['SCF']['ADDED_MOS'] = 20
parameters['FORCE_EVAL']['DFT']['SCF']['CHOLESKY'] = 'OFF'

# If insulator then employ OT.
elif electronic_type == ElectronicType.INSULATOR:
parameters['FORCE_EVAL']['DFT']['SCF']['OT'] = {
'PRECONDITIONER': 'FULL_SINGLE_INVERSE',
'MINIMIZER': 'CG',
}

# Magnetic calculation.
if spin_type == SpinType.NONE:
parameters['FORCE_EVAL']['DFT']['UKS'] = False
if magnetization_per_site is not None:
import warnings
warnings.warn('`magnetization_per_site` will be ignored as `spin_type` is set to SpinType.NONE')

elif spin_type == SpinType.COLLINEAR:
parameters['FORCE_EVAL']['DFT']['UKS'] = True
structure, magnetization_tags = tags_and_magnetization(structure, magnetization_per_site)
parameters['FORCE_EVAL']['DFT']['MULTIPLICITY'] = guess_multiplicity(structure, magnetization_per_site)
if 'sirius' not in protocol: # It is not possible to disable smearing for sirius.
# If metal then add smearing, unoccupied orbitals, and employ diagonalization.
if electronic_type == ElectronicType.METAL:
parameters['FORCE_EVAL']['DFT']['SCF']['SMEAR'] = {
'_': 'ON',
'METHOD': 'FERMI_DIRAC',
'ELECTRONIC_TEMPERATURE': '[K] 710.5',
}
parameters['FORCE_EVAL']['DFT']['SCF']['DIAGONALIZATION'] = {
'EPS_ADAPT': '1',
}
parameters['FORCE_EVAL']['DFT']['SCF']['MIXING'] = {
'METHOD': 'BROYDEN_MIXING',
'ALPHA': '0.1',
'BETA': '1.5',
}
parameters['FORCE_EVAL']['DFT']['SCF']['ADDED_MOS'] = 20
parameters['FORCE_EVAL']['DFT']['SCF']['CHOLESKY'] = 'OFF'

# If insulator then employ OT.
elif electronic_type == ElectronicType.INSULATOR:
parameters['FORCE_EVAL']['DFT']['SCF']['OT'] = {
'PRECONDITIONER': 'FULL_SINGLE_INVERSE',
'MINIMIZER': 'CG',
}

# Magnetic calculation. Changing this configuration has no effect on sirius calculations.
if spin_type == SpinType.NONE:
parameters['FORCE_EVAL']['DFT']['UKS'] = False
if magnetization_per_site is not None:
import warnings
warnings.warn('`magnetization_per_site` will be ignored as `spin_type` is set to SpinType.NONE')

elif spin_type == SpinType.COLLINEAR:
parameters['FORCE_EVAL']['DFT']['UKS'] = True
structure, magnetization_tags = tags_and_magnetization(structure, magnetization_per_site)
parameters['FORCE_EVAL']['DFT']['MULTIPLICITY'] = guess_multiplicity(structure, magnetization_per_site)

# Starting magnetization.
basis_pseudo = protocol_dict.pop('basis_pseudo')
dict_merge(parameters, get_kinds_section(structure, basis_pseudo, magnetization_tags))
if 'sirius' in protocol:
dict_merge(
parameters,
get_kinds_section(structure, basis_pseudo=None, magnetization_tags=magnetization_tags, use_sirius=True)
)
else:
dict_merge(
parameters,
get_kinds_section(
structure, basis_pseudo=basis_pseudo, magnetization_tags=magnetization_tags, use_sirius=False
)
)

# Relaxation type.
if relax_type == RelaxType.POSITIONS:
Expand Down Expand Up @@ -281,7 +305,10 @@ def _construct_builder(self, **kwargs) -> engine.ProcessBuilder:
builder.handler_overrides = orm.Dict(dict={'restart_incomplete_calculation': True})

# Files.
builder.cp2k.file = get_file_section()
if 'sirius' in protocol:
builder.cp2k.pseudos_upf = get_upf_pseudos_section(structure, basis_pseudo)
else:
builder.cp2k.file = get_file_section()

# Input structure.
builder.cp2k.structure = structure
Expand Down
112 changes: 112 additions & 0 deletions aiida_common_workflows/workflows/relax/cp2k/protocol.yml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -467,3 +467,115 @@ verification-pbe-v1:
RMS_FORCE: "[bohr^-1*hartree] 0.00030" #default: [bohr^-1*hartree] 0.00030
BFGS:
TRUST_RADIUS: "[angstrom] 0.5"
verification-pbe-v1-sirius:
description: 'Protocol to relax a structure with SIRIUS and UPF pseudopotentials'
kpoints_distance: 0.06
basis_pseudo: SSSP/1.2/PBE/precision
input:
GLOBAL:
EXTENDED_FFT_LENGTHS: True
PREFERRED_DIAG_LIBRARY: ScaLAPACK
PRINT_LEVEL: MEDIUM #default: MEDIUM. Important to have the correct parsing
FORCE_EVAL:
METHOD: SIRIUS
STRESS_TENSOR: ANALYTICAL
PRINT:
FORCES:
FILENAME: requested-forces
ADD_LAST: SYMBOLIC
PW_DFT:
CONTROL:
VERBOSITY: 2
PARAMETERS:
ELECTRONIC_STRUCTURE_METHOD: pseudopotential
USE_SYMMETRY: TRUE
GK_CUTOFF: 10 # a.u.^-1
PW_CUTOFF: 55 # a.u.^-1
ENERGY_TOL: 1e-15
NUM_DFT_ITER: 400
SMEARING: FERMI_DIRAC
SMEARING_WIDTH: 0.00225
ITERATIVE_SOLVER:
ENERGY_TOLERANCE: 0.001
NUM_STEPS: 20
SUBSPACE_SIZE: 4
CONVERGE_BY_ENERGY: 1
DFT:
XC:
XC_FUNCTIONAL:
GGA_X_PBE:
_: ''
GGA_C_PBE:
_: ''
PRINT:
MO:
_: 'OFF'
ADD_LAST: SYMBOLIC
EIGENVALUES: ON
OCCUPATION_NUMBERS: ON
NDIGITS: 12
EACH:
CELL_OPT: 0
GEO_OPT: 0
MD: 0
QS_SCF: 0
MULLIKEN:
_: 'ON'
ADD_LAST: SYMBOLIC
EACH:
CELL_OPT: 0
GEO_OPT: 0
MD: 0
LOWDIN:
_: 'OFF'
HIRSHFELD:
_: 'OFF'
SUBSYS:
MOTION:
PRINT:
TRAJECTORY:
FORMAT: DCD_ALIGNED_CELL
EACH:
CELL_OPT: 1
GEO_OPT: 1
MD: 1
RESTART:
BACKUP_COPIES: 0
EACH:
CELL_OPT: 1
GEO_OPT: 1
MD: 1
RESTART_HISTORY:
_: 'OFF'
CELL:
_: 'OFF'
VELOCITIES:
_: 'OFF'
FORCES:
_: 'ON'
STRESS:
_: 'ON'
GEO_OPT:
TYPE: MINIMIZATION #default: MINIMIZATION
OPTIMIZER: BFGS #default: BFGS
MAX_ITER: 1000 #default: 200
MAX_DR: "[bohr] 0.0030" #default: [bohr] 0.0030
RMS_DR: "[bohr] 0.0015" #default: [bohr] 0.0015
MAX_FORCE: "[bohr^-1*hartree] 0.00045" #default: [bohr^-1*hartree] 0.00045
RMS_FORCE: "[bohr^-1*hartree] 0.00030" #default: [bohr^-1*hartree] 0.00030
BFGS:
TRUST_RADIUS: "[angstrom] 0.25" #default: [angstrom] 0.25
CELL_OPT:
TYPE: DIRECT_CELL_OPT #default: DIRECT_CELL_OPT
EXTERNAL_PRESSURE: "[bar] 0.0"
PRESSURE_TOLERANCE: "[bar] 100" #default: 100
KEEP_ANGLES: False
KEEP_SYMMETRY: False
OPTIMIZER: BFGS #default: BFGS
MAX_ITER: 1000 #default: 200
MAX_DR: "[bohr] 0.0030" #default: [bohr] 0.0030
RMS_DR: "[bohr] 0.0015" #default: [bohr] 0.0015
MAX_FORCE: "[bohr^-1*hartree] 0.00045" #default: [bohr^-1*hartree] 0.00045
RMS_FORCE: "[bohr^-1*hartree] 0.00030" #default: [bohr^-1*hartree] 0.00030
BFGS:
TRUST_RADIUS: "[angstrom] 0.5"
2 changes: 1 addition & 1 deletion pyproject.toml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,7 @@ dependencies = [
'aiida-bigdft>=0.2.6',
'aiida-castep>=1.2.0a5',
'aiida-core[atomic_tools]~=1.6',
'aiida-cp2k~=1.3',
'aiida-cp2k~=1.6',
'aiida-fleur>=1.3.0',
'aiida-gaussian',
'aiida-nwchem>=2.1.0',
Expand Down

0 comments on commit 3015631

Please sign in to comment.