scnc
/
qmcpack
mirror of https://github.com/QMCPACK/qmcpack.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity
qmcpack/nexus/lib/qmcpack_input.py

8915 lines
309 KiB
Python
Raw Permalink Blame History

##################################################################
## (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
import inspect
import keyword
import numpy as np
from numpy import fromstring,empty,array,float64,\
loadtxt,ndarray,dtype,sqrt,pi,arange,exp,eye,\
ceil,mod,dot,abs,identity,floor,linalg,where,isclose
from io import StringIO
from superstring import string2val
from generic import obj,hidden
from xmlreader import XMLreader,XMLelement
from developer import DevBase,error
from periodic_table import is_element
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 debug import ci as interact
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}
def is_int(var):
try:
int(var)
return True
except ValueError:
return False
#end try
#end def is_int
def is_float(var):
try:
float(var)
return True
except ValueError:
return False
#end try
#end def is_float
def is_array(var,type):
try:
if isinstance(var,str):
array(var.split(),type)
else:
array(var,type)
#end if
return True
except ValueError:
return False
#end try
#end def is_float_array
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 = array(attr.split(),int)
if val.size==9:
val.shape = 3,3
#end if
elif is_array(attr,float):
val = array(attr.split(),float)
else:
val = attr
#end if
return val
#end def attribute_to_value
#local write types
def yesno(var):
return render_bool(var,'yes','no')
#end def yesno
def yesnostr(var):
if isinstance(var,str):
return var
else:
return yesno(var)
#end if
#end def yesnostr
def onezero(var):
return render_bool(var,'1','0')
#end def onezero
def truefalse(var):
return render_bool(var,'true','false')
#end def onezero
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])
class QIobj(DevBase):
afqmc_mode = False
# user settings
permissive_read = False
permissive_write = False
permissive_init = False
@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
class meta(obj):
None
#end class meta
class section(QIobj):
def __init__(self,*args,**kwargs):
self.args = args
self.kwargs = kwargs
#end def __init__
#end class section
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__
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 not element.tag 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
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
def get_single(self,preference=None):
if len(self)>0:
if preference!=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
def list(self):
lst = []
for key in self.hidden().order:
lst.append(self[key])
#end for
return lst
#end def list
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
def make_collection(elements):
return collection(*elements)
#end def make_collection
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
class QmcpackInputCollections(QIobj):
def add(self,element):
if element.tag in plurals_inv:
cname = plurals_inv[element.tag]
if not cname 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
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()
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)
@staticmethod
def set_expanded_names(**kwargs):
exnames = obj(**kwargs)
Names.expanded_names = exnames
Names.rsqmc_expanded_names = exnames
#end def set_expanded_names
@staticmethod
def set_afqmc_expanded_names(**kwargs):
Names.afqmc_expanded_names = obj(**kwargs)
#end def set_afqmc_expanded_names
@staticmethod
def use_rsqmc_expanded_names():
Names.expanded_names = Names.rsqmc_expanded_names
#end def use_rsqmc_expanded_names
@staticmethod
def use_afqmc_expanded_names():
Names.expanded_names = Names.afqmc_expanded_names
#end def use_afqmc_expanded_names
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
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
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
def condensed_name_report(self):
print()
print('Condensed Name Report:')
print('----------------------')
keylist = array(list(self.condensed_names.keys()))
order = 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
class QIxml(Names):
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
@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
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==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__
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!=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!=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!=None:
self[self.text] = param(xml)
#end if
#end def init_from_xml
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
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!=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
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!=None:
if overwrite or not name 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
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
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 not xname in pc:
pc[xname] = obj()
#end if
pc[xname].append(profile)
#end def collect_profile
def get_single(self,preference):
return self
#end def get_single
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!=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 = not name 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!=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 = not n 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
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
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,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
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,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==None:
self[var] = vnew
else:
not_coll = not isinstance(val,collection)
not_xml = not isinstance(val,QIxml)
not_arr = not isinstance(val,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
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!=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==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
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!=None and dest!=None and id(host)!=id(dest):
if not name in host:
name = plurals_inv[name]
#end if
dest[name] = host[name]
del host[name]
#end if
#end for
#end def move
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
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!=None:
kdiff[kk]=kkdiff
#end if
if kkd1!=None:
kd1[kk]=kkd1
#end if
if kkd2!=None:
kd2[kk]=kkd2
#end if
#end if
#end for
else:
if isinstance(value1,ndarray):
a1 = value1.ravel()
else:
a1 = array([value1])
#end if
if isinstance(value2,ndarray):
a2 = value2.ravel()
else:
a2 = 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 = abs(a1-a2).max()/max(1e-99,abs(a1).max(),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!=None:
diff[k] = kdiff
#end if
if kd1!=None:
d1[k] = kd1
#end if
if kd2!=None:
d2[k] = kd2
#end if
#end for
#end if
if root:
if diff!=None:
diff.remove_empty()
#end if
d1.remove_empty()
d2.remove_empty()
#end if
return different,diff,d1,d2
#end def difference
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
def get_host(self,names):
return self.get(names,host=True)
#end def get_host
def get_precision(self):
return self.__class__.class_get('precision')
#end def get_precision
#end class QIxml
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!=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__
def init_class(self):
None # this is for compatibility with QIxml only (do not overwrite)
#end def init_class
#end class QIxmlFactory
class Param(Names):
metadata = None
def __init__(self):
self.reset_precision()
#end def __init__
def reset_precision(self):
self.precision = None
self.prec_format = None
#end def reset_precision
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__
def read(self,xml):
val = ''
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:
token = xml.text.split('\n',1)[0].split(None,1)[0]
try:
if is_int(token):
val = loadtxt(StringIO(xml.text),int)
elif is_float(token):
val = loadtxt(StringIO(xml.text),float)
else:
val = array(xml.text.split())
#end if
except:
if is_int(token):
val = array(xml.text.split(),dtype=int)
elif is_float(token):
val = array(xml.text.split(),dtype=float)
else:
val = array(xml.text.split())
#end if
#end try
if val.size==1:
val = val.ravel()[0]
#end if
#end if
return val
#end def read
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 = array(value)
#end if
if attr_mode:
if isinstance(value,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,ndarray)
is_single = not (is_array and value.size>1)
if tag!=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!=None:
c+='\n'
#end if
ndim = len(value.shape)
if ndim==1:
line_len = 70
if tag!=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 == dtype(float):
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!=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!=None:
c+=pad+'</'+tag+'>\n'
#end if
#end if
return c
#end def write
def write_val(self,val):
if self.precision!=None and isinstance(val,float):
return self.prec_format.format(val)
else:
return str(val)
#end if
#end def write_val
def init_class(self):
None
#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'
)
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']
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','href']
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',
'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',
'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,
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'
)
def set_rsqmc_mode():
QIobj.afqmc_mode = False
Names.use_rsqmc_expanded_names()
#end def set_rsqmc_mode
def set_afqmc_mode():
QIobj.afqmc_mode = True
Names.use_afqmc_expanded_names()
#end def set_afqmc_mode
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)
)
@staticmethod
def settings(**kwargs):
QIobj.settings(**kwargs)
#end def settings
def __init__(self,arg0=None,arg1=None):
Param.metadata = None
filepath = None
metadata = None
element = None
if arg0==None and arg1==None:
None
elif isinstance(arg0,str) and arg1==None:
filepath = arg0
elif isinstance(arg0,QIxml) and arg1==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!=None:
self._metadata = metadata
else:
self._metadata = meta()
#end if
if filepath!=None:
self.read(filepath)
elif element!=None:
#simulation = arg0
#self.simulation = self.simulation_type(simulation)
elem_class = element.__class__
if elem_class.identifier!=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__
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
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
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
def read(self,filepath=None,xml=None):
if xml!=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!=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
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
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
def get(self,*names):
base = self.get_base()
return base.get(names)
#end def get
def remove(self,*names):
base = self.get_base()
base.remove(*names)
#end def remove
def assign(self,**kwargs):
base = self.get_base()
base.assign(**kwargs)
#end def assign
def replace(self,*args,**kwargs):# input is list of keyword=(oldval,newval)
base = self.get_base()
base.replace(*args,**kwargs)
#end def replace
def move(self,**elemdests):
base = self.get_base()
base.move(**elemdests)
#end def move
def get_host(self,names):
base = self.get_base()
return base.get_host(names)
#end if
def incorporate_defaults(self,elements=False,overwrite=False,propagate=False):
base = self.get_base()
base.incorporate_defaults(elements,overwrite,propagate)
#end def incorporate_defaults
def pluralize(self):
base = self.get_base()
base.pluralize()
#end def pluralize
def standard_placements(self):
self.move(particleset='qmcsystem',wavefunction='qmcsystem',hamiltonian='qmcsystem')
#end def standard_placements
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!=None:
diff.remove_empty()
#end if
d1.remove_empty()
d2.remove_empty()
return different,diff,d1,d2
#end def difference
def remove_empty(self):
base = self.get_base()
base.remove_empty()
#end def remove_empty
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
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==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!=None:
del host[section_name]
#end if
host[qname] = qxml
else:
if section_name==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.
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 not element_type 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 not placement 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 not ekey 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!=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 not pbef is None and paft is None:
cur_elems[elem,'on' ] = pbef
cur_elems[elem,'before'] = None
elif not paft is 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!=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
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
def generate_jastrows(self,size=None,j1func='bspline',j1size=8,j2func='bspline',j2size=8):
if size!=None:
j1size = size
j2size = size
#end if
#self.remove('jastrow')
lattice,particlesets,wavefunction = self.get('lattice','particleset','wavefunction')
no_lattice = lattice==None
no_particleset = particlesets==None
no_wavefunction = wavefunction==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 not 'e' 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 = sqrt(4*pi*density)
dr = rcut/size
r = .02 + dr*arange(size)
uuc = .5/(wp*r)*(1.-exp(-r*sqrt(wp/2)))*exp(-(2*r/rcut)**2)
udc = .5/(wp*r)*(1.-exp(-r*sqrt(wp)))*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 = abs(array(c.uu.coefficients.coeff)).sum() + abs(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
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!=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 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==None:
qs = qmcsystem()
qs.incorporate_defaults(elements=False,propagate=False)
self.simulation.qmcsystem = qs
#end if
if sc==None:
sc = simulationcell()
sc.incorporate_defaults(elements=False,propagate=False)
qs.simulationcell = sc
#end if
if ham==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!=None:
fs = structure.folded_structure
axes = array(pwscf_array_string(fs.axes).split(),dtype=float)
axes.shape = fs.axes.shape
axes = dot(structure.tmatrix,axes)
if 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 = array(pwscf_array_string(structure.axes).split(),dtype=float)
axes.shape = 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!=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==None:
pp = ham.get('PseudoPot')
if pp!=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!=None and not ion.name 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!=None:
self.replace(old_eps_name,'e')
#end if
if old_ips_name!=None and len(ions)>0:
self.replace(old_ips_name,'ion0')
#end if
udet,ddet = self.get('updet','downdet')
if udet!=None:
udet.size = elns.up_electron.count
#end if
if ddet!=None:
ddet.size = elns.down_electron.count
#end if
if abs(net_spin) > 1e-1:
if ddet!=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
def get_electron_particle_set(self):
input = self.copy()
input.pluralize()
return input.get('particlesets').e
#end def get_electron_particle_set
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!=None and 'groups' in elns:
simcell = input.get('simulationcell')
if simcell!=None and 'rs' in simcell:
have_ions = False
have_jellium = True
elif not 'pairpots' in H:
have_ions = False
#end if
#end if
if elns==None:
self.error('could not find electron particleset')
#end if
if ions==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 = array(elem)
pos = 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 = 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:
pos.shape = (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
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
def get_pp_files(self):
pp_files = []
h = self.get('hamiltonian')
if h != 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!=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
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
def cusp_correction(self):
cc = False
if not self.is_afqmc_input():
ds = self.get('determinantset')
cc_var = ds!=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
def get_qmc(self,series):
qmc = None
calcs = self.get('calculations')
series_start = self.get('series')
if calcs!=None:
if series_start is None:
qmc = calcs[series]
else:
qmc = calcs[series-series_start]
#end if
#end if
return qmc
#end def get_qmc
def bundle(self,inputs,filenames):
return BundledQmcpackInput(inputs,filenames)
#end def bundle
def trace(self,quantity,values):
return TracedQmcpackInput(quantity,values,self)
#end def trace
def twist_average(self,twistnums):
return self.trace('twistnum',twistnums)
#end def twist_average
#end class QmcpackInput
# base class for bundled qmcpack input
# not used on its own
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__
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
def generate_filenames(self,infile):
self.not_implemented()
#end def generate_filenames
def write(self,filepath=None):
if filepath!=None and not 'filenames' 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==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
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__
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
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
class QmcpackInputTemplate(SimulationInputTemplate):
def preprocess(self,contents,filepath=None):
if filepath!=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):
icont = open(include_path,'r').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
def get_output_info(self,*args,**kwargs):
# just pretend
return []
#end def get_output_info
#end class QmcpackInputTemplate
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!=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!=None:
fs = structure.folded_structure
axes = array(pwscf_array_string(fs.axes).split(),dtype=float)
axes.shape = fs.axes.shape
axes = dot(structure.tmatrix,axes)
if 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 = array(pwscf_array_string(structure.axes).split(),dtype=float)
axes.shape = structure.axes.shape
#end if
structure.adjust_axes(axes)
sc.lattice = axes
#end if
#end if
return sc
#end def generate_simulationcell
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
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
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
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 = 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!=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
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!=None:
hb.twist = tuple(twist)
#end if
return hb
#end def generate_heg_builder
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
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
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 = array(exc2.split(),dtype=int)
except:
format_failed = True
#end try
if not tmp is 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
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 = identity(3,dtype=int)
if system!=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 = 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 / abs(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 = 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[where(kpath_label == k_1)][0]
k_2 = kpath_rel[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 isclose(k_1, k).all():
k_1 = knum
found_k1 = True
#end if
if 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
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!=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,symbol = is_element(ion.name,symbol=True)
if label in ppfiles:
ppfile = ppfiles[label]
elif symbol in ppfiles:
ppfile = ppfiles[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!=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
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
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 not 'type' in spo:
QmcpackInput.class_error('cannot generate estimator dm1b\n basis sposet provided does not have a "type" attribute')
#end if
if not 'name' in spo:
spo.name = 'spo_dm'
#end if
builders = wf.get('sposet_builders')
if not spo.type 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
def generate_jastrows(jastrows,system=None,return_list=False,check_ions=False):
jin = []
have_ions = True
if check_ions and system!=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 not 'type' 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 not jtype 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
def generate_jastrows_alt(
J1 = False,
J2 = False,
J3 = False,
J1_size = None,
J1_rcut = None,
J1_dr = 0.5,
J2_size = None,
J2_rcut = None,
J2_dr = 0.5,
J2_init = 'zero',
J3_isize = 3,
J3_esize = 3,
J3_rcut = 5.0,
J1_rcut_open = 5.0,
J2_rcut_open = 10.0,
J1k = False,
J1k_kcut = 5.0,
J1k_symm = 'crystal',
J2k = False,
J2k_kcut = 5.0,
J2k_symm = 'crystal',
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(ceil(J1_rcut/J1_dr))
#end if
J = generate_jastrow('J1','bspline',J1_size,J1_rcut,system=system)
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(ceil(J2_rcut/J2_dr))
#end if
J = generate_jastrow('J2','bspline',J2_size,J2_rcut,init=J2_init,system=system)
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
J = generate_jastrow('J3','polynomial',J3_esize,J3_isize,J3_rcut,system=system)
jastrows.append(J)
#end if
if J1k or J2k:
Jk_inp = obj()
if J1k:
Jk_inp.kc1 = J1k_kcut
Jk_inp.symm1 = J1k_symm
#end if
if J2k:
Jk_inp.kc2 = J2k_kcut
Jk_inp.symm2 = J2k_symm
#end if
J = generate_kspace_jastrow(**Jk_inp)
jastrows.append(J)
#end if
return jastrows
#end def generate_jastrows_alt
def generate_jastrow(descriptor,*args,**kwargs):
keywords = set(['function','size','rcut','elements','coeff','cusp','ename',
'iname','spins','density','Buu','Bud','system','isize','esize','init'])
if not 'system' in kwargs:
kwargs['system'] = None
#end if
system = kwargs['system']
del kwargs['system']
if system!=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
def generate_jastrow1(function='bspline',size=8,rcut=None,coeff=None,cusp=0.,ename='e',iname='ion0',elements=None,system=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!=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 not 'coeff' 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 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
def generate_bspline_jastrow2(size=8,rcut=None,coeff=None,spins=('u','d'),density=None,system=None,init='rpa'):
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!=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 = sqrt(4*pi*density)
dr = rcut/size
r = .02 + dr*arange(size)
uuc = .5/(wp*r)*(1.-exp(-r*sqrt(wp/2)))*exp(-(2*r/rcut)**2)
udc = .5/(wp*r)*(1.-exp(-r*sqrt(wp)) )*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 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
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
def generate_jastrow2(function='bspline',*args,**kwargs):
if not 'spins' 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
def generate_jastrow3(function='polynomial',esize=3,isize=3,rcut=4.,coeff=None,iname='ion0',spins=('u','d'),elements=None,system=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!=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!=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=[]
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',optimize=True))
)
corrs.append(
correlation(
especies1=uname,especies2=dname,ispecies=element,esize=esize,
isize=isize,rcut=rcut,
coefficients=section(id=udname+element,type='Array',optimize=True))
)
#end for
jastrow = jastrow3(
name = 'J3',type='eeI',function=function,print=True,source=iname,
correlations = corrs
)
return jastrow
#end def generate_jastrow3
def generate_kspace_jastrow(
kc1 = None,
kc2 = None,
nk1 = 0,
nk2 = 0,
symm1 = 'isotropic',
symm2 = 'isotropic',
coeff1 = None,
coeff2 = None,
):
"""Generate <jastrow type="kSpace">
Parameters
----------
kc1 : float, optional
kcut for one-body Jastrow, default 0
kc2 : float, optional
kcut for two-body Jastrow, default 0
nk1 : int, optional
number of coefficients for one-body Jastrow, default 0
nk2 : int, optional
number of coefficients for two-body Jastrow, default 0
symm1 : str, optional
one of ['crystal', 'isotropic', 'none'], default 'isotropic'
symm2 : str, optional
one of ['crystal', 'isotropic', 'none'], default is 'isotropic'
coeff1 : list, optional
one-body Jastrow coefficients, default None
coeff2 : list, optional
list, optional two-body Jastrow coefficients, default None
Returns
-------
jk: QIxml
kspace_jastrow qmcpack_input element
"""
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),
)
corrs.append(corr1)
#end if
if J2k:
corr2 = correlation(
type = 'Two-Body',
symmetry = symm2,
kc = kc2,
coefficients = section(id='cG2', type='Array', coeff=coeff2),
)
corrs.append(corr2)
#end if
jk = kspace_jastrow(
type = 'kSpace',
name = 'Jk',
source = 'ion0',
correlations = collection(corrs),
)
return jk
# end def generate_kspace_jastrow
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
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)
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 not method 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!=None:
tot_walkers = max(min_walkers,tot_walkers)
walkers = int(ceil(float(tot_walkers)/processes-.001))
if threads!=None and mod(walkers,threads)!=0:
walkers = threads*int(ceil(float(walkers)/threads-.001))
#end if
#end if
#blocks = int(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!=None:
opt.energy = energy
#end if
if rw_variance!=None:
opt.reweightedvariance = rw_variance
#end if
if urw_variance!=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
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 = 'variance',
cycles = 12,
var_cycles = 0,
var_samples = None,
init_cycles = 0,
init_samples = None,
init_minwalkers = 1e-4,
)
shared_opt_batched_defaults = obj(
samples = 204800,
#nonlocalpp = True,
#use_nonlocalpp_deriv = True,
warmupsteps = 300,
blocks = 100,
#steps = 1,
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
)
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,
})
del shared_opt_batched_defaults
del linear_quartic_batched_defaults
del linear_oneshift_batched_defaults
del linear_adaptive_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
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
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
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
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
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
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
def generate_batched_opt_calculations(
method ,
cost ,
cycles ,
var_cycles ,
var_samples,
init_cycles,
init_samples,
init_minwalkers,
loc = 'generate_opt_calculations',
**opt_inputs
):
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
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 '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_batched_opt_calculations
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, # no checkpointing support yet
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
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:
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, # not supported yet
)
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,
#checkpoint = eq_checkpoint, # not supported yet
)
)
#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, # not supported yet
)
)
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,
)
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
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
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',
interactions = 'all',
corrections = 'default',
observables = None,
estimators = None,
traces = None,
calculations = None,
det_format = 'new',
J1 = False,
J2 = False,
J3 = False,
J1_size = None,
J1_rcut = None,
J1_dr = 0.5,
J2_size = None,
J2_rcut = None,
J2_dr = 0.5,
J2_init = 'zero',
J3_isize = 3,
J3_esize = 3,
J3_rcut = 5.0,
J1_rcut_open = 5.0,
J2_rcut_open = 10.0,
J1k = False,
J1k_kcut = 5.0,
J1k_symm = 'crystal',
J2k = False,
J2k_kcut = 5.0,
J2k_symm = 'crystal',
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
)
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 ,
J2_size = kw.J2_size ,
J2_rcut = kw.J2_rcut ,
J2_dr = kw.J2_dr ,
J2_init = kw.J2_init ,
J3_isize = kw.J3_isize ,
J3_esize = kw.J3_esize ,
J3_rcut = kw.J3_rcut ,
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 ,
J2k = kw.J2k ,
J2k_kcut = kw.J2k_kcut ,
J2k_symm = kw.J2k_symm ,
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
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
if spobuilders is not None:
wfn.sposet_builders = make_collection(spobuilders)
#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
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,
)
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
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 = array([[1,1,0],[1,0,1],[0,1,1]]),
reciprocal = 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 = 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 = 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 = 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 = 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 = 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
Reference in New Issue View Git Blame Copy Permalink