Source code for nexus.qmcpack_input
##################################################################
## (c) Copyright 2015- by Jaron T. Krogel ##
##################################################################
#====================================================================#
# qmcpack_input.py #
# Supports I/O and manipulation of QMCPACK's xml input file. #
# #
# Content summary: #
# QmcpackInput #
# SimulationInput class for QMCPACK. #
# Represents the QMCPACK input file as a nested collection of #
# objects mirroring the structure of the XML input. #
# XML attributes and <parameter/> elements are joined into #
# a common keyword representation. #
# Full tranlations of test QMCPACK input files into executable #
# Python code can be found at the end of this file for #
# reference. #
# #
# BundledQmcpackInput #
# Class represents QMCPACK input when provided in 'bundled' #
# form, i.e. an input file containing a list of XML input #
# files. #
# A BundledQmcpackInput object contains many QmcpackInput #
# objects. #
# #
# QmcpackInputTemplate #
# Class supports Nexus interaction with text input files #
# provided by users as 'templates'. #
# Users can mark keywords in a template input file, generate #
# variations on that input file through this class, and #
# then use the resulting input files in Nexus workflows. #
# #
# generate_qmcpack_input #
# User-facing function to compose specific QMCPACK input files. #
# Though any QMCPACK input file can be composed directly with #
# Python code, only a more limited selection can be #
# generated directly. #
# Calls many other functions to generate specific XML elements: #
# generate_simulationcell #
# generate_particlesets #
# generate_sposets #
# generate_sposet_builder #
# generate_bspline_builder #
# generate_heg_builder #
# generate_determinantset #
# generate_determinantset_old #
# generate_hamiltonian #
# generate_jastrows #
# generate_jastrow #
# generate_jastrow1 #
# generate_bspline_jastrow2 #
# generate_pade_jastrow2 #
# generate_jastrow2 #
# generate_jastrow3 #
# generate_kspace_jastrow #
# generate_opt #
# generate_opts #
# #
# QIxml, Names #
# Class represents a generic XML element. #
# Supports read/write and setting default values. #
# Derived classes represent specific elements and contain #
# a keyword specification for each element. #
# Support for new XML elements is enabled by creating #
# a corresponding class with the allowed keywords, etc., #
# and then adding it to the global 'classes' list below. #
# See classes simulation, project, application, random, include, #
# mcwalkerset, qmcsystem, simulationcell, particleset, group, #
# sposet, bspline_builder, heg_builder, composite_builder, #
# wavefunction, determinantset, basisset, grid, atomicbasisset,#
# basisgroup, radfunc, slaterdeterminant, determinant, #
# occupation, multideterminant, detlist, ci, jastrow1, #
# jastrow2, jastrow3, correlation, var, coefficients, #
# coefficient, hamiltonian, coulomb, constant, pseudopotential,#
# pseudo, mpc, localenergy, energydensity, reference_points, #
# spacegrid, origin, axis, chiesa, density, nearestneighbors, #
# neighbor_trace, dm1b, spindensity, magnetizationdensity, #
# structurefactor, init, scalar_traces, array_traces, #
# particle_traces, traces, loop, linear, cslinear, vmc, dmc. #
# #
# QIxmlFactory #
# Class supports comprehension of XML elements that share the #
# same XML tag (e.g. <pairpot/>), but have many different #
# and distinct realizations (e.g. coulomb, mpc, etc.). #
# See QIxmlFactory instances pairpot, estimator, #
# sposet_builder, jastrow, and qmc. #
# #
# collection #
# Container class representing an ordered set of plural XML #
# elements named according to an XML attribute (e.g. named #
# particlesets). #
# XML elements are listed by name in the collection to allow #
# intuitive interactive navigation by the user. #
# #
# Param (single global instance is param) #
# Handles all reads/writes of attributes and text contents #
# of XML elements. #
# Automatically distinguishes string, int, float, and array #
# data of any of these types. #
# #
# Note: those seeking to add new XML elements should be aware #
# of the global data structures described below. #
# #
# classes #
# Global list of QIxml classes. #
# Each new QIxml class must be added here. #
# #
# types #
# Global dict of elements that are actually simple types (to be #
# interpreted the same way as <parameter/> elements) and also #
# factory instances. #
# #
# plurals #
# Global obj of elements that can be plural (i.e. multiple can #
# be present with the same tag. A QmcpackInput object will #
# contain collection instances containing the related XML #
# elements. Each collection is named according to the keys in #
# plurals. #
# #
# Names.set_expanded_names #
# Function to specify mappings from all-lowercase names used #
# in QmcpackInput objects to the names expected by QMCPACK in #
# the actual XML input file. QMCPACK does not follow any #
# consistent naming convention (camel-case, hyphenation, and #
# separation via underscore are all present in names) and some #
# names are case-sensitive while others are not. This function #
# allows developers to cope with this diversity upon write, #
# while preserving a uniform representation (all lowercase) in #
# QmcpackInput objects. #
# #
#====================================================================#
import os
from pathlib import Path
import inspect
import keyword
import numpy as np
from .numpy_extensions import reshape_inplace
from .xmlreader import XMLreader, XMLelement
from .developer import DevBase, obj, hidden, error
from .periodic_table import Elements
from .structure import Structure, Jellium, get_kpath
from .physical_system import PhysicalSystem
from .simulation import SimulationInput, SimulationInputTemplate
from .pwscf_input import array_to_string as pwscf_array_string
from .utilities import path_string
from . import numpy_extensions as npe
yesno_dict = {True:'yes' ,False:'no'}
truefalse_dict = {True:'true',False:'false'}
onezero_dict = {True:'1' ,False:'0'}
boolmap={'yes':True,'no':False,'true':True,'false':False,'1':True,'0':False}
#end try
#end def is_int
#end try
#end def is_float
[docs]
def is_array(var,type):
try:
if isinstance(var,str):
np.array(var.split(),type)
else:
np.array(var,type)
#end if
return True
except ValueError:
return False
#end try
#end def is_float_array
[docs]
def attribute_to_value(attr):
if is_int(attr):
val = int(attr)
elif is_float(attr):
val = float(attr)
elif is_array(attr,int):
val = np.array(attr.split(),int)
if val.size==9:
npe.reshape_inplace(val, (3, 3))
#end if
elif is_array(attr,float):
val = np.array(attr.split(),float)
else:
val = attr
#end if
return val
#end def attribute_to_value
#local write types
#end def yesno
#end if
#end def yesnostr
#end def onezero
#end def onezero
[docs]
def render_bool(var,T,F):
if isinstance(var,bool) or var in (1,0):
if var:
return T
else:
return F
#end if
elif var in (T,F):
return var
else:
error('Invalid QMCPACK input encountered.\nUser provided an invalid value of "{}" when yes/no was expected.\nValid options are: "{}", "{}", True, False, 1, 0'.format(var,T,F))
#end if
#end def render_bool
bool_write_types = set([yesno,onezero,truefalse])
[docs]
class QIobj(DevBase):
afqmc_mode = False
# user settings
permissive_read = False
permissive_write = False
permissive_init = False
[docs]
@staticmethod
def settings(
permissive_read = False,
permissive_write = False,
permissive_init = False,
):
QIobj.permissive_read = permissive_read
QIobj.permissive_write = permissive_write
QIobj.permissive_init = permissive_init
#end def settings
#end class QIobj
#end class meta
[docs]
class section(QIobj):
def __init__(self,*args,**kwargs):
self.args = args
self.kwargs = kwargs
#end def __init__
#end class section
[docs]
class collection(hidden):
def __init__(self,*elements):
hidden.__init__(self)
if len(elements)==1 and isinstance(elements[0],list):
elements = elements[0]
elif len(elements)==1 and isinstance(elements[0],collection):
elements = elements[0].__dict__.values()
#end if
self.hidden().order = []
for element in elements:
self.add(element)
#end for
#end def __init__
def __setitem__(self,name,value):
#self.error('elements can only be set via the add function')
self.add(value,key=name)
#end def __setitem__
def __delitem__(self,name):
#self.error('elements can only be deleted via the remove function')
self.remove(name)
#end def __delitem__
def __setattr__(self,name,value):
#self.error('elements can only be set via the add function')
self.add(value,key=name)
#end def __setattr__
def __delattr__(self,name):
#self.error('elements can only be deleted via the remove function')
self.remove(name)
#end def __delattr__
[docs]
def add(self,element,strict=True,key=None):
if not isinstance(element,QIxml):
self.error('collection cannot be formed\nadd attempted for non QIxml element\ntype received: {0}'.format(element.__class__.__name__))
#end if
keyin = key
key = None
public = self.public()
identifier = element.identifier
missing_identifier = False
if element.tag not in plurals_inv and element.collection_id is None:
self.error('collection cannot be formed\n encountered non-plural element\n element class: {0}\n element tag: {1}\n tags allowed in a collection: {2}'.format(element.__class__.__name__,element.tag,sorted(plurals_inv.keys())))
elif identifier is None:
key = len(public)
elif isinstance(identifier,str):
if identifier in element:
key = element[identifier]
else:
missing_identifier = True
#end if
else:
key = ''
for ident in identifier:
if ident in element:
key+=element[ident]
#end if
#end for
missing_identifier = key==''
#end if
if missing_identifier:
key = len(public)
#end if
if keyin is not None and not isinstance(key,int) and keyin.lower()!=key.lower():
self.error('attempted to add key with incorrect name\nrequested key: {0}\ncorrect key: {1}'.format(keyin,key))
#end if
#if key in public:
# self.error('attempted to add duplicate key to collection: {0}\n keys present: {1}'.format(key,sorted(public.keys())))
##end if
public[key] = element
self.hidden().order.append(key)
return True
#end def add
[docs]
def remove(self,key):
public = self.public()
if key in public:
del public[key]
self.hidden().order.remove(key)
else:
raise KeyError
#end if
#end def remove
[docs]
def get_single(self,preference=None):
if len(self)>0:
if preference is not None and preference in self:
return self[preference]
else:
return self.list()[0]
#end if
else:
#return self
return None
#end if
#end def get_single
[docs]
def list(self):
lst = []
for key in self.hidden().order:
lst.append(self[key])
#end for
return lst
#end def list
[docs]
def pairlist(self):
pairs = []
for key in self.hidden().order:
pairs.append((key,self[key]))
#end for
return pairs
#end def pairlist
#end class collection
#end def make_collection
[docs]
class classcollection(QIobj):
def __init__(self,*classes):
if len(classes)==1 and isinstance(classes[0],list):
classes = classes[0]
#end if
self.classes = classes
#end def __init__
#end class classcollection
[docs]
class QmcpackInputCollections(QIobj):
[docs]
def add(self,element):
if element.tag in plurals_inv:
cname = plurals_inv[element.tag]
if cname not in self:
coll = collection()
success = coll.add(element,strict=False)
if success:
self[cname] = coll
#end if
else:
self[cname].add(element,strict=False)
#end if
#end if
#end def add
[docs]
def get(self,cname,label=None):
v = None
if cname in self:
if label is None:
v = self[cname]
elif label in self[cname]:
v = self[cname][label]
#end if
#end if
return v
#end def get
#end class QmcpackInputCollections
QIcollections = QmcpackInputCollections()
[docs]
class Names(QIobj):
condensed_names = obj()
expanded_names = None
rsqmc_expanded_names = None
afqmc_expanded_names = None
escape_names = set(keyword.kwlist+['write'])
escaped_names = list(escape_names)
for i in range(len(escaped_names)):
escaped_names[i]+='_'
#end for
escaped_names = set(escaped_names)
[docs]
@staticmethod
def set_expanded_names(**kwargs):
exnames = obj(**kwargs)
Names.expanded_names = exnames
Names.rsqmc_expanded_names = exnames
#end def set_expanded_names
[docs]
@staticmethod
def set_afqmc_expanded_names(**kwargs):
Names.afqmc_expanded_names = obj(**kwargs)
#end def set_afqmc_expanded_names
[docs]
@staticmethod
def use_rsqmc_expanded_names():
Names.expanded_names = Names.rsqmc_expanded_names
#end def use_rsqmc_expanded_names
[docs]
@staticmethod
def use_afqmc_expanded_names():
Names.expanded_names = Names.afqmc_expanded_names
#end def use_afqmc_expanded_names
[docs]
def expand_name(self,condensed):
expanded = condensed
cname = self.condense_name(condensed)
if cname in self.escaped_names:
cname = cname[:-1]
expanded = cname
#end if
if cname in self.expanded_names:
expanded = self.expanded_names[cname]
#end if
return expanded
#end def expand_name
[docs]
def condense_name(self,expanded):
condensed = expanded
condensed = condensed.replace('___','_').replace('__','_')
condensed = condensed.replace('-','_').replace(' ','_')
condensed = condensed.lower()
if condensed in self.escape_names:
condensed += '_'
#end if
self.condensed_names[expanded]=condensed
return condensed
#end def condense_name
[docs]
def condense_names(self,*namelists):
out = []
for namelist in namelists:
exp = obj()
for expanded in namelist:
condensed = self.condense_name(expanded)
exp[condensed]=expanded
#end for
out.append(exp)
#end for
return out
#end def condense_names
[docs]
def condensed_name_report(self):
print()
print('Condensed Name Report:')
print('----------------------')
keylist = np.array(list(self.condensed_names.keys()))
order = np.array(list(self.condensed_names.values())).argsort()
keylist = keylist[order]
for expanded in keylist:
condensed = self.condensed_names[expanded]
if expanded!=condensed:
print(" {0:15} = '{1}'".format(condensed,expanded))
#end if
#end for
print()
print()
#end def condensed_name_report
#end class Names
[docs]
class QIxml(Names):
[docs]
def init_from_args(self,args):
print()
print('In init from args (not implemented).')
print('Possible reasons for incorrect entry: ')
print(' Is xml element {0} meant to be plural?'.format(self.__class__.__name__))
print(' If so, add it to the plurals object.')
print()
print('Arguments received:')
print(args)
print()
self.not_implemented()
#end def init_from_args
[docs]
@classmethod
def init_class(cls):
cls.class_set_optional(
tag = cls.__name__,
identifier = None,
attributes = [],
elements = [],
text = None,
parameters = [],
attribs = [],
costs = [],
h5tags = [],
types = obj(),
write_types = obj(),
attr_types = None,
precision = None,
defaults = obj(),
collection_id = None,
exp_names = None,
)
for v in ['attributes','elements','parameters','attribs','costs','h5tags']:
names = cls.class_get(v)
for i in range(len(names)):
if names[i] in cls.escape_names:
names[i]+='_'
#end if
#end for
#end for
cls.params = cls.parameters + cls.attribs + cls.costs + cls.h5tags
cls.plurals_inv = obj()
for e in cls.elements:
if e in plurals_inv:
cls.plurals_inv[e] = plurals_inv[e]
#end if
#end for
cls.plurals = cls.plurals_inv.inverse()
if cls.exp_names is not None:
cls.expanded_names = obj(Names.expanded_names,cls.exp_names)
#end if
#end def init_class
[docs]
def write(self,indent_level=0,pad=' ',first=False):
param.set_precision(self.get_precision())
if not QIobj.permissive_write:
self.check_junk(exit=True)
#end if
indent = indent_level*pad
ip = indent+pad
ipp= ip+pad
expanded_tag = self.expand_name(self.tag)
c = indent+'<'+expanded_tag
for a in self.attributes:
if a in self:
val = self[a]
if isinstance(val,str):
val = self.expand_name(val)
#end if
c += ' '+self.expand_name(a)+'='
if a in self.write_types:
c += '"'+self.write_types[a](val)+'"'
else:
c += '"'+param.write(val)+'"'
#end if
#end if
#end for
#if first:
# c+=' xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xsi:noNamespaceSchemaLocation="http://www.mcc.uiuc.edu/qmc/schema/molecu.xsd"'
##end if
#no_contents = len(set(self.keys())-set(self.elements)-set(self.plurals.keys()))==0
no_contents = len(set(self.keys())-set(self.attributes))==0
if no_contents:
c += '/>\n'
else:
c += '>\n'
for v in self.h5tags:
if v in self:
if v in self.write_types:
write_type = self.write_types[v]
else:
write_type = None
#end if
c += param.write(self[v],name=self.expand_name(v),tag='h5tag',mode='elem',pad=ip,write_type=write_type)
#end if
#end for
for v in self.costs:
if v in self:
c += param.write(self[v],name=self.expand_name(v),tag='cost',mode='elem',pad=ip)
#end if
#end for
for p in self.parameters:
if p in self:
if p in self.write_types:
write_type = self.write_types[p]
else:
write_type = None
#end if
c += param.write(self[p],name=self.expand_name(p),mode='elem',pad=ip,write_type=write_type)
#end if
#end for
for a in self.attribs:
if a in self:
if a in self.write_types:
write_type = self.write_types[a]
else:
write_type = None
#end if
c += param.write(self[a],name=self.expand_name(a),tag='attrib',mode='elem',pad=ip,write_type=write_type)
#end if
#end for
elements = self.elements
if self.afqmc_mode and 'afqmc_order' in self.__class__.__dict__:
elements = self.afqmc_order
#end if
for e in elements:
if e in self:
elem = self[e]
if isinstance(elem,QIxml):
c += elem.write(indent_level+1)
else:
begin = '<'+e+'>'
contents = param.write(elem)
end = '</'+e+'>'
if contents.strip()=='':
c += ip+begin+end+'\n'
else:
c += ip+begin+'\n'
c += ipp+contents+'\n'
c += ip+end+'\n'
#end if
#end if
elif e in plurals_inv and plurals_inv[e] in self:
coll = self[plurals_inv[e]]
if not isinstance(coll,collection):
self.error('write failed\n element {0} is not a collection\n contents of element {0}:\n{1}'.format(plurals_inv[e],str(coll)))
#end if
for instance in coll.list():
c += instance.write(indent_level+1)
#end for
#end if
#end for
if self.text is not None:
c = c.rstrip('\n')
c+=param.write(self[self.text],mode='elem',pad=ip,tag=None,normal_elem=True)
#end if
c+=indent+'</'+expanded_tag+'>\n'
#end if
param.reset_precision()
return c
#end def write
def __init__(self,*args,**kwargs):
if Param.metadata is None:
Param.metadata = meta()
#end if
if len(args)==1:
a = args[0]
if isinstance(a,XMLelement):
self.init_from_xml(a)
elif isinstance(a,section):
self.init_from_inputs(a.args,a.kwargs)
elif isinstance(a,self.__class__):
self.transfer_from(a)
else:
self.init_from_inputs(args,kwargs)
#end if
else:
self.init_from_inputs(args,kwargs)
#end if
QIcollections.add(self)
#end def __init__
[docs]
def init_from_xml(self,xml):
al,el = self.condense_names(xml._attributes.keys(),xml._elements.keys())
xa,sa = set(al.keys()) , set(self.attributes)
attr = xa & sa
junk = xa-attr
junk_elem = []
for e,ecap in el.items():
value = xml._elements[ecap]
if (isinstance(value,list) or isinstance(value,tuple)) and e in self.plurals_inv.keys():
if e not in types:
self.error('input element "{}" is unknown'.format(e))
#end if
p = self.plurals_inv[e]
plist = []
for instance in value:
plist.append(types[e](instance))
#end for
self[p] = make_collection(plist)
elif e in self.elements:
if e not in types:
self.error('input element "{}" is unknown'.format(e))
#end if
self[e] = types[e](value)
elif e in ['parameter','attrib','cost','h5tag']:
if isinstance(value,XMLelement):
value = [value]
#end if
for p in value:
name = self.condense_name(p.name)
if name in self.params:
self[name] = param(p)
else:
junk_elem.append(name)
#end if
#end for
else:
junk_elem.append(e)
#end if
#end for
junk = junk | set(junk_elem)
if QmcpackInput.profile_collection is not None:
self.collect_profile(xml,al,el,junk)
#end for
if not QIobj.permissive_read:
self.check_junk(junk)
#end if
if self.attr_types is not None:
typed_attr = attr & set(self.attr_types.keys())
attr -= typed_attr
for a in typed_attr:
self[a] = self.attr_types[a](xml._attributes[al[a]])
#end for
#end if
for a in attr:
if a in self.write_types and self.write_types[a] in bool_write_types:
aval = xml._attributes[al[a]]
if aval in boolmap:
self[a] = boolmap[aval]
else:
self.error('{0} is not a valid value for boolean attribute {1}\n valid values are: {2}'.format(aval,a,boolmap.keys()))
#end if
else:
self[a] = attribute_to_value(xml._attributes[al[a]])
#end if
#end for
if self.text is not None:
self[self.text] = param(xml)
#end if
#end def init_from_xml
[docs]
def init_from_inputs(self,args,kwargs):
if len(args)>0:
if len(args)==1 and isinstance(args[0],self.__class__):
self.transfer_from(args[0])
elif len(args)==1 and isinstance(args[0],dict):
self.init_from_kwargs(args[0])
else:
self.init_from_args(args)
#end if
#end if
self.init_from_kwargs(kwargs)
#end def init_from_inputs
[docs]
def init_from_kwargs(self,kwargs):
ks=[]
kmap = dict()
for key,val in kwargs.items():
ckey = self.condense_name(key)
ks.append(ckey)
kmap[ckey] = val
#end for
ks = set(ks)
kwargs = kmap
h5tags = ks & set(self.h5tags)
costs = ks & set(self.costs)
parameters = ks & set(self.parameters)
attribs = ks & set(self.attribs)
attr = ks & set(self.attributes)
elem = ks & set(self.elements)
plur = ks & set(self.plurals.keys())
if self.text is not None:
text = ks & set([self.text])
else:
text = set()
#end if
if not QIobj.permissive_init:
junk = ks -attr -elem -plur -h5tags -costs -parameters -attribs -text
self.check_junk(junk,exit=True)
#end if
for v in h5tags:
self[v] = param(kwargs[v])
#end for
for v in costs:
self[v] = param(kwargs[v])
#end for
for v in parameters:
self[v] = param(kwargs[v])
#end for
for v in attribs:
self[v] = param(kwargs[v])
#end for
for a in attr:
self[a] = kwargs[a]
#end for
for e in elem:
self[e] = types[e](kwargs[e])
#end for
for p in plur:
e = self.plurals[p]
kwcoll = kwargs[p]
if isinstance(kwcoll,collection):
plist = kwcoll.list()
elif isinstance(kwcoll,(list,tuple)):
plist = kwcoll
else:
self.error('init failed\n encountered non-list collection')
#end if
ilist = []
for instance in plist:
ilist.append(types[e](instance))
#end for
self[p] = make_collection(ilist)
#end for
for t in text:
self[t] = kwargs[t]
#end for
#end def init_from_kwargs
[docs]
def incorporate_defaults(self,elements=False,overwrite=False,propagate=True):
for name,value in self.defaults.items():
defval=None
if isinstance(value,classcollection):
if elements:
coll=[]
for cls in value.classes:
ins = cls()
ins.incorporate_defaults()
coll.append(ins)
#end for
defval = make_collection(coll)
#end if
elif inspect.isclass(value):
if elements:
defval = value()
#end if
elif inspect.isfunction(value):
defval = value()
else:
defval = value
#end if
if defval is not None:
if overwrite or name not in self:
self[name] = defval
#end if
#end if
#end for
if propagate:
for name,value in self.items():
if isinstance(value,QIxml):
value.incorporate_defaults(elements,overwrite)
elif isinstance(value,collection):
for v in value:
if isinstance(v,QIxml):
v.incorporate_defaults(elements,overwrite)
#end if
#end for
#end if
#end for
#end if
#end def incorporate_defaults
[docs]
def check_junk(self,junk=None,exit=False):
if junk is None:
ks = set(self.keys())
h5tags = ks & set(self.h5tags)
costs = ks & set(self.costs)
parameters = ks & set(self.parameters)
attribs = ks & set(self.attribs)
attr = ks & set(self.attributes)
elem = ks & set(self.elements)
plur = ks & set(self.plurals.keys())
if self.text is not None:
text = ks & set([self.text])
else:
text = set()
#end if
junk = ks -attr -elem -plur -h5tags -costs -parameters -attribs -text
#end if
if len(junk)>0:
oname = ''
if self.tag!=self.__class__.__name__:
oname = ' ('+self.__class__.__name__+')'
#end if
msg = '{0}{1} does not have the following attributes/elements:\n'.format(self.tag,oname)
for jname in junk:
msg+=' '+jname+'\n'
#end for
#if QmcpackInput.profile_collection is None:
# self.error(msg,'QmcpackInput',exit=exit,trace=exit)
##end if
#print(obj(dict(self.__class__.__dict__)))
self.error(msg,'QmcpackInput',exit=exit,trace=exit)
#end if
#end def check_junk
[docs]
def collect_profile(self,xml,al,el,junk):
attributes = obj(**al)
parameters = obj()
elements = obj()
for e,ecap in el.items():
if e=='parameter':
parameters[e]=ecap
else:
elements[e]=ecap
#end if
#end for
profile = obj(
attributes = attributes,
parameters = parameters,
elements = elements,
xml = xml,
junk = junk
)
#xml = xml.copy()
xname = xml._name
if xname[-1:].isdigit():
xname = xname[:-1]
elif xname[-2:].isdigit():
xname = xname[:-2]
#end if
#xml._name = xname
if len(profile.junk)>0:
print(' '+xname+' (found '+str(junk)+')')
for sector in 'attributes elements'.split():
missing = []
for n in profile.junk:
if n in profile[sector]:
missing.append(profile[sector][n])
#end if
#end for
if len(missing)>0:
ms= ' '+sector+':'
for m in missing:
ms+=' '+m
#end for
print(ms)
#end if
#end for
if 'parameter' in profile.xml:
params = obj()
for p in profile.xml.parameter:
params[p.name.lower()] = p.text.strip()
#end for
missing = []
for n in profile.junk:
if n in params:
missing.append(n)
#end if
#end for
if len(missing)>0:
ms= ' parameters:'
for m in missing:
ms+=' '+m
#end for
print(ms)
#end if
#end if
if junk!=set(['analysis']) and junk!=set(['ratio']) and junk!=set(['randmo']) and junk!=set(['printeloc', 'source']) and junk!=set(['warmup_steps']) and junk!=set(['sposet_collection']) and junk!=set(['eigensolve', 'atom']) and junk!=set(['maxweight', 'reweightedvariance', 'unreweightedvariance', 'energy', 'exp0', 'stabilizerscale', 'minmethod', 'alloweddifference', 'stepsize', 'beta', 'minwalkers', 'nstabilizers', 'bigchange', 'usebuffer']) and junk!=set(['loop2']) and junk!=set(['random']) and junk!=set(['max_steps']):
exit()
#end if
#end if
pc = QmcpackInput.profile_collection
if xname not in pc:
pc[xname] = obj()
#end if
pc[xname].append(profile)
#end def collect_profile
#end def get_single
[docs]
def get(self,names,namedict=None,host=False,root=True):
if namedict is None:
namedict = {}
#end if
if isinstance(names,str):
names = [names]
#end if
if root and not host:
if self.identifier is not None and self.identifier in self:
identity = self[self.identifier]
else:
identity = None
#end if
for name in names:
if name==self.tag:
namedict[name]=self
elif name==identity:
namedict[name]=self
#end if
#end for
#end if
for name in names:
loc = None
if name in self:
loc = name
elif name in plurals_inv and plurals_inv[name] in self:
loc = plurals_inv[name]
#end if
name_absent = name not in namedict
not_element = False
if not name_absent:
not_xml = not isinstance(namedict[name],QIxml)
not_coll = not isinstance(namedict[name],collection)
not_element = not_xml and not_coll
#end if
if loc is not None and (name_absent or not_element):
if host:
namedict[name] = self
else:
namedict[name] = self[loc]
#end if
#end if
#end for
for name,value in self.items():
if isinstance(value,QIxml):
value.get(names,namedict,host,root=False)
elif isinstance(value,collection):
for n,v in value.items():
name_absent = n not in namedict
not_element = False
if not name_absent:
not_xml = not isinstance(namedict[n],QIxml)
not_coll = not isinstance(namedict[n],collection)
not_element = not_xml and not_coll
#end if
if n in names and (name_absent or not_element):
if host:
namedict[n] = value
else:
namedict[n] = v
#end if
#end if
if isinstance(v,QIxml):
v.get(names,namedict,host,root=False)
#end if
#end if
#end if
#end for
if root:
namelist = []
for name in names:
if name in namedict:
namelist.append(namedict[name])
else:
namelist.append(None)
#end if
#end for
if len(namelist)==1:
return namelist[0]
else:
return namelist
#end if
#end if
#end def get
[docs]
def remove(self,*names):
if len(names)==1 and not isinstance(names[0],str):
names = names[0]
#end if
remove = []
for name in names:
attempt = True
if name in self:
rname = name
elif name in plurals_inv and plurals_inv[name] in self:
rname = plurals_inv[name]
else:
attempt = False
#end if
if attempt:
val = self[rname]
if isinstance(val,QIxml) or isinstance(val,collection):
remove.append(rname)
#end if
#end if
#end for
for name in remove:
del self[name]
#end for
for name,value in self.items():
if isinstance(value,QIxml):
value.remove(*names)
elif isinstance(value,collection):
for element in value:
if isinstance(element,QIxml):
element.remove(*names)
#end if
#end for
#end if
#end for
#end def remove
[docs]
def assign(self,**kwargs):
for var,vnew in kwargs.items():
if var in self:
val = self[var]
not_coll = not isinstance(val,collection)
not_xml = not isinstance(val,QIxml)
not_arr = not isinstance(val,np.ndarray)
if not_coll and not_xml and not_arr:
self[var] = vnew
#end if
#end if
#end for
for vname,val in self.items():
if isinstance(val,QIxml):
val.assign(**kwargs)
elif isinstance(val,collection):
for v in val:
if isinstance(v,QIxml):
v.assign(**kwargs)
#end if
#end for
#end if
#end for
#end def assign
[docs]
def replace(self,*args,**kwargs):
if len(args)==2 and isinstance(args[0],str) and isinstance(args[1],str):
vold,vnew = args
args = [(vold,vnew)]
#end for
for valpair in args:
vold,vnew = valpair
for var,val in self.items():
not_coll = not isinstance(val,collection)
not_xml = not isinstance(val,QIxml)
not_arr = not isinstance(val,np.ndarray)
if not_coll and not_xml and not_arr and val==vold:
self[var] = vnew
#end if
#end for
#end for
for var,valpair in kwargs.items():
vold,vnew = valpair
if var in self:
val = self[var]
if vold is None:
self[var] = vnew
else:
not_coll = not isinstance(val,collection)
not_xml = not isinstance(val,QIxml)
not_arr = not isinstance(val,np.ndarray)
if not_coll and not_xml and not_arr and val==vold:
self[var] = vnew
#end if
#end if
#end if
#end for
for vname,val in self.items():
if isinstance(val,QIxml):
val.replace(*args,**kwargs)
elif isinstance(val,collection):
for v in val:
if isinstance(v,QIxml):
v.replace(*args,**kwargs)
#end if
#end for
#end if
#end for
#end def replace
[docs]
def combine(self,other):
#elemental combine only
for name,element in other.items():
plural = isinstance(element,collection)
single = isinstance(element,QIxml)
if single or plural:
elem = []
single_name = None
plural_name = None
if single:
elem.append(element)
single_name = name
if name in plurals_inv:
plural_name = plurals_inv[name]
#end if
else:
elem.extend(element.values())
plural_name = name
single_name = plurals[name]
#end if
if single_name in self:
elem.append(self[single_name])
del self[single_name]
elif plural_name is not None and plural_name in self:
elem.append(self[plural_name])
del self[plural_name]
#end if
if len(elem)==1:
self[single_name]=elem[0]
elif plural_name is None:
self.error('attempting to combine non-plural elements: '+single_name)
else:
self[plural_name] = make_collection(elem)
#end if
#end if
#end for
#end def combine
[docs]
def move(self,**elemdests):
names = list(elemdests.keys())
hosts = self.get_host(names)
dests = self.get(list(elemdests.values()))
if len(names)==1:
hosts = [hosts]
dests = [dests]
#end if
for i in range(len(names)):
name = names[i]
host = hosts[i]
dest = dests[i]
if host is not None and dest is not None and id(host)!=id(dest):
if name not in host:
name = plurals_inv[name]
#end if
dest[name] = host[name]
del host[name]
#end if
#end for
#end def move
[docs]
def pluralize(self):
make_plural = []
for name,value in self.items():
if isinstance(value,QIxml):
if name in plurals_inv:
make_plural.append(name)
#end if
value.pluralize()
elif isinstance(value,collection):
if name in plurals_inv:
make_plural.append(name)
#end if
for v in value:
if isinstance(v,QIxml):
v.pluralize()
#end if
#end for
#end if
#end for
for name in make_plural:
value = self[name]
del self[name]
plural_name = plurals_inv[name]
self[plural_name] = make_collection([value])
#end for
#end def pluralize
[docs]
def difference(self,other,root=True):
if root:
q1 = self.copy()
q2 = other.copy()
else:
q1 = self
q2 = other
#end if
if q1.__class__!=q2.__class__:
different = True
diff = None
d1 = q1
d2 = q2
else:
cls = q1.__class__
s1 = set(q1.keys())
s2 = set(q2.keys())
shared = s1 & s2
unique1 = s1 - s2
unique2 = s2 - s1
different = len(unique1)>0 or len(unique2)>0
diff = cls()
d1 = cls()
d2 = cls()
d1.transfer_from(q1,unique1)
d2.transfer_from(q2,unique2)
for k in shared:
value1 = q1[k]
value2 = q2[k]
is_coll1 = isinstance(value1,collection)
is_coll2 = isinstance(value2,collection)
is_qxml1 = isinstance(value1,QIxml)
is_qxml2 = isinstance(value2,QIxml)
if is_coll1!=is_coll2 or is_qxml1!=is_qxml2:
self.error('values for '+k+' are of inconsistent types\n difference could not be taken')
#end if
if is_qxml1 and is_qxml2:
kdifferent,kdiff,kd1,kd2 = value1.difference(value2,root=False)
elif is_coll1 and is_coll2:
ks1 = set(value1.keys())
ks2 = set(value2.keys())
kshared = ks1 & ks2
kunique1 = ks1 - ks2
kunique2 = ks2 - ks1
kdifferent = len(kunique1)>0 or len(kunique2)>0
kd1 = collection()
kd2 = collection()
kd1.transfer_from(value1,kunique1)
kd2.transfer_from(value2,kunique2)
kdiff = collection()
for kk in kshared:
v1 = value1[kk]
v2 = value2[kk]
if isinstance(v1,QIxml) and isinstance(v2,QIxml):
kkdifferent,kkdiff,kkd1,kkd2 = v1.difference(v2,root=False)
kdifferent = kdifferent or kkdifferent
if kkdiff is not None:
kdiff[kk]=kkdiff
#end if
if kkd1 is not None:
kd1[kk]=kkd1
#end if
if kkd2 is not None:
kd2[kk]=kkd2
#end if
#end if
#end for
else:
if isinstance(value1,np.ndarray):
a1 = value1.ravel()
else:
a1 = np.array([value1])
#end if
if isinstance(value2,np.ndarray):
a2 = value2.ravel()
else:
a2 = np.array([value2])
#end if
if len(a1)!=len(a2):
kdifferent = True
elif len(a1)==0:
kdifferent = False
elif (isinstance(a1[0],float) or isinstance(a2[0],float)) and not (isinstance(a1[0],str) or isinstance(a2[0],str)):
kdifferent = np.abs(a1-a2).max()/max(1e-99,np.abs(a1).max(),np.abs(a2).max()) > 1e-6
else:
kdifferent = not (a1==a2).all()
#end if
if kdifferent:
kdiff = (value1,value2)
kd1 = value1
kd2 = value2
else:
kdiff = None
kd1 = None
kd2 = None
#end if
#end if
different = different or kdifferent
if kdiff is not None:
diff[k] = kdiff
#end if
if kd1 is not None:
d1[k] = kd1
#end if
if kd2 is not None:
d2[k] = kd2
#end if
#end for
#end if
if root:
if diff is not None:
diff.remove_empty()
#end if
d1.remove_empty()
d2.remove_empty()
#end if
return different,diff,d1,d2
#end def difference
[docs]
def remove_empty(self):
names = list(self.keys())
for name in names:
value = self[name]
if isinstance(value,QIxml):
value.remove_empty()
if len(value)==0:
del self[name]
#end if
elif isinstance(value,collection):
ns = list(value.keys())
for n in ns:
v = value[n]
if isinstance(v,QIxml):
v.remove_empty()
if len(v)==0:
del value[n]
#end if
#end if
#end for
if len(value)==0:
del self[name]
#end if
#end if
#end for
#end def remove_empty
#end def get_host
#end def get_precision
#end class QIxml
[docs]
class QIxmlFactory(Names):
def __init__(self,name,types,typekey='',typeindex=-1,typekey2='',default=None):
self.name = name
self.types = types
self.typekey = typekey
self.typekey2 = typekey2
self.typeindex = typeindex
self.default = default
#end def __init__
def __call__(self,*args,**kwargs):
#emulate QIxml.__init__
#get the value of the typekey
a = args
kw = kwargs
found_type = False
if len(args)>0:
v = args[0]
if isinstance(v,XMLelement):
kw = v._attributes
elif isinstance(v,section):
a = v.args
kw = v.kwargs
elif isinstance(v,tuple(self.types.values())):
found_type = True
type = v.__class__.__name__
#end if
#end if
if not found_type:
if self.typekey in kw.keys():
type = kw[self.typekey]
elif self.typekey2 in kw.keys():
type = kw[self.typekey2]
elif self.default is not None:
type = self.default
elif self.typeindex==-1:
self.error('QMCPACK input file is misformatted\ncannot identify type for <{0}/> element\nwith contents:\n{1}\nplease find the XML element matching this description in the input file to identify the problem\nmost likely, it is missing attributes "{2}" or "{3}"'.format(self.name,str(v).rstrip(),self.typekey,self.typekey2))
else:
type = a[self.typeindex]
#end if
#end if
type = self.condense_name(type)
if type in self.types:
return self.types[type](*args,**kwargs)
else:
msg = self.name+' factory is not aware of the following subtype:\n'
msg+= ' '+type+'\n'
self.error(msg,exit=False,trace=False)
#end if
#end def __call__
#end def init_class
#end class QIxmlFactory
[docs]
class Param(Names):
metadata = None
def __init__(self):
self.reset_precision()
#end def __init__
#end def reset_precision
[docs]
def set_precision(self,precision):
if precision is None:
self.reset_precision()
elif not isinstance(precision,str):
self.error('attempted to set precision with non-string: {0}'.format(precision))
else:
self.precision = precision
self.prec_format = '{0:'+precision+'}'
#end if
#end def set_precision
def __call__(self,*args,**kwargs):
if len(args)==0:
self.error('no arguments provided, should have received one XMLelement')
elif not isinstance(args[0],XMLelement):
return args[0]
#self.error('first argument is not an XMLelement')
#end if
return self.read(args[0])
#end def __call__
[docs]
def read(self,xml):
val = None
attr = set(xml._attributes.keys())
other_attr = attr-set(['name'])
if 'name' in attr and len(other_attr)>0:
oa = obj()
for a in other_attr:
oa[a] = xml._attributes[a]
#end for
self.metadata[xml.name] = oa
#end if
if 'text' in xml:
text = xml.text.strip()
# nothing in text
if len(text)==0:
return text
# scalar value
if ' ' not in text:
try:
val = int(text)
except ValueError:
try:
val = float(text)
except ValueError:
val = text
return val
# array value
tokens = []
rowlens = []
for line in text.splitlines():
if len(line.strip())==0:
continue
t = line.split()
rowlens.append(len(t))
tokens.extend(t)
try:
val = np.array(tokens,dtype=int)
except ValueError:
try:
val = np.array(tokens,dtype=float)
except ValueError:
val = np.array(tokens,dtype=object)
if len(set(rowlens))==1 and len(rowlens)>1:
# rows have identical size: 2d
reshape_inplace(val,len(rowlens),rowlens[0])
if val.size==1:
self.error('scalar value interpreted as an array. This is a developer error.')
if val is None:
val = ''
return val
#end def read
[docs]
def write(self,value,mode='attr',tag='parameter',name=None,pad=' ',write_type=None,normal_elem=False):
c = ''
attr_mode = mode=='attr'
elem_mode = mode=='elem'
if not attr_mode and not elem_mode:
self.error(mode+' is not a valid mode. Options are attr,elem.')
#end if
if isinstance(value,list) or isinstance(value,tuple):
value = np.array(value)
#end if
if attr_mode:
if isinstance(value,np.ndarray):
arr = value.ravel()
for v in arr:
c+=self.write_val(v)+' '
#end for
c=c[:-1]
else:
c = self.write_val(value)
#end if
elif elem_mode:
c+=pad
is_array = isinstance(value,np.ndarray)
is_single = not (is_array and value.size>1)
if tag is not None:
if is_single:
max_len = 20
rem_len = max(0,max_len-len(name))
else:
rem_len = 0
#end if
other=''
if name in self.metadata:
for a,v in self.metadata[name].items():
other +=' '+self.expand_name(a)+'="'+self.write_val(v)+'"'
#end for
#end if
c+='<'+tag+' name="'+name+'"'+other+rem_len*' '+'>'
pp = pad+' '
else:
pp = pad
#end if
if is_array:
if normal_elem:
c+='\n'
#end if
if tag is not None:
c+='\n'
#end if
ndim = len(value.shape)
if ndim==1:
line_len = 70
if tag is not None:
c+=pp
#end if
line = ''
for v in value:
line+=self.write_val(v)+' '
if len(line)>line_len:
c+=line+'\n'
line = ''
#end if
#end for
if len(line)>0:
c+=line
#end if
c=c[:-1]+'\n'
elif ndim==2:
nrows,ncols = value.shape
fmt=pp
if value.dtype == np.float64:
if self.precision is None:
vfmt = ':16.8f' # must have 8 digits of post decimal accuracy to meet qmcpack tolerance standards
#vfmt = ':16.8e'
else:
vfmt = ': '+self.precision
#end if
else:
vfmt = ''
#end if
for nc in range(ncols):
fmt+='{'+str(nc)+vfmt+'} '
#end for
fmt = fmt[:-2]+'\n'
for nr in range(nrows):
c+=fmt.format(*value[nr])
#end for
else:
self.error('only 1 and 2 dimensional arrays are supported for xml formatting.\n Received '+ndim+' dimensional array.')
#end if
else:
if write_type is not None:
val = write_type(value)
else:
val = value
#end if
#c += ' '+str(val)
c += ' {0:<10}'.format(self.write_val(val))
#end if
if tag is not None:
c+=pad+'</'+tag+'>\n'
#end if
#end if
return c
#end def write
[docs]
def write_val(self,val):
if self.precision is not None and isinstance(val,float):
return self.prec_format.format(val)
else:
return str(val)
#end if
#end def write_val
#end def init_class
#end class Param
param = Param()
class simulation(QIxml):
# afqmc
attributes = ['method']
# rsqmc
elements = ['project','random','include','qmcsystem','particleset',
'wavefunction','hamiltonian','init','traces',
'qmc','loop','mcwalkerset','cmc']+\
['afqmcinfo','walkerset','propagator','execute'] # afqmc
afqmc_order = ['project','random','afqmcinfo','hamiltonian',
'wavefunction','walkerset','propagator','execute']
write_types = obj(random=yesno)
#end class simulation
class project(QIxml):
attributes = ['id','series']
parameters = ['driver_version','maxcpusecs','max_seconds']
elements = ['application','host','date','user']
#end class project
class application(QIxml):
attributes = ['name','role','class','version']
#end class application
class host(QIxml):
text = 'value'
#end class host
class date(QIxml):
text = 'value'
#end class date
class user(QIxml):
text = 'value'
#end class user
class random(QIxml):
attributes = ['seed','parallel']
write_types= obj(parallel=truefalse)
#end class random
class include(QIxml):
attributes = ['href']
#end def include
class mcwalkerset(QIxml):
attributes = ['fileroot','version','collected','node','nprocs','href','target','file','walkers']
write_types = obj(collected=yesno)
#end class mcwalkerset
class qmcsystem(QIxml):
attributes = ['dim'] #,'wavefunction','hamiltonian'] # breaks QmcpackInput
elements = ['simulationcell','particleset','wavefunction','hamiltonian','random','init','mcwalkerset','estimators']
#end class qmcsystem
class simulationcell(QIxml):
attributes = ['name','tilematrix']
parameters = ['lattice','reciprocal','bconds','lr_dim_cutoff','lr_tol','lr_handler','rs','nparticles','scale','uc_grid']
#end class simulationcell
class particleset(QIxml):
attributes = ['name','size','random','random_source','randomsrc','charge','source','spinor']
elements = ['group','simulationcell']
attribs = ['ionid','position']
write_types= obj(random=yesno,spinor=yesno)
identifier = 'name'
#end class particleset
class group(QIxml):
attributes = ['name','size','mass'] # mass attr and param, bad bad bad!!!
parameters = ['charge','valence','atomicnumber','mass','lmax',
'cutoff_radius','spline_radius','spline_npoints']
attribs = ['position']
identifier = 'name'
#end class group
class sposet(QIxml):
attributes = ['basisset','type','name','group','size',
'index_min','index_max','energy_min','energy_max',
'spindataset','cuspinfo','sort','gpu','href','twistnum',
'gs_sposet','basis_sposet','same_k','frequency','mass',
'source','version','precision','tilematrix',
'meshfactor']
elements = ['occupation','coefficient','coefficients']
text = None
identifier = 'name'
#end class sposet
class rotated_sposet(QIxml):
attributes = ['name']
elements = ['sposet']
identifier = 'name'
#end class rotated_sposet
class bspline_builder(QIxml):
tag = 'sposet_builder'
identifier = 'type'
attributes = ['type','href','sort','tilematrix','twistnum','twist','source',
'version','meshfactor','gpu','transform','precision','truncate',
'lr_dim_cutoff','shell','randomize','key','buffer','rmax_core','dilation','tag','hybridrep','gpusharing']
elements = ['sposet','rotated_sposet']
write_types = obj(gpu=yesno,sort=onezero,transform=yesno,truncate=yesno,randomize=truefalse,hybridrep=yesno,gpusharing=yesno)
#end class bspline_builder
class heg_builder(QIxml):
tag = 'sposet_builder'
identifier = 'type'
attributes = ['type','twist']
elements = ['sposet']
#end class heg_builder
class molecular_orbital_builder(QIxml):
tag = 'sposet_builder'
identifier = 'type'
attributes = ['name','type','transform','source','cuspcorrection','href']
elements = ['basisset','sposet']
elements = ['basisset','sposet','rotated_sposet']
write_types = obj(transform=yesno,cuspcorrection=yesno)
#end class molecular_orbital_builder
class composite_builder(QIxml):
tag = 'sposet_builder'
identifier = 'type'
attributes = ['type']
elements = ['sposet']
#end class composite_builder
sposet_builder = QIxmlFactory(
name = 'sposet_builder',
types = dict(bspline=bspline_builder,
einspline=bspline_builder,
heg=heg_builder,
composite=composite_builder,
molecularorbital = molecular_orbital_builder),
typekey = 'type'
)
sposet_collection = QIxmlFactory(
name = 'sposet_collection',
types = dict(bspline=bspline_builder,
einspline=bspline_builder,
heg=heg_builder,
composite=composite_builder,
molecularorbital = molecular_orbital_builder),
typekey = 'type'
)
class wavefunction(QIxml):
# rsqmc afqmc
attributes = ['name','target','id','ref']+['info','type']
# afqmc
parameters = ['filetype','filename','cutoff']
elements = ['sposet_builder','determinantset','jastrow','override_variational_parameters','sposet_collection']
identifier = 'name','id'
#end class wavefunction
class determinantset(QIxml):
attributes = ['type','href','sort','tilematrix','twistnum','twist','source','version','meshfactor','gpu','transform','precision','truncate','lr_dim_cutoff','shell','randomize','key','rmax_core','dilation','name','cuspcorrection','tiling','usegrid','meshspacing','shell2','src','buffer','bconds','keyword','hybridrep','pbcimages','gpusharing']
elements = ['basisset','sposet','slaterdeterminant','multideterminant','spline','backflow','cubicgrid']
h5tags = ['twistindex','twistangle','rcut']
write_types = obj(gpu=yesno,sort=onezero,transform=yesno,truncate=yesno,randomize=truefalse,cuspcorrection=yesno,usegrid=yesno,gpusharing=yesno)
#end class determinantset
class spline(QIxml):
attributes = ['method']
elements = ['grid']
#end class spline
class cubicgrid(QIxml):
attributes = ['method']
elements = ['grid']
#end class cubicgrid
class basisset(QIxml):
attributes = ['ecut','name','ref','type','source','transform','key']
elements = ['grid','atomicbasisset']
write_types = obj(transform=yesno)
#end class basisset
class grid(QIxml):
attributes = ['dir','npts','closed','type','ri','rf','rc','step']
#identifier = 'dir'
#end class grid
class atomicbasisset(QIxml):
attributes = ['type','elementtype','expandylm','href','normalized','name','angular']
elements = ['grid','basisgroup']
identifier = 'elementtype'
write_types= obj(#expandylm=yesno,
normalized=yesno)
#end class atomicbasisset
class basisgroup(QIxml):
attributes = ['rid','ds','n','l','m','zeta','type','s','imin','source']
parameters = ['b']
elements = ['radfunc']
#identifier = 'rid'
#end class basisgroup
class radfunc(QIxml):
attributes = ['exponent','node','contraction','id','type']
precision = '16.12e'
#end class radfunc
class slaterdeterminant(QIxml):
attributes = ['optimize','delay_rank','gpu','matrix_inverter','batch']
elements = ['determinant']
write_types = obj(optimize=yesno,gpu=yesno,batch=yesno)
#end class slaterdeterminant
class determinant(QIxml):
attributes = ['id','group','sposet','size','ref','spin','href','orbitals','spindataset','name','cuspinfo','debug']
elements = ['occupation','coefficient']
identifier = 'id'
write_types = obj(debug=yesno)
#end class determinant
class occupation(QIxml):
attributes = ['mode','spindataset','size','pairs','format']
text = 'contents'
#end class occupation
class multideterminant(QIxml):
attributes = ['optimize','spo_up','spo_dn']
elements = ['detlist']
#end class multideterminant
class detlist(QIxml):
attributes = ['size','type','nca','ncb','nea','neb','nstates','cutoff','ext_level','href','optimize']
elements = ['ci','csf']
#end class detlist
class ci(QIxml):
attributes = ['id','coeff','qc_coeff','alpha','beta']
#identifier = 'id'
attr_types = obj(alpha=str,beta=str)
precision = '16.12e'
#end class ci
class csf(QIxml):
attributes = ['id','exctlvl','coeff','coeff_real','coeff_imag','qchem_coeff','occ']
elements = ['det']
attr_types = obj(occ=str)
#end class csf
class det(QIxml):
attributes = ['id','coeff','alpha','beta']
attr_types = obj(alpha=str,beta=str)
#end class det
class backflow(QIxml):
attributes = ['optimize']
elements = ['transformation']
write_types = obj(optimize=yesno)
#end class backflow
class transformation(QIxml):
attributes = ['name','type','function','source']
elements = ['correlation']
identifier = 'name'
#end class transformation
class jastrow1(QIxml):
tag = 'jastrow'
attributes = ['type','name','function','source','print','spin','transform']
elements = ['correlation','distancetable','grid']
identifier = 'name'
write_types = obj(print=yesno,spin=yesno,transform=yesno)
#end class jastrow1
class jastrow2(QIxml):
tag = 'jastrow'
attributes = ['type','name','function','print','spin','init','kc','transform','source','optimize']
elements = ['correlation','distancetable','basisset','grid','basisgroup']
parameters = ['b','longrange']
identifier = 'name'
write_types = obj(print=yesno,transform=yesno,optimize=yesno)
#end class jastrow2
class jastrow3(QIxml):
tag = 'jastrow'
attributes = ['type','name','function','print','source']
elements = ['correlation']
identifier = 'name'
write_types = obj(print=yesno)
#end class jastrow3
class kspace_jastrow(QIxml):
tag = 'jastrow'
attributes = ['type','name','source']
elements = ['correlation']
identifier = 'name'
write_types = obj(optimize=yesno)
#end class kspace_jastrow
class rpa_jastrow(QIxml):
tag = 'jastrow'
attributes = ['type','name','source','function','kc']
parameters = ['longrange']
identifier = 'name'
write_types = obj(longrange=yesno)
#end class rpa_jastrow
class correlation(QIxml):
attributes = ['elementtype','speciesa','speciesb','size','ispecies','especies',
'especies1','especies2','isize','esize','rcut','cusp','pairtype',
'kc','type','symmetry','cutoff','spindependent','dimension','init',
'species']
parameters = ['a','b','c','d']
elements = ['coefficients','var','coefficient']
identifier = 'speciesa','speciesb','elementtype','especies1','especies2','ispecies'
write_types = obj(init=yesno)
#end class correlation
class var(QIxml):
attributes = ['id','name','optimize']
text = 'value'
identifier = 'id'
write_types=obj(optimize=yesno)
#end class var
class coefficients(QIxml):
attributes = ['id','type','optimize','state','size','cusp','rcut']
text = 'coeff'
write_types= obj(optimize=yesno)
exp_names = obj(array='Array')
#end class coefficients
class coefficient(QIxml): # this is bad!!! coefficients/coefficient
attributes = ['id','type','size','dataset','spindataset']
text = 'coeff'
precision = '16.12e'
#end class coefficient
class distancetable(QIxml):
attributes = ['source','target']
#end class distancetable
jastrow = QIxmlFactory(
name = 'jastrow',
types = dict(one_body=jastrow1,two_body=jastrow2,jastrow1=jastrow1,jastrow2=jastrow2,eei=jastrow3,jastrow3=jastrow3,kspace=kspace_jastrow,kspace_jastrow=kspace_jastrow,rpa=rpa_jastrow,rpa_jastrow=rpa_jastrow),
typekey = 'type'
)
class override_variational_parameters(QIxml):
attributes = ['href']
#end class override_variational_parameters
class estimators(QIxml):
elements = ['estimator']
#end class estimators
class hamiltonian(QIxml):
# rsqmc afqmc
attributes = ['name','type','target','default']+['info']
# afqmc
parameters = ['filetype','filename']
elements = ['pairpot','constant','estimator']
identifier = 'name'
#end class hamiltonian
class coulomb(QIxml):
tag = 'pairpot'
attributes = ['type','name','source','target','physical','forces']
write_types = obj(physical=yesno)
identifier = 'name'
#end class coulomb
class constant(QIxml):
attributes = ['type','name','source','target','forces']
write_types= obj(forces=yesno)
#end class constant
class pseudopotential(QIxml):
tag = 'pairpot'
attributes = ['type','name','source','wavefunction','format','target','forces','dla','algorithm']
elements = ['pseudo']
write_types= obj(forces=yesno,dla=yesno)
identifier = 'name'
#end class pseudopotential
class pseudo(QIxml):
attributes = ['elementtype','href','format','cutoff','lmax','nrule','l_local']
elements = ['header','local','grid']
identifier = 'elementtype'
#end class pseudo
class mpc(QIxml):
tag='pairpot'
attributes=['type','name','source','target','ecut','physical']
write_types = obj(physical=yesno)
identifier='name'
#end class mpc
class cpp(QIxml):
tag = 'pairpot'
attributes = ['type','name','source','target']
elements = ['element']
identifier = 'name'
#end class cpp
class element(QIxml):
attributes = ['name','alpha','rb']
#end class element
pairpot = QIxmlFactory(
name = 'pairpot',
types = dict(coulomb=coulomb,pseudo=pseudopotential,
pseudopotential=pseudopotential,mpc=mpc,
cpp=cpp),
typekey = 'type'
)
class header(QIxml):
attributes = ['symbol','atomic-number','zval']
#end class header
class local(QIxml):
elements = ['grid']
#end class local
class localenergy(QIxml):
tag = 'estimator'
attributes = ['name','hdf5']
write_types= obj(hdf5=yesno)
identifier = 'name'
#end class localenergy
class energydensity(QIxml):
tag = 'estimator'
attributes = ['type','name','dynamic','static','ion_points']
elements = ['reference_points','spacegrid']
identifier = 'name'
write_types = obj(ion_points=yesno)
#end class energydensity
class reference_points(QIxml):
attributes = ['coord']
text = 'points'
#end class reference_points
class spacegrid(QIxml):
attributes = ['coord','min_part','max_part']
elements = ['origin','axis']
#end class spacegrid
class origin(QIxml):
attributes = ['p1','p2']
#end class origin
class axis(QIxml):
attributes = ['p1','p2','scale','label','grid']
identifier = 'label'
#end class axis
class chiesa(QIxml):
tag = 'estimator'
attributes = ['name','type','source','psi','wavefunction','target']
identifier = 'name'
#end class chiesa
class density(QIxml):
tag = 'estimator'
attributes = ['name','type','delta','x_min','x_max','y_min','y_max','z_min','z_max']
identifier = 'type'
#end class density
class nearestneighbors(QIxml):
tag = 'estimator'
attributes = ['type']
elements = ['neighbor_trace']
identifier = 'type'
#end class nearestneighbors
class neighbor_trace(QIxml):
attributes = ['count','neighbors','centers']
identifier = 'neighbors','centers'
#end class neighbor_trace
class spindensity(QIxml):
tag = 'estimator'
attributes = ['type','name','report']
parameters = ['dr','grid','cell','center','corner','voronoi','test_moves']
write_types = obj(report=yesno)
identifier = 'name'
#end class spindensity
class magnetizationdensity(QIxml):
tag = 'estimator'
attributes = ['type','name','report']
parameters = ['dr','grid','center','corner','integrator','samples']
write_types = obj(report=yesno)
identifier = 'name'
#end class magnetizationdensity
class structurefactor(QIxml):
tag = 'estimator'
attributes = ['type','name','report']
write_types = obj(report=yesno)
identifier = 'name'
#end class structurefactor
class force(QIxml):
tag = 'estimator'
attributes = ['type','name','mode','source','species','target','addionion',
'fast_derivatives','spacewarp','epsilon']
parameters = ['rcut','nbasis','weightexp']
identifier = 'name'
write_types= obj(addionion=yesno,fast_derivatives=yesno,spacewarp=yesno)
#end class force
class forwardwalking(QIxml):
tag = 'estimator'
attributes = ['type','blocksize']
elements = ['observable']
identifier = 'name'
#end class forwardwalking
class pressure(QIxml):
tag = 'estimator'
attributes = ['type','potential','etype','function']
parameters = ['kc']
identifier = 'type'
#end class pressure
class dmccorrection(QIxml):
tag = 'estimator'
attributes = ['type','blocksize','max','frequency']
elements = ['observable']
identifier = 'type'
#end class dmccorrection
class nofk(QIxml):
tag = 'estimator'
attributes = ['type','name','wavefunction']
identifier = 'name'
#end class nofk
class mpc_est(QIxml):
tag = 'estimator'
attributes = ['type','name','physical']
write_types = obj(physical=yesno)
identifier = 'name'
#end class mpc_est
class sk(QIxml):
tag = 'estimator'
attributes = ['name','type','hdf5']
identifier = 'name'
write_types = obj(hdf5=yesno)
#end class sk
class skall(QIxml):
tag = 'estimator'
attributes = ['name','type','hdf5','source','target','writeionion']
identifier = 'name'
write_types = obj(hdf5=yesno,writeionion=yesno)
#end class skall
class gofr(QIxml):
tag = 'estimator'
attributes = ['type','name','num_bin','rmax','source']
identifier = 'name'
#end class gofr
class flux(QIxml):
tag = 'estimator'
attributes = ['type','name']
identifier = 'name'
#end class flux
class orbitalimages(QIxml):
tag = 'estimator'
attributes = ['type','name','ions']
parameters = ['sposets','grid','center_grid','value','corner','cell','center','batch_size']
write_types = obj(center_grid=yesno)
identifier = 'name'
#end class orbitalimages
class momentum(QIxml): # legacy
tag = 'estimator'
attributes = ['type','name','grid','samples','hdf5','wavefunction','kmax','kmax0','kmax1','kmax2']
identifier = 'name'
write_types = obj(hdf5=yesno)
#end class momentum
class momentumdistribution(QIxml): # batched
tag = 'estimator'
attributes = ['type','name','grid','samples','hdf5','wavefunction','kmax','kmax0','kmax1','kmax2']
identifier = 'name'
write_types = obj(hdf5=yesno)
#end class momentumdistribution
class dm1b(QIxml): # legacy
tag = 'estimator'
identifier = 'type'
attributes = ['type','name','reuse']#reuse is a temporary dummy keyword
parameters = ['energy_matrix','basis_size','integrator','points','scale','basis','evaluator','center','check_overlap','check_derivatives','acceptance_ratio','rstats','normalized','volume_normed','samples']
write_types = obj(energy_matrix=yesno,check_overlap=yesno,check_derivatives=yesno,acceptance_ratio=yesno,rstats=yesno,normalized=yesno,volume_normed=yesno)
#end class dm1b
class onebodydensitymatrices(QIxml): # batched
tag = 'estimator'
identifier = 'type'
attributes = ['type','name','reuse']#reuse is a temporary dummy keyword
parameters = ['energy_matrix','basis_size','integrator','points','scale','basis','evaluator','center','check_overlap','check_derivatives','acceptance_ratio','rstats','normalized','volume_normed','samples']
write_types = obj(energy_matrix=yesno,check_overlap=yesno,check_derivatives=yesno,acceptance_ratio=yesno,rstats=yesno,normalized=yesno,volume_normed=yesno)
#end class onebodydensitymatrices
# afqmc estimators
class back_propagation(QIxml):
tag = 'estimator'
attributes = ['name']
parameters = ['naverages','block_size','ortho','nsteps']
elements = ['onerdm']
identifier = 'name'
#end class back_propagation
estimator = QIxmlFactory(
name = 'estimator',
types = dict(localenergy = localenergy,
energydensity = energydensity,
chiesa = chiesa,
density = density,
nearestneighbors = nearestneighbors,
dm1b = dm1b,
spindensity = spindensity,
magnetizationdensity = magnetizationdensity,
structurefactor = structurefactor,
force = force,
forwardwalking = forwardwalking,
pressure = pressure,
dmccorrection = dmccorrection,
nofk = nofk,
mpc = mpc_est,
sk = sk,
skall = skall,
gofr = gofr,
flux = flux,
orbitalimages = orbitalimages,
momentum = momentum,
momentumdistribution = momentumdistribution,
onebodydensitymatrices = onebodydensitymatrices,
# afqmc estimators
back_propagation = back_propagation,
),
typekey = 'type',
typekey2 = 'name'
)
class observable(QIxml):
attributes = ['name','max','frequency']
identifier = 'name'
#end class observable
class init(QIxml):
attributes = ['source','target']
#end class
class scalar_traces(QIxml):
attributes = ['defaults']
text = 'quantities'
write_types = obj(defaults=yesno)
#end class scalar_traces
class array_traces(QIxml):
attributes = ['defaults']
text = 'quantities'
write_types = obj(defaults=yesno)
#end class array_traces
class particle_traces(QIxml): # legacy
attributes = ['defaults']
text = 'quantities'
write_types = obj(defaults=yesno)
#end class particle_traces
class traces(QIxml):
attributes = ['write','throttle','format','verbose','scalar','array',
'scalar_defaults','array_defaults',
'particle','particle_defaults']
elements = ['scalar_traces','array_traces','particle_traces']
write_types = obj(write_=yesno,verbose=yesno,scalar=yesno,array=yesno,
scalar_defaults=yesno,array_defaults=yesno,
particle=yesno,particle_defaults=yesno)
#end class traces
class record(QIxml):
attributes = ['name','stride']
#end class record
class loop(QIxml):
collection_id = 'qmc'
attributes = ['max']
elements = ['qmc','init']
def unroll(self):
calculations=[]
calcs = []
if 'qmc' in self:
calcs = [self.qmc]
elif 'calculations' in self:
calcs = self.calculations
#end if
for n in range(self.max):
for i in range(len(calcs)):
calculations.append(calcs[i].copy())
#end for
#end for
return make_collection(calculations)
#end def unroll
#end class loop
class optimize(QIxml):
text = 'parameters'
#end class optimize
class cg_optimizer(QIxml):
tag = 'optimizer'
attributes = ['method']
parameters = ['max_steps','tolerance','stepsize','friction','epsilon',
'xybisect','verbose','max_linemin','tolerance_g','length_cg',
'rich','xypolish','gfactor']
#end class cg_optimizer
class flex_optimizer(QIxml):
tag = 'optimizer'
attributes = ['method']
parameters = ['max_steps','tolerance','stepsize','epsilon',
'xybisect','verbose','max_linemin','tolerance_g','length_cg',
'rich','xypolish','gfactor']
#end class flex_optimizer
optimizer = QIxmlFactory(
name = 'optimizer',
types = dict(cg=cg_optimizer,flexopt=flex_optimizer),
typekey = 'method',
)
class optimize_qmc(QIxml):
collection_id = 'qmc'
tag = 'qmc'
attributes = ['method','move','renew','completed','checkpoint','gpu']
parameters = ['blocks','steps','timestep','walkers','minwalkers','useweight',
'power','correlation','maxweight','usedrift','min_walkers',
'minke','samples','warmupsteps','minweight','warmupblocks',
'maxdispl','tau','tolerance','stepsize','epsilon',
'en_ref','usebuffer','substeps','stepsbetweensamples',
'samplesperthread','max_steps','nonlocalpp']
elements = ['optimize','optimizer','estimator']
costs = ['energy','variance','difference','weight','unreweightedvariance','reweightedvariance']
write_types = obj(renew=yesno,completed=yesno)
#end class optimize_qmc
class linear(QIxml):
collection_id = 'qmc'
tag = 'qmc'
attributes = ['method','move','profiling','kdelay', # batched
'checkpoint','gpu','trace'] # legacy - batched
elements = ['estimator']
parameters = ['total_walkers','walkers_per_rank','crowds','opt_num_crowds', # batched
'walkers','warmupsteps','blocks','steps','substeps','timestep', # who knows
'usedrift','stepsbetweensamples','samples','minmethod',
'minwalkers','maxweight','nonlocalpp','use_nonlocalpp_deriv',
'usebuffer','alloweddifference','gevmethod','beta','exp0',
'bigchange','stepsize','stabilizerscale','nstabilizers',
'max_its','cgsteps','eigcg','stabilizermethod',
'rnwarmupsteps','walkersperthread','minke','gradtol','alpha',
'tries','min_walkers','samplesperthread',
'shift_i','shift_s','max_relative_change','max_param_change',
'chase_lowest','chase_closest','block_lm','nblocks','nolds',
'nkept','max_seconds','spin_mass',
'sr_tau','sr_tolerance','sr_regularization','line_search',
]
costs = ['energy','unreweightedvariance','reweightedvariance','variance','difference']
write_types = obj(gpu=yesno,usedrift=yesno,nonlocalpp=yesno,usebuffer=yesno,
use_nonlocalpp_deriv=yesno,chase_lowest=yesno,
chase_closest=yesno,block_lm=yesno,line_search=yesno)
#end class linear
class cslinear(QIxml):
collection_id = 'qmc'
tag = 'qmc'
attributes = ['method','move','checkpoint','gpu','trace']
elements = ['estimator']
parameters = ['walkers','warmupsteps','blocks','steps','substeps','timestep',
'usedrift','stepsbetweensamples','samples','minmethod',
'minwalkers','maxweight','nonlocalpp','usebuffer',
'alloweddifference','gevmethod','beta','exp0','bigchange',
'stepsize','stabilizerscale','nstabilizers','max_its',
'stabilizermethod','cswarmupsteps','alpha_error','gevsplit',
'beta_error','use_nonlocalpp_deriv']
costs = ['energy','unreweightedvariance','reweightedvariance']
write_types = obj(gpu=yesno,usedrift=yesno,nonlocalpp=yesno,use_nonlocalpp_deriv=yesno,usebuffer=yesno)
#end class cslinear
class vmc(QIxml):
collection_id = 'qmc'
tag = 'qmc'
attributes = ['method','move','profiling','kdelay', # batched
'multiple','warp','gpu','checkpoint','trace', # legacy - batched
'target','completed','id']
elements = ['estimator', # batched
'record'] # legacy - batched
parameters = ['total_walkers','walkers_per_rank','crowds','warmupsteps', # batched
'blocks','steps','substeps','timestep','maxcpusecs','rewind',
'storeconfigs','checkproperties','recordconfigs','current',
'stepsbetweensamples','samplesperthread','samples','usedrift',
'spin_mass','estimator_period',
'walkers','nonlocalpp','tau','walkersperthread','reconfiguration', # legacy - batched
'dmcwalkersperthread','current','ratio','firststep',
'minimumtargetwalkers','max_seconds']
write_types = obj(usedrift=yesno,profiling=yesno, # batched
gpu=yesno,nonlocalpp=yesno,reconfiguration=yesno, # legacy - batched
ratio=yesno,completed=yesno)
#end class vmc
class dmc(QIxml):
collection_id = 'qmc'
tag = 'qmc'
attributes = ['method','move','profiling','kdelay', # batched
'gpu','multiple','warp','checkpoint','trace', # legacy - batched
'target','completed','id','continue']
elements = ['estimator']
parameters = ['total_walkers','walkers_per_rank','crowds','warmupsteps',
'crowd_serialize_walkers', # batched
'blocks','steps','substeps','timestep','maxcpusecs','rewind',
'storeconfigs','checkproperties','recordconfigs','current',
'stepsbetweensamples','samplesperthread','samples','reconfiguration',
'nonlocalmoves','maxage','alpha','gamma','reserve','use_nonblocking',
'branching_cutoff_scheme','feedback','sigmabound',
'spin_mass','estimator_period',
'walkers','nonlocalmove','pop_control','targetwalkers', # legacy - batched
'minimumtargetwalkers','energybound','feedback','recordwalkers',
'fastgrad','popcontrol','branchinterval','usedrift','storeconfigs',
'en_ref','tau','alpha','gamma','max_branch','killnode','swap_walkers',
'swap_trigger','branching_cutoff_scheme','l2_diffusion','maxage',
'max_seconds']
write_types = obj(usedrift=yesno,profiling=yesno,reconfiguration=yesno,
crowd_serialize_walkers=yesno, # batched
nonlocalmoves=yesnostr,use_nonblocking=yesno,
gpu=yesno,fastgrad=yesno,completed=yesno,killnode=yesno, # legacy - batched
swap_walkers=yesno,l2_diffusion=yesno)
#end class dmc
class rmc(QIxml):
collection_id = 'qmc'
tag = 'qmc'
attributes = ['method','multiple','target','observables','target','warp']
parameters = ['blocks','steps','chains','cuts','bounce','clone','walkers','timestep','trunclength','maxtouch','mass','collect']
elements = ['qmcsystem']
write_types = obj(collect=yesno)
#end class rmc
class vmc_batch(QIxml):
# Do not assume all of the parameters below are supported.
# These were simply copied over from legacy drivers because the
# batched driver compatible inputs have yet not been listed anywhere.
collection_id = 'qmc'
tag = 'qmc'
attributes = ['method','move','profiling','kdelay','checkpoint']
elements = ['estimator']
parameters = ['total_walkers','walkers_per_rank','crowds','warmupsteps','blocks','steps','substeps','timestep','maxcpusecs','rewind','storeconfigs','checkproperties','recordconfigs','current','stepsbetweensamples','samplesperthread','samples','usedrift']
write_types = obj(usedrift=yesno,profiling=yesno)
#end class vmc_batch
class dmc_batch(QIxml):
# Do not assume all of the parameters below are supported.
# These were simply copied over from legacy drivers because the
# batched driver compatible inputs have yet not been listed anywhere.
collection_id = 'qmc'
tag = 'qmc'
attributes = ['method','move','profiling','kdelay','checkpoint']
elements = ['estimator']
parameters = ['total_walkers','walkers_per_rank','crowd_serialize_walkers','crowds','warmupsteps','blocks','steps','substeps','timestep','maxcpusecs','rewind','storeconfigs','checkproperties','recordconfigs','current','stepsbetweensamples','samplesperthread','samples','reconfiguration','nonlocalmoves','maxage','alpha','gamma','reserve','use_nonblocking','branching_cutoff_scheme','feedback','sigmabound']
write_types = obj(usedrift=yesno,profiling=yesno,reconfiguration=yesno,nonlocalmoves=yesnostr,use_nonblocking=yesno, crowd_serialize_walkers=yesno)
#end class dmc_batch
class linear_batch(QIxml):
# Do not assume all of the parameters below are supported.
# These were simply copied over from legacy drivers because the
# batched driver compatible inputs have yet not been listed anywhere.
collection_id = 'qmc'
tag = 'qmc'
attributes = ['method','move','profiling','kdelay']
elements = ['estimator']
parameters = ['walkers','warmupsteps','blocks','steps','substeps','timestep',
'usedrift','stepsbetweensamples','samples','minmethod',
'minwalkers','maxweight','nonlocalpp','use_nonlocalpp_deriv',
'usebuffer','alloweddifference','gevmethod','beta','exp0',
'bigchange','stepsize','stabilizerscale','nstabilizers',
'max_its','cgsteps','eigcg','stabilizermethod',
'rnwarmupsteps','walkersperthread','minke','gradtol','alpha',
'tries','min_walkers','samplesperthread',
'shift_i','shift_s','max_relative_change','max_param_change',
'chase_lowest','chase_closest','block_lm','nblocks','nolds',
'nkept',
'crowds','opt_num_crowds'
]
costs = ['energy','unreweightedvariance','reweightedvariance','variance','difference']
write_types = obj(usedrift=yesno,nonlocalpp=yesno,usebuffer=yesno,use_nonlocalpp_deriv=yesno,chase_lowest=yesno,chase_closest=yesno,block_lm=yesno)
#end class linear_batch
class wftest(QIxml):
collection_id = 'qmc'
tag = 'qmc'
attributes = ['method','checkpoint', 'gpu', 'move', 'multiple', 'warp']
parameters = ['ratio','walkers','clone','source','hamiltonianpbyp','orbitalutility','printeloc','basic','virtual_move']
#elements = ['printeloc','source']
write_types = obj(ratio=yesno,clone=yesno,hamiltonianpbyp=yesno,orbitalutility=yesno,printeloc=yesno,basic=yesno,virtual_move=yesno)
#end class wftest
class setparams(QIxml):
collection_id = 'qmc'
tag = 'qmc'
attributes = ['method','move','checkpoint','gpu']
parameters = ['alpha','blocks','warmupsteps','stepsbetweensamples','timestep','samples','usedrift']
elements = ['estimator']
#end class setparams
qmc = QIxmlFactory(
name = 'qmc',
types = dict(linear=linear,cslinear=cslinear,vmc=vmc,dmc=dmc,loop=loop,optimize=optimize_qmc,wftest=wftest,rmc=rmc,setparams=setparams,vmc_batch=vmc_batch,dmc_batch=dmc_batch,linear_batch=linear_batch),
typekey = 'method',
default = 'loop'
)
class cmc(QIxml):
attributes = ['method','target']
#end class cmc
# afqmc elements
class afqmcinfo(QIxml):
attributes = ['name']
parameters = ['nmo','naea','naeb']
#end class afqmcinfo
class walkerset(QIxml):
attributes = ['name','type']
parameters = ['walker_type']
#end class walkerset
class propagator(QIxml):
attributes = ['name','info']
parameters = ['hybrid']
write_types = obj(hybrid=yesno)
#end class propagator
class execute(QIxml):
attributes = ['info','ham','wfn','wset','prop']
parameters = ['ncores','nwalkers','blocks','steps','timestep']
elements = ['estimator']
#end class execute
class onerdm(QIxml):
None
#end class onerdm
class gen(QIxml):
attributes = []
elements = []
#end class gen
classes = [ #standard classes
simulation,project,application,random,qmcsystem,simulationcell,particleset,
group,hamiltonian,constant,pseudopotential,coulomb,pseudo,mpc,chiesa,density,
localenergy,energydensity,spacegrid,origin,axis,wavefunction,
determinantset,slaterdeterminant,basisset,grid,determinant,occupation,
jastrow1,jastrow2,jastrow3,
correlation,coefficients,loop,linear,cslinear,vmc,dmc,vmc_batch,dmc_batch,linear_batch,
atomicbasisset,basisgroup,init,var,traces,scalar_traces,particle_traces,array_traces,
reference_points,nearestneighbors,neighbor_trace,dm1b,
coefficient,radfunc,spindensity,structurefactor,magnetizationdensity,
sposet,bspline_builder,composite_builder,heg_builder,include,
multideterminant,detlist,ci,mcwalkerset,csf,det,
optimize,cg_optimizer,flex_optimizer,optimize_qmc,wftest,kspace_jastrow,
header,local,force,forwardwalking,observable,record,rmc,pressure,dmccorrection,
nofk,mpc_est,flux,orbitalimages,distancetable,cpp,element,spline,setparams,
backflow,transformation,cubicgrid,molecular_orbital_builder,cmc,sk,skall,gofr,
host,date,user,rpa_jastrow,momentum,override_variational_parameters,
momentumdistribution,onebodydensitymatrices,estimators,rotated_sposet,
# afqmc classes
afqmcinfo,walkerset,propagator,execute,back_propagation,onerdm
]
types = dict( #simple types and factories
#host = param,
#date = param,
#user = param,
pairpot = pairpot,
estimator = estimator,
sposet_builder = sposet_builder,
sposet_collection = sposet_collection,
jastrow = jastrow,
qmc = qmc,
optimizer = optimizer,
)
plurals = obj(
particlesets = 'particleset',
groups = 'group',
hamiltonians = 'hamiltonian',
pairpots = 'pairpot',
pseudos = 'pseudo',
estimators = 'estimator',
spacegrids = 'spacegrid',
axes = 'axis',
wavefunctions = 'wavefunction',
grids = 'grid',
determinants = 'determinant',
correlations = 'correlation',
jastrows = 'jastrow',
basisgroups = 'basisgroup',
calculations = 'qmc',
vars = 'var',
neighbor_traces = 'neighbor_trace',
sposet_builders = 'sposet_builder',
sposets = 'sposet',
radfuncs = 'radfunc',
#qmcsystems = 'qmcsystem', # not a good idea
atomicbasissets = 'atomicbasisset',
cis = 'ci',
csfs = 'csf',
dets = 'det',
observables = 'observable',
optimizes = 'optimize',
#coefficientss = 'coefficients', # bad plurality of qmcpack
constants = 'constant',
mcwalkersets = 'mcwalkerset',
transformations = 'transformation',
rotated_sposets = 'rotated_sposet',
)
plurals_inv = plurals.inverse()
plural_names = set(plurals.keys())
single_names = set(plurals.values())
Names.set_expanded_names(
elementtype = 'elementType',
energydensity = 'EnergyDensity',
gevmethod = 'GEVMethod',
localenergy = 'LocalEnergy',
lr_dim_cutoff = 'LR_dim_cutoff',
lr_tol = 'LR_tol',
lr_handler = 'LR_handler',
minmethod = 'MinMethod',
one_body = 'One-Body',
speciesa = 'speciesA',
speciesb = 'speciesB',
substeps = 'subSteps',
two_body = 'Two-Body',
usedrift = 'useDrift',
maxweight = 'maxWeight',
warmupsteps = 'warmupSteps',
twistindex = 'twistIndex',
twistangle = 'twistAngle',
usebuffer = 'useBuffer',
mpc = 'MPC',
kecorr = 'KEcorr',
ionion = 'IonIon',
elecelec = 'ElecElec',
pseudopot = 'PseudoPot',
posarray = 'posArray',
#array = 'Array', # handle separately, namespace collision
atomicbasisset = 'atomicBasisSet',
basisgroup = 'basisGroup',
expandylm = 'expandYlm',
mo = 'MO',
numerical = 'Numerical',
nearestneighbors = 'NearestNeighbors',
cuspcorrection = 'cuspCorrection',
cuspinfo = 'cuspInfo',
exctlvl = 'exctLvl',
pairtype = 'pairType',
printeloc = 'printEloc',
spindependent = 'spinDependent',
l_local = 'l-local',
pbcimages = 'PBCimages',
dla = 'DLA',
l2_diffusion = 'L2_diffusion',
maxage = 'MaxAge',
sigmabound = 'sigmaBound',
spin_mass = 'spin_mass',
)
# afqmc names
Names.set_afqmc_expanded_names(
afqmcinfo = 'AFQMCInfo',
nmo = 'NMO',
naea = 'NAEA',
naeb = 'NAEB',
hamiltonian = 'Hamiltonian',
wavefunction = 'Wavefunction',
walkerset = 'WalkerSet',
propagator = 'Propagator',
onerdm = 'OneRDM',
nwalkers = 'nWalkers',
estimator = 'Estimator',
)
for c in classes:
c.init_class()
types[c.__name__] = c
#end for
#set default values
simulation.defaults.set(
project = project,
qmcsystem = qmcsystem,
calculations = lambda:list()
)
project.defaults.set(
series=0,
application = application
)
application.defaults.set(
name='qmcpack',role='molecu',class_='serial',version='1.0'
)
#simulationcell.defaults.set(
# bconds = 'p p p',lr_dim_cutoff=15
# )
wavefunction.defaults.set(
name='psi0'
)
#determinantset.defaults.set(
# type='einspline',tilematrix=lambda:eye(3,dtype=int),meshfactor=1.,gpu=False,precision='double'
# )
#occupation.defaults.set(
# mode='ground',spindataset=0
# )
jastrow1.defaults.set(
name='J1',type='one-body',function='bspline',print=True,source='ion0',
correlation=correlation
)
jastrow2.defaults.set(
name='J2',type='two-body',function='bspline',print=True,
correlation=correlation
)
jastrow3.defaults.set(
name='J3',type='eeI',function='polynomial',print=True,source='ion0',
correlation=correlation
)
correlation.defaults.set(
coefficients=coefficients
)
coefficients.defaults.set(
type='Array'
)
#hamiltonian.defaults.set(
# name='h0',type='generic',target='e',
# constant = constant,
# pairpots = classcollection(coulomb,pseudopotential,mpc),
# estimators = classcollection(chiesa),
# )
#coulomb.defaults.set(
# name='ElecElec',type='coulomb',source='e',target='e'
# )
#constant.defaults.set(
# name='IonIon',type='coulomb',source='ion0',target='ion0'
# )
#pseudopotential.defaults.set(
# name='PseudoPot',type='pseudo',source='ion0',wavefunction='psi0',format='xml'
# )
#mpc.defaults.set(
# name='MPC',type='MPC',ecut=60.0,source='e',target='e',physical=False
# )
localenergy.defaults.set(
name='LocalEnergy',hdf5=True
)
#chiesa.defaults.set(
# name='KEcorr',type='chiesa',source='e',psi='psi0'
# )
#energydensity.defaults.set(
# type='EnergyDensity',name='EDvoronoi',dynamic='e',static='ion0',
# spacegrid = spacegrid
# )
#spacegrid.defaults.set(
# coord='voronoi'
# )
density.defaults.set(
type='density',name='Density'
)
spindensity.defaults.set(
type='spindensity',name='SpinDensity'
)
magnetizationdensity.defaults.set(
type='magnetizationdensity',name='MagnetizationDensity'
)
skall.defaults.set(
type='skall',name='skall',source='ion0',target='e',hdf5=True
)
force.defaults.set(
type='Force',name='force'
)
pressure.defaults.set(
type='Pressure'
)
momentum.defaults.set(
type='momentum'
)
momentumdistribution.defaults.set(
type='MomentumDistribution',name='nofk'
)
dm1b.defaults.set(
type = 'dm1b',name='DensityMatrices',energy_matrix=False,
evaluator='matrix',
)
onebodydensitymatrices.defaults.set(
type = 'OneBodyDensityMatrices',name='DensityMatrices',energy_matrix=False,
evaluator='matrix',
)
linear.defaults.set(
method = 'linear',move='pbyp',checkpoint=-1,
#estimators = classcollection(localenergy)
# #jtk
# method='linear',move='pbyp',checkpoint=-1,gpu=True,
# energy=0, reweightedvariance=0, unreweightedvariance=0,
# warmupsteps = 20,
# usedrift = True,
# timestep = .5,
# minmethod ='rescale',
# stepsize = .5,
# beta = 0.05,
# alloweddifference = 1e-8,
# bigchange = 1.1,
# cgsteps = 3,
# eigcg = 1,
# exp0 = -6,
# maxweight = 1e9,
# minwalkers = .5,
# nstabilizers = 10,
# stabilizerscale = .5,
# usebuffer = True,
)
cslinear.defaults.set(
method='cslinear', move='pbyp', checkpoint=-1,
#estimators = classcollection(localenergy)
#jtk
#method='cslinear',move='pbyp',checkpoint=-1,gpu=True,
#energy=0,reweightedvariance=0,unreweightedvariance=1.,
#warmupsteps=5,steps=2,usedrift=True,timestep=.5,
#minmethod='quartic',gevmethod='mixed',exp0=-15,
#nstabilizers=5,stabilizerscale=3,stepsize=.35,
#alloweddifference=1e-5,beta=.05,bigchange=5.,
#estimators=classcollection(localenergy)
#lschulen
#method='cslinear', move='pbyp', checkpoint=-1, gpu=True,
#energy=0, reweightedvariance=0, unreweightedvariance=0,
#warmupsteps = 20,
##steps = 5,
#usedrift = True,
#timestep = .8,
#nonlocalpp = False,
#minmethod = 'rescale',
#stepsize = .4,
#beta = .05,
#gevmethod = 'mixed',
#alloweddifference = 1e-4,
#bigchange = 9.,
#exp0 = -16,
#max_its = 1,
#maxweight = 1e9,
#minwalkers = .5,
#nstabilizers = 3,
#stabilizerscale = 1,
#usebuffer = False,
#estimators = classcollection(localenergy)
#jmm
#method='cslinear', move='pbyp', checkpoint=-1, gpu=True,
#energy=0, reweightedvariance=0, unreweightedvariance=0,
#warmupsteps = 20,
#usedrift = True,
#timestep = .5,
#nonlocalpp = True,
#minmethod = 'quartic',
#stepsize = .4,
#beta = 0.0,
#gevmethod = 'mixed',
#alloweddifference = 1.0e-4,
#bigchange = 9.0,
#exp0 = -16,
#max_its = 1,
#maxweight = 1e9,
#minwalkers = 0.5,
#nstabilizers = 3,
#stabilizerscale = 1.0,
#usebuffer = True,
#estimators = classcollection(localenergy)
)
vmc.defaults.set(
method='vmc',move='pbyp',
#walkers = 1,
#warmupsteps = 50,
#substeps = 3,
#usedrift = True,
#timestep = .5,
#estimators = classcollection(localenergy)
)
dmc.defaults.set(
method='dmc',move='pbyp',
#warmupsteps = 20,
#timestep = .01,
#nonlocalmoves = True,
#estimators = classcollection(localenergy)
)
vmc_batch.defaults.set(
method='vmc_batch',move='pbyp',
)
dmc_batch.defaults.set(
method='dmc_batch',move='pbyp',
)
linear_batch.defaults.set(
method='linear_batch',move='pbyp',
)
# afqmc defaults
afqmcinfo.defaults.set(
name = 'info0',
)
walkerset.defaults.set(
name = 'wset0',
)
propagator.defaults.set(
name = 'prop0',
info = 'info0',
)
execute.defaults.set(
info = 'info0',
ham = 'ham0',
wfn = 'wfn0',
wset = 'wset0',
prop = 'prop0',
)
back_propagation.defaults.set(
name='back_propagation'
)
#end def set_rsqmc_mode
#end def set_afqmc_mode
[docs]
class QmcpackInput(SimulationInput,Names):
profile_collection = None
opt_methods = set(['opt','linear','cslinear','linear_batch'])
simulation_type = simulation
default_metadata = meta(
lattice = dict(units='bohr'),
reciprocal = dict(units='2pi/bohr'),
ionid = dict(datatype='stringArray'),
position = dict(datatype='posArray', condition=0)
)
#end def settings
def __init__(self,arg0=None,arg1=None):
Param.metadata = None
filepath = None
metadata = None
element = None
if arg0 is None and arg1 is None:
None
elif isinstance(arg0,(str, Path)) and arg1 is None:
filepath = path_string(arg0)
elif isinstance(arg0,QIxml) and arg1 is None:
element = arg0
elif isinstance(arg0,meta) and isinstance(arg1,QIxml):
metadata = arg0
element = arg1
else:
self.error('input arguments of types '+arg0.__class__.__name__+' and '+arg0.__class__.__name__+' cannot be used to initialize QmcpackInput')
#end if
if metadata is not None:
self._metadata = metadata
else:
self._metadata = meta()
#end if
if filepath is not None:
self.read(filepath)
elif element is not None:
#simulation = arg0
#self.simulation = self.simulation_type(simulation)
elem_class = element.__class__
if elem_class.identifier is not None:
name = elem_class.identifier
else:
name = elem_class.__name__
#end if
self[name] = elem_class(element)
#end if
Param.metadata = None
QIcollections.clear()
#end def __init__
[docs]
def is_afqmc_input(self):
is_afqmc = False
if 'simulation' in self:
sim = self.simulation
is_afqmc = 'method' in sim and sim.method.lower()=='afqmc'
#end if
return is_afqmc
#end def is_afqmc_input
[docs]
def get_base(self):
elem_names = list(self.keys())
elem_names.remove('_metadata')
if len(elem_names)>1:
self.error('qmcpack input cannot have more than one base element\n You have provided '+str(len(elem_names))+': '+str(elem_names))
#end if
return self[elem_names[0]]
#end def get_base
[docs]
def get_basename(self):
elem_names = list(self.keys())
elem_names.remove('_metadata')
if len(elem_names)>1:
self.error('qmcpack input cannot have more than one base element\n You have provided '+str(len(elem_names))+': '+str(elem_names))
#end if
return elem_names[0]
#end def get_basename
[docs]
def read(self,filepath=None,xml=None):
if xml is not None or os.path.exists(filepath):
element_joins=['qmcsystem']
element_aliases=dict(loop='qmc')
xml = XMLreader(filepath,element_joins,element_aliases,warn=False,xml=xml).obj
xml.condense()
self._metadata = meta() #store parameter/attrib attribute metadata
Param.metadata = self._metadata
if 'simulation' in xml:
self.simulation = simulation(xml.simulation)
else:
#try to determine the type
elements = []
keys = []
error = False
for key,value in xml.items():
if isinstance(key,str) and key[0]!='_':
if key in types:
elements.append(types[key](value))
keys.append(key)
else:
self.error('element '+key+' is not a recognized type',exit=False)
error = True
#end if
#end if
#end for
if error:
self.error('cannot read input xml file')
#end if
if len(elements)==0:
self.error('no valid elements were found for input xml file')
#end if
for i in range(len(elements)):
elem = elements[i]
key = keys[i]
if isinstance(elem,QIxml):
if elem.identifier is not None:
name = elem.identifier
else:
name = elem.tag
#end if
else:
name = key
#end if
self[name] = elem
#end for
#end if
Param.metadata = None
else:
self.error('the filepath you provided does not exist.\n Input filepath: '+filepath)
#end if
return self
#end def read
[docs]
def write_text(self,filepath=None):
set_rsqmc_mode()
if self.is_afqmc_input():
set_afqmc_mode()
#end if
c = ''
header = '''<?xml version="1.0"?>
'''
c+= header
if len(self._metadata)==0:
Param.metadata = self.default_metadata
else:
Param.metadata = self._metadata
#end if
base = self.get_base()
c+=base.write(first=True)
Param.metadata = None
set_rsqmc_mode()
return c
#end def write_text
[docs]
def unroll_calculations(self,modify=True):
qmc = []
sim = self.simulation
if 'calculations' in sim:
calcs = sim.calculations
elif 'qmc' in sim:
calcs = [sim.qmc]
else:
calcs = []
#end if
for i in range(len(calcs)):
c = calcs[i]
if isinstance(c,loop):
qmc.extend(c.unroll())
else:
qmc.append(c)
#end if
#end for
qmc = make_collection(qmc)
if modify:
self.simulation.calculations = qmc
#end if
return qmc
#end def unroll_calculations
#end def get
#end def remove
#end def assign
[docs]
def replace(self,*args,**kwargs):# input is list of keyword=(oldval,newval)
base = self.get_base()
base.replace(*args,**kwargs)
#end def replace
#end def move
#end def get_host
[docs]
def incorporate_defaults(self,elements=False,overwrite=False,propagate=False):
base = self.get_base()
base.incorporate_defaults(elements,overwrite,propagate)
#end def incorporate_defaults
#end def pluralize
[docs]
def standard_placements(self):
self.move(particleset='qmcsystem',wavefunction='qmcsystem',hamiltonian='qmcsystem')
#end def standard_placements
[docs]
def difference(self,other):
s1 = self.copy()
s2 = other.copy()
b1 = s1.get_basename()
b2 = s2.get_basename()
q1 = s1[b1]
q2 = s2[b2]
if b1!=b2:
different = True
d1 = q1
d2 = q2
diff = None
else:
s1.standard_placements()
s2.standard_placements()
s1.pluralize()
s2.pluralize()
different,diff,d1,d2 = q1.difference(q2,root=False)
#end if
if diff is not None:
diff.remove_empty()
#end if
d1.remove_empty()
d2.remove_empty()
return different,diff,d1,d2
#end def difference
#end def remove_empty
[docs]
def read_xml(self,filepath=None,xml=None):
if os.path.exists(filepath):
element_joins=['qmcsystem']
element_aliases=dict(loop='qmc')
if xml is None:
xml = XMLreader(filepath,element_joins,element_aliases,warn=False).obj
else:
xml = XMLreader(None,element_joins,element_aliases,warn=False,xml=xml).obj
#end if
xml.condense()
else:
self.error('the filepath you provided does not exist.\n Input filepath: '+filepath)
#end if
return xml
#end def read_xml
[docs]
def include_xml(self,xmlfile,replace=True,exists=True):
xml = self.read_xml(xmlfile)
Param.metadata = self._metadata
for name,exml in xml.items():
if not name.startswith('_'):
qxml = types[name](exml)
qname = qxml.tag
host = self.get_host(qname)
if host is None and exists:
self.error('host xml section for '+qname+' not found','QmcpackInput')
#end if
if qname in host:
section_name = qname
elif qname in plurals_inv and plurals_inv[qname] in host:
section_name = plurals_inv[qname]
else:
section_name = None
#end if
if replace:
if section_name is not None:
del host[section_name]
#end if
host[qname] = qxml
else:
if section_name is None:
host[qname] = qxml
else:
section = host[section_name]
if isinstance(section,collection):
section[qxml.identifier] = qxml
elif section_name in plurals_inv:
coll = collection()
coll[section.identifier] = section
coll[qxml.identifier] = qxml
del host[section_name]
host[plurals_inv[section_name]] = coll
else:
section.combine(qxml)
#end if
#end if
#end if
#end if
#end for
Param.metadata = None
#end def include_xml
# This include functionality is currently not being used
# The rationale is essentially this:
# -Having includes explicitly represented in the input file object
# makes it very difficult to search for various components
# i.e. where is the particleset? the wavefunction? a particular determinant?
# -Difficulty in locating components makes it difficult to modify them
# -Includes necessarily introduce greater variability in input file structure
# and it is difficult to ensure every possible form is preserved each and
# every time a modification is made
# -The only time it is undesirable to incorporate the contents of an
# include directly into the input file object is if the data is large
# e.g. for an xml wavefunction or pseudopotential.
# In these cases, an external file should be provided that contains
# only the large object in question (pseudo or wavefunction).
# This is already done for pseudopotentials and should be done for
# wavefunctions, e.g. multideterminants.
# Until that time, wavefunctions will be explicitly read into the full
# input file.
[docs]
def add_include(self,element_type,href,placement='on'):
# check the element type
elems = ['cell','ptcl','wfs','ham']
emap = obj(
simulationcell = 'cell',
particleset = 'ptcl',
wavefunction = 'wfs',
hamiltonian = 'ham'
)
if element_type not in elems:
self.error('cannot add include for element of type {0}\n valid element types are {1}'.format(element_type,elems))
#end if
# check the requested placement
placements = ('before','on','after')
if placement not in placements:
self.error('cannot add include for element with placement {0}\n valid placements are {1}'.format(placement,list(placements)))
#end if
# check that the base element is a simulation
base = self.get_base()
if not isinstance(base,simulation):
self.error('an include can only be added to simulation\n attempted to add to {0}'.format(base.__class__.__name__))
#end if
# gather a list of current qmcsystems
if 'qmcsystem' in base:
qslist = [(0,base.qmcsystem)]
del base.qmcsystem
elif 'qmcsystems' in base:
qslist = base.qmcsystems.pairlist()
del base.qmcsystems
else:
qslist = []
#end if
# organize the elements of the qmcsystems
cur_elems = obj()
for elem in elems:
for place in placements:
cur_elems[elem,place] = None
#end for
#end for
for qskey,qs in qslist:
if isinstance(qs,include):
inc = qs
ekey = qskey.split('_')[1]
if ekey not in elems:
self.error('encountered invalid element key: {0}\n valid keys are: {1}'.format(ekey,elems))
#end if
if cur_elems[ekey,'on'] is None:
cur_elems[ekey,'before'] = ekey,inc
else:
cur_elems[ekey,'after' ] = ekey,inc
#end if
elif not isinstance(qs,qmcsystem):
self.error('expected qmcsystem element, got {0}'.format(qs.__class__.__name__))
else:
for elem in qmcsystem.elements:
elem_plural = elem+'s'
name = None
if elem in qs:
name = elem
elif elem_plural in qs:
name = elem_plural
#end if
if name is not None:
cur_elems[emap[elem],'on'] = name,qs[name]
del qs[name]
#end if
#end for
residue = list(qs.keys())
if len(residue)>0:
self.error('extra keys found in qmcsystem: {0}'.format(sorted(residue)))
#end if
#end if
#end for
for elem in elems:
pbef = cur_elems[elem,'before']
pon = cur_elems[elem,'on' ]
paft = cur_elems[elem,'after' ]
if pon is None:
if pbef is not None and paft is None:
cur_elems[elem,'on' ] = pbef
cur_elems[elem,'before'] = None
elif paft is not None and pbef is None:
cur_elems[elem,'on' ] = paft
cur_elems[elem,'after' ] = None
#end if
#end if
#end for
# insert the new include
inc_name = 'include_'+element_type
inc_value = include(href=href)
cur_elems[element_type,placement] = inc_name,inc_value
# create a collection of qmcsystems
qmcsystems = collection()
qskey = ''
qs = qmcsystem()
for elem in elems:
for place in placements:
cur_elem = cur_elems[elem,place]
if cur_elem is not None:
name,value = cur_elem
if isinstance(value,include):
if len(qskey)>0:
qmcsystems.add(qs,key=qskey)
qskey = ''
qs = qmcsystem()
#end if
qmcsystems.add(value,key=name)
else:
qskey += elem[0]
qs[name] = value
#end if
#end if
#end for
#end for
if len(qskey)>0:
qmcsystems.add(qs,key=qskey)
#end if
# attach the collection to the input file
base.qmcsystems = qmcsystems
#end def add_include
[docs]
def get_output_info(self,*requests):
project = self.simulation.project
prefix = project.id
series = project.series
qmc = []
calctypes = set()
outfiles = []
if not self.is_afqmc_input():
qmc_ur = self.unroll_calculations(modify=False)
n=0
for qo in qmc_ur:
q = obj()
q.prefix = prefix
q.series = series+n
n+=1
method = qo.method
if method in self.opt_methods:
q.type = 'opt'
else:
q.type = method
#end if
calctypes.add(q.type)
q.method = method
fprefix = prefix+'.s'+str(q.series).zfill(3)+'.'
files = obj()
files.scalar = fprefix+'scalar.dat'
files.stat = fprefix+'stat.h5'
# apparently this one is no longer generated by default as of r5756
#files.config = fprefix+'storeConfig.h5'
if q.type=='opt':
files.opt = fprefix+'opt.xml'
elif q.type=='dmc':
files.dmc = fprefix+'dmc.dat'
#end if
outfiles.extend(files.values())
q.files = files
qmc.append(q)
#end for
else:
q = obj()
q.prefix = prefix
q.series = series
q.type = 'afqmc'
q.method = 'afqmc'
calctypes.add(q.type)
fprefix = prefix+'.s'+str(q.series).zfill(3)+'.'
files = obj()
files.scalar = fprefix+'scalar.dat'
outfiles.extend(files.values())
q.files = files
qmc.append(q)
#end if
res = dict(qmc=qmc,calctypes=calctypes,outfiles=outfiles)
values = []
for req in requests:
if req in res:
values.append(res[req])
else:
self.error(req+' is not a valid output info request')
#end if
#end for
if len(values)==1:
return values[0]
else:
return values
#end if
#end def get_output_info
[docs]
def generate_jastrows(self,size=None,j1func='bspline',j1size=8,j2func='bspline',j2size=8):
if size is not None:
j1size = size
j2size = size
#end if
#self.remove('jastrow')
lattice,particlesets,wavefunction = self.get('lattice','particleset','wavefunction')
no_lattice = lattice is None
no_particleset = particlesets is None
no_wavefunction = wavefunction is None
if no_lattice:
self.error('a simulationcell lattice must be present to generate jastrows',exit=False)
#end if
if no_particleset:
self.error('a particleset must be present to generate jastrows',exit=False)
#end if
if no_wavefunction:
self.error('a wavefunction must be present to generate jastrows',exit=False)
#end if
if no_lattice or no_particleset or no_wavefunction:
self.error('jastrows cannot be generated')
#end if
if isinstance(particlesets,QIxml):
particlesets = make_collection([particlesets])
#end if
if 'e' not in particlesets:
self.error('electron particleset (e) not found\n particlesets: '+str(particlesets.keys()))
#end if
jastrows = collection()
cell = Structure(lattice)
volume = cell.volume()
rcut = cell.rmin()
#use the rpa jastrow for electrons (modeled after Luke's tool)
size = j2size
e = particlesets.e
nelectrons = 0
for g in e.groups:
nelectrons += g.size
#end for
density = nelectrons/volume
wp = np.sqrt(4*np.pi*density)
dr = rcut/size
r = .02 + dr*np.arange(size)
uuc = .5/(wp*r)*(1.-np.exp(-r*np.sqrt(wp/2)))*np.exp(-(2*r/rcut)**2)
udc = .5/(wp*r)*(1.-np.exp(-r*np.sqrt(wp)))*np.exp(-(2*r/rcut)**2)
jastrows.J2 = jastrow2(
name = 'J2',type='Two-Body',function=j2func,print='yes',
correlations = collection(
uu = correlation(speciesA='u',speciesB='u',size=size,
coefficients=section(id='uu',type='Array',coeff=uuc)),
ud = correlation(speciesA='u',speciesB='d',size=size,
coefficients=section(id='ud',type='Array',coeff=udc))
)
)
#generate electron-ion jastrows, if ions present
ions = []
for name in particlesets.keys():
if name=='i' or name.startswith('ion'):
ions.append(name)
#end if
#end for
if len(ions)>0:
size = j1size
j1 = []
for ion in ions:
i = particlesets[ion]
if 'group' in i:
groups = [i.group]
else:
groups = i.groups
#end if
corr = []
for g in groups:
elem = g.name
c=correlation(
elementtype=elem,
cusp=0.,
size=size,
coefficients=section(
id='e'+elem,
type='Array',
coeff=size*[0]
)
)
corr.append(c)
#end for
j=jastrow1(
name='J1_'+ion,
type='One-Body',
function=j1func,
source=ion,
print='yes',
correlations = corr
)
j1.append(j)
#end for
if len(j1)==1:
j1[0].name='J1'
#end if
for j in j1:
jastrows[j.name]=j
#end for
#end if
if 'J2' in wavefunction.jastrows:
J2 = wavefunction.jastrows.J2
if 'function' in J2 and J2.function.lower()=='bspline':
c = wavefunction.jastrows.J2.correlations
ctot = np.abs(np.array(c.uu.coefficients.coeff)).sum() + np.abs(np.array(c.ud.coefficients.coeff)).sum()
if ctot < 1e-3:
wavefunction.jastrows.J2 = jastrows.J2
#end if
#end if
#end if
#only add the jastrows if ones of the same type
# (one-body,two-body,etc) are not already present
for jastrow in jastrows:
jtype = jastrow.type.lower().replace('-','_')
has_jtype = False
for wjastrow in wavefunction.jastrows:
wjtype = wjastrow.type.lower().replace('-','_')
has_jtype = has_jtype or wjtype==jtype
#end for
if not has_jtype:
wavefunction.jastrows[jastrow.name] = jastrow
#end if
#end for
#end def generate_jastrows
[docs]
def incorporate_system(self,system):
self.warn('incorporate_system may or may not work\n please check the qmcpack input produced\n if it is wrong, please contact the developer')
system = system.copy()
system.check_folded_system()
system.change_units('B')
#system.structure.group_atoms()
system.structure.order_by_species()
particles = system.particles
structure = system.structure
net_charge = system.net_charge
net_spin = system.net_spin
qs,sc,ham,ps = self.get('qmcsystem','simulationcell','hamiltonian','particleset')
old_eps_name = None
old_ips_name = None
if ps is not None:
if isinstance(ps,particleset):
ps = make_collection([ps])
#end if
for pname,pset in ps.items():
g0name = list(pset.groups.keys())[0]
g0 = pset.groups[g0name]
if np.abs(-1-g0.charge)<1e-2:
old_eps_name = pname
elif 'ionid' in pset:
old_ips_name = pname
#end if
#end for
#end if
del ps
self.remove('particleset')
if qs is None:
qs = qmcsystem()
qs.incorporate_defaults(elements=False,propagate=False)
self.simulation.qmcsystem = qs
#end if
if sc is None:
sc = simulationcell()
sc.incorporate_defaults(elements=False,propagate=False)
qs.simulationcell = sc
#end if
if ham is None:
ham = hamiltonian()
ham.incorporate_defaults(elements=False,propagate=False)
qs.hamiltonian = ham
elif isinstance(ham,collection):
if 'h0' in ham:
ham = ham.h0
elif len(ham)==1:
ham = ham.list()[0]
else:
self.error('cannot find hamiltonian for system incorporation')
#end if
#end if
elem = structure.elem
pos = structure.pos
if len(structure.axes)>0: #exclude systems with open boundaries
#setting the 'lattice' (cell axes) requires some delicate care
# qmcpack will fail if this is even 1e-10 off of what is in
# the wavefunction hdf5 file from pwscf
if structure.folded_structure is not None:
fs = structure.folded_structure
axes = np.array(pwscf_array_string(fs.axes).split(),dtype=float)
npe.reshape_inplace(axes, fs.axes.shape)
axes = np.dot(structure.tmatrix,axes)
if np.abs(axes-structure.axes).sum()>1e-5:
self.error('supercell axes do not match tiled version of folded cell axes\n you may have changed one set of axes (super/folded) and not the other\n folded cell axes:\n'+str(fs.axes)+'\n supercell axes:\n'+str(structure.axes)+'\n folded axes tiled:\n'+str(axes))
#end if
else:
axes = np.array(pwscf_array_string(structure.axes).split(),dtype=float)
npe.reshape_inplace(axes, structure.axes.shape)
#end if
structure.adjust_axes(axes)
sc.lattice = axes
#end if
elns = particles.get_electrons()
ions = particles.get_ions()
eup = elns.up_electron
edn = elns.down_electron
particlesets = []
eps = particleset(
name='e',random=True,
groups = [
group(name='u',charge=-1,mass=eup.mass,size=eup.count),
group(name='d',charge=-1,mass=edn.mass,size=edn.count)
]
)
particlesets.append(eps)
if len(ions)>0:
if sc is not None and 'bconds' in sc and tuple(sc.bconds)!=('p','p','p'):
eps.randomsrc = 'ion0'
#end if
ips = particleset(
name='ion0',
)
groups = []
ham.pluralize()
pseudos = ham.get('pseudo')
if pseudos is None:
pp = ham.get('PseudoPot')
if pp is not None:
pseudos = collection()
pp.pseudos = pseudos
#end if
#end if
for ion in ions:
gpos = pos[elem==ion.name]
g = group(
name = ion.name,
charge = ion.charge,
valence = ion.charge,
atomicnumber = ion.protons,
mass = ion.mass,
position = gpos,
size = len(gpos)
)
groups.append(g)
if pseudos is not None and ion.name not in pseudos:
pseudos[ion.name] = pseudo(elementtype=ion.name,href='MISSING.xml')
#end if
#end for
ips.groups = make_collection(groups)
particlesets.append(ips)
#end if
qs.particlesets = make_collection(particlesets)
if old_eps_name is not None:
self.replace(old_eps_name,'e')
#end if
if old_ips_name is not None and len(ions)>0:
self.replace(old_ips_name,'ion0')
#end if
udet,ddet = self.get('updet','downdet')
if udet is not None:
udet.size = elns.up_electron.count
#end if
if ddet is not None:
ddet.size = elns.down_electron.count
#end if
if np.abs(net_spin) > 1e-1:
if ddet is not None:
if 'occupation' in ddet:
ddet.occupation.spindataset = 1
else:
ss = self.get('sposets')
ss[ddet.sposet].spindataset = 1
#end if
#end if
#end if
#end def incorporate_system
[docs]
def get_electron_particle_set(self):
input = self.copy()
input.pluralize()
return input.get('particlesets').e
#end def get_electron_particle_set
[docs]
def return_system(self,structure_only=False):
input = self.copy()
input.pluralize()
axes,ps,H = input.get('lattice','particlesets','hamiltonian')
if ps is None:
return None
#end if
# find electrons and ions
have_ions = True
have_jellium = False
ions = None
elns = None
ion_list = []
for name,p in ps.items():
if 'ionid' in p:
ion_list.append(p)
elif name.startswith('e'):
elns = p
#end if
#end for
if len(ion_list)==0: #try to identify ions by positive charged groups
for name,p in ps.items():
if 'groups' in p:
for g in p.groups:
if 'charge' in g and g.charge>0:
ion_list.append(p)
break
#end if
#end for
#end if
#end for
#end if
if len(ion_list)==1:
ions = ion_list[0]
elif len(ion_list)>1:
self.error('ability to handle multiple ion particlesets has not been implemented')
#end if
if ions is None and elns is not None and 'groups' in elns:
simcell = input.get('simulationcell')
if simcell is not None and 'rs' in simcell:
have_ions = False
have_jellium = True
elif 'pairpots' not in H:
have_ions = False
#end if
#end if
if elns is None:
self.error('could not find electron particleset')
#end if
if ions is None and have_ions:
self.error('could not find ion particleset')
#end if
#compute spin and electron charge
net_spin = 0
eln_charge = 0
for spin,eln in elns.groups.items():
if spin[0]=='u':
net_spin+=eln.size
elif spin[0]=='d':
net_spin-=eln.size
#end if
eln_charge += eln.charge*eln.size
#end if
#get structure and ion charge
structure = None
ion_charge = 0
valency = dict()
if have_ions:
elem = None
if 'ionid' in ions:
if isinstance(ions.ionid,str):
elem = [ions.ionid]
else:
elem = list(ions.ionid)
#end if
pos = ions.position
elif 'size' in ions and ions.size==1:
elem = [ions.groups.list()[0].name]
pos = [[0,0,0]]
elif 'groups' in ions:
elem = []
pos = []
for group in ions.groups:
if 'position' in group:
nions = group.size
elem.extend(nions*[group.name])
if group.size==1:
pos.extend([list(group.position)])
else:
pos.extend(list(group.position))
#end if
#end if
#end for
if len(elem)==0:
elem = None
pos = None
else:
elem = np.array(elem)
pos = np.array(pos)
order = elem.argsort()
elem = elem[order]
pos = pos[order]
#end if
#end if
if elem is None:
self.error('could not read ions from ion particleset')
#end if
if axes is None:
center = (0,0,0)
else:
md = input._metadata
if 'position' in md and 'condition' in md['position'] and md['position']['condition']==1:
pos = np.dot(pos,axes)
#end if
center = axes.sum(0)/2
#end if
# pos must be a 2D array, shape (N,3)
# reshape single atom case, shape (3,) as shape (1,3)
pos = np.asarray(pos)
if len(pos.flatten())==3:
npe.reshape_inplace(pos, (1, 3))
#end if
structure = Structure(axes=axes,elem=elem,pos=pos,center=center,units='B')
for name,element in ions.groups.items():
if 'charge' in element:
valence = element.charge
elif 'valence' in element:
valence = element.valence
elif 'atomic_number' in element:
valence = element.atomic_number
else:
self.error('could not identify valency of '+name)
#end if
valency[name] = valence
count = list(elem).count(name)
ion_charge += valence*count
#end for
elif have_jellium:
structure = Jellium(rs=simcell.rs,background_charge=-eln_charge)
ion_charge = structure.background_charge
#end if
net_charge = ion_charge + eln_charge
system = PhysicalSystem(structure,net_charge,net_spin,**valency)
if structure_only:
return structure
else:
return system
#end if
#end def return_system
[docs]
def get_ion_particlesets(self):
ions = obj()
ps = self.get('particlesets')
#try to identify ions by positive charged groups
for name,p in ps.items():
if name.startswith('ion') or name.startswith('atom'):
ions[name] = p
elif 'groups' in p:
for g in p.groups:
if 'charge' in g and g.charge>0:
ions[name] = p
break
#end if
#end for
#end if
#end for
return ions
#end def get_ion_particlesets
[docs]
def get_pp_files(self):
pp_files = []
h = self.get('hamiltonian')
if h is not None:
pp = None
if 'pairpots' in h:
for pairpot in h.pairpots:
if 'type' in pairpot and pairpot.type=='pseudo':
pp = pairpot
#end if
#end for
elif 'pairpot' in h and 'type' in h.pairpot and h.pairpot.type=='pseudo':
pp = h.pairpot
#end if
if pp is not None:
if 'pseudo' in pp and 'href' in pp.pseudo:
pp_files.append(pp.pseudo.href)
elif 'pseudos' in pp:
for pseudo in pp.pseudos:
if 'href' in pseudo:
pp_files.append(pseudo.href)
#end if
#end for
#end if
#end if
#end if
return pp_files
#end def get_pp_files
[docs]
def remove_physical_system(self):
qs = self.simulation.qmcsystem
if 'simulationcell' in qs:
del qs.simulationcell
#end if
if 'particlesets' in qs:
del qs.particlesets
#end if
for name in qs.keys():
if isinstance(qs[name],particleset):
del qs[name]
#end if
#end for
self.replace('ion0','i')
#end def remove_physical_system
[docs]
def cusp_correction(self):
cc = False
if not self.is_afqmc_input():
ds = self.get('determinantset')
cc_var = ds is not None and 'cuspcorrection' in ds and ds.cuspcorrection == True
cc_run = len(self.simulation.calculations)==0
cc = cc_var and cc_run
#end if
return cc
#end def cusp_correction
[docs]
def get_driver(self):
driver = self.get('driver_version')
if driver is None or driver.startswith('batch'):
driver = 'batched'
assert driver in ('batched','legacy')
return driver
#end def get_driver()
[docs]
def set_driver(self,driver):
if driver.startswith('batch'):
driver = 'batched'
assert driver in ('batched','legacy')
proj = self.get('project')
proj.driver_version = driver
#end set_driver
#end def has_jastrows
#end def get_jastrows
#end def remove_jastrows
[docs]
def remove_J1(self):
J = self.get('jastrow')
if J is not None:
Jrem = []
for name,Jn in J.items():
if Jn.type=='One-Body':
Jrem.append(name)
for name in Jrem:
del J[name]
#end def remove_J1
[docs]
def remove_J2(self):
J = self.get('jastrow')
if J is not None:
Jrem = []
for name,Jn in J.items():
if Jn.type=='Two-Body':
Jrem.append(name)
for name in Jrem:
del J[name]
#end def remove_J2
[docs]
def remove_J3(self):
J = self.get('jastrow')
if J is not None:
Jrem = []
for name,Jn in J.items():
if Jn.type=='eeI':
Jrem.append(name)
for name in Jrem:
del J[name]
#end def remove_J3
[docs]
def gen_jastrows(self,**kwargs):
self.remove_jastrows()
system = kwargs.pop('system',None)
if system is None:
system = self.return_system()
else:
assert isinstance(system,PhysicalSystem)
jastrows = generate_jastrows_alt(system=system,**kwargs)
wfn = self.get('wavefunction')
if wfn is None:
self.error('cannot set jastrows.\nWavefunction is missing.')
wfn.jastrows = make_collection(jastrows)
#end def gen_jastrows
[docs]
def optimize_jastrows(self,opt=True):
opt = bool(opt)
jastrows = self.get_jastrows()
if jastrows is not None:
jastrow_classes = tuple(jastrow.types.values())
for n,Jn in jastrows.items():
if not isinstance(Jn,QIxml):
continue
assert isinstance(Jn,jastrow_classes)
for corr in Jn.correlations.values():
corr.coefficients.optimize = opt
#end def optimize_jastrows
[docs]
def set_orbitals_h5(self,orbitals_h5):
assert isinstance(orbitals_h5,str)
assert ' ' not in orbitals_h5
assert orbitals_h5.endswith('.h5')
wfn = self.get('wavefunction')
assert wfn is not None
assert 'determinantset' in wfn
dset = wfn.determinantset
spob = self.get('sposet_builder')
if spob is None:
spob = self.get('sposet_collection')
if spob is None:
dset.href = orbitals_h5
elif 'bspline' in spob:
spob.bspline.href = orbitals_h5
else:
self.error('orbital file assignment is only supported for sposet_builder with B-spline orbitals')
#end def set_orbitals_h5
[docs]
def has_lcao_orbitals(self):
dset = self.get('determinantset')
assert dset is not None
lcao = 'type' in dset and dset.type=='MolecularOrbital'
return lcao
#end def has_lcao_orbitals
[docs]
def set_lcao_orbital_file(self,filepath):
assert filepath.endswith('.h5')
if not self.has_lcao_orbitals():
self.error('calculation type is not LCAO. Cannot assign LCAO orbital file.')
dset = self.get('determinantset')
assert dset is not None
dset.href = filepath
#end def set_lcao_orbital_file
#end def has_multidet
#end def get_multidet
[docs]
def optimize_multidet(self,opt=True):
opt = bool(opt)
md = self.get_multidet()
if md is None:
self.error('input file has no multideterminant')
assert isinstance(md,multideterminant)
assert 'detlist' in md
md.detlist.optimize = opt
#end def optimize_multidet
[docs]
def set_multidet_params(self,**kwargs):
md = self.get_multidet()
if md is None:
self.error('input file has no multideterminant')
dl = md.detlist
names = set(list(kwargs.keys()))
mdc = multideterminant
md_names = set(mdc.attributes)|set(mdc.parameters)
dl_names = set(detlist.attributes)|set(detlist.parameters)
allowed_names = md_names|dl_names
invalid = names - allowed_names
if len(invalid)>0:
self.error('unrecognized multideterminant parameters encountered.\n Allowed params are: {}\nYou provided:{}'.format(list(sorted(allowed_names)),list(sorted(invalid))))
for name in md_names:
if name in kwargs:
md[name] = kwargs[name]
for name in dl_names:
if name in kwargs:
dl[name] = kwargs[name]
#end def set_mutidet_params
[docs]
def set_multidet_h5(self,filepath):
assert isinstance(filepath,str)
assert ' ' not in filepath
assert filepath.endswith('.h5')
self.set_multidet_params(href=filepath)
#end def set_multidet_h5
[docs]
def set_pseudo_files(self,**pseudo_files):
pps = self.get('pseudo')
assert pps is not None
species = list(pps.keys())
for spec1,filepath in pseudo_files.items():
for spec2 in species:
if spec1.lower()==spec2.lower():
pp = pps[spec2]
pp.href = filepath
del pps[spec2]
pps[spec1] = pp
#end def set_pseudo_files
#end def has_qmc
[docs]
def get_qmc(self,series):
series = int(series)
qmc = None
calcs = self.get('calculations')
series_start = self.get('series')
if calcs is not None:
if series_start is not None:
series -= series_start
if series in calcs:
qmc = calcs[series]
return qmc
#end def get_qmc
[docs]
def remove_qmc(self,series):
series = int(series)
qmc = None
calcs = self.get('calculations')
series_start = self.get('series')
if calcs is not None:
if series_start is not None:
series -= series_start
if series in calcs:
qmc = calcs[series]
del calcs[series]
while series+1 in calcs:
calcs[series] = calcs[series+1]
del calcs[series+1]
series += 1
else:
self.error('qmc method with series {} not found'.format(series))
#return qmc
#end def remove_qmc
#end def has_calculations
#end def get_calculations
#end def remove_calculations
[docs]
def gen_calculations(self,qmc,**kw):
allowed_qmc = ('opt','vmc','vmc_test','vmc_noJ',
'dmc','dmc_test','dmc_noJ')
if qmc not in allowed_qmc:
self.error('calculation type "{}" is unrecognized.\nValid options are: {}'.format(qmc,allowed_qmc))
kw = obj(**kw)
driver = self.get_driver()
kw.set_optional(**qmc_defaults[driver][qmc])
kw.driver = driver
#self.remove_calculations()
if qmc=='opt':
calcs = generate_opt_calculations(**kw)
elif 'vmc' in qmc:
calcs = generate_vmc_calculations(**kw)
elif 'dmc' in qmc:
calcs = generate_dmc_calculations(**kw)
assert isinstance(calcs,list)
for calc in calcs:
calc.incorporate_defaults(elements=False,overwrite=False,propagate=True)
calcs = make_collection(calcs)
self.simulation.calculations = calcs
#end def gen_calculations
[docs]
def modify(self,
driver = None,
remove_system = False,
change_system = False,
remove_jastrows = False,
remove_J1 = False,
remove_J2 = False,
remove_J3 = False,
remove_determinants = False,
remove_multidet = False,
remove_calculations = False,
# generate_jastrow_alt inputs
J1 = False,
J2 = False,
J3 = False,
J1_size = None,
J1_rcut = None,
J1_dr = 0.5,
J1_opt = True,
J2_size = None,
J2_rcut = None,
J2_dr = 0.5,
J2_init = 'zero',
J2_opt = True,
J3_isize = 3,
J3_esize = 3,
J3_rcut = 5.0,
J3_opt = None,
J1_rcut_open = 5.0,
J2_rcut_open = 10.0,
J1k = False,
J1k_kcut = 5.0,
J1k_symm = 'crystal',
J1k_opt = None,
J2k = False,
J2k_kcut = 5.0,
J2k_symm = 'crystal',
J2k_opt = None,
system = None,
# other jastrow
jastrow_opt = None,
# determinant inputs
orbitals_h5 = None,
# multidet inputs
multidet_h5 = None,
multidet_cutoff = None,
multidet_opt = None,
# other wavefunction
optimize = None,
# hamiltonian
pseudo_files = None,
# calculations input
calculations = None,
qmc = None,
**gen_calcs
):
"""Modify the parameters and xml elements of a QMCPACK input file.
Parameters
----------
driver : {'batched', 'legacy'} or None, default=None
Sets ``driver_version`` in QMCPACK input. If ``None``, ``'batched'``
is assumed.
remove_system : bool, default=False
Removes ``<simulationcell/>`` and ``<particleset/>``.
change_system : PhysicalSystem (such as from `generate_physical_system`)
Updates physical system information in ``<simulationcell/>`` and
``<particleset/>`` to match the contents of the ``PhysicalSystem``
object.
remove_jastrows : bool, default=False
Removes all ``<jastrow/>`` elements.
remove_J1 : bool, default=False
Remove only the one-body jastrow, ``<jastrow type="One-Body"/>``
remove_J2 : bool, default=False
Remove only the two-body jastrow, ``<jastrow type="Two-Body"/>``
remove_J3 : bool, default=False
Remove only the three-body jastrow, ``<jastrow type="eeI"/>``
remove_determinants : bool, default=False
Removes ``<determinantset/>``
remove_multidet : bool, default=False
Removes ``<multideterminant/>``
remove_calculations : bool, default=False
Removes all ``<qmc/>`` and ``<loop/>`` elements.
optimize : bool or None, default=None
Sets ``optimize`` parameter in all wavefunction components.
If ``True`` or ``False`` ``optimize`` is set accordingly.
If ``None`` no changes are made.
jastrow_opt : bool or None, default=None
Sets ``optimize`` parameters in all ``<jastrow/>`` elements.
Logic is identical to ``optimize``.
orbitals_h5 : str or None
Sets path to an HDF5 file containing single particle orbitals.
If type ``str``, ``href`` is set in ``<sposet_builder/>`` or
``<sposet_collection/>`` if present and in ``<determinantset/>``
otherwise.
If ``None``, no action is taken.
multidet_h5 : bool or None, default=None
Set path to an HDF5 file containing multideterminat coefficents.
If type ``str``, ``href`` is set in ``<multideterminant/>``.
If ``None``, no action is taken.
multidet_cutoff : float or None, default=None
Sets the multideterminant coefficient cutoff, which itself
determints to include based on their magnitude relative to
the cutoff.
If type ``float``, ``cutoff`` in ``<detlist/>`` is set.
If ``None``, no action is taken.
pseudo_files : dict of str:str, default={}
Sets paths to pseudopotential files.
Any atomic species as keywords and pseudopotential filepaths
as values.
For example:
.. code-block:: python
qi.modify(
pseudo_files = dict(
Mo = 'Mo.ccECP.xml',
S = 'S.ccECP.xml'))
Atomic species matching is case insensitive.
calculations : None or list of qmc or loop objects
Overwrite all ``<qmc/>`` or ``<loop/>`` elements with those provided.
If ``None``, no action is taken.
Notes
-----
The remaining input parameters are broken into sections based on
what part of the input file they modify.
Jastrow Generation Parameters
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
Generate Jastrow factors based on the parameters given. Existing
Jastrows are overwritten. The parameter signature is identical
to ``generate_jastrows_alt`` called both here and by
``generate_qmcpack_input`` or ``generate_qmcpack``.
J1 : bool, default=False
Creates a one-body B-spline Jastrow if ``True``.
If no other ``J1_`` parameters are given, sensible defaults are set:
cutoff set to the Wigner-Seitz radius for periodic systems or to
5 Bohr for open boundary conditions.
By default, one knot is placed every 0.5 Bohr up to the cutoff.
J2 : bool, default=False.
Creates both one-body and two-body B-spline Jastrows if ``True``.
If no other ``J2_`` parameters are given, sensible defaults are set:
cutoff set to the Wigner-Seitz radius for periodic systems or to
10 Bohr for open boundary conditions.
By default, one knot is placed every 0.5 Bohr up to the cutoff.
J3 : bool, default=False.
Creates one-, two-, and three-body B-spline Jastrows if ``True``.
If no other ``J3_`` parameters are given, sensible defaults are set:
cutoff set to 5 Bohr, ``isize=3``, ``esize=3``.
J1_rcut : float or None, default=None
Sets the cutoff (``rcut``) in the one-body Jastrow.
If ``None``, the Wigner-Seitz radius is used for periodic systems,
or the value of ``J1_rcut_open`` for open boundary conditions.
J1_rcut_open : float, default=5.0
Sets the cutoff (``rcut``) in the one-body Jastrow for open systems.
J1_size : int or None, default=None
Sets the number of knots in the one-body B-spline Jastrow.
If ``int``, knots are placed up to the cutoff.
If ``None``, ``J1_dr`` is used instead.
J1_dr : float, default=0.5
Sets B-spline knots every ``J1_dr`` up to the cutoff.
J1_opt : bool or None, default=None
If ``bool``, sets ``optimize`` flag in the one-body Jastrow.
If ``None``, no action is taken.
J2_rcut : float or None, default=None
Sets the cutoff (``rcut``) in the two-body Jastrow.
If ``None``, the Wigner-Seitz radius is used for periodic systems,
or the value of ``J2_rcut_open`` for open boundary conditions.
J2_rcut_open : float, default=10.0
Sets the cutoff (``rcut``) in the two-body Jastrow for open systems.
J2_size : int or None, default=None
Sets the number of knots in the one-body B-spline Jastrow.
If ``int``, knots are placed up to the cutoff.
If ``None``, ``J2_dr`` is used instead.
J2_dr : float, default=0.5
Sets B-spline knots every ``J2_dr`` up to the cutoff.
J2_opt : bool or None, default=None
If ``bool``, sets ``optimize`` flag in the two-body Jastrow.
If ``None``, no action is taken.
J2_init : {'zero', 'rpa'}, default='zero'
If ``zero``, set all B-spline coefficients to 0.0.
If ``rpa``, set B-spline coefficents based on the RPA Jastrow for
a homogeneous electron gas with the same electron density as the
current atomic system.
For an open system only ``'zero'`` is allowed.
J1k : bool, default=False
Creates a one-body k-space Jastrow with defaults below if ``True``.
J1k_kcut : float, default=5.0
Sets the k-space cutoff which determines how many plane-waves
and coefficients are used.
J1k_symm : {'crystal', 'isotropic', 'none'}
Whether to use symmetries to constrain the plane-wave coefficients.
If ``'crystal'``, enforce translation symmetries.
If ``'isotropic'``, impose symmetry based on identical :math:`|k|`.
If ``'none'``, the coefficients are fully unconstrained.
J1k_opt : bool or None, default=None
If ``bool``, sets ``optimize`` flag in the one-body k-space Jastrow.
If ``None``, no action is taken.
J2k : bool, default=False
Creates a two-body k-space Jastrow with defaults below if ``True``.
J2k_kcut : float, default=5.0
Sets the k-space cutoff which determines how many plane-waves
and coefficients are used.
J2k_symm : {'crystal', 'isotropic', 'none'}
Whether to use symmetries to constrain the plane-wave coefficients.
If ``'crystal'``, enforce translation symmetries.
If ``'isotropic'``, impose symmetry based on identical :math:`|k-k'|`.
If ``'none'``, the coefficients are fully unconstrained.
J2k_opt : bool or None, default=None
If ``bool``, sets ``optimize`` flag in the two-body k-space Jastrow.
If ``None``, no action is taken.
QMC Calculation Generation Parameters
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
Generate ``<qmc/>`` and/or ``<loop/>`` elements. Any existing elements
are overwritten.
The number of input parameters depends on the value of ``qmc`` and
are given as keyword inputs represented by ``gen_calcs``.
Only inputs for batched drivers are described below.
Parameters at the top are shared by nearly all ``qmc`` methods.
qmc : {'vmc', 'vmc_test', 'vmc_noJ', 'dmc', 'dmc_test', 'dmc_noJ', 'opt'} or None
If ``None``, no action is taken.
Otherwise calculations are generated as detailed below.
Shared Parameters
^^^^^^^^^^^^^^^^^
total_walkers : int or None, default=None
If not ``None``, set the ``total_walkers`` parameter, which is the
number of independent walker configuration trajectories across
within each VMC sampling. If using MPI or threads, the walkers
will be divided roughly evenly between each MPI rank/thread.
walkers_per_rank : int or None, default=None
If not ``None``, set the ``walkers_per_rank`` parameter, which is
the number of independent walker configuration trajectories
within each MPI rank. In this case, the total number of
walkers is ``#MPI_ranks*walkers_per_rank``.
Only one of {``walkers_per_rank``, ``total_walkers``} should be
provided.
warmupsteps : int
Number of VMC steps used to move the walker population toward
the equilibrium distribution before sampling estimators
such as the total energy.
blocks : int
Sets ``blocks`` parameter, the outer loop in the VMC/DMC
sampling process.
steps : int or None, default=None
If not None, set the ``steps`` parameter, the inner loop in
the VMC/DMC sampling. The resulting number of samples per
walker is ``blocks*steps``.
Only one of {``samples``, ``steps``} should be provided.
substeps : int
Sets the ``substeps`` parameter, which is the number of VMC
steps in between the generation of each sample. Used to
decorrelate walker configurations between the collection of
each sample (energy evaluation). Does not apply to DMC
calculations.
timestep : float
Sets the ``timestep`` parameter, which is the width of the
gaussian used to generate the next configuration in each
walker's configuration trajectory. Affects the acceptance
ratio, and hence the efficiency of the sampling. In production
DMC a small value should be used (e.g. 0.01) to prioritize the
accuracy of the solution (minimize timestep) over apparent
gains in efficiency.
usedrift : bool
Sets the ``usedrift`` parameter.
If ``True``, use the logarithmic gradient to shift the gaussian
center for a more efficient sampling (higher acceptance ratio).
Used only in VMC.
checkpoint : int or None, default=None
If not ``None``, set the ``checkpoint`` parameter. A checkpoint
HDF5 file will be written every ``checkpoint`` blocks.
maxcpusecs : float or None
If not ``None``, set the ``maxcpusecs`` parameter. QMCPACK will
terminate gracefully if the walltime exceeds this value.
crowds : int or None, default=None
If not ``None``, set the ``crowds`` parameter, which controls
the partitioning of walkers for parallel (thread/gpu)
execution.
spinmass : float or None, default=None
If not ``None``, set the ``spinmass`` parameter. Generally only
used in calculations including spin-orbit coupling.
Case ``qmc='vmc'``
^^^^^^^^^^^^^^^^^^
As in "Shared Parameters" above, but with the defaults below.
- ``warmupsteps : int, default=50``
- ``blocks : int, default=800``
- ``steps : int, default=10``
- ``substeps : int, default=3``
- ``timestep : float, default=0.3``
- ``usedrift : bool, default=False``
Case ``qmc='vmc_test'``
^^^^^^^^^^^^^^^^^^^^^^^
As in ``qmc='vmc'``, but with the defaults below. Intended to
make a quick test run to check for successful execution or to
obtain timing estimates to design production runs.
Case ``qmc='vmc_noJ'``
^^^^^^^^^^^^^^^^^^^^^^
As in ``qmc='vmc'``, but with the defaults below. Uses increased
sampling intended to better deal with the increased variance
present in Jastrow-free runs.
- ``warmupsteps : int, default=200``
- ``blocks : int, default=800``
- ``steps : int, default=100``
Case ``qmc='dmc'``
^^^^^^^^^^^^^^^^^^
As in "Shared Parameters" above, but with the defaults below.
These parameter names and defaults refer to the DMC sections.
- ``warmupsteps : int, default=20``
- ``blocks : int, default=200``
- ``steps : int, default=10``
- ``timestep : float, default=0.01``
nonlocalmoves : {'v0', 'v1', 'v3'} or bool or None, default=None
Perform T-moves or the locality approximation.
If ``None``, use QMCPACK's default (locality approx)
If ``False``, use the locality approximation.
If ``True`` or ``'v0'``, use the first developed T-moves algorithm.
If ``'v1'``, use the second developed T-moves algorithm.
If ``'v3'``, use a modified T-moves algorithm, courtesy Ye Luo.
branching_cutoff_scheme
See QMCPACK manual.
crowd_serialize_walkers
See QMCPACK manual.
reconfiguration
See QMCPACK manual.
maxage
See QMCPACK manual.
feedback
See QMCPACK manual.
sigmabound
See QMCPACK manual.
vmc_warmupsteps : int, default=30
Set ``warmupsteps`` in the VMC block executed prior to DMC.
The parameters below set the respective params in VMC.
vmc_blocks : int, default=40
vmc_steps : int, default=10
vmc_substeps : int, default=3
vmc_timestep : float, default=0.3
vmc_usedrift : bool, default=False
vmc_checkpoint : int or None, default=None
vmc_spin_mass : float or None, default=None
eq_dmc : bool, default=False
Insert a DMC block following VMC for the purpose of
rapid equilibration prior to the subsequent production
DMC sections.
eq_warmupsteps : int, default=20
eq_blocks : int, default=20
eq_steps : int, default=5
eq_timestep : float, default=0.02
The timestep should be greater than or equal to the ones used
in the subsequent DMC sections.
eq_checkpoint : int or None, default=None
ntimesteps : int, default=1
If greater than one, create a sequence of ``ntimesteps`` DMC
sections with successively smaller timesteps. Intended for
DMC timestep extrapolation.
timestep_factor : float, default=0.5
The first timestep is given by ``timestep``, the following ones
are reduced by successive multiplication of ``timestep_factor``.
Case ``qmc='dmc_test'``
^^^^^^^^^^^^^^^^^^^^^^^
As in ``qmc='dmc'``, but with the defaults below. Intended to
make a quick test run to check for successful execution or to
obtain timing estimates to design production runs.
- ``warmupsteps : int, default=2``
- ``blocks : int, default=10``
- ``steps : int, default=2``
- ``vmc_warmupsteps : int, default=10``
- ``vmc_blocks : int, default=4``
- ``eq_dmc : bool, default=False``
- ``eq_warmupsteps : int, default=2``
- ``eq_blocks : int, default=5``
- ``eq_steps : int, default=2``
Case ``qmc='dmc_noJ'``
^^^^^^^^^^^^^^^^^^^^^^
As in ``qmc='dmc'``, but with the defaults below. Uses increased
sampling intended to better deal with the increased variance
present in Jastrow-free runs. Note that Jastrow-free runs are
much more likely to be unstable due to large fluctations in the
branching weights.
- ``warmupsteps : int, default=40``
- ``blocks : int, default=400``
- ``steps : int, default=20``
Case ``qmc='opt'``
^^^^^^^^^^^^^^^^^^
Generate calculation elements for wavefunction optimization.
The parameter signature is identical to ``generate_opt_calculations``,
which depends on the value of ``method`` and ``minmethod``.
method : {'linear', 'cslinear'}, default='linear'
If ``'linear'``, use one of the versions of the linear method.
If ``'cslinear'``, use the correlated sampling linear method.
minmethod : {'quartic' , 'rescale' , 'linemin', 'adaptive', 'oneshift', 'sr_cg'}, default='quartic'
minwalkers : float, default=0.3
Minimum threshold to accept a parameter update based on the
ratio of wavefunction values between internal sub-iterations.
The value of ``minwalkers`` should be given in the range (0,1].
A small value of ``minwalkers`` will easily accept parameter
updates, likely resulting in an unstable run.
cost : {'energy', 'variance'} or tuple
If ``'energy'``, energy minization is performed.
If ``'variance'``, variance minimization is performed.
If length 2 tuple of floats ``(we, wv)``,
``cost = we*energy + wv*variance``
If length 3 tuple of floats ``(we, wv, wuv)``,
``cost = we*energy + wv*variance + wuv*unreweightedvariance``
When ``minmethod='oneshift'``, no cost function is being
minimized, but instead the parameter updates are determined
solely by ``minwalkers``.
cycles : int, default=12
Number of top level optimization iterations to perform.
Sets ``<loop max="cycles"/>``.
samples : int or None, default=None
If not ``None`` set the ``samples`` parameter, i.e. the total
number of VMC walker configurations to use in each optimization
cycle.
init_cycles : int, default=0
If ``init_cycles>0``, introduce a preceding optimization loop of
the same type (same ``minmethod``, ``cost`` and most other
parameters).
Sets ``<loop max="init_cycles"/>`` in this prior loop.
A few parameters can be set to different values from the
subsequent/main loop as listed below.
init_samples : int or None, default=None
If not ``None`` set the ``samples`` parameter, i.e. the total
number of VMC walker configurations to use in the preceding
optimization loop.
init_steps : int or None
If not ``None`` set the ``steps`` parameter in the preceding
optimization loop.
init_minwalkers : float, default=0.1
If not ``None`` set the ``minwalkers`` parameter in the preceding
optimization loop. Often set to a smaller value than in the
subsequent/main loop to allow more aggressive parameter updates
in hopes of a faster convergence to the general vicinity
of the cost minimum.
init_line_search : bool, default=False
Only applicable to ``minmethod='sr_cg'``, see below.
If ``True``, perform a linesearch along the direction of the
parameter gradient using the minimum cost to determine the
parameter stepsize.
init_sr_tau : float, default=0.1
Only applicable to ``minmethod='opt_sr'``, see below.
Set the ``sr_tau`` parameter appearing in the stochastic
reconfiguration projector.
Case ``qmc='opt'`` ``method={'linear', 'cslinear'}`` ``minmethod={'quartic', 'rescale', 'linemin'}``
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
minmethod : {'quartic', 'rescale', 'linemin'}, default='quartic'
Sets ``minmethod`` parameter. See QMCPACK manual.
usebuffer : bool, default=True
Sets ``usebuffer`` parameter. See QMCPACK manual.
exp0 : float, default=-6
Sets ``exp0`` parameter. See QMCPACK manual.
bigchange : float, default=10.0
Sets ``bigchange`` parameter. See QMCPACK manual.
alloweddifference : float, default=1e-4
Sets ``alloweddifference`` parameter. See QMCPACK manual.
stepsize : float, default=0.15
Sets ``stepsize`` parameter. See QMCPACK manual.
nstabilizers : int, default=1
Sets ``nstabilizers`` parameter. See QMCPACK manual.
var_cycles : int, default=0
If ``var_cycles>0``, introduce a preceding loop of variance
minmization to obtain a preconditioned starting point, e.g.
to stabilize subsequent energy minimization.
Sets ``<loop max="var_cycles"/>`` in this prior loop.
Uses all other parameters as set for the subsequent/main loop,
perhaps excepting ``samples``.
var_samples : int or None, default=None
If not ``None`` set the ``samples`` parameter, i.e. the total
number of VMC walker configurations to use in the preceding
variance minimization cycle.
Case ``qmc='opt'`` ``method='linear'`` ``minmethod='oneshift'``
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
Use the ``"oneshift"`` variant of the linear method, courtesy Ye Luo.
shift_i : float
Set the ``shift_i`` parameter. See QMCPACK manual.
shift_s : float
Set the ``shift_s`` parameter. See QMCPACK manual.
Case ``qmc='opt'`` ``method='linear'`` ``minmethod='adaptive'``
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
max_relative_change : float, default=10.0
Sets ``max_relative_change`` parameter. See QMCPACK manual.
max_param_change : float, default=0.3
Sets ``max_param_change`` parameter. See QMCPACK manual.
shift_i : float, default=0.01
Set the ``shift_i`` parameter. See QMCPACK manual.
shift_s : float, default=1.0
Set the ``shift_s`` parameter. See QMCPACK manual.
Case ``qmc='opt'`` ``method='linear'`` ``minmethod='sr_cg'``
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
Use a preliminary implementation of stochastic reconfiguration,
courtesy Cody Melton.
sr_tau : float, default=0.01
Set the ``sr_tau`` parameter, which is the timestep in the
stochastic reconfiguration projector.
sr_tolerance : float, default=0.001.
Set the ``sr_tolerance`` parameter. See QMCPACK manual.
sr_regularization : float, default=0.01.
Set the ``sr_regularization`` parameter. See QMCPACK manual.
linesearch : bool, default=False
Perform a correlated sampling linesearch to determine tau
automatically for each iteration.
If ``True``, the default for ``sr_tau`` is 0.1 instead.
"""
if optimize is not None:
optimize = bool(optimize)
J1_opt = optimize
J2_opt = optimize
J3_opt = optimize
J1k_opt = optimize
J2k_opt = optimize
jastrow_opt = optimize
multidet_opt = optimize
if jastrow_opt is not None:
jastrow_opt = bool(jastrow_opt)
if multidet_opt is not None:
multidet_opt = bool(multidet_opt)
if calculations is not None and qmc is not None:
self.error('cannot both provide calculation list ("calculations keyword") and request calculation generation ("qmc" keyword)')
# set driver version
if driver is not None:
self.set_driver(driver)
# remove system
remove_system = bool(remove_system)
if remove_system:
self.remove_physical_system()
# change system
change_system = bool(change_system)
if change_system:
assert system is not None
self.remove_physical_system()
self.incorporate_system(system)
# remove jastrows
remove_jastrows = bool(remove_jastrows)
if remove_jastrows:
self.remove_jastrows()
remove_J1 = bool(remove_J1)
if remove_J1:
self.remove_J1()
remove_J2 = bool(remove_J2)
if remove_J2:
self.remove_J2()
remove_J3 = bool(remove_J3)
if remove_J3:
self.remove_J3()
# generate jastrows
if J1 or J2 or J3 or J1k or J2k:
self.gen_jastrows(
J1 = J1 ,
J2 = J2 ,
J3 = J3 ,
J1_size = J1_size ,
J1_rcut = J1_rcut ,
J1_dr = J1_dr ,
J1_opt = J1_opt ,
J2_size = J2_size ,
J2_rcut = J2_rcut ,
J2_dr = J2_dr ,
J2_init = J2_init ,
J2_opt = J2_opt ,
J3_isize = J3_isize ,
J3_esize = J3_esize ,
J3_rcut = J3_rcut ,
J3_opt = J3_opt ,
J1_rcut_open = J1_rcut_open,
J2_rcut_open = J2_rcut_open,
J1k = J1k ,
J1k_kcut = J1k_kcut ,
J1k_symm = J1k_symm ,
J1k_opt = J1k_opt ,
J2k = J2k ,
J2k_kcut = J2k_kcut ,
J2k_symm = J2k_symm ,
J2k_opt = J2k_opt ,
system = system ,
)
# remove deteriminants
remove_determinants = bool(remove_determinants)
if remove_determinants:
self.remove('determinantset')
# remove multidet
remove_multidet = bool(remove_multidet)
if remove_multidet:
self.remove('multideterminant')
# set orbital file
if orbitals_h5 is not None:
self.set_orbitals_h5(orbitals_h5)
# set jastrow params
if self.has_jastrows():
if jastrow_opt is not None:
self.optimize_jastrows(jastrow_opt)
# set multidet params
if self.has_multidet():
if multidet_h5 is not None:
self.set_multidet_h5(multidet_h5)
if multidet_cutoff is not None:
self.set_multidet_params(cutoff=multidet_cutoff)
if multidet_opt is not None:
self.optimize_multidet(multidet_opt)
# set hamiltonian params
if pseudo_files is not None:
assert isinstance(pseudo_files,(dict,obj))
self.set_pseudo_files(**pseudo_files)
# remove calculations
remove_calculations = bool(remove_calculations)
if remove_calculations:
if 'simulation' in self and 'calculations' in self.simulation:
del self.simulation.calculations
else:
self.remove('calculation')
# set calculations
if qmc is not None:
self.gen_calculations(qmc,**gen_calcs)
elif len(gen_calcs)>0:
self.error('invalid keywords provided to the modify function:\n{}\n'.format(sorted(gen_calcs.keys()))+' Please see the documentation. If you are trying to generate qmc calculation sections, please provide the "qmc" keyword.')
elif calculations is not None:
self.simulation.calculations = make_collection(calculations).copy()
#end def modify
#end def bundle
#end def trace
#end def twist_average
#end class QmcpackInput
# base class for bundled qmcpack input
# not used on its own
[docs]
class BundledQmcpackInput(SimulationInput):
def __init__(self,inputs,filenames):
self.inputs = obj()
for input in inputs:
self.inputs.append(input)
#end for
self.filenames = filenames
#end def __init__
[docs]
def get_output_info(self,*requests):
outfiles = []
for index,input in self.inputs.items():
outfs = input.get_output_info('outfiles')
infile = self.filenames[index]
outfile= infile.rsplit('.',1)[0]+'.g'+str(index).zfill(3)+'.qmc'
outfiles.append(infile)
outfiles.append(outfile)
twfile = infile.rsplit('.',3)[0]+'.twist_info.dat'
outfiles.append(twfile)
for outf in outfs:
prefix,rest = outf.split('.',1)
outfiles.append(prefix+'.g'+str(index).zfill(3)+'.'+rest)
#end for
#end for
values = []
for req in requests:
if req=='outfiles':
values.append(outfiles)
else:
values.append(None)
#end if
#end for
if len(values)==1:
return values[0]
else:
return values
#end if
#end def get_output_info
#end def generate_filenames
[docs]
def write(self,filepath=None):
if filepath is not None and 'filenames' not in self:
infile = os.path.split(filepath)[1]
if not infile.endswith('.xml'):
infile+='.xml'
#end if
self.generate_filenames(infile)
#end if
if filepath is None:
c = ''
for i in range(len(self.inputs)):
c += self.filenames[i]+'\n'
#end for
return c
else:
path,file = os.path.split(filepath)
#if file!=self.filenames[-1]:
# self.error('main filenames do not match\n internal: '+self.filenames[-1]+'\n inputted: '+file)
##end if
c = ''
for i in range(len(self.inputs)):
input = self.inputs[i]
bfile = self.filenames[i]
c += bfile+'\n'
bfilepath = os.path.join(path,bfile)
input.write(bfilepath)
#end for
fobj = open(filepath,'w')
fobj.write(c)
fobj.close()
#end if
#end def write
#end class BundledQmcpackInput
[docs]
class TracedQmcpackInput(BundledQmcpackInput):
def __init__(self,quantity=None,values=None,input=None):
self.quantities = obj()
self.variables = obj()
self.inputs = obj()
self.filenames = None
if quantity is not None and values is not None and input is not None:
self.bundle_inputs(quantity,values,input)
#end if
#end def __init__
[docs]
def bundle_inputs(self,quantity,values,input):
range = len(self.inputs),len(self.inputs)+len(values)
self.quantities.append(obj(quantity=quantity,range=range))
for value in values:
inp = input.copy()
qhost = inp.get_host(quantity)
#print(qhost)
if qhost is not None:
if not isinstance(value,obj):
qhost[quantity] = value
else:
for k,v in value.items():
qhost[k] = v
#end for
#end if
else:
self.error('quantity '+quantity+' was not found in '+input.__class__.__name__)
#end if
self.variables.append(obj(quantity=quantity,value=value))
self.inputs.append(inp)
#end for
#end def bundle_inputs
[docs]
def generate_filenames(self,infile):
prefix,ext = infile.split('.',1)
if not ext.endswith('xml'):
ext+='.xml'
#end if
self.filenames = []
for i in range(len(self.variables)):
var = self.variables[i]
q = var.quantity
v = var.value
if isinstance(v,obj):
v = i
#end if
bfile = prefix+'.g'+str(i).zfill(3)+'.'+q+'_'+str(v)+'.'+ext
self.filenames.append(bfile)
#end if
self.filenames.append(prefix+'.in')
#end def generate_filenames
#end class TracedQmcpackInput
[docs]
class QmcpackInputTemplate(SimulationInputTemplate):
[docs]
def preprocess(self,contents,filepath=None):
if filepath is not None:
basepath,filename = os.path.split(filepath)
c = contents
contents=''
for line in c.splitlines():
if '<include' in line and '/>' in line:
tokens = line.replace('<include','').replace('/>','').split()
for token in tokens:
if token.startswith('href'):
include_file = token.replace('href','').replace('=','').replace('"','').strip()
include_path = os.path.join(basepath,include_file)
if os.path.exists(include_path):
with open(include_path, "r") as f:
icont = f.read() + "\n"
line = ''
for iline in icont.splitlines():
if '<?' not in iline:
line+=iline+'\n'
#end if
#end for
#end if
#end if
#end for
#end if
contents+=line+'\n'
#end for
#end if
return contents
#end def preprocess
#end def get_output_info
#end class QmcpackInputTemplate
[docs]
def generate_simulationcell(bconds='ppp',lr_dim_cutoff=15,lr_tol=None,lr_handler=None,system=None):
bconds = tuple(bconds)
sc = simulationcell(bconds=bconds)
periodic = 'p' in bconds
axes_valid = system is not None and len(system.structure.axes)>0
if periodic:
sc.lr_dim_cutoff = lr_dim_cutoff
if lr_tol is not None:
sc.lr_tol = lr_tol
#end if
if lr_handler is not None:
sc.lr_handler = lr_handler
#end if
if not axes_valid:
QmcpackInput.class_error('invalid axes in generate_simulationcell\nargument system must be provided\naxes of the structure must have non-zero dimension')
#end if
#end if
if axes_valid:
system.check_folded_system()
system.change_units('B')
structure = system.structure
if isinstance(structure,Jellium):
sc.rs = structure.rs()
sc.nparticles = system.particles.count_electrons()
else:
#setting the 'lattice' (cell axes) requires some delicate care
# qmcpack will fail if this is even 1e-10 off of what is in
# the wavefunction hdf5 file from pwscf
if structure.folded_structure is not None:
fs = structure.folded_structure
axes = np.array(pwscf_array_string(fs.axes).split(),dtype=float)
npe.reshape_inplace(axes, fs.axes.shape)
axes = np.dot(structure.tmatrix,axes)
if np.abs(axes-structure.axes).sum()>1e-5:
QmcpackInput.class_error('in generate_simulationcell\nsupercell axes do not match tiled version of folded cell axes\nyou may have changed one set of axes (super/folded) and not the other\nfolded cell axes:\n'+str(fs.axes)+'\nsupercell axes:\n'+str(structure.axes)+'\nfolded axes tiled:\n'+str(axes))
#end if
else:
axes = np.array(pwscf_array_string(structure.axes).split(),dtype=float)
npe.reshape_inplace(axes, structure.axes.shape)
#end if
structure.adjust_axes(axes)
sc.lattice = axes
#end if
#end if
return sc
#end def generate_simulationcell
[docs]
def generate_particlesets(electrons = 'e',
ions = 'ion0',
up = 'u',
down = 'd',
spinor = None,
system = None,
randomsrc = True,
hybrid_rcut = None,
hybrid_lmax = None,
):
if system is None:
QmcpackInput.class_error('generate_particlesets argument system must not be None')
#end if
ename = electrons
iname = ions
uname = up
dname = down
del electrons
del ions
del up
del down
system.check_folded_system()
system.change_units('B')
particles = system.particles
structure = system.structure
net_charge = system.net_charge
net_spin = system.net_spin
elns = particles.get_electrons()
ions = particles.get_ions()
eup = elns.up_electron
edn = elns.down_electron
use_spinor = spinor is not None and spinor
particleset_groups = []
if not use_spinor:
if eup.count > 0:
particleset_groups.append(group(name=uname,charge=-1,mass=eup.mass,size=eup.count))
#end if
if edn.count > 0:
particleset_groups.append(group(name=dname,charge=-1,mass=edn.mass,size=edn.count))
#end if
else:
ecount = eup.count+edn.count
if ecount>0:
particleset_groups.append(group(name=uname,charge=-1,mass=eup.mass,size=ecount))
#end if
#end if
particlesets = []
eps = particleset(
name = ename,
random = True,
groups = particleset_groups,
)
if use_spinor:
eps.spinor = True
#end if
particlesets.append(eps)
if len(ions)>0:
# maintain consistent order
ion_species,ion_counts = structure.order_by_species()
elem = structure.elem
pos = structure.pos
if randomsrc:
eps.randomsrc = iname
#end if
ips = particleset(name=iname)
# handle hybrid rep
hybridrep = hybrid_rcut is not None or hybrid_lmax is not None
if hybridrep:
hybrid_vars = (
('hybrid_rcut',hybrid_rcut),
('hybrid_lmax',hybrid_lmax),
)
for hvar,hval in hybrid_vars:
if not isinstance(hval,obj):
QmcpackInput.class_error('generate_particlesets argument "{0}" must be of type obj\nyou provided type: {1}\nwith value: {2}'.format(hvar,hval.__class__.__name__,hval))
#end if
if set(hval.keys())!=set(ion_species):
QmcpackInput.class_error('generate_particsets argument "{0}" is incorrect\none entry must be present for each atomic species\natomic species present in the simulation: {1}\nvalues provided for the following species: {2}'.format(hvar,sorted(ion_species),sorted(hval.keys())))
#end if
#end for
#end if
# make groups
groups = []
for ion_spec in ion_species:
ion = ions[ion_spec]
gpos = pos[elem==ion.name]
g = group(
name = ion.name,
charge = ion.charge,
valence = ion.charge,
atomicnumber = ion.protons,
mass = ion.mass,
position = gpos,
size = len(gpos)
)
if hybridrep:
g.lmax = hybrid_lmax[ion_spec]
g.cutoff_radius = hybrid_rcut[ion_spec]
#end if
groups.append(g)
#end for
ips.groups = make_collection(groups)
particlesets.append(ips)
#end if
particlesets = make_collection(particlesets)
return particlesets
#end def generate_particlesets
[docs]
def generate_sposets(type = None,
occupation = None,
spin_polarized = False,
nup = None,
ndown = None,
spo_up = 'spo_u',
spo_down = 'spo_d',
system = None,
sposets = None,
spindatasets = False,
spinor = None,
rotate = False,
):
ndn = ndown
if type is None:
QmcpackInput.class_error('cannot generate sposets\n type of sposet not specified')
#end if
if sposets is not None:
for spo in sposets:
spo.type = type
#end for
elif occupation=='slater_ground':
have_counts = not (nup is None or ndown is None)
if system is None and not have_counts:
QmcpackInput.class_error('cannot generate sposets in occupation mode {0}\n arguments nup & ndown or system must be given to generate_sposets'.format(occupation))
elif not have_counts:
elns = system.particles.get_electrons()
nup = elns.up_electron.count
ndn = elns.down_electron.count
else:
ndn = ndown
#end if
use_spinor = spinor is not None and spinor
if not use_spinor:
if not spin_polarized:
if nup==ndn:
sposets = [sposet(type=type,name='spo_ud',spindataset=0,size=nup)]
else:
sposets = [sposet(type=type,name=spo_up, spindataset=0,size=nup),
sposet(type=type,name=spo_down,spindataset=0,size=ndn)]
#end if
else:
sposets_list = []
if nup > 0:
sposets_list.append(sposet(type=type,name=spo_up, spindataset=0,size=nup))
#end if
if ndn > 0:
sposets_list.append(sposet(type=type,name=spo_down,spindataset=1,size=ndn))
#end if
sposets = sposets_list
#end if
else:
sposets = [sposet(type=type,name='spo_u',spindataset=0,size=nup+ndn)]
#end if
if not spindatasets:
for spo in sposets:
del spo.spindataset
#end for
#end if
else:
QmcpackInput.class_error('cannot generate sposets in occupation mode {0}\n generate_sposets currently supports the following occupation modes:\n slater_ground'.format(occupation))
#end if
if rotate:
rotated_sposets = []
for spo in sposets:
rot_spo = rotated_sposet(name='rot_'+spo.name, sposet=spo)
rotated_sposets.append(rot_spo)
#end for
sposets = rotated_sposets
#end if
return make_collection(sposets)
#end def generate_sposets
[docs]
def generate_sposet_builder(type,*args,**kwargs):
if type=='bspline' or type=='einspline':
return generate_bspline_builder(type,*args,**kwargs)
elif type=='heg':
return generate_heg_builder(*args,**kwargs)
else:
QmcpackInput.class_error('cannot generate sposet_builder\n sposet_builder of type {0} is unrecognized'.format(type))
#end if
#end def generate_sposet_builder
[docs]
def generate_bspline_builder(type = 'bspline',
meshfactor = 1.0,
precision = 'float',
twistnum = None,
twist = None,
sort = None,
version = '0.10',
truncate = False,
buffer = None,
spin_polarized = False,
hybridrep = None,
href = 'MISSING.h5',
rotate = False,
ions = 'ion0',
spo_up = 'spo_u',
spo_down = 'spo_d',
sposets = None,
system = None,
orbitals_cpu = None,
gpusharing = None,
spinor = None,
):
tilematrix = np.identity(3,dtype=int)
if system is not None:
tilematrix = system.structure.tilematrix()
#end if
bsb = bspline_builder(
type = type,
meshfactor = meshfactor,
precision = precision,
tilematrix = tilematrix,
href = href,
version = version,
truncate = truncate,
source = ions,
)
sposets = generate_sposets(
type = type,
occupation = 'slater_ground',
spin_polarized = spin_polarized,
system = system,
sposets = sposets,
spindatasets = True,
spinor = spinor,
rotate = rotate,
)
if not rotate:
bsb.sposets = sposets
else:
bsb.rotated_sposets = sposets
#end if
if sort is not None:
bsb.sort = sort
#end if
if truncate and buffer is not None:
bsb.buffer = buffer
#end if
if hybridrep is not None:
bsb.hybridrep = hybridrep
#end if
if twist is not None:
bsb.twistnum = system.structure.select_twist(twist)
elif twistnum is not None:
bsb.twistnum = twistnum
elif len(system.structure.kpoints)==1:
bsb.twistnum = 0
else:
bsb.twistnum = None
#end if
if orbitals_cpu is not None and orbitals_cpu:
bsb.gpu = False
#end if
if gpusharing is not None:
bsb.gpusharing = gpusharing
#end if
return bsb
#end def generate_bspline_builder
[docs]
def generate_heg_builder(twist = None,
spin_polarized = False,
spo_up = 'spo_u',
spo_down = 'spo_d',
sposets = None,
system = None
):
type = 'heg'
hb = heg_builder(
type = type,
sposets = generate_sposets(
type = type,
occupation = 'slater_ground',
spin_polarized = spin_polarized,
system = system,
sposets = sposets
)
)
if twist is not None:
hb.twist = tuple(twist)
#end if
return hb
#end def generate_heg_builder
[docs]
def partition_sposets(sposet_builder,partition,partition_meshfactors=None):
ssb = sposet_builder
spos_in =ssb.sposets
del ssb.sposets
if isinstance(partition,(dict,obj)):
partition_indices = sorted(partition.keys())
partition_contents = partition
else:
partition_indices = list(partition)
partition_contents = None
#end if
if partition_meshfactors is not None:
if partition_contents is None:
partition_contents = obj()
for p in partition_indices:
partition_contents[p] = obj()
#end for
#end if
for p,mf in zip(partition_indices,partition_meshfactors):
partition_contents[p].meshfactor = mf
#end for
#end if
# partition each spo in the builder and create a corresponding composite spo
comp_spos = []
part_spos = []
for spo in spos_in.list():
part_spo_names = []
part_ranges = partition_indices+[spo.size]
for i in range(len(partition_indices)):
index_min = part_ranges[i]
index_max = part_ranges[i+1]
if index_min>spo.size:
break
elif index_max>spo.size:
index_max = spo.size
#end if
part_spo_name = spo.name+'_'+str(index_min)
part_spo = sposet(**spo)
part_spo.name = part_spo_name
if index_min==0:
part_spo.size = index_max
else:
part_spo.index_min = index_min
part_spo.index_max = index_max
del part_spo.size
#end if
if partition_contents is not None:
part_spo.set(**partition_contents[index_min])
#end if
part_spos.append(part_spo)
part_spo_names.append(part_spo_name)
#end for
comp_spo = sposet(
name = spo.name,
size = spo.size,
spos = part_spo_names,
)
comp_spos.append(comp_spo)
#end for
ssb.sposets = make_collection(part_spos)
cssb = composite_builder(
type = 'composite',
sposets = make_collection(comp_spos),
)
return [ssb,cssb]
#end def partition_sposets
[docs]
def generate_determinantset(up = 'u',
down = 'd',
spo_up = 'spo_u',
spo_down = 'spo_d',
spin_polarized = False,
delay_rank = None,
det_batch = None,
matrix_inv_cpu = None,
system = None,
spinor = None,
rotate = False,
):
if system is None:
QmcpackInput.class_error('generate_determinantset argument system must not be None')
#end if
elns = system.particles.get_electrons()
nup = elns.up_electron.count
ndn = elns.down_electron.count
use_spinor = spinor is not None and spinor
if not spin_polarized and nup==ndn and not use_spinor:
spo_u = 'spo_ud'
spo_d = 'spo_ud'
else:
spo_u = spo_up
spo_d = spo_down
#end if
if rotate:
spo_u = 'rot_'+spo_u
spo_d = 'rot_'+spo_d
#end if
determinants_list = []
if not use_spinor:
if nup > 0:
determinants_list.append(
determinant(
id = 'updet',
group = up,
sposet = spo_u,
size = nup
)
)
#end if
if ndn > 0:
determinants_list.append(
determinant(
id = 'downdet',
group = down,
sposet = spo_d,
size = ndn
)
)
#end if
else:
if nup+ndn > 0:
determinants_list.append(
determinant(
id = 'updet',
group = up,
sposet = spo_u,
size = nup+ndn,
)
)
#end if
#end if
dset = determinantset(
slaterdeterminant = slaterdeterminant(
determinants = collection(*determinants_list)
)
)
if delay_rank is not None:
dset.slaterdeterminant.delay_rank = delay_rank
#end if
if det_batch is not None:
dset.slaterdeterminant.batch = det_batch
#end if
if matrix_inv_cpu is not None and matrix_inv_cpu:
dset.slaterdeterminant.matrix_inverter = 'host'
#end if
return dset
#end def generate_determinantset
[docs]
def check_excitation_type(excitation):
# Possible spin channels or spin states
exc_spins = obj(
up = 1, # 'up'
down = 2, # 'down'
singlet = 3, # 'singlet'
triplet = 4, # 'triplet'
)
# Possible orbital excitation types
exc_types = obj(
band = 1, # '0 45 3 46' # Type 1
energy = 2, # '-215 +216' # Type 2
kpoint = 3, # 'L vb F cb' # Type 3
lowest = 4, # 'lowest' # Type 4
)
exc_spin = None
exc_type = None
# Check that 'excitation' is correctly formated
format_failed = False
# Extract elements form excitation
if not isinstance(excitation,(tuple,list)) or len(excitation) != 2:
format_failed = True
else:
exc1,exc2 = excitation
if not isinstance(exc1,str) or not isinstance(exc2,str):
format_failed = True
#end if
#end if
# Check first element
if not format_failed:
if exc1.lower() not in ('up','down','singlet','triplet'):
format_failed = True
else:
exc_spin = exc_spins[exc1.lower()]
#end if
#end if
# Check second element
if not format_failed:
if any(substr in exc2.lower() for substr in ('vb','cb','lowest')):
if exc2.lower()=='lowest':
exc_type = exc_types.lowest
elif len(exc2.split())!=4:
format_failed = True
else:
exc_type = exc_types.kpoint
#end if
else:
tmp = None
try:
tmp = np.array(exc2.split(),dtype=int)
except:
format_failed = True
#end try
if tmp is not None:
if len(tmp)==4:
# '0 45 3 46'
if not tmp[0]>=0 or not tmp[1]>=0 or not tmp[2]>=0 or not tmp[3]>=0:
format_failed = True
#end if
exc_type = exc_types.band
elif len(tmp)==2:
# '-215 +216'
if not tmp[0]<0 or not tmp[1]>0:
format_failed = True
#end if
exc_type = exc_types.energy
else:
format_failed = True
#end if
#end if
#end if
#end if
if format_failed:
msg = 'excitation must be a tuple or list with with two elements.\n'
msg += 'The first element must be either "up", "down", "singlet", or "triplet"\n'
msg += 'and the second element must be a band format (e.g. "0 45 3 46"),\n'
msg += 'energy format (e.g. "-215 +216"), kpoint format (e.g. "L vb F cb"),\n'
msg += 'or lowest format (e.g. "lowest").\n'
msg += 'You Provided: {0}'
msg = msg.format(excitation)
QmcpackInput.class_error(msg)
#end if
return exc_spin,exc_type,exc_spins,exc_types,exc1,exc2
#end def check_excitation_type
[docs]
def generate_determinantset_old(type = 'bspline',
meshfactor = 1.0,
precision = 'float',
twistnum = None,
twist = None,
spin_polarized = False,
hybridrep = None,
source = 'ion0',
href = 'MISSING.h5',
excitation = None,
delay_rank = None,
gpusharing = None,
system = None,
spinor = None,
):
if system is None:
QmcpackInput.class_error('generate_determinantset argument system must not be None')
#end if
elns = system.particles.get_electrons()
down_spin = 0
if spin_polarized:
down_spin=1
#end if
tilematrix = np.identity(3,dtype=int)
if system is not None:
tilematrix = system.structure.tilematrix()
#end if
use_spinor = spinor is not None and spinor
nup = elns.up_electron.count
ndn = elns.down_electron.count
determinants_list = []
if not use_spinor:
if nup > 0:
determinants_list.append(
determinant(
id = 'updet',
size = nup,
occupation=section(mode='ground',spindataset=0)
),
)
#end if
if ndn > 0:
determinants_list.append(
determinant(
id = 'downdet',
size = ndn,
occupation=section(mode='ground',spindataset=down_spin)
)
)
#end if
else:
if nup+ndn > 0:
determinants_list.append(
determinant(
id = 'updet',
size = nup+ndn,
occupation=section(mode='ground',spindataset=0)
),
)
#end if
#end if
dset = determinantset(
type = type,
meshfactor = meshfactor,
precision = precision,
tilematrix = tilematrix,
href = href,
source = source,
slaterdeterminant = slaterdeterminant(
determinants = collection(*determinants_list)
)
)
if twist is not None:
dset.twistnum = system.structure.select_twist(twist)
elif twistnum is not None:
dset.twistnum = twistnum
elif len(system.structure.kpoints)==1:
dset.twistnum = 0
else:
dset.twistnum = None
#end if
if hybridrep is not None:
if hybridrep=='yes' or hybridrep=='no':
dset.hybridrep = hybridrep
else:
dset.hybridrep = yesno_dict[hybridrep]
#end if
#end if
if delay_rank is not None:
dset.slaterdeterminant.delay_rank = delay_rank
#end if
if gpusharing is not None:
dset.gpusharing = gpusharing
#end if
if excitation is not None:
exc_spin,exc_type,exc_spins,exc_types,exc1,exc2 = check_excitation_type(excitation)
if exc_spin==exc_spins.up:
sdet = dset.get('updet')
elif exc_spin==exc_spins.down:
sdet = dset.get('downdet')
elif exc_spin in (exc_spins.singlet,exc_spins.triplet):
# Are there an equal number of up and down electrons?
# If no, then exit. Currently, singlet and triplet
# excitations are assumed to have ms = 0.
if elns.down_electron.count != elns.up_electron.count:
QmcpackInput.class_error('The \'singlet\' and \'triplet\' excitation types currently assume number of up and down electrons is the same for the reference ground state. Otherwise, one should use \'up\' or \'down\' types.\nFor your system: Nup={} and Ndown={}.\nWe plan to expand to additional cases in the future.'.format(elns.up_electron.count,elns.down_electron.count))
#end if
coeff_sign = ''
if exc_spin==exc_spins.triplet:
coeff_sign = '-'
#end if
if down_spin:
sposet_list = [sposet(name = 'spo_u',
spindataset = 0,
size = elns.up_electron.count+1,
occupation = section(mode='ground'),
coefficient = section(size=90,spindataset=0),
spos = ''
),
sposet(name = 'spo_d',
spindataset = 1,
size = elns.up_electron.count+1,
occupation = section(mode='ground'),
coefficient = section(spindataset=1),
spos = ''
)]
else:
sposet_list = [sposet(name = 'spo_ud',
spindataset = 0,
size = elns.up_electron.count+1,
occupation = section(mode='ground'),
coefficient = section(spindataset=0),
spos = ''
)]
#end if
dset = determinantset(
type = type,
meshfactor = meshfactor,
precision = precision,
tilematrix = tilematrix,
twistnum = twistnum,
href = href,
source = source,
sposets = sposet_list,
multideterminant = multideterminant(
optimize = 'no',
spo_up='spo_u' if down_spin else 'spo_ud',
spo_dn='spo_d' if down_spin else 'spo_ud',
detlist = detlist(
size = '1',
type = 'CSF',
nca = '0',
ncb = '0',
nea = elns.up_electron.count,
neb = elns.down_electron.count,
cutoff = '0.001',
csf = csf(
id = 'CSF_0',
exctLvl = '1',
coeff = '1.0',
coeff_real = '1.0',
qchem_coeff = '1.0',
dets = collection(
det(
id='csf_00',
coeff='0.70710678118654752440',
),
det(
id='csf_01',
coeff=coeff_sign+'0.70710678118654752440',
),
)
)
),
)
)
if exc_type in (exc_types.energy,exc_types.lowest):
nup = elns.up_electron.count
if exc_type==exc_types.lowest:
exc_orbs = [nup,nup+1]
else:
# assume excitation of form '-216 +217' or '-216 217'
exc_orbs = np.array(exc2.split(),dtype=int)
exc_orbs[0] *= -1
#end if
for sp in dset.sposets:
sp.size=exc_orbs[1]
#end for
dset.multideterminant.detlist.nstates = exc_orbs[1]
dset.multideterminant.detlist.csf.occ = '2'*nup+'0'*(exc_orbs[1]-nup-1)+'1'
dset.multideterminant.detlist.csf.occ = dset.multideterminant.detlist.csf.occ[:exc_orbs[0]-1]+'1'+dset.multideterminant.detlist.csf.occ[exc_orbs[0]:]
dset.multideterminant.detlist.csf.dets[0].alpha = '1'*(exc_orbs[0]-1)+'0'+'1'*(nup-exc_orbs[0])+'0'*(exc_orbs[1]-nup-1)+'1'
dset.multideterminant.detlist.csf.dets[0].beta = '1'*nup+'0'*(exc_orbs[1]-nup)
dset.multideterminant.detlist.csf.dets[1].alpha = '1'*nup+'0'*(exc_orbs[1]-nup)
dset.multideterminant.detlist.csf.dets[1].beta = '1'*(exc_orbs[0]-1)+'0'+'1'*(nup-exc_orbs[0])+'0'*(exc_orbs[1]-nup-1)+'1'
elif exc_type == exc_types.kpoint:
QmcpackInput.class_error('{} excitation is not yet available for kpoint type'.format(exc1))
else:
QmcpackInput.class_error('{} excitation is not yet available for band type'.format(exc1))
#end if
return dset
#end if
occ = sdet.occupation
occ.pairs = 1
occ.mode = 'excited'
occ.contents = '\n'+exc2+'\n'
# add new input format
if exc_type == exc_types.kpoint:
# assume excitation of form 'gamma vb k cb' or 'gamma vb-1 k cb+1'
excitation = exc2.upper().split(' ')
if len(excitation) == 4:
k_1, band_1, k_2, band_2 = excitation
else:
QmcpackInput.class_error('excitation with vb-cb band format works only with special k-points')
#end if
vb = int(sdet.size / np.abs(np.linalg.det(tilematrix))) -1 # Separate for each spin channel
cb = vb+1
# Convert band_1, band_2 to band indexes
bands = [band_1, band_2]
for bnum, b in enumerate(bands):
b = b.lower()
if 'cb' in b:
if '-' in b:
b = b.split('-')
bands[bnum] = cb - int(b[1])
elif '+' in b:
b = b.split('+')
bands[bnum] = cb + int(b[1])
else:
bands[bnum] = cb
#end if
elif 'vb' in b:
if '-' in b:
b = b.split('-')
bands[bnum] = vb - int(b[1])
elif '+' in b:
b = b.split('+')
bands[bnum] = vb + int(b[1])
else:
bands[bnum] = vb
#end if
else:
QmcpackInput.class_error('{0} in excitation has the wrong formatting'.format(b))
#end if
#end for
band_1, band_2 = bands
# Convert k_1 k_2 to wavevector indexes
structure = system.structure.get_smallest().copy()
structure.change_units('A')
kpath = get_kpath(structure=structure)
kpath_label = np.array(kpath['explicit_kpoints_labels'])
kpath_rel = kpath['explicit_kpoints_rel']
k1_in = k_1
k2_in = k_2
if k_1 in kpath_label and k_2 in kpath_label:
k_1 = kpath_rel[np.where(kpath_label == k_1)][0]
k_2 = kpath_rel[np.where(kpath_label == k_2)][0]
#kpts = nscf.input.k_points.kpoints
kpts = structure.kpoints_unit()
found_k1 = False
found_k2 = False
for knum, k in enumerate(kpts):
if np.isclose(k_1, k).all():
k_1 = knum
found_k1 = True
#end if
if np.isclose(k_2, k).all():
k_2 = knum
found_k2 = True
#end if
#end for
if not found_k1 or not found_k2:
QmcpackInput.class_error('Requested special kpoint is not in the tiled cell\nRequested "{}", present={}\nRequested "{}", present={}\nAvailable kpoints: {}'.format(k1_in,found_k1,k2_in,found_k2,sorted(set(kpath_label))))
#end if
else:
QmcpackInput.class_error('Excitation wavevectors are not found in the kpath\nlabels requested: {} {}\nlabels present: {}'.format(k_1,k_2,sorted(set(kpath_label))))
#end if
#Write everything in band (ti,bi) format
occ.contents = '\n'+str(k_1)+' '+str(band_1)+' '+str(k_2)+' '+str(band_2)+'\n'
occ.format = 'band'
elif exc_type == exc_types.energy:
# assume excitation of form '-216 +217'
occ.format = 'energy'
elif exc_type == exc_types.lowest: # Type 4
occ.format = 'energy'
if exc_spin == exc_spins.up:
nel = elns.up_electron.count
else:
nel = elns.down_electron.count
#end if
excitation = '-{} +{}'.format(nel,nel+1)
occ.contents = '\n'+excitation+'\n'
else: #Type 1
# assume excitation of form '6 36 6 37'
occ.format = 'band'
#end if
#end if
return dset
#end def generate_determinantset_old
[docs]
def generate_hamiltonian(name = 'h0',
type = 'generic',
electrons = 'e',
ions = 'ion0',
wavefunction = 'psi0',
pseudos = None,
algorithm = None,
dla = None,
format = 'xml',
estimators = None,
system = None,
wf_elem = None,
interactions = 'default',
):
if system is None:
QmcpackInput.class_error('generate_hamiltonian argument system must not be None')
#end if
ename = electrons
iname = ions
wfname = wavefunction
ppfiles = pseudos
del electrons
del ions
del pseudos
del wavefunction
particles = system.particles
if particles.count_electrons()==0:
QmcpackInput.class_error('cannot generate hamiltonian, no electrons present')
#end if
pairpots = []
if interactions is not None:
pairpots.append(coulomb(name='ElecElec',type='coulomb',source=ename,target=ename))
if particles.count_ions()>0:
pairpots.append(coulomb(name='IonIon',type='coulomb',source=iname,target=iname))
ions = particles.get_ions()
if not system.pseudized:
pairpots.append(coulomb(name='ElecIon',type='coulomb',source=iname,target=ename))
else:
if ppfiles is None or len(ppfiles)==0:
QmcpackInput.class_error('cannot generate hamiltonian\n system is pseudized, but no pseudopotentials have been provided\n please provide pseudopotential files via the pseudos keyword')
#end if
if isinstance(ppfiles,list):
pplist = ppfiles
ppfiles = obj()
for pppath in pplist:
if '/' in pppath:
ppfile = pppath.split('/')[-1]
else:
ppfile = pppath
#end if
element = ppfile.split('.')[0]
if len(element)>2:
element = element[0:2]
#end if
ppfiles[element] = pppath
#end for
#end if
pseudos = collection()
for ion in ions:
label = ion.name
iselem, element = Elements.is_element(ion.name, return_element=True)
if label in ppfiles:
ppfile = ppfiles[label]
elif element.symbol in ppfiles:
ppfile = ppfiles[element.symbol]
else:
QmcpackInput.class_error('pseudos provided to generate_hamiltonian are incomplete\n a pseudopotential for ion of type {0} is missing\n pseudos provided:\n{1}'.format(ion.name,str(ppfiles)))
#end if
pseudos.add(pseudo(elementtype=label,href=ppfile))
#end for
pp = pseudopotential(name='PseudoPot',type='pseudo',source=iname,wavefunction=wfname,format=format,pseudos=pseudos)
if algorithm is not None:
pp.algorithm = algorithm
#end if
if dla is not None:
pp.dla = dla
#end if
pairpots.append(pp)
#end if
#end if
#end if
ests = []
if estimators is not None:
for estimator in estimators:
if isinstance(estimator,QIxml):
estimator = estimator.copy()
#end if
est=estimator
if isinstance(estimator,str):
estname = estimator.lower().replace(' ','_').replace('-','_').replace('__','_')
if estname=='mpc':
pairpots.append(mpc(name='MPC',type='MPC',ecut=60.0,source=ename,target=ename,physical=False))
est = None
elif estname=='chiesa':
est = chiesa(name='KEcorr',type='chiesa',source=ename,psi=wfname)
elif estname=='localenergy':
est = localenergy(name='LocalEnergy')
elif estname=='skall':
est = skall(name='SkAll',type='skall',source=iname,target=ename,hdf5=True)
elif estname=='energydensity':
est = energydensity(
type='EnergyDensity',name='EDvoronoi',dynamic=ename,static=iname,
spacegrid = spacegrid(coord='voronoi')
)
elif estname=='pressure':
est = pressure(type='Pressure')
else:
QmcpackInput.class_error('estimator '+estimator+' has not yet been enabled in generate_basic_input')
#end if
elif not isinstance(estimator,QIxml):
QmcpackInput.class_error('generate_hamiltonian received an invalid estimator\n an estimator must either be a name or a QIxml object\n inputted estimator type: {0}\n inputted estimator contents: {1}'.format(estimator.__class__.__name__,estimator))
elif isinstance(estimator,energydensity):
est.set_optional(
type = 'EnergyDensity',
dynamic = ename,
static = iname,
)
elif isinstance(estimator,dm1b):
est = process_dm1b_estimator(estimator,wfname,wf_elem=wf_elem)
#end if
if est is not None:
ests.append(est)
#end if
#end for
#end if
estimators = ests
hmltn = hamiltonian(
name = name,
type = type,
target = ename
)
if len(pairpots)>0:
hmltn.pairpots = make_collection(pairpots)
#end if
if len(estimators)>0:
hmltn.estimators = make_collection(estimators)
#end if
return hmltn
#end def generate_hamiltonian
[docs]
def generate_estimators_batched(estimators,
electrons = 'e',
ions = 'ion0',
wavefunction = 'psi0',
wf_elem = None,
):
assert len(estimators)>0
ename = electrons
iname = ions
wfname = wavefunction
del electrons
del ions
ests = []
for estimator in estimators:
if isinstance(estimator,QIxml):
estimator = estimator.copy()
#end if
est = estimator
if isinstance(estimator,str):
estname = estimator.lower().replace(' ','_').replace('-','_').replace('__','_')
#if estname=='chiesa':
# est = chiesa(name='KEcorr',type='chiesa',source=ename,psi=wfname)
#else:
QmcpackInput.class_error('estimator '+estimator+' has not yet been enabled in generate_estimators')
##end if
elif not isinstance(estimator,QIxml):
QmcpackInput.class_error('generate_estimators received an invalid estimator\n an estimator must either be a name or a QIxml object\n inputted estimator type: {0}\n inputted estimator contents: {1}'.format(estimator.__class__.__name__,estimator))
elif isinstance(estimator,momentum):
estimator.type = 'MomentumDistribution'
elif isinstance(estimator,onebodydensitymatrices):
est = process_dm1b_estimator(estimator,wfname,wf_elem)
#end if
if est is not None:
ests.append(est)
#end if
#end for
estimators = make_collection(ests)
return estimators
#end def generate_estimators_batched
[docs]
def process_dm1b_estimator(dm,wfname,wf_elem):
reuse = False
if 'reuse' in dm:
reuse = bool(dm.reuse)
del dm.reuse
#end if
basis = []
builder = None
maxed = False
wf = wf_elem
if reuse and 'basis' in dm and isinstance(dm.basis,sposet):
spo = dm.basis
# get sposet size
if 'size' in dm.basis:
size = spo.size
del spo.size
elif 'index_max' in dm.basis:
size = spo.index_max
del spo.index_max
else:
QmcpackInput.class_error('cannot generate estimator dm1b\n basis sposet provided does not have a "size" attribute')
#end if
try:
# get sposet from wavefunction
dets = wf.get('determinant')
det = dets.get_single()
if 'sposet' in det:
rsponame = det.sposet
else:
rsponame = det.id
#end if
builders = wf.get('sposet_builders')
if builders is None:
builders = [wf.sposet_builders.bspline]
#end if
rspo = None
for bld in builders:
if rsponame in bld.sposets:
builder = bld
rspo = bld.sposets[rsponame]
break
#end if
#end for
basis.append(rsponame)
# adjust current sposet
spo.index_min = rspo.size
spo.index_max = size
maxed = rspo.size>=size
except Exception as e:
msg = 'cannot generate estimator dm1b\n '
if wf is None:
QmcpackInput.class_error(msg+'wavefunction {0} not found'.format(wfname))
elif dets is None or det is None:
QmcpackInput.class_error(msg+'determinant not found')
elif builders is None:
QmcpackInput.class_error(msg+'sposet_builders not found')
elif rspo is None:
QmcpackInput.class_error(msg+'sposet {0} not found'.format(rsponame))
else:
QmcpackInput.class_error(msg+'cause of failure could not be determined\n see the following error message:\n{0}'.format(e))
#end if
#end if
#end if
# put the basis sposet in the appropriate builder
if isinstance(dm.basis,sposet) and not maxed:
spo = dm.basis
del dm.basis
if 'type' not in spo:
QmcpackInput.class_error('cannot generate estimator dm1b\n basis sposet provided does not have a "type" attribute')
#end if
if 'name' not in spo:
spo.name = 'spo_dm'
#end if
builders = wf.get('sposet_builders')
if spo.type not in builders:
bld = generate_sposet_builder(spo.type,sposets=[spo])
builders.add(bld)
else:
bld = builders[spo.type]
bld.sposets.add(spo)
#end if
basis.append(spo.name)
#end if
dm.basis = basis
dm.incorporate_defaults(elements=False,overwrite=False,propagate=False)
return dm
#end def process_dm1b_estimator
[docs]
def generate_jastrows(jastrows,system=None,return_list=False,check_ions=False):
jin = []
have_ions = True
if check_ions and system is not None:
have_ions = system.particles.count_ions()>0
#end if
if isinstance(jastrows,str):
jorders = set(jastrows.replace('generate',''))
if '1' in jorders and have_ions:
jterm = generate_jastrow('J1','bspline',8,system=system)
#end if
if '2' in jorders:
jterm = generate_jastrow('J2','bspline',8,system=system)
#end if
if '3' in jorders and have_ions:
jterm = generate_jastrow('J3','polynomial',3,3,4.0,system=system)
#end if
if 'k' in jorders:
kcut = max(system.rpa_kf())
nksh = system.structure.count_kshells(kcut)
jterm = generate_kspace_jastrow(kc1=0, kc2=kcut, nk1=0, nk2=nksh)
#end if
jin.append(jterm)
if len(jin)==0:
QmcpackInput.class_error('jastrow generation requested but no orders specified (1,2,and/or 3)')
#end if
else:
jset = set(['J1','J2','J3'])
for jastrow in jastrows:
if isinstance(jastrow,QIxml):
jin.append(jastrow)
elif isinstance(jastrow,dict) or isinstance(jastrow,obj):
jdict = dict(**jastrow)
if 'type' not in jastrow:
QmcpackInput.class_error("could not determine jastrow type from input\n field 'type' must be 'J1', 'J2', or 'J3'\n object you provided: "+str(jastrow))
#end if
jtype = jdict['type']
if jtype not in jset:
QmcpackInput.class_error("invalid jastrow type provided\n field 'type' must be 'J1', 'J2', or 'J3'\n object you provided: "+str(jdict))
#end if
del jdict['type']
if 'system' in jdict:
jsys = jdict['system']
del jdict['system']
else:
jsys = system
#end if
jterm = generate_jastrow(jtype,system=jsys,**jdict)
if jterm is not None:
jin.append(jterm)
#end if
del jtype
del jsys
elif jastrow[0] in jset:
jin.append(generate_jastrow(jastrow,system=system))
else:
QmcpackInput.class_error('starting jastrow unrecognized:\n '+str(jastrow))
#end if
#end for
#end if
if return_list:
return jin
else:
wf = wavefunction(jastrows=jin)
wf.pluralize()
return wf.jastrows
#end if
#end def generate_jastrows
[docs]
def generate_jastrows_alt(
J1 = False,
J2 = False,
J3 = False,
J1_size = None,
J1_rcut = None,
J1_dr = 0.5,
J1_opt = None,
J2_size = None,
J2_rcut = None,
J2_dr = 0.5,
J2_init = 'zero',
J2_opt = True,
J3_isize = 3,
J3_esize = 3,
J3_rcut = 5.0,
J3_opt = None,
J1_rcut_open = 5.0,
J2_rcut_open = 10.0,
J1k = False,
J1k_kcut = 5.0,
J1k_symm = 'crystal',
J1k_opt = None,
J2k = False,
J2k_kcut = 5.0,
J2k_symm = 'crystal',
J2k_opt = None,
system = None,
):
if system is None:
QmcpackInput.class_error('input variable "system" is required to generate jastrows','generate_jastrows_alt')
elif system.structure.units!='B':
system = system.copy()
system.structure.change_units('B')
#end if
openbc = system.structure.is_open()
natoms = system.particles.count_ions()
nelec = system.particles.count_electrons()
jastrows = []
J2 |= J3
J1 |= J2
if not J1 and not J2 and not J3:
J1 = True
J2 = True
#end if
rwigner = None
if J1:
if natoms<1:
QmcpackInput.class_error('One-body Jastrow (J1) requested, but no atoms are present','generate_jastrows_alt')
#end if
if J1_rcut is None:
if openbc:
J1_rcut = J1_rcut_open
else:
if rwigner is None:
rwigner = system.structure.rwigner(1)
#end if
J1_rcut = rwigner
#end if
#end if
if J1_size is None:
J1_size = int(np.ceil(J1_rcut/J1_dr))
#end if
kw = obj(system=system)
if J1_opt is not None:
kw.opt = bool(J1_opt)
J = generate_jastrow('J1','bspline',J1_size,J1_rcut,**kw)
jastrows.append(J)
#end if
if J2:
if nelec<2:
QmcpackInput.class_error('Two-body Jastrow (J2) requested, but not enough electrons are present.\nElectrons required: 2 or more\nElectrons present: {}'.format(nelec),'generate_jastrows_alt')
#end if
if J2_rcut is None:
if openbc:
J2_rcut = J2_rcut_open
else:
if rwigner is None:
rwigner = system.structure.rwigner(1)
#end if
J2_rcut = rwigner
#end if
#end if
if J2_size is None:
J2_size = int(np.ceil(J2_rcut/J2_dr))
#end if
kw = obj(system=system)
if J2_opt is not None:
kw.opt = bool(J2_opt)
J = generate_jastrow('J2','bspline',J2_size,J2_rcut,init=J2_init,**kw)
jastrows.append(J)
#end if
if J3:
if natoms<1 or nelec<2:
QmcpackInput.class_error('Three-body Jastrow (J3) requested, but not enough particles are present.\nAtoms required: 1 or more\nElectrons required: 2 or more\nAtoms present: {}\nElectrons present: {}'.format(natoms,nelec),'generate_jastrows_alt')
#end if
if not openbc:
if rwigner is None:
rwigner = system.structure.rwigner(1)
#end if
J3_rcut = min(J3_rcut,rwigner)
#end if
kw = obj(system=system)
if J3_opt is not None:
kw.opt = bool(J3_opt)
J = generate_jastrow('J3','polynomial',J3_esize,J3_isize,J3_rcut,**kw)
jastrows.append(J)
#end if
if J1k or J2k:
Jk_inp = obj()
if J1k:
Jk_inp.kc1 = J1k_kcut
Jk_inp.symm1 = J1k_symm
if J1k_opt is not None:
Jk_inp.opt1 = bool(J1k_opt)
#end if
if J2k:
Jk_inp.kc2 = J2k_kcut
Jk_inp.symm2 = J2k_symm
if J2k_opt is not None:
Jk_inp.opt2 = bool(J2k_opt)
#end if
J = generate_kspace_jastrow(**Jk_inp)
jastrows.append(J)
#end if
return jastrows
#end def generate_jastrows_alt
[docs]
def generate_jastrow(descriptor,*args,**kwargs):
keywords = set(['function','size','rcut','elements','coeff','cusp','ename',
'iname','spins','density','Buu','Bud','opt','system','isize','esize','init'])
if not 'system' in kwargs:
kwargs['system'] = None
#end if
system = kwargs['system']
del kwargs['system']
if system is not None:
system.change_units('B')
#end if
if isinstance(descriptor,str):
descriptor = [descriptor]
#end if
ikw=0
for i in range(len(descriptor)):
if descriptor[i] in keywords:
break
#end if
ikw += 1
#end for
dargs = descriptor[1:ikw]
if len(dargs)>0:
args = dargs
#end if
for i in range(ikw,len(descriptor),2):
d = descriptor[i]
if isinstance(d,str):
if d in keywords:
kwargs[d] = descriptor[i+1]
else:
QmcpackInput.class_error('keyword {0} is unrecognized\n valid options are: {1}'.format(d,str(keywords)),'generate_jastrow')
#end if
#end if
#end for
kwargs['system'] = system
jtype = descriptor[0]
if jtype=='J1':
jastrow = generate_jastrow1(*args,**kwargs)
elif jtype=='J2':
jastrow = generate_jastrow2(*args,**kwargs)
elif jtype=='J3':
jastrow = generate_jastrow3(*args,**kwargs)
else:
QmcpackInput.class_error('jastrow type unrecognized: '+jtype)
#end if
return jastrow
#end def generate_jastrow
[docs]
def generate_jastrow1(function='bspline',size=8,rcut=None,coeff=None,cusp=0.,ename='e',iname='ion0',elements=None,system=None,opt=None,**elemargs):
noelements = elements is None
nosystem = system is None
noelemargs = len(elemargs)==0
isopen = False
isperiodic = False
rwigner = 1e99
if noelements and nosystem and noelemargs:
QmcpackInput.class_error('must specify elements or system','generate_jastrow1')
#end if
if noelements:
elements = []
#end if
if not nosystem:
elements.extend(list(set(system.structure.elem)))
isopen = system.structure.is_open()
isperiodic = system.structure.is_periodic()
if not isopen and isperiodic:
rwigner = system.structure.rwigner()
#end if
#end if
if not noelemargs:
elements.extend(elemargs.keys())
#end if
# remove duplicate elements
eset = set()
elements = [ e for e in elements if e not in eset and not eset.add(e) ]
corrs = []
for i in range(len(elements)):
element = elements[i]
if cusp == 'Z':
QmcpackInput.class_error('need to implement Z cusp','generate_jastrow1')
else:
lcusp = cusp
#end if
lrcut = rcut
lcoeff = size*[0]
if coeff is not None:
if element in coeff:
lcoeff = coeff[element]
else:
lcoeff = coeff[i]
#end if
#end if
if element in elemargs:
v = elemargs[element]
if 'cusp' in v:
lcusp = v['cusp']
#end if
if 'rcut' in v:
lrcut = v['rcut']
#end if
if 'size' in v and 'coeff' not in v:
lcoeff = v['size']*[0]
#end if
if 'coeff' in v:
lcoeff = v['coeff']
#end if
#end if
corr = correlation(
elementtype = element,
size = len(lcoeff),
cusp = cusp,
coefficients=section(
id = ename+element,
type = 'Array',
coeff = lcoeff,
)
)
if opt is not None:
corr.coefficients.optimize = bool(opt)
if lrcut!=None:
if isperiodic and lrcut>rwigner:
QmcpackInput.class_error('rcut must not be greater than the simulation cell wigner radius\nyou provided: {0}\nwigner radius: {1}'.format(lrcut,rwigner),'generate_jastrow1')
corr.rcut = lrcut
elif isopen:
QmcpackInput.class_error('rcut must be provided for an open system','generate_jastrow1')
elif isperiodic:
corr.rcut = rwigner
#end if
corrs.append(corr)
#end for
j1 = jastrow1(
name = 'J1',
type = 'One-Body',
function = function,
source = iname,
print = True,
correlations = corrs
)
return j1
#end def generate_jastrow1
[docs]
def generate_bspline_jastrow2(size=8,rcut=None,coeff=None,spins=('u','d'),density=None,system=None,init='rpa',opt=None):
if coeff is None and system is None and (init=='rpa' and density is None or rcut is None):
QmcpackInput.class_error('rcut and density or system must be specified','generate_bspline_jastrow2')
#end if
isopen = False
isperiodic = False
allperiodic = False
rwigner = 1e99
if system is not None:
isopen = system.structure.is_open()
isperiodic = system.structure.is_periodic()
allperiodic = system.structure.all_periodic()
if not isopen and isperiodic:
rwigner = system.structure.rwigner()
#end if
volume = system.structure.volume()
if isopen:
if rcut is None:
QmcpackInput.class_error('rcut must be provided for an open system','generate_bspline_jastrow2')
#end if
if init=='rpa':
init = 'zero'
#end if
else:
if rcut is None and isperiodic:
rcut = rwigner
#end if
nelectrons = system.particles.count_electrons()
density = nelectrons/volume
#end if
elif init=='rpa':
init = 'zero'
#end if
if coeff is None:
if init=='rpa':
if not allperiodic:
QmcpackInput.class_error('rpa initialization can only be used for fully periodic systems','generate_bspline_jastrow2')
#end if
wp = np.sqrt(4*np.pi*density)
dr = rcut/size
r = .02 + dr*np.arange(size)
uuc = .5/(wp*r)*(1.-np.exp(-r*np.sqrt(wp/2)))*np.exp(-(2*r/rcut)**2)
udc = .5/(wp*r)*(1.-np.exp(-r*np.sqrt(wp)) )*np.exp(-(2*r/rcut)**2)
coeff = [uuc,udc]
elif init=='zero' or init==0:
coeff = [size*[0],size*[0]]
else:
QmcpackInput.class_error(str(init)+' is not a valid value for parameter init\n valid options are: rpa, zero','generate_bspline_jastrow2')
#end if
elif len(coeff)!=2:
QmcpackInput.class_error('must provide 2 sets of coefficients (uu,ud)','generate_bspline_jastrow2')
#end if
size = len(coeff[0])
uname,dname = spins
uuname = uname+uname
udname = uname+dname
corrs = [
correlation(speciesA=uname,speciesB=uname,size=size,
coefficients=section(id=uuname,type='Array',coeff=coeff[0])),
correlation(speciesA=uname,speciesB=dname,size=size,
coefficients=section(id=udname,type='Array',coeff=coeff[1]))
]
if opt is not None:
for corr in corrs:
corr.coefficients.optimize = bool(opt)
if rcut!=None:
if isperiodic and rcut>rwigner:
QmcpackInput.class_error('rcut must not be greater than the simulation cell wigner radius\nyou provided: {0}\nwigner radius: {1}'.format(rcut,rwigner),'generate_jastrow2')
#end if
for corr in corrs:
corr.rcut=rcut
#end for
#end if
j2 = jastrow2(
name = 'J2',type='Two-Body',function='bspline',print=True,
correlations = corrs
)
return j2
#end def generate_bspline_jastrow2
[docs]
def generate_pade_jastrow2(Buu=None,Bud=None,spins=('u','d'),system=None):
if Buu is None:
Buu = 2.0
#end if
if Bud is None:
Bud = float(Buu)
#end if
uname,dname = spins
uuname = uname+uname
udname = uname+dname
cuu = var(id=uuname+'_b',name='B',value=Buu)
cud = var(id=udname+'_b',name='B',value=Bud)
corrs = [
correlation(speciesA=uname,speciesB=uname,
vars=[cuu]),
correlation(speciesA=uname,speciesB=dname,
vars=[cud])
]
j2 = jastrow2(
name = 'J2',type='Two-Body',function='pade',
correlations = corrs
)
return j2
#end def generate_pade_jastrow2
[docs]
def generate_jastrow2(function='bspline',*args,**kwargs):
if 'spins' not in kwargs:
kwargs['spins'] = ('u','d')
#end if
spins = kwargs['spins']
if not isinstance(spins,tuple) and not isinstance(spins,list):
QmcpackInput.class_error('spins must be a list or tuple of u/d spin names\n you provided: '+str(spins))
#end if
if len(spins)!=2:
QmcpackInput.class_error('name for up and down spins must be specified\n you provided: '+str(spins))
#end if
if not isinstance(function,str):
QmcpackInput.class_error('function must be a string\n you provided: '+str(function),'generate_jastrow2')
#end if
if function=='bspline':
j2 = generate_bspline_jastrow2(*args,**kwargs)
elif function=='pade':
j2 = generate_pade_jastrow2(*args,**kwargs)
else:
QmcpackInput.class_error('function is invalid\n you provided: {0}\n valid options are: bspline or pade'.format(function),'generate_jastrow2')
#end if
if 'system' in kwargs and kwargs['system'] is not None:
nup,ndn = kwargs['system'].particles.electron_counts()
if nup<2:
del j2.correlations.uu
#end if
if nup<1 or ndn<1:
del j2.correlations.ud
#end if
if len(j2.correlations)==0:
j2=None
#end if
#end if
return j2
#end def generate_jastrow2
[docs]
def generate_jastrow3(function='polynomial',esize=3,isize=3,rcut=4.,coeff=None,iname='ion0',spins=('u','d'),elements=None,system=None,opt=None):
if elements is None and system is None:
QmcpackInput.class_error('must specify elements or system','generate_jastrow3')
elif elements is None:
elements = list(set(system.structure.elem))
#end if
if coeff is not None:
QmcpackInput.class_error('handling coeff is not yet implemented for generate jastrow3')
#end if
if len(spins)!=2:
QmcpackInput.class_error('must specify name for up and down spins\n provided: '+str(spins),'generate_jastrow3')
#end if
if rcut is None:
QmcpackInput.class_error('must specify rcut','generate_jastrow3')
#end if
if system is not None and system.structure.is_periodic():
rwigner = system.structure.rwigner()
if rcut>rwigner:
QmcpackInput.class_error('rcut must not be greater than the simulation cell wigner radius\nyou provided: {0}\nwigner radius: {1}'.format(rcut,rwigner),'generate_jastrow3')
#end if
#end if
uname,dname = spins
uuname = uname+uname
udname = uname+dname
corrs=[]
kw = obj()
if opt is not None:
kw.opt = bool(opt)
for element in elements:
corrs.append(
correlation(
especies1=uname,especies2=uname,ispecies=element,esize=esize,
isize=isize,rcut=rcut,
coefficients=section(id=uuname+element,type='Array',**kw))
)
corrs.append(
correlation(
especies1=uname,especies2=dname,ispecies=element,esize=esize,
isize=isize,rcut=rcut,
coefficients=section(id=udname+element,type='Array',**kw))
)
#end for
jastrow = jastrow3(
name = 'J3',type='eeI',function=function,print=True,source=iname,
correlations = corrs
)
return jastrow
#end def generate_jastrow3
[docs]
def generate_kspace_jastrow(
kc1: float | None = None,
kc2: float | None = None,
nk1: int = 0,
nk2: int = 0,
symm1: str = 'isotropic',
symm2: str = 'isotropic',
coeff1: list = None,
coeff2: list = None,
opt1: bool | None = None,
opt2: bool | None = None,
):
"""Generate ``<jastrow type="kSpace">``
Parameters
----------
kc1 : float, optional
kcut for one-body Jastrow. Must provide this and/or ``kc2``.
kc2 : float, optional
kcut for two-body Jastrow. Must provide this and/or ``kc1``.
nk1 : int, default=0
Number of coefficients for one-body Jastrow.
nk2 : int, default=0
Number of coefficients for two-body Jastrow.
symm1 : {'crystal', 'isotropic', 'none'}, default='isotropic'
Impose specified symmetry on 1-body Jastrow coefficients.
See Notes for description.
symm2 : {'crystal', 'isotropic', 'none'}, default='isotropic'
Impose specified symmetry on 2-body Jastrow coefficients.
See Notes for description.
coeff1 : list, optional
One-body Jastrow coefficients, optional.
coeff2 : list, optional
Two-body Jastrow coefficients, optional.
opt1 : bool or None, default=None
Set whether or not the one-body Jastrow coefficients are optimizable.
See Notes for more information.
opt2 : bool or None, default=None
Set whether or not the two-body Jastrow coefficients are optimizable.
See Notes for more information.
Returns
-------
jk : QIxml
``kspace_jastrow`` qmcpack_input element
Notes
-----
The ``symm1`` and ``symm2`` parameters yield the following behavior:
``"crystal"``
Impose crystal symmetry on coefficients according to the
structure factor.
``"isotropic"``
Impose spherical symmetry on coefficients according to G-vector
magnitude.
``"none"``
Impose no symmetry on the coefficients.
The parameters ``opt1`` and ``opt2`` are not necessarily guaranteed
to control the behavior of QMCPACK. See QMCPACK's documentation on
`k-space Jastrow`_.
.. _k-space Jastrow: https://qmcpack.readthedocs.io/en/develop/intro_wavefunction.html#long-ranged-jastrow-k-space-jastrow
"""
J1k = kc1 is not None
J2k = kc2 is not None
if not J1k and not J2k:
QmcpackInput.class_error('must have at least one term', 'generate_kspace_jastrow')
#end if
if coeff1 is None:
coeff1 = [0]*nk1
if coeff2 is None:
coeff2 = [0]*nk2
if len(coeff1) != nk1:
QmcpackInput.class_error('coeff1 mismatch', 'generate_kspace_jastrow')
#end if
if len(coeff2) != nk2:
QmcpackInput.class_error('coeff2 mismatch', 'generate_kspace_jastrow')
#end if
corrs = []
if J1k:
corr1 = correlation(
type = 'One-Body',
symmetry = symm1,
kc = kc1,
coefficients = section(id='cG1', type='Array', coeff=coeff1),
)
if opt1 is not None:
corr1.coefficients.optimize = bool(opt1)
corrs.append(corr1)
#end if
if J2k:
corr2 = correlation(
type = 'Two-Body',
symmetry = symm2,
kc = kc2,
coefficients = section(id='cG2', type='Array', coeff=coeff2),
)
if opt2 is not None:
corr2.coefficients.optimize = bool(opt2)
corrs.append(corr2)
#end if
jk = kspace_jastrow(
type = 'kSpace',
name = 'Jk',
source = 'ion0',
correlations = collection(corrs),
)
return jk
# end def generate_kspace_jastrow
[docs]
def count_jastrow_params(jastrows):
if isinstance(jastrows,QIxml):
jastrows = [jastrows]
#end if
params = 0
for jastrow in jastrows:
name = jastrow.name
if 'type' in jastrow:
type = jastrow.type.lower()
else:
type = ''
#end if
jastrow.pluralize()
if name=='J1' or type=='one-body':
for correlation in jastrow.correlations:
params += correlation.size
#end for
elif name=='J2' or type=='two-body':
for correlation in jastrow.correlations:
params += correlation.size
#end for
elif name=='J3' or type=='eeI':
for correlation in jastrow.correlations:
params += correlation.esize
params += correlation.isize
#end for
#end if
#end for
return params
#end def count_jastrow_params
[docs]
def generate_energydensity(
name = None,
dynamic = None,
static = None,
coord = None,
grid = None,
scale = None,
ion_grids = None,
system = None,
):
if dynamic is None:
dynamic = 'e'
#end if
if static is None:
static = 'ion0'
#end if
refp = None
sg = []
if coord is None:
QmcpackInput.class_error('coord must be provided','generate_energydensity')
elif coord=='voronoi':
if name is None:
name = 'EDvoronoi'
#end if
sg.append(spacegrid(coord=coord))
elif coord=='cartesian':
if name is None:
name = 'EDcell'
#end if
if grid is None:
QmcpackInput.class_error('grid must be provided for cartesian coordinates','generate_energydensity')
#end if
axes = [
axis(p1='a1',scale='.5',label='x'),
axis(p1='a2',scale='.5',label='y'),
axis(p1='a3',scale='.5',label='z'),
]
n=0
for ax in axes:
ax.grid = '-1 ({0}) 1'.format(grid[n])
n+=1
#end for
sg.append(spacegrid(coord=coord,origin=origin(p1='zero'),axes=axes))
elif coord=='spherical':
if name is None:
name = 'EDatom'
#end if
if ion_grids is None:
QmcpackInput.class_error('ion_grids must be provided for spherical coordinates','generate_energydensity')
#end if
refp = reference_points(coord='cartesian',points='\nr1 1 0 0\nr2 0 1 0\nr3 0 0 1\n')
if system is None:
i=1
for scale,g1,g2,g3 in ion_grids:
grid = g1,g2,g3
axes = [
axis(p1='r1',scale=scale,label='r'),
axis(p1='r2',scale=scale,label='phi'),
axis(p1='r3',scale=scale,label='theta'),
]
n=0
for ax in axes:
ax.grid = '0 ({0}) 1'.format(grid[n])
n+=1
#end for
sg.append(spacegrid(coord=coord,origin=origin(p1=static+str(i)),axes=axes))
i+=1
#end for
else:
ig = ion_grids
ion_grids = obj()
for e,s,g1,g2,g3 in ig:
ion_grids[e] = s,(g1,g2,g3)
#end for
missing = set(ion_grids.keys())
i=1
for e in system.structure.elem:
if e in ion_grids:
scale,grid = ion_grids[e]
axes = [
axis(p1='r1',scale=scale,label='r'),
axis(p1='r2',scale=scale,label='phi'),
axis(p1='r3',scale=scale,label='theta'),
]
n=0
for ax in axes:
ax.grid = '0 ({0}) 1'.format(grid[n])
n+=1
#end for
sg.append(spacegrid(coord=coord,origin=origin(p1=static+str(i)),axes=axes))
if e in missing:
missing.remove(e)
#end if
#end if
i+=1
#end for
if len(missing)>0:
QmcpackInput.class_error('ion species not found for spherical grid\nspecies not found: {0}\nspecies present: {1}'.format(sorted(missing),sorted(set(list(system.structure.elem)))),'generate_energydensity')
#end if
#end if
else:
QmcpackInput.class_error('unsupported coord type\ncoord type provided: {0}\nsupported coord types: voronoi, cartesian, spherical'.format(coord),'generate_energydensity')
#end if
ed = energydensity(
type = 'EnergyDensity',
name = name,
dynamic = dynamic,
static = static,
spacegrids = sg,
)
if refp is not None:
ed.reference_points = refp
#end if
return ed
#end def generate_energydensity
opt_map = dict(linear=linear,cslinear=cslinear,linear_batch=linear_batch)
[docs]
def generate_opt(method,
repeat = 1,
energy = None,
rw_variance = None,
urw_variance = None,
params = None,
jastrows = None,
processes = None,
walkers_per_proc = None,
threads = None,
blocks = 2000,
#steps = 5,
decorr = 10,
min_walkers = None, #use e.g. 128 for gpu's
timestep = .5,
nonlocalpp = False,
sample_factor = 1.0):
if method not in opt_map:
QmcpackInput.class_error('section cannot be generated for optimization method '+method)
#end if
if energy is None and rw_variance is None and urw_variance is None:
QmcpackInput.class_error('at least one cost parameter must be specified\n options are: energy, rw_variance, urw_variance')
#end if
if params is None and jastrows is None:
QmcpackInput.class_error('must provide either number of opt parameters (params) or a list of jastrow objects (jastrows)')
#end if
if processes is None:
QmcpackInput.class_error('must specify total number of processes')
elif walkers_per_proc is None and threads is None:
QmcpackInput.class_error('must specify walkers_per_proc or threads')
#end if
if params is None:
params = count_jastrow_params(jastrows)
#end if
samples = max(100000,100*params**2)
samples = int(round(sample_factor*samples))
samples_per_proc = int(round(float(samples)/processes))
if walkers_per_proc is None:
walkers = 1
walkers_per_proc = threads
else:
walkers = walkers_per_proc
#end if
tot_walkers = processes*walkers_per_proc
if min_walkers is not None:
tot_walkers = max(min_walkers,tot_walkers)
walkers = int(np.ceil(float(tot_walkers)/processes-.001))
if threads is not None and np.mod(walkers,threads)!=0:
walkers = threads*int(np.ceil(float(walkers)/threads-.001))
#end if
#end if
#blocks = int(np.ceil(float(decorr*samples)/(steps*tot_walkers)))
blocks = min(blocks,samples_per_proc*decorr)
opt = opt_map[method]()
opt.set(
walkers = walkers,
blocks = blocks,
#steps = steps,
samples = samples,
substeps = decorr,
timestep = timestep,
nonlocalpp = nonlocalpp,
stepsbetweensamples = 1
)
if energy is not None:
opt.energy = energy
#end if
if rw_variance is not None:
opt.reweightedvariance = rw_variance
#end if
if urw_variance is not None:
opt.unreweightedvariance = urw_variance
#end if
opt.incorporate_defaults(elements=True)
if repeat>1:
opt = loop(max=repeat,qmc=opt)
#end if
return opt
#end def generate_opt
[docs]
def generate_opts(opt_reqs,**kwargs):
opts = []
for opt_req in opt_reqs:
opts.append(generate_opt(*opt_req,**kwargs))
#end for
return opts
#end def generate_opts
# legacy driver defaults
opt_legacy_defaults = obj(
method = 'linear',
minmethod = 'quartic',
cost = 'variance',
cycles = 12,
var_cycles = 0,
var_samples = None,
init_cycles = 0,
init_samples = None,
init_minwalkers = 1e-4,
)
shared_opt_legacy_defaults = obj(
samples = 204800,
nonlocalpp = True,
use_nonlocalpp_deriv = True,
warmupsteps = 300,
blocks = 100,
steps = 1,
substeps = 10,
timestep = 0.3,
usedrift = False,
max_seconds = None,
spin_mass = None,
)
linear_quartic_legacy_defaults = obj(
minwalkers = 0.3,
usebuffer = True,
exp0 = -6,
bigchange = 10.0,
alloweddifference = 1e-04,
stepsize = 0.15,
nstabilizers = 1,
**shared_opt_legacy_defaults
)
linear_oneshift_legacy_defaults = obj(
minwalkers = 0.5,
#shift_i = 0.01,
#shift_s = 1.00,
**shared_opt_legacy_defaults
)
linear_adaptive_legacy_defaults = obj(
minwalkers = 0.3,
max_relative_change = 10.0,
max_param_change = 0.3,
shift_i = 0.01,
shift_s = 1.00,
**shared_opt_legacy_defaults
)
opt_method_legacy_defaults = obj({
('linear' ,'quartic' ) : linear_quartic_legacy_defaults,
('linear' ,'rescale' ) : linear_quartic_legacy_defaults,
('linear' ,'linemin' ) : linear_quartic_legacy_defaults,
('cslinear','quartic' ) : linear_quartic_legacy_defaults,
('cslinear','rescale' ) : linear_quartic_legacy_defaults,
('cslinear','linemin' ) : linear_quartic_legacy_defaults,
('linear' ,'adaptive') : linear_adaptive_legacy_defaults,
('linear' ,'oneshift') : linear_oneshift_legacy_defaults,
('linear' ,'oneshiftonly') : linear_oneshift_legacy_defaults,
})
del shared_opt_legacy_defaults
del linear_quartic_legacy_defaults
del linear_oneshift_legacy_defaults
del linear_adaptive_legacy_defaults
allowed_opt_method_legacy_inputs = set(linear.attributes+linear.parameters
+cslinear.attributes+cslinear.parameters)
vmc_legacy_defaults = obj(
walkers = 1,
warmupsteps = 50,
blocks = 800,
steps = 10,
substeps = 3,
timestep = 0.3,
checkpoint = -1,
usedrift = None,
max_seconds = None,
spin_mass = None,
)
vmc_test_legacy_defaults = obj(
warmupsteps = 10,
blocks = 20,
steps = 4,
).set_optional(**vmc_legacy_defaults)
vmc_noJ_legacy_defaults = obj(
warmupsteps = 200,
blocks = 800,
steps = 100,
).set_optional(**vmc_legacy_defaults)
dmc_legacy_defaults = obj(
warmupsteps = 20,
blocks = 200,
steps = 10,
timestep = 0.01,
checkpoint = -1,
vmc_samples = 2048,
vmc_samplesperthread = None,
vmc_walkers = None,
vmc_warmupsteps = 30,
vmc_blocks = 40,
vmc_steps = 10,
vmc_substeps = 3,
vmc_timestep = 0.3,
vmc_usedrift = None,
vmc_checkpoint = -1,
vmc_spin_mass = None,
eq_dmc = False,
eq_warmupsteps = 20,
eq_blocks = 20,
eq_steps = 5,
eq_timestep = 0.02,
eq_checkpoint = -1,
ntimesteps = 1,
timestep_factor = 0.5,
nonlocalmoves = None,
branching_cutoff_scheme = None,
maxage = None,
feedback = None,
sigmabound = None,
max_seconds = None,
spin_mass = None,
)
dmc_test_legacy_defaults = obj(
vmc_warmupsteps = 10,
vmc_blocks = 20,
vmc_steps = 4,
eq_warmupsteps = 2,
eq_blocks = 5,
eq_steps = 2,
warmupsteps = 2,
blocks = 10,
steps = 2,
).set_optional(**dmc_legacy_defaults)
dmc_noJ_legacy_defaults = obj(
warmupsteps = 40,
blocks = 400,
steps = 20,
).set_optional(**dmc_legacy_defaults)
# batched driver defaults
opt_batched_defaults = obj(
method = 'linear',
minmethod = 'quartic',
cost = None, # default: energy if sr_cg else variance
cycles = 12,
var_cycles = 0,
var_samples = None,
init_cycles = 0,
init_samples = None,
init_steps = None,
init_minwalkers = 0.1,
init_line_search = True,
init_sr_tau = 0.1,
)
shared_opt_batched_defaults = obj(
samples = None, # 204800 if steps is None
#nonlocalpp = True,
#use_nonlocalpp_deriv = True,
warmupsteps = 300,
blocks = 100,
steps = None,
substeps = 10,
timestep = 0.3,
usedrift = False,
spin_mass = None,
walkers_per_rank = None,
total_walkers = None,
)
linear_quartic_batched_defaults = obj(
minwalkers = 0.3,
usebuffer = True,
exp0 = -6,
bigchange = 10.0,
alloweddifference = 1e-04,
stepsize = 0.15,
nstabilizers = 1,
**shared_opt_batched_defaults
)
linear_oneshift_batched_defaults = obj(
minwalkers = 0.5,
#shift_i = 0.01,
#shift_s = 1.00,
**shared_opt_batched_defaults
)
linear_adaptive_batched_defaults = obj(
minwalkers = 0.3,
max_relative_change = 10.0,
max_param_change = 0.3,
shift_i = 0.01,
shift_s = 1.00,
**shared_opt_batched_defaults
)
linear_sr_cg_batched_defaults = obj(
sr_tau = None, # projector: 1-tau*H (0.01/0.1 if line_search=no/yes)
sr_tolerance = 0.001, # conjugate gradient convergence tolerance
sr_regularization = 0.01, # ~diagonal shift to overlap matrix
line_search = False, # corr samp line search on cost along sr param direction
**shared_opt_batched_defaults
)
opt_method_batched_defaults = obj({
('linear' ,'quartic' ) : linear_quartic_batched_defaults,
('linear' ,'rescale' ) : linear_quartic_batched_defaults,
('linear' ,'linemin' ) : linear_quartic_batched_defaults,
('cslinear','quartic' ) : linear_quartic_batched_defaults,
('cslinear','rescale' ) : linear_quartic_batched_defaults,
('cslinear','linemin' ) : linear_quartic_batched_defaults,
('linear' ,'adaptive') : linear_adaptive_batched_defaults,
('linear' ,'oneshift') : linear_oneshift_batched_defaults,
('linear' ,'oneshiftonly') : linear_oneshift_batched_defaults,
('linear' ,'sr_cg') : linear_sr_cg_batched_defaults,
})
del shared_opt_batched_defaults
del linear_quartic_batched_defaults
del linear_oneshift_batched_defaults
del linear_adaptive_batched_defaults
del linear_sr_cg_batched_defaults
allowed_opt_method_batched_inputs = set(linear.attributes+linear.parameters
+cslinear.attributes+cslinear.parameters)
vmc_batched_defaults = obj(
total_walkers = None,
walkers_per_rank = None,
warmupsteps = 50,
blocks = 800,
steps = 10,
substeps = 3,
timestep = 0.3,
usedrift = False,
checkpoint = None,
maxcpusecs = None,
crowds = None,
spin_mass = None,
)
vmc_test_batched_defaults = obj(
warmupsteps = 10,
blocks = 20,
steps = 4,
).set_optional(**vmc_batched_defaults)
vmc_noJ_batched_defaults = obj(
warmupsteps = 200,
blocks = 800,
steps = 100,
).set_optional(**vmc_batched_defaults)
dmc_batched_defaults = obj(
total_walkers = None,
walkers_per_rank = None,
warmupsteps = 20,
blocks = 200,
steps = 10,
substeps = None,
timestep = 0.01,
checkpoint = None,
vmc_warmupsteps = 30,
vmc_blocks = 40,
vmc_steps = 10,
vmc_substeps = 3,
vmc_timestep = 0.3,
vmc_usedrift = False,
vmc_checkpoint = None,
vmc_spin_mass = None,
eq_dmc = False,
eq_warmupsteps = 20,
eq_blocks = 20,
eq_steps = 5,
eq_timestep = 0.02,
eq_checkpoint = None,
ntimesteps = 1,
timestep_factor = 0.5,
nonlocalmoves = None,
branching_cutoff_scheme = None,
crowd_serialize_walkers = None,
crowds = None,
reconfiguration = None,
maxage = None,
feedback = None,
sigmabound = None,
spin_mass = None,
)
dmc_test_batched_defaults = obj(
vmc_warmupsteps = 10,
vmc_blocks = 20,
vmc_steps = 4,
eq_warmupsteps = 2,
eq_blocks = 5,
eq_steps = 2,
warmupsteps = 2,
blocks = 10,
steps = 2,
).set_optional(**dmc_batched_defaults)
dmc_noJ_batched_defaults = obj(
warmupsteps = 40,
blocks = 400,
steps = 20,
).set_optional(**dmc_batched_defaults)
# collected defaults for opt methods
opt_method_defaults = obj(
legacy = opt_method_legacy_defaults,
batched = opt_method_batched_defaults,
)
# collected defaults for all drivers
qmc_defaults = obj(
legacy = obj(
opt = opt_legacy_defaults,
vmc = vmc_legacy_defaults,
vmc_test = vmc_test_legacy_defaults,
vmc_noJ = vmc_noJ_legacy_defaults,
dmc = dmc_legacy_defaults,
dmc_test = dmc_test_legacy_defaults,
dmc_noJ = dmc_noJ_legacy_defaults,
),
batched = obj(
opt = opt_batched_defaults,
vmc = vmc_batched_defaults,
vmc_test = vmc_test_batched_defaults,
vmc_noJ = vmc_noJ_batched_defaults,
dmc = dmc_batched_defaults,
dmc_test = dmc_test_batched_defaults,
dmc_noJ = dmc_noJ_batched_defaults,
),
)
del opt_legacy_defaults
del vmc_legacy_defaults
del vmc_test_legacy_defaults
del vmc_noJ_legacy_defaults
del dmc_legacy_defaults
del dmc_test_legacy_defaults
del dmc_noJ_legacy_defaults
del opt_batched_defaults
del vmc_batched_defaults
del vmc_test_batched_defaults
del vmc_noJ_batched_defaults
del dmc_batched_defaults
del dmc_test_batched_defaults
del dmc_noJ_batched_defaults
[docs]
def generate_opt_calculations(driver,**kwargs):
if driver=='legacy':
calcs = generate_legacy_opt_calculations(**kwargs)
elif driver=='batched':
calcs = generate_batched_opt_calculations(**kwargs)
else:
error('Cannot generate calculations for unrecognized driver.\nUnrecognized driver: {}'.format(driver))
#end if
return calcs
#end def generate_opt_calculations
[docs]
def generate_vmc_calculations(driver,**kwargs):
if driver=='legacy':
calcs = generate_legacy_vmc_calculations(**kwargs)
elif driver=='batched':
calcs = generate_batched_vmc_calculations(**kwargs)
else:
error('Cannot generate calculations for unrecognized driver.\nUnrecognized driver: {}'.format(driver))
#end if
return calcs
#end def generate_vmc_calculations
[docs]
def generate_dmc_calculations(driver,**kwargs):
if driver=='legacy':
calcs = generate_legacy_dmc_calculations(**kwargs)
elif driver=='batched':
calcs = generate_batched_dmc_calculations(**kwargs)
else:
error('Cannot generate calculations for unrecognized driver.\nUnrecognized driver: {}'.format(driver))
#end if
return calcs
#end def generate_dmc_calculations
[docs]
def generate_legacy_opt_calculations(
method ,
cost ,
cycles ,
var_cycles ,
var_samples,
init_cycles,
init_samples,
init_minwalkers,
loc = 'generate_opt_calculations',
**opt_inputs
):
methods = obj(linear=linear,cslinear=cslinear)
if method not in methods:
error('invalid optimization method requested\ninvalid method: {0}\nvalid options are: {1}'.format(method,sorted(methods.keys())),loc)
#end if
opt = methods[method]
opt_inputs = obj(opt_inputs)
invalid = set(opt_inputs.keys())-allowed_opt_method_legacy_inputs
oneshift = False
if len(invalid)>0:
error('invalid optimization inputs provided\ninvalid inputs: {}\nvalid options are: {}'.format(sorted(invalid),sorted(allowed_opt_method_legacy_inputs)))
#end if
for k in list(opt_inputs.keys()):
if opt_inputs[k] is None:
del opt_inputs[k]
#end if
#end for
if 'minmethod' in opt_inputs and opt_inputs.minmethod.lower().startswith('oneshift'):
opt_inputs.minmethod = 'OneShiftOnly'
oneshift = True
#end if
if cost=='variance':
cost = (0.0,1.0,0.0)
elif cost=='energy':
cost = (1.0,0.0,0.0)
elif isinstance(cost,(tuple,list)) and (len(cost)==2 or len(cost)==3):
if len(cost)==2:
cost = (cost[0],0.0,cost[1])
#end if
else:
error('invalid optimization cost function encountered\ninvalid cost fuction: {0}\nvalid options are: variance, energy, (0.95,0.05), etc'.format(cost),loc)
#end if
opt_calcs = []
if var_cycles>0:
vmin_opt = opt(
energy = 0.0,
unreweightedvariance = 1.0,
reweightedvariance = 0.0,
**opt_inputs
)
if var_samples is not None:
vmin_opt.samples = var_samples
#end if
opt_calcs.append(loop(max=var_cycles,qmc=vmin_opt))
#end if
if init_cycles>0:
init_opt = opt(**opt_inputs)
if init_samples is not None:
init_opt.samples = init_samples
#end if
init_opt.minwalkers = init_minwalkers
if not oneshift:
init_opt.energy = cost[0]
init_opt.unreweightedvariance = cost[1]
init_opt.reweightedvariance = cost[2]
#end if
opt_calcs.append(loop(max=init_cycles,qmc=init_opt))
#end if
cost_opt = opt(**opt_inputs)
if not oneshift:
cost_opt.energy = cost[0]
cost_opt.unreweightedvariance = cost[1]
cost_opt.reweightedvariance = cost[2]
#end if
opt_calcs.append(loop(max=cycles,qmc=cost_opt))
return opt_calcs
#end def generate_legacy_opt_calculations
[docs]
def generate_legacy_vmc_calculations(
walkers ,
warmupsteps,
blocks ,
steps ,
substeps ,
timestep ,
checkpoint ,
usedrift ,
max_seconds,
spin_mass,
loc = 'generate_vmc_calculations',
):
vmc_calc = vmc(
walkers = walkers,
warmupsteps = warmupsteps,
blocks = blocks,
steps = steps,
substeps = substeps,
timestep = timestep,
checkpoint = checkpoint,
)
if usedrift is not None:
vmc_calc.usedrift = usedrift
#end if
if max_seconds is not None:
vmc_calc.max_seconds = max_seconds
#end if
if spin_mass is not None:
vmc_calc.spin_mass = spin_mass
#end if
vmc_calcs = [vmc_calc]
return vmc_calcs
#end def generate_legacy_vmc_calculations
[docs]
def generate_legacy_dmc_calculations(
warmupsteps ,
blocks ,
steps ,
timestep ,
checkpoint ,
vmc_samples ,
vmc_samplesperthread ,
vmc_walkers ,
vmc_warmupsteps ,
vmc_blocks ,
vmc_steps ,
vmc_substeps ,
vmc_timestep ,
vmc_usedrift ,
vmc_checkpoint ,
vmc_spin_mass ,
eq_dmc ,
eq_warmupsteps ,
eq_blocks ,
eq_steps ,
eq_timestep ,
eq_checkpoint ,
ntimesteps ,
timestep_factor ,
nonlocalmoves ,
branching_cutoff_scheme,
maxage ,
feedback ,
sigmabound ,
max_seconds ,
spin_mass ,
loc = 'generate_dmc_calculations',
):
if vmc_samples is None and vmc_samplesperthread is None and vmc_walkers is None:
error('vmc samples (dmc walkers) not specified\nplease provide one of the following keywords: vmc_samples, vmc_samplesperthread, vmc_walkers',loc)
#end if
if vmc_walkers is None:
vmc_walkers = 1
#end if
vmc_calc = vmc(
walkers = vmc_walkers,
warmupsteps = vmc_warmupsteps,
blocks = vmc_blocks,
steps = vmc_steps,
substeps = vmc_substeps,
timestep = vmc_timestep,
checkpoint = vmc_checkpoint,
)
if vmc_samplesperthread is not None:
vmc_calc.samplesperthread = vmc_samplesperthread
elif vmc_samples is not None:
vmc_calc.samples = vmc_samples
#end if
if vmc_usedrift is not None:
vmc_calc.usedrift = vmc_usedrift
#end if
if max_seconds is not None:
vmc_calc.max_seconds = max_seconds
#end if
if vmc_spin_mass is not None:
vmc_calc.spin_mass = vmc_spin_mass
#end if
dmc_calcs = [vmc_calc]
if eq_dmc:
dmc_calcs.append(
dmc(
warmupsteps = eq_warmupsteps,
blocks = eq_blocks,
steps = eq_steps,
timestep = eq_timestep,
checkpoint = eq_checkpoint,
)
)
#end if
tfac = 1.0
for n in range(ntimesteps):
sfac = 1.0/tfac
dmc_calcs.append(
dmc(
warmupsteps = int(sfac*warmupsteps),
blocks = blocks,
steps = int(sfac*steps),
timestep = tfac*timestep,
checkpoint = checkpoint,
)
)
tfac *= timestep_factor
#end for
optional_dmc_inputs = obj(
nonlocalmoves = nonlocalmoves,
branching_cutoff_scheme = branching_cutoff_scheme,
maxage = maxage ,
feedback = feedback ,
sigmabound = sigmabound,
max_seconds = max_seconds,
spin_mass = spin_mass,
)
for calc in dmc_calcs:
if isinstance(calc,dmc):
for k,v in optional_dmc_inputs.items():
if v is not None:
calc[k] = v
#end if
#end for
#end if
#end for
return dmc_calcs
#end def generate_legacy_dmc_calculations
[docs]
def generate_batched_opt_calculations(
method ,
cost ,
cycles ,
var_cycles ,
var_samples,
init_cycles,
init_samples,
init_steps,
init_minwalkers,
init_line_search,
init_sr_tau,
loc = 'generate_opt_calculations',
**opt_inputs
):
opt_inputs = obj(opt_inputs)
has = obj(
minmethod = 'minmethod' in opt_inputs,
samples = 'samples' in opt_inputs,
steps = 'steps' in opt_inputs,
sr_tau = 'sr_tau' in opt_inputs,
linesearch = 'linesearch' in opt_inputs,
)
for k,v in has.items():
has[k] = v and opt_inputs[k] is not None
#end for
if cost is None:
if has.minmethod and opt_inputs.minmethod=='sr_cg':
cost = 'energy'
else:
cost = 'variance'
#end if
#end if
if not has.samples and not has.steps:
opt_inputs.samples = 204800
#end if
if init_samples is None and init_steps is None and not has.steps:
init_samples = 204800
#end if
if not has.sr_tau:
if has.linesearch and opt_inputs.line_search:
opt_inputs.sr_tau = 0.1
else:
opt_inputs.sr_tau = 0.01
#end if
#end if
for k in list(opt_inputs.keys()):
if opt_inputs[k] is None:
del opt_inputs[k]
#end if
#end for
methods = obj(linear=linear)
if method not in methods:
error('invalid optimization method requested\ninvalid method: {0}\nvalid options are: {1}'.format(method,sorted(methods.keys())),loc)
#end if
opt = methods[method]
opt_inputs = obj(opt_inputs)
invalid = set(opt_inputs.keys())-allowed_opt_method_batched_inputs
oneshift = False
sr_cg = False
if len(invalid)>0:
error('invalid optimization inputs provided\ninvalid inputs: {}\nvalid options are: {}'.format(sorted(invalid),sorted(allowed_opt_method_batched_inputs)))
#end if
if has.minmethod:
if opt_inputs.minmethod.lower().startswith('oneshift'):
opt_inputs.minmethod = 'OneShiftOnly'
oneshift = True
elif opt_inputs.minmethod=='sr_cg':
sr_cg = True
#end if
#end if
if cost=='variance':
cost = (0.0,1.0,0.0)
elif cost=='energy':
cost = (1.0,0.0,0.0)
elif isinstance(cost,(tuple,list)) and (len(cost)==2 or len(cost)==3):
if len(cost)==2:
cost = (cost[0],0.0,cost[1])
#end if
else:
error('invalid optimization cost function encountered\ninvalid cost fuction: {0}\nvalid options are: variance, energy, (0.95,0.05), etc'.format(cost),loc)
#end if
opt_calcs = []
if var_cycles>0:
vmin_opt = opt(
energy = 0.0,
unreweightedvariance = 1.0,
reweightedvariance = 0.0,
**opt_inputs
)
if var_samples is not None:
vmin_opt.samples = var_samples
#end if
opt_calcs.append(loop(max=var_cycles,qmc=vmin_opt))
#end if
if init_cycles>0:
init_opt = opt(**opt_inputs)
if init_samples is not None:
init_opt.samples = init_samples
#end if
if init_steps is not None:
init_opt.steps = init_steps
#end if
if init_sr_tau is not None:
init_opt.sr_tau = init_sr_tau
#end if
if init_line_search is not None:
init_opt.line_search = init_line_search
if not sr_cg:
init_opt.minwalkers = init_minwalkers
#end if
if not oneshift:
init_opt.energy = cost[0]
init_opt.unreweightedvariance = cost[1]
init_opt.reweightedvariance = cost[2]
#end if
opt_calcs.append(loop(max=init_cycles,qmc=init_opt))
#end if
cost_opt = opt(**opt_inputs)
if not oneshift:
cost_opt.energy = cost[0]
cost_opt.unreweightedvariance = cost[1]
cost_opt.reweightedvariance = cost[2]
#end if
opt_calcs.append(loop(max=cycles,qmc=cost_opt))
return opt_calcs
#end def generate_batched_opt_calculations
[docs]
def generate_batched_vmc_calculations(
total_walkers ,
walkers_per_rank ,
warmupsteps ,
blocks ,
steps ,
substeps ,
timestep ,
usedrift ,
checkpoint ,
maxcpusecs ,
crowds ,
spin_mass ,
loc = 'generate_vmc_calculations',
):
if total_walkers is not None and walkers_per_rank is not None:
error('Only one of "total_walkers" and "walkers_per_rank" may be provided.',loc)
#end if
vmc_inputs = obj(
warmupsteps = warmupsteps,
blocks = blocks,
steps = steps,
substeps = substeps,
timestep = timestep,
usedrift = usedrift,
)
optional_vmc_inputs = obj(
total_walkers = total_walkers,
walkers_per_rank = walkers_per_rank,
checkpoint = checkpoint,
maxcpusecs = maxcpusecs,
crowds = crowds,
spin_mass = spin_mass,
)
for name,value in optional_vmc_inputs.items():
if value is not None:
vmc_inputs[name] = value
#end if
#end for
vmc_calcs = [vmc(**vmc_inputs)]
return vmc_calcs
#end def generate_batched_vmc_calculations
[docs]
def generate_batched_dmc_calculations(
total_walkers ,
walkers_per_rank ,
warmupsteps ,
blocks ,
steps ,
substeps ,
timestep ,
checkpoint ,
vmc_warmupsteps ,
vmc_blocks ,
vmc_steps ,
vmc_substeps ,
vmc_timestep ,
vmc_usedrift ,
vmc_checkpoint ,
vmc_spin_mass ,
eq_dmc ,
eq_warmupsteps ,
eq_blocks ,
eq_steps ,
eq_timestep ,
eq_checkpoint ,
ntimesteps ,
timestep_factor ,
nonlocalmoves ,
branching_cutoff_scheme,
crowd_serialize_walkers,
crowds ,
reconfiguration ,
maxage ,
feedback ,
sigmabound ,
spin_mass ,
loc = 'generate_dmc_calculations',
):
if total_walkers is None and walkers_per_rank is None:
total_walkers = 2048
#error('DMC walker count not specified via "total_walkers" or "walkers_per_rank".\nPlease provide at least one of these.\n\nWarning: use care in the selection of these parameters.\nPerformance critically depends on the walker count and the batched QMCPACK \ndrivers make no effort to prevent substantial under-utilization.',loc)
elif total_walkers is not None and walkers_per_rank is not None:
error('Only one of "total_walkers" and "walkers_per_rank" may be provided.',loc)
#end if
vmc_inputs = obj(
warmupsteps = vmc_warmupsteps,
blocks = vmc_blocks,
steps = vmc_steps,
substeps = vmc_substeps,
timestep = vmc_timestep,
usedrift = vmc_usedrift,
)
optional_vmc_inputs = obj(
total_walkers = total_walkers,
walkers_per_rank = walkers_per_rank,
crowds = crowds,
spin_mass = vmc_spin_mass,
checkpoint = vmc_checkpoint,
)
for name,value in optional_vmc_inputs.items():
if value is not None:
vmc_inputs[name] = value
#end if
#end for
vmc_calc = vmc(**vmc_inputs)
dmc_calcs = [vmc_calc]
if eq_dmc:
dmc_calcs.append(
dmc(
warmupsteps = eq_warmupsteps,
blocks = eq_blocks,
steps = eq_steps,
timestep = eq_timestep,
)
)
#end if
tfac = 1.0
for n in range(ntimesteps):
sfac = 1.0/tfac
dmc_calcs.append(
dmc(
warmupsteps = int(sfac*warmupsteps),
blocks = blocks,
steps = int(sfac*steps),
timestep = tfac*timestep,
)
)
tfac *= timestep_factor
#end for
optional_dmc_inputs = obj(
total_walkers = total_walkers,
walkers_per_rank = walkers_per_rank,
substeps = substeps,
nonlocalmoves = nonlocalmoves,
branching_cutoff_scheme = branching_cutoff_scheme,
crowd_serialize_walkers = crowd_serialize_walkers,
crowds = crowds,
reconfiguration = reconfiguration,
maxage = maxage,
feedback = feedback,
sigmabound = sigmabound,
spin_mass = spin_mass,
checkpoint = checkpoint,
)
for calc in dmc_calcs:
if isinstance(calc,dmc):
for name,value in optional_dmc_inputs.items():
if value is not None:
calc[name] = value
#end if
#end for
#end if
#end for
return dmc_calcs
#end def generate_batched_dmc_calculations
[docs]
def generate_qmcpack_input(**kwargs):
QIcollections.clear()
system = kwargs.get('system',None)
if isinstance(system,PhysicalSystem):
system.update_particles()
#end if
selector = kwargs.pop('input_type','basic')
if selector=='basic':
inp = generate_basic_input(**kwargs)
elif selector=='basic_afqmc':
inp = generate_basic_afqmc_input(**kwargs)
elif selector=='opt_jastrow':
inp = generate_opt_jastrow_input(**kwargs)
else:
QmcpackInput.class_error('selection '+str(selector)+' has not been implemented for qmcpack input generation')
#end if
return inp
#end def generate_qmcpack_input
[docs]
def read_jastrows(filepath):
qi = QmcpackInput(filepath)
qi.pluralize()
jastrows = qi.get('jastrows')
return jastrows
#end def read_jastrows
gen_basic_input_defaults = obj(
id = 'qmc',
series = 0,
purpose = '',
maxcpusecs = None,
max_seconds = None,
seed = None,
bconds = None,
truncate = False,
buffer = None,
lr_dim_cutoff = 15,
lr_tol = None,
lr_handler = None,
remove_cell = False,
randomsrc = True,
meshfactor = 1.0,
orbspline = None,
precision = 'float',
twistnum = None,
twist = None,
gcta = None,
spin_polarized = None,
partition = None,
partition_mf = None,
hybridrep = None,
hybrid_rcut = None,
hybrid_lmax = None,
orbitals_h5 = 'MISSING.h5',
rotated_orbitals = False,
run_path = None,
check_paths = True,
excitation = None,
system = 'missing',
pseudos = None,
pseudo_algorithm = None,
spinor = None,
dla = None,
delay_rank = None,
det_batch = None,
jastrows = 'generateJ12',
opt_params = None,
interactions = 'all',
corrections = 'default',
observables = None,
estimators = None,
estimator_period = None,
traces = None,
calculations = None,
det_format = 'new',
J1 = False,
J2 = False,
J3 = False,
J1_size = None,
J1_rcut = None,
J1_dr = 0.5,
J1_opt = None,
J2_size = None,
J2_rcut = None,
J2_dr = 0.5,
J2_init = 'zero',
J2_opt = None,
J3_isize = 3,
J3_esize = 3,
J3_rcut = 5.0,
J3_opt = None,
J1_rcut_open = 5.0,
J2_rcut_open = 10.0,
J1k = False,
J1k_kcut = 5.0,
J1k_symm = 'crystal',
J1k_opt = None,
J2k = False,
J2k_kcut = 5.0,
J2k_symm = 'crystal',
J2k_opt = None,
driver = 'batched', # legacy,batched
# batched driver inputs
orbitals_cpu = None, # place/evaluate orbitals on cpu if on gpu
matrix_inv_cpu = None, # evaluate matrix inverse on cpu if on gpu
# legacy cuda inputs
gpusharing = None,
qmc = None, # opt,vmc,vmc_test,dmc,dmc_test
)
[docs]
def generate_basic_input(**kwargs):
# capture inputs
kw = obj(kwargs)
# apply general defaults
kw.set_optional(**gen_basic_input_defaults)
valid = set(gen_basic_input_defaults.keys())
# apply method specific defaults
if kw.qmc is not None:
if kw.driver not in qmc_defaults:
QmcpackInput.class_error('Invalid input for argument "driver".\nInvalid input: {}\nValid options are: {}'.format(kw.driver,sorted(qmc_defaults.keys())),'generate_qmcpack_input')
#end if
qmc_driver_defaults = qmc_defaults[kw.driver]
if kw.qmc not in qmc_driver_defaults:
QmcpackInput.class_error('Invalid input for argument "qmc".\nInvalid input: {}\nValid options are: {}'.format(kw.qmc,sorted(qmc_driver_defaults.keys())),'generate_qmcpack_input')
#end if
qmc_keys = ['driver']
kw.set_optional(**qmc_driver_defaults[kw.qmc])
qmc_keys += list(qmc_driver_defaults[kw.qmc].keys())
if kw.qmc=='opt':
opt_method_driver_defaults = opt_method_defaults[kw.driver]
key = (kw.method,kw.minmethod.lower())
if key not in opt_method_driver_defaults:
QmcpackInput.class_error('invalid input for arguments "method,minmethod".\nInvalid input: {}\nValid options are: {}'.format(key,sorted(opt_method_driver_defaults.keys())),'generate_qmcpack_input')
#end if
kw.set_optional(**opt_method_driver_defaults[key])
qmc_keys += list(opt_method_driver_defaults[key].keys())
del key
#end if
valid |= set(qmc_keys)
#end if
# screen for invalid keywords
invalid_kwargs = set(kw.keys())-valid
if len(invalid_kwargs)>0:
QmcpackInput.class_error('invalid input parameters encountered.\nInvalid input parameters: {0}\nValid options are: {1}'.format(sorted(invalid_kwargs),sorted(valid)),'generate_qmcpack_input')
#end if
batched = kw.driver=='batched'
legacy = kw.driver=='legacy'
if kw.system=='missing':
QmcpackInput.class_error('argument "system" is missing.\nIf you really do not want particlesets to be generated, set system to None.','generate_qmcpack_input')
#end if
if kw.bconds is None:
if kw.system is not None:
s = kw.system.structure
kw.bconds = s.bconds
if len(kw.bconds)==0 or not s.has_axes():
kw.bconds = 'nnn'
#end if
else:
kw.bconds = 'ppp'
#end if
#end if
if kw.corrections=='default' and tuple(kw.bconds)==tuple('ppp') and not kw.spinor:
if not batched:
kw.corrections = ['mpc','chiesa']
else:
kw.corrections = ['mpc']
#end if
elif isinstance(kw.corrections,(list,tuple)):
None
else:
kw.corrections = []
#end if
if kw.observables is None:
#observables = ['localenergy']
kw.observables = []
#end if
if kw.estimators is None:
kw.estimators = []
#end if
kw.estimators = kw.estimators + kw.observables + kw.corrections
if kw.calculations is None:
kw.calculations = []
#end if
if kw.spin_polarized is None:
kw.spin_polarized = kw.system.net_spin>0
#end if
if kw.partition is not None:
kw.det_format = 'new'
#end if
if kw.hybrid_rcut is not None or kw.hybrid_lmax is not None:
kw.hybridrep = True
#end if
metadata = QmcpackInput.default_metadata.copy()
proj = project(
id = kw.id,
series = kw.series,
application = application(),
driver_version = kw.driver,
)
if batched:
if kw.maxcpusecs is not None:
proj.maxcpusecs = kw.maxcpusecs
#end if
if kw.max_seconds is not None:
proj.max_seconds = kw.max_seconds
#end if
#end if
simcell = generate_simulationcell(
bconds = kw.bconds,
lr_dim_cutoff = kw.lr_dim_cutoff,
lr_tol = kw.lr_tol,
lr_handler = kw.lr_handler,
system = kw.system,
)
if kw.system is not None:
kw.system.structure.set_bconds(kw.bconds)
particlesets = generate_particlesets(
system = kw.system,
randomsrc = kw.randomsrc,
hybrid_rcut = kw.hybrid_rcut,
hybrid_lmax = kw.hybrid_lmax,
spinor = kw.spinor,
)
#end if
if kw.det_format=='new':
if kw.excitation is not None:
QmcpackInput.class_error('user provided "excitation" input argument with new style determinant format.\nPlease add det_format="old" and try again','generate_qmcpack_input')
#end if
if kw.system is not None and isinstance(kw.system.structure,Jellium):
ssb = generate_sposet_builder(
type = 'heg',
twist = kw.twist,
spin_polarized = kw.spin_polarized,
system = kw.system,
)
else:
if kw.orbspline is None:
kw.orbspline = 'bspline'
#end if
if kw.orbitals_h5!='MISSING.h5':
orbfile_exists = os.path.exists(kw.orbitals_h5)
if kw.check_paths and not orbfile_exists:
QmcpackInput.class_error('user provided "orbitals_h5" path does not exist\nPath provided: {}\nTo disable this check, set check_paths=False'.format(kw.orbitals_h5),'generate_qmcpack_input')
#end if
if kw.run_path is not None:
kw.orbitals_h5 = os.path.relpath(kw.orbitals_h5,kw.run_path)
#end if
#end if
ssb = generate_sposet_builder(
type = kw.orbspline,
twist = kw.twist,
twistnum = kw.twistnum,
meshfactor = kw.meshfactor,
precision = kw.precision,
truncate = kw.truncate,
buffer = kw.buffer,
hybridrep = kw.hybridrep,
href = kw.orbitals_h5,
rotate = kw.rotated_orbitals,
spin_polarized = kw.spin_polarized,
system = kw.system,
orbitals_cpu = kw.orbitals_cpu,
gpusharing = kw.gpusharing,
spinor = kw.spinor,
)
#end if
if kw.partition is None:
spobuilders = [ssb]
else:
spobuilders = partition_sposets(
sposet_builder = ssb,
partition = kw.partition,
partition_meshfactors = kw.partition_mf,
)
#end if
dset = generate_determinantset(
spin_polarized = kw.spin_polarized,
delay_rank = kw.delay_rank,
det_batch = kw.det_batch,
matrix_inv_cpu = kw.matrix_inv_cpu,
system = kw.system,
spinor = kw.spinor,
rotate = kw.rotated_orbitals,
)
elif kw.det_format=='old':
spobuilders = None
if kw.orbspline is None:
kw.orbspline = 'einspline'
#end if
dset = generate_determinantset_old(
type = kw.orbspline,
twistnum = kw.twistnum,
meshfactor = kw.meshfactor,
precision = kw.precision,
hybridrep = kw.hybridrep,
href = kw.orbitals_h5,
spin_polarized = kw.spin_polarized,
excitation = kw.excitation,
delay_rank = kw.delay_rank,
gpusharing = kw.gpusharing,
system = kw.system,
spinor = kw.spinor,
)
else:
QmcpackInput.class_error('argument "det_format" is invalid.\nReceived: {0}\nValid options are: new, old'.format(det_format),'generate_qmcpack_input')
#end if
wfn = wavefunction(
name = 'psi0',
target = 'e',
determinantset = dset,
)
if isinstance(kw.jastrows,str) and kw.jastrows.endswith('.xml'):
if not os.path.exists(kw.jastrows):
QmcpackInput.class_error('user provided "jastrows" file path does not exist\nFile path provided: {}'.format(kw.jastrows),'generate_qmcpack_input')
#end if
jastrows = read_jastrows(kw.jastrows)
if jastrows is None:
QmcpackInput.class_error('no jastrows found at user provided "jastrows" file.\nFile path provided: {}'.format(kw.jastrows),'generate_qmcpack_input')
#end if
kw.jastrows = jastrows
elif kw.J1 or kw.J2 or kw.J3:
kw.jastrows = generate_jastrows_alt(
J1 = kw.J1 ,
J2 = kw.J2 ,
J3 = kw.J3 ,
J1_size = kw.J1_size ,
J1_rcut = kw.J1_rcut ,
J1_dr = kw.J1_dr ,
J1_opt = kw.J1_opt ,
J2_size = kw.J2_size ,
J2_rcut = kw.J2_rcut ,
J2_dr = kw.J2_dr ,
J2_init = kw.J2_init ,
J2_opt = kw.J2_opt ,
J3_isize = kw.J3_isize ,
J3_esize = kw.J3_esize ,
J3_rcut = kw.J3_rcut ,
J3_opt = kw.J3_opt ,
J1_rcut_open = kw.J1_rcut_open,
J2_rcut_open = kw.J2_rcut_open,
J1k = kw.J1k ,
J1k_kcut = kw.J1k_kcut ,
J1k_symm = kw.J1k_symm ,
J1k_opt = kw.J1k_opt ,
J2k = kw.J2k ,
J2k_kcut = kw.J2k_kcut ,
J2k_symm = kw.J2k_symm ,
J2k_opt = kw.J2k_opt ,
system = kw.system ,
)
#end if
if kw.jastrows is not None:
wfn.jastrows = generate_jastrows(kw.jastrows,kw.system,check_ions=True)
#end if
if kw.spinor is not None and kw.spinor:
# remove u-d
# also set correct cusp
J2 = wfn.jastrows.get('J2')
if J2 is not None:
corr = J2.get('correlation')
if 'ud' in corr:
del corr.ud
if 'uu' in corr:
corr.uu.cusp = -0.5
#end if
#end if
#end if
J3 = wfn.jastrows.get('J3')
if J3 is not None:
corr = J3.get('correlation')
j3_ids = []
for j3_term in corr:
j3_id = j3_term.coefficients.id
j3_ids.append(j3_id)
#end for
for j3_id in j3_ids:
if 'ud' in j3_id:
delattr(corr, j3_id)
#end if
#end for
#end if
#end if
if spobuilders is not None:
wfn.sposet_builders = make_collection(spobuilders)
#end if
if kw.opt_params is not None:
if not isinstance(kw.opt_params,str):
QmcpackInput.class_error('opt_params must be a file path.\nYou provided: {}'.format(kw.opt_params),'generate_qmcpack_input')
#end if
if not kw.opt_params.endswith('vp.h5'):
QmcpackInput.class_error('opt_params must a vp.h5 file.\nYou provided: {}'.format(kw.opt_params),'generate_qmcpack_input')
#end if
if kw.check_paths and not os.path.exists(kw.opt_params):
QmcpackInput.class_error('opt_params file does not exist.\nFile path provided: {}\nTo disable this check, set check_paths=False'.format(kw.opt_params),'generate_qmcpack_input')
#end if
wfn.override_variational_parameters = override_variational_parameters(
href = os.path.abspath(kw.opt_params)
)
#end if
h_estimators = kw.estimators
d_estimators = None
if batched:
h_estimators = []
d_estimators = []
if isinstance(kw.estimators,list):
for est in kw.estimators:
if isinstance(est,str) and est.lower()=='mpc':
h_estimators.append(est)
else:
d_estimators.append(est)
#end if
#end for
#end if
#end if
hmltn = generate_hamiltonian(
system = kw.system,
pseudos = kw.pseudos,
algorithm = kw.pseudo_algorithm,
dla = kw.dla,
interactions = kw.interactions,
estimators = h_estimators,
wf_elem = wfn,
)
qmcsys = qmcsystem(
simulationcell = simcell,
wavefunction = wfn,
hamiltonian = hmltn,
)
if kw.system is not None:
qmcsys.particlesets = particlesets
#end if
sim = simulation(
project = proj,
qmcsystem = qmcsys,
)
if kw.seed is not None:
sim.random = random(seed=kw.seed)
#end if
if kw.traces is not None:
sim.traces = kw.traces
#end if
if len(kw.calculations)==0 and kw.qmc is not None:
qmc_inputs = kw.obj(*qmc_keys)
if kw.qmc=='opt':
kw.calculations = generate_opt_calculations(**qmc_inputs)
elif 'vmc' in kw.qmc:
kw.calculations = generate_vmc_calculations(**qmc_inputs)
elif 'dmc' in kw.qmc:
kw.calculations = generate_dmc_calculations(**qmc_inputs)
#end if
#end if
if batched and d_estimators is not None and len(d_estimators)>0:
ests = generate_estimators_batched(d_estimators,wf_elem=wfn)
ests_elem = estimators()
ests_elem.estimators = ests
sim.qmcsystem.estimators = ests_elem
#end if
if kw.estimator_period is not None:
for c in kw.calculations:
if isinstance(c,(vmc,dmc)):
c.estimator_period = kw.estimator_period
#end if
#end for
#end if
sim.calculations = make_collection(kw.calculations).copy()
qi = QmcpackInput(metadata,sim)
qi.incorporate_defaults(elements=False,overwrite=False,propagate=True)
if kw.remove_cell:
qi.remove_physical_system()
#end if
for calc in sim.calculations:
if isinstance(calc,loop):
calc = calc.qmc
#end if
if isinstance(calc,(linear,cslinear,linear_batch)) and 'nonlocalpp' not in calc and not batched:
calc.nonlocalpp = True
calc.use_nonlocalpp_deriv = True
#end if
#end for
return qi
#end def generate_basic_input
gen_basic_afqmc_input_defaults = obj(
id = 'qmc',
series = 0,
seed = None,
nmo = None,
naea = None,
naeb = None,
ham_file = None,
wfn_file = None,
wfn_type = 'NOMSD',
cutoff = 1e-8,
wset_type = 'shared',
walker_type = 'CLOSED',
hybrid = True,
ncores = 1,
nwalkers = 10,
blocks = 10000,
steps = 10,
timestep = 0.005,
estimators = None,
info_name = 'info0',
ham_name = 'ham0',
wfn_name = 'wfn0',
wset_name = 'wset0',
prop_name = 'prop0',
system = None,
run_path = None,
)
[docs]
def generate_basic_afqmc_input(**kwargs):
# capture inputs
kw = obj(kwargs)
gen_info = obj()
for k,v in kw.items():
if not isinstance(v,obj):
gen_info[k] = v
#end if
#end for
# apply general defaults
kw.set_optional(**gen_basic_afqmc_input_defaults)
valid = set(gen_basic_afqmc_input_defaults.keys())
# screen for invalid keywords
invalid_kwargs = set(kw.keys())-valid
if len(invalid_kwargs)>0:
QmcpackInput.class_error('invalid input parameters encountered\ninvalid input parameters: {0}\nvalid options are: {1}'.format(sorted(invalid_kwargs),sorted(valid)),'generate_qmcpack_input')
#end if
metadata = meta(
generation_info = gen_info.copy(),
)
sim = simulation(
method = 'afqmc',
)
sim.project = project(
id = kw.id,
series = kw.series,
)
if kw.seed is not None:
sim.random = random(seed=kw.seed)
#end if
info = afqmcinfo(
name = kw.info_name,
)
if kw.nmo is not None:
info.nmo = kw.nmo
#end if
if kw.naea is not None:
info.naea = kw.naea
#end if
if kw.naeb is not None:
info.naeb = kw.naeb
#end if
sim.afqmcinfo = info
if kw.ham_file is None and kw.wfn_file is not None:
kw.ham_file = kw.wfn_file
elif kw.ham_file is not None and kw.wfn_file is None:
kw.wfn_file = kw.ham_file
elif kw.ham_file is None and kw.wfn_file is None:
kw.ham_file = 'MISSING.h5'
kw.wfn_file = 'MISSING.h5'
#end if
def get_filetype(filename,loc):
if filename.endswith('.h5'):
filetype = 'hdf5'
else:
QmcpackInput.class_error('Type of {} file "{}" is unrecognized.\n The following file extensions are allowed: .h5'.format(loc,filename))
#end if
return filetype
#end def get_filetype
ham = hamiltonian(
name = kw.ham_name,
info = info.name,
filetype = get_filetype(kw.ham_file,'hamiltonian'),
filename = kw.ham_file,
)
sim.hamiltonian = ham
wfn = wavefunction(
name = kw.wfn_name,
info = info.name,
filetype = get_filetype(kw.wfn_file,'wavefunction'),
filename = kw.wfn_file,
)
if kw.wfn_type is not None:
wfn.type = kw.wfn_type
#end if
if kw.cutoff is not None:
wfn.cutoff = kw.cutoff
#end if
sim.wavefunction = wfn
wset = walkerset(
name = kw.wset_name,
)
if kw.wset_type is not None:
wset.type = kw.wset_type
#end if
if kw.walker_type is not None:
wset.walker_type = kw.walker_type
#end if
sim.walkerset = wset
prop = propagator(
name = kw.prop_name,
info = info.name,
)
if kw.hybrid is not None:
prop.hybrid = kw.hybrid
#end if
sim.propagator = prop
exe = execute(
info = info.name,
ham = ham.name,
wfn = wfn.name,
wset = wset.name,
prop = prop.name,
)
for k in execute.parameters:
if k in kw and kw[k] is not None:
exe[k] = kw[k]
#end if
#end for
estimators = []
valid_estimators = (back_propagation,)
if kw.estimators is not None:
for est in kw.estimators:
invalid = False
if isinstance(est,QIxml):
est = est.copy()
else:
invalid = True
#end if
invalid |= not isinstance(est,valid_estimators)
if invalid:
valid_names = [e.__class__.__name__ for e in valid_estimators]
QmcpackInput.class_error('invalid estimator input encountered\nexpected one of the following: {}\ninputted type: {}\ninputted value: {}'.format(valid_names,est.__class__.__name__,est))
#end if
est.incorporate_defaults()
estimators.append(est)
#end for
#end if
if len(estimators)>0:
exe.estimators = make_collection(estimators)
#end if
sim.execute = exe
qi = QmcpackInput(metadata,sim)
return qi
#end def generate_basic_afqmc_input
[docs]
def generate_opt_jastrow_input(id = 'qmc',
series = 0,
purpose = '',
seed = None,
bconds = None,
remove_cell = False,
meshfactor = 1.0,
precision = 'float',
twistnum = None,
twist = None,
spin_polarized = False,
orbitals_h5 = 'MISSING.h5',
system = None,
pseudos = None,
jastrows = 'generateJ12',
corrections = None,
observables = None,
processes = None,
walkers_per_proc = None,
threads = None,
decorr = 10,
min_walkers = None, #use e.g. 128 for gpu's
timestep = 0.5,
nonlocalpp = False,
sample_factor = 1.0,
opt_calcs = None,
det_format = 'new'):
jastrows = generate_jastrows(jastrows,system)
if opt_calcs is None:
opt_calcs = [
('linear', 4, 0, 0, 1.0),
('linear', 4, .8, .2, 0)
]
#end if
opts = []
for opt_calc in opt_calcs:
if isinstance(opt_calc,QIxml):
opts.append(opt_calc)
elif len(opt_calc)==5:
if opt_calc[0] in opt_map:
opts.append(
generate_opt(
*opt_calc,
jastrows = jastrows,
processes = processes,
walkers_per_proc = walkers_per_proc,
threads = threads,
decorr = decorr,
min_walkers = min_walkers,
timestep = timestep,
nonlocalpp = nonlocalpp,
sample_factor = sample_factor
)
)
else:
QmcpackInput.class_error('optimization method '+opt_calc[0]+' has not yet been implemented')
#end if
else:
QmcpackInput.class_error('optimization calculation is ill formatted\n opt calc provided: \n'+str(opt_calc))
#end if
#end if
input = generate_basic_input(
id = id ,
series = series ,
purpose = purpose ,
seed = seed ,
bconds = bconds ,
remove_cell = remove_cell ,
meshfactor = meshfactor ,
precision = precision ,
twistnum = twistnum ,
twist = twist ,
spin_polarized = spin_polarized ,
orbitals_h5 = orbitals_h5 ,
system = system ,
pseudos = pseudos ,
jastrows = jastrows ,
corrections = corrections ,
observables = observables ,
calculations = opts ,
det_format = det_format ,
)
return input
#end def generate_opt_jastrow_input
if __name__=='__main__':
filepath = './example_input_files/c_boron/qmcpack.in.xml'
element_joins=['qmcsystem']
element_aliases=dict(loop='qmc')
xml = XMLreader(filepath,element_joins,element_aliases,warn=False).obj
xml.condense()
qi = QmcpackInput()
qi.read(filepath)
s = qi.simulation
q = s.qmcsystem
c = s.calculations
h = q.hamiltonian
p = q.particlesets
w = q.wavefunction
j = w.jastrows
co= j.J1.correlations.B.coefficients
qi.write('./output/qmcpack.in.xml')
#qi.condensed_name_report()
#exit()
test_ret_system = 1
test_gen_input = 0
test_difference = 0
test_moves = 0
test_defaults = 0
test_substitution = 0
test_generation = 0
if test_ret_system:
from .structure import generate_structure
from .physical_system import PhysicalSystem
system = PhysicalSystem(
structure = generate_structure('diamond','fcc','Ge',(2,2,2),scale=5.639,units='A'),
net_charge = 1,
net_spin = 1,
Ge = 4
)
gi = generate_qmcpack_input('basic',system=system)
rsys = gi.return_system()
print(rsys)
#end if
if test_gen_input:
from .structure import generate_structure
from .physical_system import PhysicalSystem
system = PhysicalSystem(
structure = generate_structure('diamond','fcc','Ge',(2,2,2),scale=5.639,units='A'),
net_charge = 1,
net_spin = 1,
Ge = 4
)
gi = generate_qmcpack_input('basic',system=system)
print(gi)
print(gi.write())
#end if
if test_difference:
tstep = QmcpackInput('./example_input_files/luke_tutorial/diamond-dmcTsteps.xml')
opt = QmcpackInput('./example_input_files/luke_tutorial/opt-diamond.xml')
different,diff,d1,d2 = tstep.difference(tstep)
different,diff,d1,d2 = tstep.difference(opt)
#end if
if test_moves:
print(50*'=')
sim = qi.simulation
print(repr(sim))
print(repr(sim.qmcsystem))
print(50*'=')
qi.move(particleset='simulation')
print(repr(sim))
print(repr(sim.qmcsystem))
print(50*'=')
qi.standard_placements()
print(repr(sim))
print(repr(sim.qmcsystem))
qi.pluralize()
#end if
if test_defaults:
q=QmcpackInput(
simulation(
qmcsystem=section(
simulationcell = section(),
wavefunction = section(),
hamiltonian = section()
),
calculations = [
cslinear(),
vmc(),
dmc()
]
)
)
#q.simulation = simulation()
q.incorporate_defaults(elements=True)
print(q)
#end if
if test_substitution:
q = qi.copy()
q.remove('simulationcell','particleset','wavefunction')
q.write('./output/qmcpack.remove.xml')
q.include_xml('./example_input_files/energy_density/Si.ptcl.xml',replace=False)
q.include_xml('./example_input_files/energy_density/Si.wfs.xml',replace=False)
q.write('./output/qmcpack.replace.xml')
qnj = QmcpackInput()
qnj.read('./example_input_files/jastrowless/opt_jastrow.in.xml')
qnj.generate_jastrows(size=6)
qnj.write('./output/jastrow_gen.in.xml')
#end if
if test_generation:
q=QmcpackInput(
meta(
lattice = {'units':'bohr'},
reciprocal = {'units':'2pi/bohr'},
ionid = {'datatype':'stringArray'},
position = {'datatype':'posArray', 'condition':0}
),
simulation(
project = section(
id='C16B',
series = 0,
application = section(
name = 'qmcpack',
role = 'molecu',
class_ = 'serial',
version = .2
),
host = 'kraken',
date = '3 May 2012',
user = 'jtkrogel'
),
random = section(seed=13),
qmcsystem = section(
simulationcell = section(
name = 'global',
lattice = np.array([[1,1,0],[1,0,1],[0,1,1]]),
reciprocal = np.array([[1,1,-1],[1,-1,1],[-1,1,1]]),
bconds = 'p p p',
LR_dim_cutoff = 15
),
particlesets = [
particleset(
name = 'ion0',
size = 32,
groups=[
group(
name='C',
charge=4.
),
group(
name='B',
charge = 3.
)
],
ionid = ['B','C','C','C','C','C','C','C','C','C','C','C','C','C','C','C',
'B','C','C','C','C','C','C','C','C','C','C','C','C','C','C','C'],
position = np.array([
[ 0.00, 0.00, 0.00],[ 1.68, 1.68, 1.68],[ 3.37, 3.37, 0.00],
[ 5.05, 5.05, 1.68],[ 3.37, 0.00, 3.37],[ 5.05, 1.68, 5.05],
[ 6.74, 3.37, 3.37],[ 8.42, 5.05, 5.05],[ 0.00, 3.37, 3.37],
[ 1.68, 5.05, 5.05],[ 3.37, 6.74, 3.37],[ 5.05, 8.42, 5.05],
[ 3.37, 3.37, 6.74],[ 5.05, 5.05, 8.42],[ 6.74, 6.74, 6.74],
[ 8.42, 8.42, 8.42],[ 6.74, 6.74, 0.00],[ 8.42, 8.42, 1.68],
[10.11,10.11, 0.00],[11.79,11.79, 1.68],[10.11, 6.74, 3.37],
[11.79, 8.42, 5.05],[13.48,10.11, 3.37],[15.16,11.79, 5.05],
[ 6.74,10.11, 3.37],[ 8.42,11.79, 5.05],[10.11,13.48, 3.37],
[11.79,15.16, 5.05],[10.11,10.11, 6.74],[11.79,11.79, 8.42],
[13.48,13.48, 6.74],[15.16,15.16, 8.42]])
),
particleset(
name='e',
random = 'yes',
random_source = 'ion0',
groups=[
group(
name='u',
size=64,
charge=-1
),
group(
name='d',
size=63,
charge=-1
)
]
),
],
hamiltonians = [
hamiltonian(
name='h0',
type='generic',
target='e',
pairpots=[
pairpot(
type = 'coulomb',
name = 'ElecElec',
source = 'e',
target = 'e'
),
pairpot(
type = 'pseudo',
name = 'PseudoPot',
source = 'ion0',
wavefunction='psi0',
format='xml',
pseudos = [
pseudo(
elementtype='B',
href='B.pp.xml'
),
pseudo(
elementtype='C',
href='C.pp.xml'
)
]
)
],
constant = section(
type='coulomb',
name='IonIon',
source='ion0',
target='ion0'
),
estimators = [
estimator(
type='energydensity',
name='edvoronoi',
dynamic='e',
static='ion0',
spacegrid = section(
coord = 'voronoi'
)
),
energydensity(
name='edchempot',
dynamic='e',
static='ion0',
spacegrid=spacegrid(
coord='voronoi',
min_part=-4,
max_part=5
)
),
estimator(
type='energydensity',
name='edcell',
dynamic='e',
static='ion0',
spacegrid = section(
coord = 'cartesian',
origin = section(p1='zero'),
axes = (
axis(label='x',p1='a1',scale=.5,grid='-1 (192) 1'),
axis(label='y',p1='a2',scale=.5,grid='-1 (1) 1'),
axis(label='z',p1='a3',scale=.5,grid='-1 (1) 1')
)
# axes = collection(
# x = section(p1='a1',scale=.5,grid='-1 (192) 1'),
# y = section(p1='a2',scale=.5,grid='-1 (1) 1'),
# z = section(p1='a3',scale=.5,grid='-1 (1) 1')
# )
)
)
]
)
],
wavefunction = section(
name = 'psi0',
target = 'e',
determinantset = section(
type='bspline',
href='Si.pwscf.h5',
sort = 1,
tilematrix = np.array([[1,0,0],[0,1,0],[0,0,1]]),
twistnum = 0,
source = 'ion0',
slaterdeterminant = section(
determinants=[
determinant(
id='updet',
size=64,
occupation = section(
mode='ground',
spindataset=0
)
),
determinant(
id='downdet',
size=63,
occupation = section(
mode='ground',
spindataset=1
)
)
]
),
),
jastrows = [
jastrow(
type='two-body',
name='J2',
function='bspline',
print='yes',
correlations = [
correlation(
speciesA='u',
speciesB='u',
size=6,
rcut=3.9,
coefficients = section(
id='uu',
type='Array',
coeff=[0,0,0,0,0,0]
)
),
correlation(
speciesA='u',
speciesB='d',
size=6,
rcut=3.9,
coefficients = section(
id='ud',
type='Array',
coeff=[0,0,0,0,0,0]
)
)
]
),
jastrow(
type='one-body',
name='J1',
function='bspline',
source='ion0',
print='yes',
correlations = [
correlation(
elementtype='C',
size=6,
rcut=3.9,
coefficients = section(
id='eC',
type='Array',
coeff=[0,0,0,0,0,0]
)
),
correlation(
elementtype='B',
size=6,
rcut=3.9,
coefficients = section(
id='eB',
type='Array',
coeff=[0,0,0,0,0,0]
)
)
]
)
]
),
),
calculations=[
loop(max=4,
qmc=qmc(
method='cslinear',
move='pbyp',
checkpoint=-1,
gpu='no',
blocks = 3125,
warmupsteps = 5,
steps = 2,
samples = 80000,
timestep = .5,
usedrift = 'yes',
minmethod = 'rescale',
gevmethod = 'mixed',
exp0=-15,
nstabilizers = 5,
stabilizerscale = 3,
stepsize=.35,
alloweddifference=1e-5,
beta = .05,
bigchange = 5.,
energy = 0.,
unreweightedvariance = 0.,
reweightedvariance = 0.,
estimators=[
estimator(
name='LocalEnergy',
hdf5='no'
)
]
)
),
qmc(
method = 'vmc',
multiple = 'no',
warp = 'no',
move = 'pbyp',
walkers = 1,
blocks = 2,
steps = 500,
substeps = 3,
timestep = .5,
usedrift = 'yes',
estimators=[
estimator(
name='LocalEnergy',
hdf5='yes'
)
]
),
qmc(
method='dmc',
move='pbyp',
walkers = 72,
blocks = 2,
steps = 50,
timestep = .01,
nonlocalmove = 'yes',
estimators=[
estimator(
name='LocalEnergy',
hdf5='no'
)
]
)
]
)
)
q.write('./output/gen.in.xml')
#something broke this, check later
exit()
qs=QmcpackInput(
simulation = section(
project = section(
id='C16B',series = 0,
application = section(name='qmcpack',role='molecu',class_='serial',version=.2),
host='kraken',date='3 May 2012',user='jtkrogel'
),
random = section(seed=13),
qmcsystem = section(
simulationcell = section(
name='global',bconds='p p p',lr_dim_cutoff=15,
lattice = [[1,1,0] ,[1,0,1] ,[0,1,1]],
reciprocal = [[1,1,-1],[1,-1,1],[-1,1,1]],
),
particlesets = collection(
ion0=particleset(
size=32,
groups=collection(
C = group(charge=4.),
B = group(charge=3.)),
ionid = ('B','C','C','C','C','C','C','C','C','C','C','C','C','C','C','C',
'B','C','C','C','C','C','C','C','C','C','C','C','C','C','C','C'),
position = [[ 0.00, 0.00, 0.00],[ 1.68, 1.68, 1.68],[ 3.37, 3.37, 0.00],
[ 5.05, 5.05, 1.68],[ 3.37, 0.00, 3.37],[ 5.05, 1.68, 5.05],
[ 6.74, 3.37, 3.37],[ 8.42, 5.05, 5.05],[ 0.00, 3.37, 3.37],
[ 1.68, 5.05, 5.05],[ 3.37, 6.74, 3.37],[ 5.05, 8.42, 5.05],
[ 3.37, 3.37, 6.74],[ 5.05, 5.05, 8.42],[ 6.74, 6.74, 6.74],
[ 8.42, 8.42, 8.42],[ 6.74, 6.74, 0.00],[ 8.42, 8.42, 1.68],
[10.11,10.11, 0.00],[11.79,11.79, 1.68],[10.11, 6.74, 3.37],
[11.79, 8.42, 5.05],[13.48,10.11, 3.37],[15.16,11.79, 5.05],
[ 6.74,10.11, 3.37],[ 8.42,11.79, 5.05],[10.11,13.48, 3.37],
[11.79,15.16, 5.05],[10.11,10.11, 6.74],[11.79,11.79, 8.42],
[13.48,13.48, 6.74],[15.16,15.16, 8.42]]
),
e=particleset(
random='yes',random_source='ion0',
groups = collection(
u=group(size=64,charge=-1),
d=group(size=63,charge=-1))
),
),
hamiltonian = section(
name='h0',type='generic',target='e',
pairpots=collection(
ElecElec = coulomb(name='ElecElec',source='e',target='e'),
PseudoPot = pseudopotential(
source='ion0',wavefunction='psi0',format='xml',
pseudos = collection(
B = pseudo(href='B.pp.xml'),
C = pseudo(href='C.pp.xml'))
)
),
constant = section(type='coulomb',name='IonIon',source='ion0',target='ion0'),
estimators = collection(
edvoronoi = energydensity(
dynamic='e',static='ion0',spacegrid=section(coord ='voronoi')
),
edchempot = energydensity(
dynamic='e',static='ion0',
spacegrid=section(coord='voronoi',min_part=-4,max_part=5)
),
edcell = energydensity(
dynamic='e',static='ion0',
spacegrid = section(
coord = 'cartesian',
origin = section(p1='zero'),
axes = collection(
x = axis(p1='a1',scale=.5,grid='-1 (192) 1'),
y = axis(p1='a2',scale=.5,grid='-1 (1) 1'),
z = axis(p1='a3',scale=.5,grid='-1 (1) 1'))
)
)
)
),
wavefunction = section(
name = 'psi0',target = 'e',
determinantset = section(
type='bspline',href='Si.pwscf.h5',sort=1,twistnum=0,source='ion0',
tilematrix=(1,0,0,0,1,0,0,0,1),
slaterdeterminant = section(
determinants=collection(
updet = determinant(
size=64,
occupation=section(mode='ground',spindataset=0)
),
downdet = determinant(
size=63,
occupation = section(mode='ground',spindataset=1))
)
),
),
jastrows = collection(
J2=jastrow2(
function='bspline',print='yes',
correlations = collection(
uu=correlation(
speciesA='u',speciesB='u',size=6,rcut=3.9,
coefficients = section(id='uu',type='Array',coeff=[0,0,0,0,0,0])
),
ud=correlation(
speciesA='u',speciesB='d',size=6,rcut=3.9,
coefficients = section(id='ud',type='Array',coeff=[0,0,0,0,0,0])
)
)
),
J1=jastrow1(
function='bspline',source='ion0',print='yes',
correlations = collection(
C=correlation(
size=6,rcut=3.9,
coefficients = section(
id='eC',type='Array',coeff=[0,0,0,0,0,0])
),
B=correlation(
size=6,rcut=3.9,
coefficients = section(id='eB',type='Array',coeff=[0,0,0,0,0,0])
)
)
)
)
),
),
calculations=(
loop(max=4,
qmc=cslinear(
move='pbyp',checkpoint=-1,gpu='no',
blocks = 3125,
warmupsteps = 5,
steps = 2,
samples = 80000,
timestep = .5,
usedrift = 'yes',
minmethod = 'rescale',
gevmethod = 'mixed',
exp0 = -15,
nstabilizers = 5,
stabilizerscale = 3,
stepsize = .35,
alloweddifference = 1e-5,
beta = .05,
bigchange = 5.,
energy = 0.,
unreweightedvariance = 0.,
reweightedvariance = 0.,
estimator = localenergy(hdf5='no')
)
),
vmc(multiple='no',warp='no',move='pbyp',
walkers = 1,
blocks = 2,
steps = 500,
substeps = 3,
timestep = .5,
usedrift = 'yes',
estimator = localenergy(hdf5='no')
),
dmc(move='pbyp',
walkers = 72,
blocks = 2,
steps = 50,
timestep = .01,
nonlocalmove = 'yes',
estimator = localenergy(hdf5='yes')
)
)
)
)
qs.write('./output/simple.in.xml')
est = qs.simulation.qmcsystem.hamiltonian.estimators
sg = est.edcell.spacegrid
print(repr(est))
exit()
q=QmcpackInput()
q.simulation = section(
project = section('C16B',0,
application = section('qmcpack','molecu','serial',.2),
host = 'kraken',
date = '3 May 2012',
user = 'jtkrogel'
),
random = (13),
qmcsystem = section(
simulationcell = section(
units = 'bohr',
lattice = np.array([[1,1,0],[1,0,1],[0,1,1]]),
bconds = 'p p p',
LR_dim_cutoff = 15
),
particlesets = [
particleset('ion0', ('C',4), ('B',3),
ionid = ['B','C','C','C','C','C','C','C','C','C','C','C','C','C','C','C',
'B','C','C','C','C','C','C','C','C','C','C','C','C','C','C','C'],
position = np.array([
[ 0.00, 0.00, 0.00],[ 1.68, 1.68, 1.68],[ 3.37, 3.37, 0.00],
[ 5.05, 5.05, 1.68],[ 3.37, 0.00, 3.37],[ 5.05, 1.68, 5.05],
[ 6.74, 3.37, 3.37],[ 8.42, 5.05, 5.05],[ 0.00, 3.37, 3.37],
[ 1.68, 5.05, 5.05],[ 3.37, 6.74, 3.37],[ 5.05, 8.42, 5.05],
[ 3.37, 3.37, 6.74],[ 5.05, 5.05, 8.42],[ 6.74, 6.74, 6.74],
[ 8.42, 8.42, 8.42],[ 6.74, 6.74, 0.00],[ 8.42, 8.42, 1.68],
[10.11,10.11, 0.00],[11.79,11.79, 1.68],[10.11, 6.74, 3.37],
[11.79, 8.42, 5.05],[13.48,10.11, 3.37],[15.16,11.79, 5.05],
[ 6.74,10.11, 3.37],[ 8.42,11.79, 5.05],[10.11,13.48, 3.37],
[11.79,15.16, 5.05],[10.11,10.11, 6.74],[11.79,11.79, 8.42],
[13.48,13.48, 6.74],[15.16,15.16, 8.42]])
),
particleset('e', ('u',-1,64), ('d',-1,63), random_source = 'ion0'),
],
hamiltonian = section('h0','e',
pairpots=[
coulomb('ElecElec','e','e'),
pseudopotential('PseudoPot','ion0','psi0',('B','B.pp.xml'),('C','C.pp.xml')),
coulomb('IonIon','ion0','ion0'),
],
estimators = [
energydensity('edvoronoi','e','ion0','voronoi',-4,5),
energydensity('edcell','e','ion0',
spacegrid('cartesian',
origin = 'zero',
x = ('a1',.5,'-1 (192) 1'),
y = ('a2',.5,'-1 (1) 1'),
z = ('a3',.5,'-1 (1) 1')
)
)
]
),
wavefunction = section('psi0','e',
determinantset = section('bspline','Si.pwscf.h5','ion0',
sort = 1,
tilematrix = np.array([[1,0,0],[0,1,0],[0,0,1]]),
twistnum = 0,
slaterdeterminant = [
determinant('updet',64,'ground',0),
determinant('downdet',63,'ground',1)
],
jastrows = [
twobody('J2','bspline',
('u','u',3.9,[0,0,0,0,0,0]),
('u','d',3.9,[0,0,0,0,0,0])),
onebody('J1','bspline','ion0',
('C',3.9,[0,0,0,0,0,0]),
('B',3.9,[0,0,0,0,0,0]))
]
)
)
),
calculations=[
loop(4,
cslinear(
blocks = 3125,
warmupsteps = 5,
steps = 2,
samples = 80000,
timestep = .5,
minmethod = 'rescale',
gevmethod = 'mixed',
exp0=-15,
nstabilizers = 5,
stabilizerscale = 3,
stepsize=.35,
alloweddifference=1e-5,
beta = .05,
bigchange = 5.,
energy = 0.,
unreweightedvariance = 0.,
reweightedvariance = 0.,
estimator = localenergy(hdf5='no')
)
),
vmc(
blocks = 2,
steps = 500,
substeps = 3,
timestep = .5,
estimator = localenergy(hdf5='yes')
),
dmc(
walkers = 72,
blocks = 2,
steps = 50,
timestep = .01,
nonlocalmove = 'yes',
estimator = localenergy(hdf5='no')
)
]
)
#end if
#end if