#!/usr/bin/env python
# -*- coding: utf-8 -*-
"""
In this module you find the kkrparams class that helps defining the KKR input parameters
Also some defaults for the parameters are defined
"""
#use print('message') instead of print 'message' in python 2.7 as well:
from __future__ import print_function
# redefine raw_input for python 3/2.7 compatilbility
from sys import version_info
if version_info[0] >= 3:
def raw_input(msg):
return input(msg)
__copyright__ = (u"Copyright (c), 2017, Forschungszentrum Jülich GmbH,"
"IAS-1/PGI-1, Germany. All rights reserved.")
__license__ = "MIT license, see LICENSE.txt file"
__version__ = "0.3"
__contributors__ = u"Philipp Rüßmann"
# This defines the default parameters for KKR used in the aiida plugin:
__kkr_default_params__ = {"LMAX": 3, # lmax-cutoff
"INS": 1, # use shape corrections (full potential)
"NSPIN": 2, # spin-polarized calculation (but by default not automatically initialized with external field)
"RMAX": 10., # Madelung sum real-space cutoff
"GMAX": 100. # Madelung sum reciprocal-space cutoff
}
[docs]class kkrparams(object):
"""
Class for creating and handling the parameter input for a KKR calculation
Optional keyword arguments are passed to init and stored in values dictionary.
Example usage:
params = kkrparams(LMAX=3, BRAVAIS=array([[1,0,0], [0,1,0], [0,0,1]]))
Alternatively values can be set afterwards either individually with
params.set_value('LMAX', 3)
or multiple keys at once with
params.set_multiple_values(EMIN=-0.5, EMAX=1)
Other useful functions:
- print the description of a keyword: params.get_description([key]) where [key] is a string for a keyword in params.values
- print a list of mandatory keywords: params.get_all_mandatory()
- print a list of keywords that are set including their value: params.get_set_values()
Note: KKR-units (e.g. atomic units with energy in Ry, length in a_Bohr) are assumed
except for the keys'<RBLEFT>', '<RBRIGHT>', 'ZPERIODL', and 'ZPERIODR' which should be given in Ang. units!
"""
[docs] def __init__(self, **kwargs):
"""
Initialize class instance with containing the attribute values that also have
a format, mandatory flags (defaults for KKRcode, changed for example via params_type='voronoi' keyword) and a description.
"""
if 'params_type' in kwargs:
self.__params_type = kwargs.pop('params_type')
else:
#parameter are set for kkr or voronoi code? (changes mandatory flags)
self.__params_type = 'kkr' #default value, also possible: 'voronoi', 'kkrimp'
#TODO: implement update_mandatory for self.__params_type=='kkrimp'
keyw = self._create_keyword_default_values(**kwargs)
#values of keywords:
self.values = {}
#formatting info
self.__format = {}
#mandatory flag
self._mandatory = {}
# description of each key
self.__description = {}
for key in keyw:
self.values[key] = keyw[key][0]
self.__format[key] = keyw[key][1]
self._mandatory[key] = keyw[key][2]
self.__description[key] = keyw[key][3]
# overwrite mandatory list for voronoi case abd update keyset
if self.__params_type == 'voronoi':
self._update_mandatory_voronoi()
[docs] @classmethod
def get_KKRcalc_parameter_defaults(self, silent=False):
"""
set defaults (defined in header of this file) and returns dict, kkrparams_version
"""
p = kkrparams()
for key, val in __kkr_default_params__.iteritems():
p.set_value(key,val,silent=silent)
return dict(p.get_set_values()), __version__
[docs] def get_dict(self, group=None, subgroup=None):
"""
Returns values dictionary.
Prints values belonging to a certain group only if the 'group' argument
is one of the following: 'lattice', 'chemistry', 'accuracy',
'external fields', 'scf cycle', 'other'
Additionally the subgroups argument allows to print only a subset of
all keys in a certain group. The following subgroups are available:
in 'lattice' group: '2D mode', 'shape functions'
in 'chemistry' group: 'Atom types', 'Exchange-correlation', 'CPA mode',
'2D mode'
in 'accuracy' group: 'Valence energy contour', 'Semicore energy contour',
'CPA mode', 'Screening clusters', 'Radial solver',
'Ewald summation', 'LLoyd'
"""
out_dict = self.values
#check for grouping
group_searchstrings = {'lattice':'Description of lattice',
'chemistry':'Chemistry',
'external fields':'External fields:',
'accuracy':'Accuracy',
'scf cycle':'Self-consistency control:',
'other':['Running and test options', 'Name of potential and shapefun file']}
subgroups_all = {'lattice':['2D mode', 'shape functions'],
'chemistry':['Atom types', 'Exchange-correlation', 'CPA mode', '2D mode'],
'accuracy':['Valence energy contour', 'Semicore energy contour',
'CPA mode', 'Screening clusters', 'Radial solver',
'Ewald summation', 'LLoyd']}
if group in ['lattice', 'chemistry', 'accuracy', 'external fields', 'scf cycle', 'other']:
print('Returning only values belonging to group %s'%group)
tmp_dict = {}
for key in out_dict.keys():
desc = self.__description[key]
key_in_group = False
if group_searchstrings[group] != 'other':
if group_searchstrings[group] in desc:
key_in_group = True
else:
if group_searchstrings[group][0] in desc or group_searchstrings[group][1] in desc:
key_in_group = True
if key_in_group:
tmp_dict[key] = self.values[key]
#check for subgrouping and overwrite tmp_dict accordingly
if group in ['lattice', 'chemistry', 'accuracy']:
if subgroup in subgroups_all[group]:
print('Restrict keys additionally to subgroup %s'%subgroup)
tmp_dict2 = {}
for key in tmp_dict.keys():
desc = self.__description[key]
key_in_group = False
if subgroup in desc:
key_in_group = True
if key_in_group:
tmp_dict2[key] = self.values[key]
tmp_dict = tmp_dict2
# overwrite out_dict with tmp_dict
out_dict = tmp_dict
return out_dict
[docs] def _get_type_from_string(self, fmtstr):
"""Helper function of get_type"""
if 'f' in fmtstr or 'e' in fmtstr:
keytype = float
elif 'i' in fmtstr:
keytype = int
elif 'l' in fmtstr:
keytype = bool
elif 's' in fmtstr:
keytype = str
else:
print('Error: type of keyvalue not found:', fmtstr)
raise TypeError
return keytype
[docs] def get_type(self, key):
"""Extract expected type of 'key' from format info"""
fmtstr = self.__format[key]
# simple format or complex pattern
simplefmt = True
if fmtstr.count('%') > 1:
simplefmt = False
if simplefmt:
keytype = self._get_type_from_string(fmtstr)
else:
fmtlist = fmtstr.replace('\n','').replace(' ','').split('%')[1:]
keytype = []
for fmtstr in fmtlist:
keytype.append(self._get_type_from_string(fmtstr))
return keytype
[docs] def _check_valuetype(self, key):
"""Consistency check if type of value matches expected type from format info"""
# check if entry is numpy array and change to list automatically:
try:
tmpval = self.values[key].flatten().tolist()
except:
tmpval = self.values[key]
tmptype = type(tmpval)
# get type of value
if tmptype == list:
valtype = []
for val in range(len(tmpval)):
valtype.append(type(tmpval[val]))
else:
valtype = tmptype
#print(key, valtype, self.get_type(key))
# check if type matches format info
cmptypes = self.get_type(key)
success = True
#print(key, type(valtype), valtype, cmptypes)
changed_type_automatically = False
if valtype == int and cmptypes == float:
changed_type_automatically = True
self.values[key] = float(self.values[key])
elif type(valtype) == list:
for ival in range(len(valtype)):
if valtype[ival] == int and cmptypes == float:
changed_type_automatically = True
self.values[key][ival] = float(self.values[key][ival])
elif valtype != cmptypes and tmpval is not None:
success = False
print('Error: type of value does not match expected type for ', key, self.values[key], cmptypes)
raise TypeError
if changed_type_automatically:
print('Warning: filling value of "%s" with integer but expects float. Converting automatically and continue'%key)
return success
[docs] def get_value(self, key):
"""Sets value of keyword 'key'"""
if key not in self.values.keys():
print('Error key not found in values dict!')
raise KeyError
else:
# deal with special cases of runopt and testopt (lists of codewords)
if key in ['RUNOPT', 'TESTOPT'] and self.values[key] is None:
return []
else:
return self.values[key]
[docs] def set_value(self, key, value, silent=False):
"""Sets value of keyword 'key'"""
if value is None:
if not silent:
print('Warning setting value None is not permitted!')
print('Use remove_value funciton instead! Ignore keyword {}'.format(key))
else:
self.values[key] = value
self._check_valuetype(key)
[docs] def remove_value(self, key):
"""Removes value of keyword 'key', i.e. resets to None"""
self.values[key] = None
[docs] def set_multiple_values(self, **kwargs):
"""Set multiple values (in example value1 and value2 of keywords 'key1' and 'key2') given as key1=value1, key2=value2"""
for key in kwargs:
key2 = key
if key not in self.values.keys():
key2 = '<'+key+'>'
#print('setting', key2, kwargs[key])
self.set_value(key2, kwargs[key])
[docs] def get_set_values(self):
"""Return a list of all keys/values that are set (i.e. not None)"""
set_values = []
added = 0
for key in self.values.keys():
if self.values[key] is not None:
set_values.append([key, self.values[key]])
added += 1
if added == 0:
print('No values set')
return set_values
[docs] def get_all_mandatory(self):
"""Return a list of mandatory keys"""
self._update_mandatory()
mandatory_list = []
for key in self.values.keys():
if self.is_mandatory(key):
mandatory_list.append(key)
return mandatory_list
[docs] def is_mandatory(self, key):
"""Returns mandatory flag (True/False) for keyword 'key'"""
return self._mandatory[key]
[docs] def get_description(self, key):
"""Returns description of keyword 'key'"""
return self.__description[key]
[docs] def _create_keyword_default_values(self, **kwargs):
"""
Creates KKR inputcard keywords dictionary and fills entry if value is given in **kwargs
entries of keyword dictionary are: 'keyword', [value, format, keyword_mandatory, description]
where
- 'value' can be a single entry or a list of entries
- 'format' contains formatting info
- 'keyword_mandatory' is a logical stating if keyword needs to be defined to run a calculation
- 'description' is a string containgin human redable info about the keyword
"""
default_keywords = dict([# complete list of keywords, detault all that are not mandatory to None
# lattice
('ALATBASIS', [None, '%f', True, 'Description of lattice: Length unit in Bohr radii usually conventional lattice parameter']),
('BRAVAIS', [None, '%f %f %f\n%f %f %f\n%f %f %f', True, 'Description of lattice: Bravais vectors in units of [ALATBASIS]']),
('NAEZ', [None, '%i', True, 'Description of lattice: Number of sites in unit cell']),
('<RBASIS>', [None, '%f %f %f', True, 'Description of lattice: Positions of sites in unit cell']),
('CARTESIAN', [None, '%l', False, 'Description of lattice: Interpret the basis vector coordinates as reduced (w. respect to bravais) or as cartesian (in lattice constant units)']),
('INTERFACE', [None, '%l', False, 'Description of lattice, 2D mode: needs to be TRUE for 2D calculation']),
('<NLBASIS>', [None, '%i', False, 'Description of lattice, 2D mode: Number of basis sites forming the half-infinite lattice to the lower (=left) part of the slab.']),
('<RBLEFT>', [None, '%f %f %f', False, 'Description of lattice, 2D mode: Positions of sites forming the basis sites of the half-infinite lattice to the lower (=left) part of the slab.']),
('ZPERIODL', [None, '%f %f %f', False, 'Description of lattice, 2D mode: Lattice vector describing the periodicity perpendicular to the slab-plane for the half-infinite lattice to the lower (=left) part of the slab (plays the role of the 3rd Bravais vector for this half-infinite lattice). The <RBLEFT> vectors are periodically repeated by the ZPERIODL vector.']),
('<NRBASIS>', [None, '%i', False, 'Description of lattice, 2D mode: Number of basis sites forming the half-infinite lattice to the upper (=right) part of the slab.']),
('<RBRIGHT>', [None, '%f %f %f', False, 'Description of lattice, 2D mode: Positions of sites forming the basis sites of the half-infinite lattice to the upper (=right) part of the slab.']),
('ZPERIODR', [None, '%f %f %f', False, 'Description of lattice, 2D mode: Lattice vector describing the periodicity perpendicular to the slab-plane for the half-infinite lattice to the upper (=right) part of the slab (plays the role of the 3rd Bravais vector for this half-infinite lattice). The <RBRIGHT> vectors are periodically repeated by the ZPERIODR vector.']),
('KSHAPE', [None, '%i', False, 'Description of lattice, shape functions: 0 for ASA ([INS]=0), 2 for full potential ([INS]=1)']),
('<SHAPE>', [None, '%i', False, 'Description of lattice, shape functions: Indexes which shape function from the shape-function file to use in which atom. Default is that each atom has its own shape function.']),
# chemistry
('<ZATOM>', [None, '%f', True, 'Chemistry, Atom types: Nuclear charge per atom. Negative value signals to use value read in from the potential file.']),
('NSPIN', [None, '%i', True, 'Chemistry, Atom types: Number of spin directions in potential. Values 1 or 2']),
('KVREL', [None, '%i', False, 'Chemistry, Atom types: Relativistic treatment of valence electrons. Takes values 0 (Schroedinger), 1 (Scalar relativistic), 2 (Dirac ; works only in ASA mode)']),
('<SOCSCL>', [None, '%f', False, 'Chemistry, Atom types: Spin-orbit coupling scaling per atom. Takes values between 0. (no spin-orbit) and 1. (full spin-orbit). Works only in combination with the Juelich spin orbit solver (runoption NEWSOSOL)']),
('KEXCOR', [None, '%i', False, 'Chemistry, Exchange-correlation: Type of exchange correlation potential. Takes values 0 (LDA, Moruzzi-Janak-Williams), 1 (LDA, von Barth-Hedin), 2 (LDA, Vosko-Wilk-Nussair), 3 (GGA, Perdew-Wang 91), 4 (GGA, PBE), 5 (GGA, PBEsol)']),
('LAMBDA_XC', [None, '%f', False, 'Chemistry, Exchange-correlation: Scale the magnetic part of the xc-potential and energy. Takes values between 0. (fully suppressed magnetisc potential) and 1. (normal magnetic potential).']),
('NAT_LDAU', [None, '%i', False, 'Chemistry, Exchange-correlation: Numer of atoms where LDA+U will be used']),
('LDAU_PARA', [None, '%i %i %f %f %f', False, 'Chemistry, Exchange-correlation: For each atom where LDA+U should be used, the entries are: [atom type] [angular mom. to apply LDA+U] [Ueff] [Jeff] [Eref] where [atom type] is between 1...[NATYP].']),
('KREADLDAU', [None, '%i', False, "Chemistry, Exchange-correlation: Takes values 0 or 1; if [KREADLDAU]=1 then read previously calculated LDA+U matrix elements from file 'ldaupot'."]),
('NATYP', [None, '%i', False, 'Chemistry, CPA mode: Number of atom types; CPA is triggered by setting [NATYP]>[NAEZ].']),
('<SITE>', [None, '%i', False, 'Chemistry, CPA mode: Takes values 1 < [<SITE>] < [NAEZ] Assigns the position (given by [<RBASIS>]) where the atom-dependent read-in potential is situated. E.g., if the 3rd-in-the-row potential should be positioned at the 2nd <RBASIS> vector, then the 3rd entry of the <SITE> list should have the value 2.']),
('<CPA-CONC>', [None, '%f', False, 'Chemistry, CPA mode: Takes values 0. < [<CPA-CONC>] < 1. Assigns the alloy-concentration corresponding to the atom-dependent read-in potential. Together with the variable <SITE>, <CPA-CONC> assigns the number and concentration of the atom-dependent potentials residing at each site form 1 to [NAEZ]. The sum of concentrations at each site should equal 1.']),
('<KAOEZL>', [None, '%i', False, 'Chemistry, 2D mode: Controls the type of t-matrix at the lower (=left) half-crystal sites in case of embedding as these are given in the left-decimation file (i.e., changes the order compared to the one in the left-decimation file).']),
('<KAOEZR>', [None, '%i', False, 'Chemistry, 2D mode: Controls the type of t-matrix at the upper (=right) half-crystal sites in case of embedding as these are given in the right-decimation file (i.e., changes the order compared to the one in the right-decimation file).']),
# external fields
('LINIPOL', [None, '%l', False, 'External fields: If TRUE, triggers an external magn. field per atom in the first iteration.']),
('HFIELD', [None, '%f', False, 'External fields: Value of an external magnetic field in the first iteration. Works only with LINIPOL, XINIPOL']),
('XINIPOL', [None, '%i', False, 'External fields: Integer multiplying the HFIELD per atom']),
('VCONST', [None, '%f', False, 'External fields: Constant potential shift in the first iteration.']),
# accuracy
('LMAX', [None, '%i', True, 'Accuracy: Angular momentum cutoff']),
('BZDIVIDE', [None, '%i %i %i', False, 'Accuracy: Maximal Brillouin zone mesh. Should not violate symmetry (e.g cubic symmetry implies i1=i2=i3; terragonal symmetry in xy implies i1=i2; i1=i2=i3 is always safe.)']),
('EMIN', [None, '%f', False, 'Accuracy, Valence energy contour: Lower value (in Ryd) for the energy contour']),
('EMAX', [None, '%f', False, 'Accuracy, Valence energy contour: Maximum value (in Ryd) for the DOS calculation Controls also [NPT2] in some cases']),
('TEMPR', [None, '%f', False, 'Accuracy, Valence energy contour: Electronic temperature in K.']),
('NPT1', [None, '%i', False, 'Accuracy, Valence energy contour: Number of energies in the 1st part of the rectangular contour ("going up").']),
('NPT2', [None, '%i', False, 'Accuracy, Valence energy contour: Number of energies in the 2nd part of the rectangular contour ("going right").']),
('NPT3', [None, '%i', False, 'Accuracy, Valence energy contour: Number of energies in the 3rd part of the rectangular contour (Fermi smearing part).']),
('NPOL', [None, '%i', False, 'Accuracy, Valence energy contour: Number of Matsubara poles For DOS calculations, set [NPOL]=0']),
('EBOTSEMI', [None, '%f', False, 'Accuracy, Semicore energy contour: Bottom of semicore contour in Ryd.']),
('EMUSEMI', [None, '%f', False, 'Accuracy, Semicore energy contour: Top of semicore contour in Ryd.']),
('TKSEMI', [None, '%f', False, 'Accuracy, Semicore energy contour: "Temperature" in K controlling height of semicore contour.']),
('NPOLSEMI', [None, '%i', False, 'Accuracy, Semicore energy contour: Control of height of semicore contour: Im z = (2 * [NPOLSEMI] * pi * kB * [TKSEMI] ) with kB=0.6333659E-5']),
('N1SEMI', [None, '%i', False, 'Accuracy, Semicore energy contour: Number of energies in first part of semicore contour ("going up").']),
('N2SEMI', [None, '%i', False, 'Accuracy, Semicore energy contour: Number of energies in second part of semicore contour ("going right").']),
('N3SEMI', [None, '%i', False, 'Accuracy, Semicore energy contour: Number of energies in third part of semicore contour ("going down").']),
('FSEMICORE', [None, '%f', False, 'Accuracy, Semicore energy contour: Initial normalization factor for semicore states (approx. 1.).']),
('CPAINFO', [None, '%f %i', False, 'Accuracy, CPA mode: CPA-error max. tolerance and max. number of CPA-cycle iterations.']),
('RCLUSTZ', [None, '%f', False, 'Accuracy, Screening clusters: Radius of screening clusters in units of [ALATBASIS], default is 11 Bohr radii.']),
('RCLUSTXY', [None, '%f', False, 'Accuracy, Screening clusters: If [RCLUSTXY] does not equal [RCLUSTZ] then cylindrical clusters are created with radius [RCLUSTXY] and height [RCLUSTZ].']),
('<RMTREF>', [None, '%f', False, 'Accuracy, Screening clusters: Muffin tin radius in Bohr radii for each site forming screening clusters. Negative value signals automatic calculation by the code.']),
('NLEFTHOS', [None, '%i', False, 'Accuracy, Screening clusters 2D mode: The vectors [<RBLEFT>] are repeated i=1,...,[NLEFTHOS] times, shifted by i*[ZPERIODL], for the later formation of screening clusters.']),
('<RMTREFL>', [None, '%f', False, 'Accuracy, Screening clusters 2D mode: Muffin-tin radius in Bohr radii for each site forming screening clusters in the lower (=left) half-crystal. Negative value signals automatic calculation by the code.']),
('NRIGHTHO', [None, '%i', False, 'Accuracy, Screening clusters 2D mode: The vectors [<RBRIGHT>] are repeated i=1,...,[NRIGHTHO] times, shifted by i*[ZPERIODR], for the later formation of screening clusters.']),
('<RMTREFR>', [None, '%f', False, 'Accuracy, Screening clusters 2D mode: Muffin-tin radius in Bohr radii for each site forming screening clusters in the upper (=right) half-crystal. Negative value signals automatic calculation by the code.']),
('INS', [None, '%i', False, 'Accuracy, Radial solver: Takes values 0 for ASA and 1 for full potential Must be 0 for Munich Dirac solver ([KREL]=2)']),
('ICST', [None, '%i', False, 'Accuracy, Radial solver: Number of iterations in the radial solver']),
('R_LOG', [None, '%f', False, 'Accuracy, Radial solver: Radius up to which log-rule is used for interval width. Used in conjunction with runopt NEWSOSOL']),
('NPAN_LOG', [None, '%i', False, 'Accuracy, Radial solver: Number of intervals from nucleus to [R_LOG] Used in conjunction with runopt NEWSOSOL']),
('NPAN_EQ', [None, '%i', False, 'Accuracy, Radial solver: Number of intervals from [R_LOG] to muffin-tin radius Used in conjunction with runopt NEWSOSOL']),
('NCHEB', [None, '%i', False, 'Accuracy, Radial solver: Number of Chebyshev polynomials per interval Used in conjunction with runopt NEWSOSOL']),
('<FPRADIUS>', [None, '%f', False, 'Accuracy, Radial solver: Full potential limit per atom (in Bohr radii); at points closer to the nucleus, the potential is assumed spherical. Negative values indicate to use values from potential file. Values larger than the muffin tin indicate to use the muffin tin radius.']),
('RMAX', [None, '%f', True, 'Accuracy, Ewald summation for Madelung potential: Max. radius in [ALATBASIS] for real space Ewald sum']),
('GMAX', [None, '%f', True, 'Accuracy, Ewald summation for Madelung potential: Max. radius in 2*pi/[ALATBASIS] for reciprocal space Ewald sum']),
('<LLOYD>', [None, '%i', False, "Accuracy, LLoyd's formula: Set to 1 in order to use Lloyd's formula"]),
('<DELTAE>', [None, '(%f, %f)', False, "Accuracy, LLoyd's formula: Energy difference for derivative calculation in Lloyd's formula"]),
('<TOLRDIF>', [None, '%e', False, 'Accuracy, Virtual atoms: For distance between scattering-centers smaller than [<TOLRDIF>], free GF is set to zero. Units are Bohr radii.']),
('<RMTCORE>', [None, '%f', False, 'Accuracy: Muffin tin radium in Bohr radii for each atom site. This sets the value of RMT used internally in the KKRcode. Needs to be smaller than the touching RMT of the cells. In particular for structure relaxations this should be kept constant.']),
# scf cycle
('NSTEPS', [None, '%i', False, 'Self-consistency control: Max. number of self-consistency iterations. Is reset to 1 in several cases that require only 1 iteration (DOS, Jij, write out GF).']),
('IMIX', [None, '%i', False, "Self-consistency control: Mixing scheme for potential. 0 means straignt (linear) mixing, 3 means Broyden's 1st method, 4 means Broyden's 2nd method, 5 means Anderson's method"]),
('STRMIX', [None, '%f', False, 'Self-consistency control: Linear mixing parameter Set to 0. if [NPOL]=0']),
('ITDBRY', [None, '%i', False, 'Self-consistency control: how many iterations to keep in the Broyden/Anderson mixing scheme.']),
('FCM', [None, '%f', False, 'Self-consistency control: Factor for increased linear mixing of magnetic part of potential compared to non-magnetic part.']),
('BRYMIX', [None, '%f', False, 'Self-consistency control: Parameter for Broyden mixing.']),
('QBOUND', [None, '%e', False, 'Self-consistency control: Lower limit of rms-error in potential to stop iterations.']),
#code options
('RUNOPT', [None, '%s%s%s%s%s%s%s%s', False, 'Running and test options: 8-character keywords in a row without spaces between them']),
('TESTOPT', [None, '%s%s%s%s%s%s%s%s\n%s%s%s%s%s%s%s%s', False, 'Running and test options: optional 8-character keywords in a row without spaces between them plus a secod row of the same.']),
#file names
('FILES', [None, '%s', False, 'Name of potential and shapefun file (list of two strings, empty string will set back to default of the one file that is supposed to be changed)'])
])
for key in kwargs:
key2 = key
if key not in default_keywords.keys():
key2 = '<'+key+'>'
default_keywords[key2][0] = kwargs[key]
return default_keywords
[docs] def _update_mandatory(self):
"""Check if mandatory flags need to be updated if certain keywords are set"""
# initialize all mandatory flags to False and update list afterwards
for key in self.values.keys():
self._mandatory[key] = False
runopts = []
if self.values['RUNOPT'] is not None:
for runopt in self.values['RUNOPT']:
runopts.append(runopt.strip())
#For a KKR calculation these keywords are always mandatory:
mandatory_list = ['ALATBASIS', 'BRAVAIS', 'NAEZ', '<RBASIS>', 'NSPIN', 'LMAX', 'RMAX', 'GMAX', '<ZATOM>']
if self.values['NPOL'] is not None and self.values['NPOL'] != 0:
mandatory_list += ['EMIN']
#Mandatory in 2D
if self.values['INTERFACE']:
mandatory_list += ['<NLBASIS>', '<RBLEFT>', 'ZPERIODL', '<NRBASIS>', '<RBRIGHT>', 'ZPERIODR']
#Mandatory in LDA+U
if 'NAT_LDAU' in self.values.keys() and 'LDAU' in runopts:
mandatory_list += ['NAT_LDAU', 'LDAU_PARA']
#Mandatory in CPA
if self.values['NATYP'] is not None and self.values['NATYP'] > self.values['NAEZ']:
mandatory_list += ['NATYP', '<SITE>', '<CPA-CONC>']
#Mandatory in SEMICORE
if 'EBOTSEMI' in self.values.keys() and 'SEMICORE' in runopts:
mandatory_list += ['EBOTSEMI', 'EMUSEMI', 'TKSEMI', 'NPOLSEMI', 'N1SEMI', 'N2SEMI', 'N3SEMI', 'FSEMICORE']
if self.values['INS'] == 1 and 'WRITEALL' not in runopts:
mandatory_list += ['<SHAPE>']
for key in mandatory_list:
self._mandatory[key] = True
# overwrite if mandatory list needs to be changed (determinded from value of self.__params_type):
if self.__params_type == 'voronoi':
self._update_mandatory_voronoi()
[docs] def _check_mandatory(self):
"""Check if all mandatory keywords are set"""
self._update_mandatory()
if self.__params_type == 'voronoi':
self._update_mandatory_voronoi()
for key in self.values.keys():
if self._mandatory[key] and self.values[key] is None:
print('Error not all mandatory keys are set!')
set_of_mandatory = set(self.get_all_mandatory())
set_of_keys = set([key[0] for key in self.get_set_values()])
print(set_of_mandatory-set_of_keys, 'missing')
raise ValueError("Missing mandatory key(s): {}".format(set_of_mandatory-set_of_keys))
[docs] def _check_array_consistency(self):
"""Check all keys in __listargs if they match their specification (mostly 1D array, except for special cases e.g. <RBASIS>)"""
from numpy import array, ndarray
vec3_entries = ['<RBASIS>', '<RBLEFT>', '<RBRIGHT>', 'ZPERIODL', 'ZPERIODR']
#success = [True]
for key in self.__listargs.keys():
if self.values[key] is not None:
tmpsuccess = True
#print('checking', key, self.values[key], self.__listargs[key])
if type(self.values[key]) not in [list, ndarray]:
self.values[key] = array([self.values[key]])
cmpdims = (self.__listargs[key], )
if key in vec3_entries:
cmpdims = (self.__listargs[key], 3)
# automatically convert if naez==1 and only 1D array is given
if self.__listargs[key] == 1 and len(array(self.values[key]).shape) == 1 and key not in ['ZPERIODL', 'ZPERIODR']:
print('Warning: expected 2D array for %s but got 1D array, converting automatically'%key)
self.values[key] = array([self.values[key]])
tmpdims = array(self.values[key]).shape
if tmpdims[0] != cmpdims[0]:
tmpsuccess = False
if len(tmpdims)==2:
if tmpdims[1] != cmpdims[1]:
tmpsuccess = False
#success.append(tmpsuccess)
if not tmpsuccess:
print('Error: array input not consistent')
print('check consistency:', key, self.values[key], cmpdims, tmpdims, tmpsuccess)
raise TypeError
[docs] def _find_value(self, charkey, txt, line=1, item=1, num=1):
"""
Search charkey in txt and return value string
parameter, input :: charkey string that is search in txt
parameter, input :: txt text that is searched (output of readlines)
parameter, input, optional :: line index in which line to start reading after key was found
parameter, input, optional :: item index which column is read
parameter, input, optional :: num number of column that are read
returns :: valtxt string or list of strings depending on num setting
"""
try:
iline = [ii for ii in range(len(txt)) if charkey in txt[ii]][0]
except IndexError:
iline = None
if iline is not None:
txtline = txt[iline]
chkeq = charkey+'='
if chkeq in txtline:
valtxt = txtline.split(chkeq)[1].split()[item-1:item-1+num]
else:
nextline = txt[iline+line]
startpos = txtline.index(charkey)
valtxt = nextline[startpos:].split()[item-1:item-1+num]
if num == 1:
return valtxt[0]
else:
return valtxt
else:
return None
# redefine _update_mandatory for voronoi code
[docs] def _update_mandatory_voronoi(self):
"""Change mandatory flags to match requirements of voronoi code"""
# initialize all mandatory flags to False and update list afterwards
for key in self.values.keys():
self._mandatory[key] = False
runopts = []
if self.values['RUNOPT'] is not None:
for runopt in self.values['RUNOPT']:
runopts.append(runopt.strip())
#For a KKR calculation these keywords are always mandatory:
mandatory_list = ['ALATBASIS', 'BRAVAIS', 'NAEZ', '<RBASIS>', 'NSPIN', 'LMAX', 'RCLUSTZ', '<ZATOM>']
#Mandatory in 2D
if self.values['INTERFACE']:
mandatory_list += ['<NLBASIS>', '<RBLEFT>', 'ZPERIODL', '<NRBASIS>', '<RBRIGHT>', 'ZPERIODR']
#Mandatory in CPA
if self.values['NATYP'] is not None and self.values['NATYP'] > self.values['NAEZ']:
mandatory_list += ['NATYP', '<SITE>', '<CPA-CONC>']
for key in mandatory_list:
self._mandatory[key] = True
[docs] def get_missing_keys(self, use_aiida=False):
"""Find list of mandatory keys that are not yet set"""
setlist = dict(self.get_set_values()).keys()
manlist = self.get_all_mandatory()
missing = []
autoset_list = ['BRAVAIS', '<RBASIS>', '<ZATOM>', 'ALATBASIS', 'NAEZ', '<SHAPE>', 'EMIN', 'RCLUSTZ']
if self.__params_type == 'voronoi':
autoset_list = ['BRAVAIS', '<RBASIS>', '<ZATOM>', 'ALATBASIS', 'NAEZ']
for key in manlist:
if key not in setlist:
if not use_aiida:
missing.append(key)
else:
if key not in autoset_list:
missing.append(key)
return missing
[docs] def update_to_voronoi(self):
"""
Update parameter settings to match voronoi specification.
Sets self.__params_type and calls _update_mandatory_voronoi()
"""
self.__params_type = 'voronoi'
self._update_mandatory_voronoi()
"""
# tests
if __name__=='__main__':
from numpy import ndarray, array
from aiida_kkr.tools.common_functions import get_Ang2aBohr
p = kkrparams(params_type='kkr')
# automatically read keywords from inpucard
p.read_keywords_from_inputcard(inputcard='/Users/ruess/sourcecodes/aiida/development/calc_import_test/inputcard')
# convert some read-in stuff back from Ang. units to alat units
rbl = p.get_value('<RBLEFT>')
rbr = p.get_value('<RBRIGHT>')
zper_l = p.get_value('ZPERIODL')
zper_r = p.get_value('ZPERIODR')
ang2alat = get_Ang2aBohr()/p.get_value('ALATBASIS')
if rbl is not None: p.set_value('<RBLEFT>', array(rbl)*ang2alat)
if rbr is not None: p.set_value('<RBRIGHT>', array(rbr)*ang2alat)
if zper_l is not None: p.set_value('ZPERIODL', array(zper_l)*ang2alat)
if zper_r is not None: p.set_value('ZPERIODR', array(zper_r)*ang2alat)
# set parameters of expected values manually
p0 = kkrparams(RUNOPT=['xigid-ef','LLOYD', 'ewald2d', 'NEWSOSOL', 'DOS'], TESTOPT=['ie','RMESH','clusters','MPIenerg','fullBZ','DOS'], LMAX=3, NSPIN=2, NATYP=80, NAEZ=80, CARTESIAN=True, ALATBASIS=20.156973053, BRAVAIS=[[0.38437499, 0., 0.], [0.19218749, -0.33287851, 0.], [0.19218749, -0.11095950, 1.]], INTERFACE=True, NRIGHTHO=10, NLEFTHOS=10, NLBASIS=10, NRBASIS=10, ZPERIODL=[-1.92187500000000e-01, 1.10959504859881e-01, -1.00000000000000e+00], ZPERIODR=[1.92187500000000e-01, -1.10959504859881e-01, 1.00000000000000e+00], RCLUSTZ=0.65, RCLUSTXY=0.65, EMIN=-1.2, EMAX=1.2, TEMPR=473., NPOL=7, NPT1=7, NPT2=40, NPT3=6, KSHAPE=2, INS=1, ICST=2, KEXCOR=2, HFIELD=0, VCONST=0, NPAN_LOG=17, NPAN_EQ=7, NCHEB=12, R_LOG=0.8, BZDIVIDE=[40, 40, 1], NSTEPS=500, IMIX=5, STRMIX=0.02, FCM=20., QBOUND=10**-7, BRYMIX=0.02, ITDBRY=30, LINIPOL=False, FILES=['potential', 'shapefun'], RMAX=15., GMAX=900.)
p0.set_value('<ZATOM>', [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 52.0, 0.0, 51.0, 0.0, 52.0, 0.0, 51.0, 0.0, 52.0, 0.0, 52.0, 0.0, 51.0, 0.0, 52.0, 0.0, 51.0, 0.0, 52.0, 0.0, 52.0, 0.0, 51.0, 0.0, 52.0, 0.0, 51.0, 0.0, 52.0, 0.0, 52.0, 0.0, 51.0, 0.0, 52.0, 0.0, 51.0, 0.0, 52.0, 0.0, 52.0, 0.0, 51.0, 0.0, 52.0, 0.0, 51.0, 0.0, 52.0, 0.0, 52.0, 0.0, 51.0, 0.0, 52.0, 0.0, 51.0, 0.0, 52.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0])
p0.set_value('<SHAPE>', [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10])
p0.set_multiple_values(KAOEZR=[i for i in range(1,11)], KAOEZL=[i for i in range(1,11)], KVREL=1, RMTREFL=[2.2671000, 2.2671000, 2.4948000, 2.3562000, 2.3562000, 2.3562000, 2.4948000, 2.2671000, 2.2671000, 2.5740000], RMTREFR=[2.2671000, 2.2671000, 2.4948000, 2.3562000, 2.3562000, 2.3562000, 2.4948000, 2.2671000, 2.2671000, 2.5740000])
p0.set_multiple_values(RMTREF=[2.2671000, 2.2671000, 2.4948000, 2.3562000, 2.3562000, 2.3562000, 2.4948000, 2.2671000, 2.2671000, 2.5740000, 2.2671000, 2.2671000, 2.4948000, 2.3562000, 2.3562000, 2.3562000, 2.4948000, 2.2671000, 2.2671000, 2.5740000, 2.2671000, 2.2671000, 2.4948000, 2.3562000, 2.3562000, 2.3562000, 2.4948000, 2.2671000, 2.2671000, 2.5740000, 2.2671000, 2.2671000, 2.4948000, 2.3562000, 2.3562000, 2.3562000, 2.4948000, 2.2671000, 2.2671000, 2.5740000, 2.2671000, 2.2671000, 2.4948000, 2.3562000, 2.3562000, 2.3562000, 2.4948000, 2.2671000, 2.2671000, 2.5740000, 2.2671000, 2.2671000, 2.4948000, 2.3562000, 2.3562000, 2.3562000, 2.4948000, 2.2671000, 2.2671000, 2.5740000, 2.2671000, 2.2671000, 2.4948000, 2.3562000, 2.3562000, 2.3562000, 2.4948000, 2.2671000, 2.2671000, 2.5740000, 2.2671000, 2.2671000, 2.4948000, 2.3562000, 2.3562000, 2.3562000, 2.4948000, 2.2671000, 2.2671000, 2.5740000])
p0.set_multiple_values(RBLEFT=[[-1.92187500000000e-01, 1.10959504859881e-01, -1.00000000000000e+00], [ 8.32667268468867e-17, 2.77555756156289e-17, -9.49500000000000e-01], [ 1.92187500000000e-01, -1.10959504859881e-01, -8.33000000000000e-01], [ 3.84375000000000e-01, -2.21919009719762e-01, -7.16500000000000e-01], [ 8.32667268468867e-17, 0.00000000000000e+00, -6.33000000000000e-01], [ 1.92187500000000e-01, -1.10959504859881e-01, -5.49500000000000e-01], [ 3.84375000000000e-01, -2.21919009719762e-01, -4.33000000000000e-01], [ 2.77555756156289e-17, 1.38777878078145e-17, -3.16500000000000e-01], [ 1.92187500000000e-01, -1.10959504859881e-01, -2.66000000000000e-01], [ 3.84375000000000e-01, -2.21919009719762e-01, -1.33000000000000e-01]],
RBRIGHT=[[1.53750000000000e+00, -8.87676038879049e-01, 8.00000000000000e+00], [1.72968750000000e+00, -9.98635543738930e-01, 8.05050000000000e+00], [1.92187500000000e+00, -1.10959504859881e+00, 8.16700000000000e+00], [2.11406250000000e+00, -1.22055455345869e+00, 8.28350000000000e+00], [1.72968750000000e+00, -9.98635543738930e-01, 8.36700000000000e+00], [1.92187500000000e+00, -1.10959504859881e+00, 8.45050000000000e+00], [2.11406250000000e+00, -1.22055455345869e+00, 8.56700000000000e+00], [1.72968750000000e+00, -9.98635543738930e-01, 8.68350000000000e+00], [1.92187500000000e+00, -1.10959504859881e+00, 8.73400000000000e+00], [2.11406250000000e+00, -1.22055455345869e+00, 8.86700000000000e+00]],
RBASIS=[[0.0, 0.0, 0.0], [0.1921875, -0.110959504859881, 0.0505000000000001], [0.384375, -0.221919009719762, 0.167], [0.5765625, -0.332878514579644, 0.2835], [0.1921875, -0.110959504859881, 0.367], [0.384375, -0.221919009719762, 0.4505], [0.5765625, -0.332878514579644, 0.567], [0.1921875, -0.110959504859881, 0.6835], [0.384375, -0.221919009719762, 0.734], [0.5765625, -0.332878514579644, 0.867], [0.1921875, -0.110959504859881, 1.0], [0.384375, -0.221919009719762, 1.0505], [0.5765625, -0.332878514579643, 1.167], [0.76875, -0.443838019439525, 1.2835], [0.384375, -0.221919009719762, 1.367], [0.5765625, -0.332878514579643, 1.4505], [0.76875, -0.443838019439525, 1.567], [0.384375, -0.221919009719762, 1.6835], [0.5765625, -0.332878514579643, 1.734], [0.76875, -0.443838019439525, 1.867], [0.384375, -0.221919009719762, 2.0], [0.5765625, -0.332878514579643, 2.0505], [0.76875, -0.443838019439525, 2.167], [0.9609375, -0.554797524299406, 2.2835], [0.5765625, -0.332878514579643, 2.367], [0.76875, -0.443838019439525, 2.4505], [0.9609375, -0.554797524299406, 2.567], [0.5765625, -0.332878514579643, 2.6835], [0.76875, -0.443838019439525, 2.734], [0.9609375, -0.554797524299406, 2.867], [0.5765625, -0.332878514579643, 3.0], [0.76875, -0.443838019439525, 3.0505], [0.9609375, -0.554797524299406, 3.167], [1.153125, -0.665757029159287, 3.2835], [0.76875, -0.443838019439525, 3.367], [0.9609375, -0.554797524299406, 3.4505], [1.153125, -0.665757029159287, 3.567], [0.76875, -0.443838019439525, 3.6835], [0.9609375, -0.554797524299406, 3.734], [1.153125, -0.665757029159287, 3.867], [0.76875, -0.443838019439525, 4.0], [0.9609375, -0.554797524299406, 4.0505], [1.153125, -0.665757029159287, 4.167], [1.3453125, -0.776716534019168, 4.2835], [0.9609375, -0.554797524299406, 4.367], [1.153125, -0.665757029159287, 4.4505], [1.3453125, -0.776716534019168, 4.567], [0.9609375, -0.554797524299406, 4.6835], [1.153125, -0.665757029159287, 4.734], [1.3453125, -0.776716534019168, 4.867], [0.9609375, -0.554797524299406, 5.0], [1.153125, -0.665757029159287, 5.0505], [1.3453125, -0.776716534019168, 5.167], [1.5375, -0.887676038879049, 5.2835], [1.153125, -0.665757029159287, 5.367], [1.3453125, -0.776716534019168, 5.4505], [1.5375, -0.887676038879049, 5.567], [1.153125, -0.665757029159287, 5.6835], [1.3453125, -0.776716534019168, 5.734], [1.5375, -0.887676038879049, 5.867], [1.153125, -0.665757029159287, 6.0], [1.3453125, -0.776716534019168, 6.0505], [1.5375, -0.887676038879049, 6.167], [1.7296875, -0.99863554373893, 6.2835], [1.3453125, -0.776716534019168, 6.367], [1.5375, -0.887676038879049, 6.4505], [1.7296875, -0.99863554373893, 6.567], [1.3453125, -0.776716534019168, 6.6835], [1.5375, -0.887676038879049, 6.734], [1.7296875, -0.99863554373893, 6.867], [1.3453125, -0.776716534019168, 7.0], [1.5375, -0.887676038879049, 7.0505], [1.7296875, -0.99863554373893, 7.167], [1.921875, -1.10959504859881, 7.2835], [1.5375, -0.887676038879049, 7.367], [1.7296875, -0.99863554373893, 7.4505], [1.921875, -1.10959504859881, 7.567], [1.5375, -0.887676038879049, 7.6835], [1.7296875, -0.99863554373893, 7.734], [1.921875, -1.10959504859881, 7.867]])
# check all values
for key in [i[0] for i in p.get_set_values()]:
v = p.get_value(key)
v0 = p0.get_value(key)
if type(v) != list and type(v) != ndarray:
if v!=v0:
print(key, v, v0)
assert v==v0
elif type(v[0]) != str:
if abs(array(v)-array(v0)).max()>=10**-14:
print(key, abs(array(v)-array(v0)).max())
assert abs(array(v)-array(v0)).max()<10**-14
else:
if set(v)-set(v0)!=set():
print(key, set(v)-set(v0))
assert set(v)-set(v0)==set()
"""