Commit 00698fbe authored by Lauri Himanen's avatar Lauri Himanen

Enabled parsing of cell information from input file.

parent bc9a0527
......@@ -9,7 +9,6 @@ metainfo_section_prefix = "x_cp2k_section_input_"
metainfo_data_prefix = "x_cp2k_input_"
#===============================================================================
class CP2KInput(object):
"""The contents of a CP2K simulation including default values and default
units from the version-specific xml file.
......@@ -154,7 +153,6 @@ class CP2KInput(object):
return parameter.lone_value
#===============================================================================
class Section(object):
"""An input section in a CP2K calculation.
"""
......@@ -216,7 +214,6 @@ class Section(object):
return value.upper()
#===============================================================================
class InputObject(object):
"""Base class for all kind of data elements in the CP2K input.
"""
......@@ -231,7 +228,6 @@ class InputObject(object):
self.default_value = None
#===============================================================================
class Keyword(InputObject):
"""Information about a keyword in a CP2K calculation.
"""
......@@ -265,13 +261,14 @@ class Keyword(InputObject):
"""
# Decode the unit and the value if not done before
proper_value = None
if self.default_unit:
if not self.value_no_unit:
if self.default_unit is not None:
if self.value_no_unit is None:
self.decode_cp2k_unit_and_value()
if self.value_no_unit is not None:
proper_value = self.value_no_unit
else:
proper_value = self.value
# if allow_default:
if proper_value is None:
proper_value = self.default_value
......@@ -331,8 +328,8 @@ class Keyword(InputObject):
"""If the units of this value can be changed, return a value and the
unit separately.
"""
if self.default_unit:
if not self.value_no_unit:
if self.default_unit is not None:
if self.value_no_unit is None:
self.decode_cp2k_unit_and_value()
return self.value_no_unit
else:
......@@ -352,22 +349,22 @@ class Keyword(InputObject):
def decode_cp2k_unit_and_value(self):
"""Given a CP2K unit name, decode it as Pint unit definition.
"""
splitted = self.value.split(None, 1)
unit_definition = splitted[0]
if unit_definition.startswith('[') and unit_definition.endswith(']'):
unit_definition = unit_definition[1:-1]
self.unit = CP2KInput.decode_cp2k_unit(self.default_unit)
self.value_no_unit = splitted[1]
elif self.default_unit:
logger.debug("No special unit definition found, returning default unit.")
self.unit = CP2KInput.decode_cp2k_unit(self.default_unit)
self.value_no_unit = self.value
else:
logger.debug("The value has no unit, returning bare value.")
self.value_no_unit = self.value
if self.value is not None:
splitted = self.value.split(None, 1)
unit_definition = splitted[0]
if unit_definition.startswith('[') and unit_definition.endswith(']'):
unit_definition = unit_definition[1:-1]
self.unit = CP2KInput.decode_cp2k_unit(unit_definition)
self.value_no_unit = splitted[1]
elif self.default_unit:
logger.debug("No special unit definition found, returning default unit.")
self.unit = CP2KInput.decode_cp2k_unit(self.default_unit)
self.value_no_unit = self.value
else:
logger.debug("The value has no unit, returning bare value.")
self.value_no_unit = self.value
#===============================================================================
class SectionParameters(InputObject):
"""Section parameters in a CP2K calculation.
......@@ -382,7 +379,6 @@ class SectionParameters(InputObject):
self.lone_keyword_value = lone_keyword_value
#===============================================================================
class DefaultKeyword(InputObject):
"""Default keyword in the CP2K input.
"""
......
......@@ -380,18 +380,19 @@ class CP2KCommonParser(CommonParser):
# correspondent, and push directly to the superBackend to avoid caching
try:
sic_cp2k = section.get_latest_value("self_interaction_correction_method")
sic_map = {
"NO": "",
"AD SIC": "SIC_AD",
"Explicit Orbital SIC": "SIC_EXPLICIT_ORBITALS",
"SPZ/MAURI SIC": "SIC_MAURI_SPZ",
"US/MAURI SIC": "SIC_MAURI_US",
}
sic_nomad = sic_map.get(sic_cp2k)
if sic_nomad is not None:
backend.superBackend.addValue('self_interaction_correction_method', sic_nomad)
else:
logger.warning("Unknown self-interaction correction method used.")
if sic_cp2k is not None:
sic_map = {
"NO": "",
"AD SIC": "SIC_AD",
"Explicit Orbital SIC": "SIC_EXPLICIT_ORBITALS",
"SPZ/MAURI SIC": "SIC_MAURI_SPZ",
"US/MAURI SIC": "SIC_MAURI_US",
}
sic_nomad = sic_map.get(sic_cp2k)
if sic_nomad is not None:
backend.superBackend.addValue('self_interaction_correction_method', sic_nomad)
else:
logger.warning("Unknown self-interaction correction method used.")
except:
pass
......@@ -515,7 +516,6 @@ class CP2KCommonParser(CommonParser):
"""
self.section_system_index = gIndex
self.cache_service.addValue("number_of_atoms")
# self.cache_service.push_array_values("simulation_cell", unit="angstrom")
self.cache_service.addArrayValues("configuration_periodic_dimensions")
self.cache_service.addArrayValues("atom_labels")
......
......@@ -8,6 +8,9 @@ import pickle
import numpy as np
from nomadcore.baseclasses import AbstractBaseParser
from cp2kparser.generic.inputparsing import metainfo_data_prefix, metainfo_section_prefix
from pint import UnitRegistry
import ase
ureg = UnitRegistry()
logger = logging.getLogger("nomad")
......@@ -66,7 +69,6 @@ class CP2KInputParser(AbstractBaseParser):
self.cache_service.add("vel_add_last")
self.cache_service.add("each_geo_opt")
self.cache_service.add("traj_add_last")
# self.cache_service.add("electronic_structure_method")
def parse(self, filepath):
......@@ -80,6 +82,91 @@ class CP2KInputParser(AbstractBaseParser):
# Parse everything in the input to cp2k specific metadata
self.fill_metadata()
# If the print level is low, parse the system geometry from the input
# file. If keywords that manipulate the system are enabled, the
# geometry is not parsed because the manipulations can be hard to
# reproduce.
print_level = self.input_tree.get_keyword("GLOBAL/PRINT_LEVEL")
multiple_unit_cell = self.input_tree.get_keyword("FORCE_EVAL/SUBSYS/MULTIPLE_UNIT_CELL")
top_multiple_unit_cell = self.input_tree.get_keyword("FORCE_EVAL/SUBSYS/TOPOLOGY/MULTIPLE_UNIT_CELL")
centering = self.input_tree.get_section("FORCE_EVAL/SUBSYS/TOPOLOGY/CENTER_COORDINATES")
if multiple_unit_cell is None:
multiple_unit_cell = np.array([1, 1, 1])
if print_level == "LOW" \
and np.array_equal(multiple_unit_cell, [1, 1, 1]) \
and np.array_equal(top_multiple_unit_cell, [1, 1, 1]) \
and not centering.accessed:
cell_section = self.input_tree.get_section("FORCE_EVAL/SUBSYS/CELL")
a_obj = cell_section.get_keyword_object("A")
b_obj = cell_section.get_keyword_object("B")
c_obj = cell_section.get_keyword_object("C")
a_val = a_obj.get_value()
b_val = b_obj.get_value()
c_val = c_obj.get_value()
abc_obj = cell_section.get_keyword_object("ABC")
angles_obj = cell_section.get_keyword_object("ALPHA_BETA_GAMMA")
abc_val = abc_obj.get_value()
angles_val = angles_obj.get_value()
# Cell given as three vectors
if a_val is not None and b_val is not None and c_val is not None:
a_unit = a_obj.get_unit()
b_unit = b_obj.get_unit()
c_unit = c_obj.get_unit()
a = (a_val * ureg(a_unit)).to("angstrom").magnitude
b = (b_val * ureg(b_unit)).to("angstrom").magnitude
c = (c_val * ureg(c_unit)).to("angstrom").magnitude
cell = np.stack((a, b, c), axis=0)
self.cache_service["simulation_cell"] = cell
self.cache_service["lattice_vectors"] = cell
# Cell given as lengths and angles
elif abc_val is not None and angles_val is not None:
abc_unit = abc_obj.get_unit()
angles_unit = angles_obj.get_unit()
abc = (abc_val * ureg(abc_unit)).to("angstrom").magnitude
angles = (angles_val * ureg(angles_unit)).to("degree").magnitude
parameters = np.hstack((abc, angles))
cell = ase.geometry.cellpar_to_cell(parameters)
self.cache_service["simulation_cell"] = cell
self.cache_service["lattice_vectors"] = cell
coord = self.input_tree.get_section("FORCE_EVAL/SUBSYS/COORD")
scaled = coord.get_keyword("SCALED")
coords = coord.default_keyword.value
print(coords)
coord_unit = coord.get_keyword("UNIT")
if len(coords) != 0:
lines = coords.split("\n")
pos = []
labels = []
for line in lines:
if len(line) != 0:
parts = line.split()
i_lab = parts[0]
i_x = float(parts[1])
i_y = float(parts[2])
i_z = float(parts[3])
labels.append(i_lab)
pos.append([i_x, i_y, i_z])
pos = np.array(pos)
coord_unit = self.get_pint_unit_string(coord_unit)
pos = (pos * ureg(coord_unit)).to("angstrom").magnitude
labels = np.array(labels)
if scaled == "F" or scaled == "FALSE":
self.cache_service["atom_positions"] = pos
self.cache_service["atom_labels"] = labels
# self.cache_service["simulation_cell"] = cell
# self.cache_service["lattice_vectors"] = cell
# coordinates = self.input_tree.get_section("FORCE_EVAL/SUBSYS/COORD")
# print(coordinates.get_section_parameter())
# self.cache_service["atom_positions"] = coordinates
# self.cache_service["atom_labels"] = labels
# Parse the used XC_functionals and their parameters
xc = self.input_tree.get_section("FORCE_EVAL/DFT/XC/XC_FUNCTIONAL")
if xc is not None:
......@@ -383,11 +470,13 @@ class CP2KInputParser(AbstractBaseParser):
# Contents (keywords, default keywords)
else:
split = line.split(None, 1)
print(split)
if len(split) <= 1:
keyword_value = ""
else:
keyword_value = split[1]
keyword_name = split[0].upper()
# print(keyword_value)
self.input_tree.set_keyword(path + "/" + keyword_name, keyword_value)
def fill_metadata(self):
......
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
&GLOBAL
PROJECT Si_bulk8
RUN_TYPE ENERGY_FORCE
PRINT_LEVEL LOW
&END GLOBAL
&FORCE_EVAL
METHOD Quickstep
&SUBSYS
&KIND Si
ELEMENT Si
BASIS_SET DZVP-GTH-PADE
POTENTIAL GTH-PADE-q4
&END KIND
&CELL
ABC 5.430697500 5.430697500 5.430697500
ALPHA_BETA_GAMMA [rad] 1.5707963268 1.5707963268 1.5707963268
&END CELL
&COORD
Si 0.000000000 0.000000000 0.000000000
Si 0.000000000 2.715348700 2.715348700
Si 2.715348700 2.715348700 0.000000000
Si 2.715348700 0.000000000 2.715348700
Si 4.073023100 1.357674400 4.073023100
Si 1.357674400 1.357674400 1.357674400
Si 1.357674400 4.073023100 4.073023100
Si 4.073023100 4.073023100 1.357674400
&END COORD
&END SUBSYS
&DFT
BASIS_SET_FILE_NAME BASIS_SET
POTENTIAL_FILE_NAME GTH_POTENTIALS
&QS
EPS_DEFAULT 1.0E-10
&END QS
&MGRID
NGRIDS 4
CUTOFF 300
REL_CUTOFF 60
&END MGRID
&XC
&XC_FUNCTIONAL PADE
&END XC_FUNCTIONAL
&END XC
&SCF
SCF_GUESS ATOMIC
EPS_SCF 1.0E-7
MAX_SCF 300
&DIAGONALIZATION ON
ALGORITHM STANDARD
&END DIAGONALIZATION
&MIXING T
METHOD BROYDEN_MIXING
ALPHA 0.4
NBROYDEN 8
&END MIXING
&END SCF
&END DFT
&PRINT
&FORCES ON
&END FORCES
&END PRINT
&END FORCE_EVAL
This diff is collapsed.
This diff is collapsed.
This source diff could not be displayed because it is too large. You can view the blob instead.
&GLOBAL
PROJECT Si_bulk8
RUN_TYPE ENERGY_FORCE
PRINT_LEVEL LOW
&END GLOBAL
&FORCE_EVAL
METHOD Quickstep
&SUBSYS
&KIND Si
ELEMENT Si
BASIS_SET DZVP-GTH-PADE
POTENTIAL GTH-PADE-q4
&END KIND
&CELL
A 5.430697500 0.000000000 0.000000000
B 0.000000000 5.430697500 0.000000000
C 0.000000000 0.000000000 5.430697500
&END CELL
&COORD
Si 0.000000000 0.000000000 0.000000000
Si 0.000000000 2.715348700 2.715348700
Si 2.715348700 2.715348700 0.000000000
Si 2.715348700 0.000000000 2.715348700
Si 4.073023100 1.357674400 4.073023100
Si 1.357674400 1.357674400 1.357674400
Si 1.357674400 4.073023100 4.073023100
Si 4.073023100 4.073023100 1.357674400
&END COORD
&END SUBSYS
&DFT
BASIS_SET_FILE_NAME BASIS_SET
POTENTIAL_FILE_NAME GTH_POTENTIALS
&QS
EPS_DEFAULT 1.0E-10
&END QS
&MGRID
NGRIDS 4
CUTOFF 300
REL_CUTOFF 60
&END MGRID
&XC
&XC_FUNCTIONAL PADE
&END XC_FUNCTIONAL
&END XC
&SCF
SCF_GUESS ATOMIC
EPS_SCF 1.0E-7
MAX_SCF 300
&DIAGONALIZATION ON
ALGORITHM STANDARD
&END DIAGONALIZATION
&MIXING T
METHOD BROYDEN_MIXING
ALPHA 0.4
NBROYDEN 8
&END MIXING
&END SCF
&END DFT
&PRINT
&FORCES ON
&END FORCES
&END PRINT
&END FORCE_EVAL
This diff is collapsed.
This source diff could not be displayed because it is too large. You can view the blob instead.
This diff is collapsed.
#!/bin/bash
#SBATCH -p test #partition
#SBATCH -t 10:00
#SBATCH --mem-per-cpu=5000
#SBATCH -o sbatch-hsw-%j.out
#SBATCH -n 24
##SBATCH -c 4
##SBATCH -N 1
#SBATCH --constraint=hsw
##SBATCH --array=50,100,150,200,250,300,350,400,450,500
module load cp2k-env/4.1
## use first input you get
input_file="$(ls -1 | grep '\.inp' | head -1)"
output_file="${input_file%inp}out"
## use first xyz-file you get
xyz_file="$(ls -1 | grep '\.xyz' | head -1)"
## define projectname in cp2k inputfile
sed -i "s/myprojectnamemy/${input_file%.inp}/g" $input_file
sed -i "s/myxyzfilemy/${xyz_file%.xyz}/g" $input_file
cp2k_bin=cp2k.popt
work_dir=.
echo "start"
echo $work_dir
echo "arguments to srun" $cp2k_bin $output_file $input_file
## runs the actual calculation
srun $cp2k_bin -o $output_file -i $input_file # option 1
############################################################################
echo "Name of the partition in which the job is running: $SLURM_JOB_PARTITION"
echo "Name of the node running the job script: $SLURMD_NODENAME"
echo " The ID of the job allocation: $SLURM_JOB_ID"
echo " List of nodes allocated to the job: $SLURM_JOB_NODELIST"
echo "finnish"
This diff is collapsed.
2
H 3.16900000 1.83000000 1.51600000
H 3.16900000 1.83000000 0.51600000
----------------------------------------------------------------------------
Setting up new environment, removing all currently loaded modules
----------------------------------------------------------------------------
Loading new modules:
binutils intelmpi/5.1.1 fftw/3.3.4 vmd/1.9.1
gcc/4.9.3 mkl/11.3.0 cp2k/4.1
CP2K version 4.1/taito_gcc revision 17462 loaded
Passes all regression tests
VMD version 1.9.1 is now in use
Due to MODULEPATH changes the following have been reloaded:
1) intelmpi/5.1.1 2) mkl/11.3.0
start
.
arguments to srun cp2k.popt sp.out sp.inp
Name of the partition in which the job is running: test
Name of the node running the job script: c984
The ID of the job allocation: 17813337
List of nodes allocated to the job: c984
finnish
# Projected DOS for atomic kind H at iteration step i = 0, E(Fermi) = 0.032245 a.u.
# MO Eigenvalue [a.u.] Occupation s p
1 -0.281587 1.000000 0.92982909 0.07017091
2 0.060055 0.000000 0.35175960 0.64824040
3 0.257068 0.000000 0.00000000 1.00000000
4 0.257068 0.000000 0.00000000 1.00000000
5 0.724722 0.000000 0.65594992 0.34405008
6 0.821978 0.000000 0.07562262 0.92437738
7 0.966780 0.000000 0.00000000 1.00000000
8 0.966780 0.000000 0.00000000 1.00000000
9 2.763651 0.000000 0.99454829 0.00545171
10 4.930257 0.000000 0.99229048 0.00770952
# Projected DOS for atomic kind H at iteration step i = 0, E(Fermi) = 0.032245 a.u.
# MO Eigenvalue [a.u.] Occupation s p
1 -0.281587 1.000000 0.92982909 0.07017091
2 0.060055 0.000000 0.35175960 0.64824040
3 0.257068 0.000000 0.00000000 1.00000000
4 0.257068 0.000000 0.00000000 1.00000000
5 0.724722 0.000000 0.65594992 0.34405008
6 0.821978 0.000000 0.07562262 0.92437738
7 0.966780 0.000000 0.00000000 1.00000000
8 0.966780 0.000000 0.00000000 1.00000000
9 2.763651 0.000000 0.99454829 0.00545171
10 4.930257 0.000000 0.99229048 0.00770952
ATOMIC FORCES in [a.u.]
# Atom Kind Element X Y Z
1 1 H 0.00000000 0.00000000 -0.11814428
2 1 H 0.00000000 0.00000000 0.11814428
SUM OF ATOMIC FORCES 0.00000000 0.00000000 0.00000000 0.00000000
LOWDIN POPULATION ANALYSIS
# Atom Element Kind Atomic population (alpha,beta) Net charge Spin moment
1 H 1 0.500000 0.500000 0.000000 0.000000
2 H 1 0.500000 0.500000 -0.000000 0.000000
# Total charge and spin 1.000000 1.000000 -0.000000 0.000000
############ General ############################
&GLOBAL
RUN_TYPE ENERGY_FORCE # ENERGY_FORCE, GEO_OPT
PREFERRED_FFT_LIBRARY FFTSG
EXTENDED_FFT_LENGTHS
PRINT_LEVEL LOW
PROJECT sp # Name of calculation
&END GLOBAL
##################################################
&FORCE_EVAL
METHOD QuickStep
&PRINT
&FORCES
FILENAME forces
&EACH
MD 1
&END EACH
&END FORCES
&END PRINT
&DFT
BASIS_SET_FILE_NAME ./BASIS_SET
POTENTIAL_FILE_NAME ./POTENTIAL
UKS
# CHARGE 0
&MGRID
CUTOFF 550 # plane-wave cutoff for the charge density [Rydbergs]
NGRIDS 5
REL_CUTOFF 60
&END MGRID
# &POISSON
# PERIODIC XY
# POISSON_SOLVER PERIODIC
# &END POISSON
&QS
METHOD GPW
EPS_DEFAULT 1.0E-12
MAP_CONSISTENT
EXTRAPOLATION ASPC
EXTRAPOLATION_ORDER 4
&END QS
### Property Section ###
&PRINT
&LOWDIN ON
FILENAME lowdin
&EACH
MD 1
&END EACH
&END LOWDIN
### CUBE FILES ###
# &MO_CUBES ! Controls which MOs are written to cube-files.
# NHOMO 1
# NLUMO 1
# &EACH MD 1
# &END EACH
# &END MO_CUBES
# &TOT_DENSITY_CUBE ON
# &EACH MD 1
# &END EACH
# &END TOT_DENSITY_CUBE
# &E_DENSITY_CUBE
# &EACH MD 1
# &END EACH
# &END E_DENSITY_CUBE
##################
&PDOS
APPEND ! Controls which MOs are included in the pdos-files.
NLUMO -1
COMPONENTS .FALSE.
FILENAME dosfile
LOG_PRINT_KEY TRUE
&END PDOS
&END PRINT
### ############### ###
&SCF
SCF_GUESS RESTART
EPS_LUMOS 0.000001 # for MO_CUBES
MAX_ITER_LUMOS 10000 # for MO_CUBES
EPS_SCF 1.0E-6 # convergence threshold for total energy
CHOLESKY INVERSE_DBCSR
### OT METHOD ###
# MAX_SCF 20
# &OT ON
# MINIMIZER DIIS
# PRECONDITIONER FULL_ALL # FULL_SINGLE_INVERSE
# #ENERGY_GAP 0.1
# ALGORITHM IRAC
# &END OT
# &OUTER_SCF ON
# EPS_SCF 1.0E-6
# MAX_SCF 100