Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

SIRIUS-CP2K protocol and workchain #312

Merged
merged 11 commits into from
Aug 30, 2023
119 changes: 73 additions & 46 deletions aiida_common_workflows/workflows/relax/cp2k/generator.py
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
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
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
Loading