mirror of https://github.com/QMCPACK/qmcpack.git
fixed twist averaged energy density
git-svn-id: https://subversion.assembla.com/svn/qmcdev/trunk@6266 e5b18d87-469d-4833-9cc0-8cdfa06e9491
This commit is contained in:
parent
69a8b77541
commit
43935b5fc5
|
@ -12,6 +12,8 @@ try:
|
|||
except ImportError:
|
||||
h5py = unavailable('h5py')
|
||||
#end try
|
||||
from debug import *
|
||||
|
||||
|
||||
|
||||
class HDFgroup(DevBase):
|
||||
|
@ -84,15 +86,17 @@ class HDFgroup(DevBase):
|
|||
#end def __init__
|
||||
|
||||
|
||||
def _remove_hidden(self):
|
||||
def _remove_hidden(self,deep=True):
|
||||
if '_parent' in self:
|
||||
del self._parent
|
||||
#end if
|
||||
for name,value in self.iteritems():
|
||||
if isinstance(value,HDFgroup):
|
||||
value._remove_hidden()
|
||||
#end if
|
||||
#end for
|
||||
if deep:
|
||||
for name,value in self.iteritems():
|
||||
if isinstance(value,HDFgroup):
|
||||
value._remove_hidden()
|
||||
#end if
|
||||
#end for
|
||||
#end if
|
||||
for name in list(self.keys()):
|
||||
if name[0]=='_':
|
||||
del self[name]
|
||||
|
@ -109,9 +113,13 @@ class HDFgroup(DevBase):
|
|||
self[name][:] = 0
|
||||
#end if
|
||||
#end for
|
||||
for value in self._groups:
|
||||
for name,value in self._groups.iteritems():
|
||||
if isinstance(value,str):
|
||||
ci(ls(),gs())
|
||||
#end if
|
||||
value.zero(*names)
|
||||
#end for
|
||||
#self.sum(*names)
|
||||
#end def zero
|
||||
|
||||
|
||||
|
@ -142,6 +150,7 @@ class HDFgroup(DevBase):
|
|||
self.error(name+' not found in minsize partner')
|
||||
#end if
|
||||
#end for
|
||||
#self.sum(*names)
|
||||
#end def minsize
|
||||
|
||||
|
||||
|
@ -180,6 +189,7 @@ class HDFgroup(DevBase):
|
|||
self.error(name+' not found in accumulate partner')
|
||||
#end if
|
||||
#end for
|
||||
#self.sum(*names)
|
||||
#end def accumulate
|
||||
|
||||
|
||||
|
@ -189,11 +199,21 @@ class HDFgroup(DevBase):
|
|||
self[name] /= normalization
|
||||
#end if
|
||||
#end for
|
||||
for value in self._groups:
|
||||
for name,value in self._groups.iteritems():
|
||||
value.normalize(normalization,*names)
|
||||
#end for
|
||||
#self.sum(*names)
|
||||
#end def normalize
|
||||
|
||||
|
||||
def sum(self,*names):
|
||||
for name in names:
|
||||
if name in self and isinstance(self[name],ndarray) and name=='value':
|
||||
s = self[name].mean(0).sum()
|
||||
print ' sum = {0}'.format(s)
|
||||
#end if
|
||||
#end for
|
||||
#end def sum
|
||||
|
||||
#end class HDFgroup
|
||||
|
||||
|
|
|
@ -9,3 +9,5 @@ except (ImportError,RuntimeError):
|
|||
figure,plot,xlabel,ylabel,title,show,ylim,legend,xlim,rcParams,savefig,bar,xticks,subplot,grid,setp,errorbar,loglog,semilogx,semilogy,text = unavailable('matplotlib.pyplot','figure','plot','xlabel','ylabel','title','show','ylim','legend','xlim','rcParams','savefig','bar','xticks','subplot','grid','setp','errorbar','loglog','semilogx','semilogy','text')
|
||||
#end try
|
||||
|
||||
|
||||
# savefig(savefile,format='png',bbox_inches ='tight',pad_inches=.3)
|
||||
|
|
|
@ -23,14 +23,15 @@ class Pwscf(Simulation):
|
|||
|
||||
def __init__(self,**sim_args):
|
||||
has_group_atoms = 'group_atoms' in sim_args
|
||||
group_atoms = True
|
||||
group_atoms = False
|
||||
if has_group_atoms:
|
||||
group_atoms = sim_args['group_atoms']
|
||||
del sim_args['group_atoms']
|
||||
#end if
|
||||
Simulation.__init__(self,**sim_args)
|
||||
if group_atoms and isinstance(self.system,PhysicalSystem):
|
||||
self.system.structure.group_atoms()
|
||||
self.warn('requested grouping by atomic species, but pwscf does not group atoms anymore!')
|
||||
#self.system.structure.group_atoms()
|
||||
#end if
|
||||
#end def post_init
|
||||
|
||||
|
|
|
@ -957,8 +957,8 @@ class PwscfInput(SimulationInput):
|
|||
self.k_points.nkpoints = nkpoints
|
||||
self.k_points.kpoints = kpoints
|
||||
self.k_points.weights = s.kweights.copy()
|
||||
self.k_points.change_specifier('crystal',self) #added to make debugging easier
|
||||
#end if
|
||||
self.k_points.change_specifier('crystal',self) #added to make debugging easier
|
||||
|
||||
atoms = p.get_ions()
|
||||
masses = obj()
|
||||
|
@ -1126,7 +1126,7 @@ def generate_scf_input(prefix = 'pwscf',
|
|||
pseudos = None,
|
||||
system = None,
|
||||
use_folded = True,
|
||||
group_atoms = True):
|
||||
group_atoms = False):
|
||||
if pseudos is None:
|
||||
pseudos = []
|
||||
#end if
|
||||
|
@ -1255,8 +1255,11 @@ def generate_scf_input(prefix = 'pwscf',
|
|||
if system is not None:
|
||||
structure = system.structure
|
||||
if group_atoms:
|
||||
structure.group_atoms()
|
||||
self.warn('requested grouping by atomic species, but pwscf does not group atoms anymore!')
|
||||
#end if
|
||||
#if group_atoms: # disabled, hopefully not needed for qmcpack
|
||||
# structure.group_atoms()
|
||||
##end if
|
||||
if structure.at_Gpoint():
|
||||
pw.k_points.clear()
|
||||
pw.k_points.set(
|
||||
|
@ -1331,7 +1334,7 @@ def generate_relax_input(prefix = 'pwscf',
|
|||
pseudos = None,
|
||||
system = None,
|
||||
use_folded = False,
|
||||
group_atoms = True):
|
||||
group_atoms = False):
|
||||
if pseudos is None:
|
||||
pseudos = []
|
||||
#end if
|
||||
|
@ -1419,8 +1422,11 @@ def generate_relax_input(prefix = 'pwscf',
|
|||
if system is not None:
|
||||
structure = system.structure
|
||||
if group_atoms:
|
||||
structure.group_atoms()
|
||||
self.warn('requested grouping by atomic species, but pwscf does not group atoms anymore!')
|
||||
#end if
|
||||
#if group_atoms: #don't group atoms for pwscf, any downstream consequences?
|
||||
# structure.group_atoms()
|
||||
##end if
|
||||
if structure.at_Gpoint():
|
||||
pw.k_points.clear()
|
||||
pw.k_points.set(
|
||||
|
|
|
@ -3,6 +3,7 @@ from numpy import minimum,resize
|
|||
from generic import obj
|
||||
from hdfreader import HDFgroup
|
||||
from qaobject import QAobject
|
||||
from debug import *
|
||||
|
||||
|
||||
class QAinformation(obj):
|
||||
|
@ -16,6 +17,7 @@ class QAdata(QAobject):
|
|||
for value in self:
|
||||
value[:] = 0
|
||||
#end for
|
||||
#self.sum()
|
||||
#end def zero
|
||||
|
||||
def minsize(self,other):
|
||||
|
@ -26,6 +28,7 @@ class QAdata(QAobject):
|
|||
self.error(name+' not found in minsize partner')
|
||||
#end if
|
||||
#end for
|
||||
#self.sum()
|
||||
#end def minsize
|
||||
|
||||
def accumulate(self,other):
|
||||
|
@ -36,20 +39,31 @@ class QAdata(QAobject):
|
|||
self.error(name+' not found in accumulate partner')
|
||||
#end if
|
||||
#end for
|
||||
#self.sum()
|
||||
#end def accumulate
|
||||
|
||||
def normalize(self,normalization):
|
||||
for value in self:
|
||||
value/=normalization
|
||||
#end for
|
||||
#self.sum()
|
||||
#end def normalize
|
||||
|
||||
|
||||
def sum(self):
|
||||
s = 0
|
||||
for value in self:
|
||||
s+=value.sum()
|
||||
#end for
|
||||
print ' sum = {0}'.format(s)
|
||||
#end def sum
|
||||
#end class QAdata
|
||||
|
||||
|
||||
|
||||
class QAHDFdata(QAdata):
|
||||
def zero(self):
|
||||
for value in self:
|
||||
for name,value in self.iteritems():
|
||||
if isinstance(value,HDFgroup):
|
||||
value.zero('value','value_squared')
|
||||
#end if
|
||||
|
|
|
@ -2573,6 +2573,7 @@ class QmcpackInput(SimulationInput,Names):
|
|||
|
||||
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()
|
||||
|
@ -2664,7 +2665,7 @@ class QmcpackInput(SimulationInput,Names):
|
|||
)
|
||||
particlesets.append(eps)
|
||||
if len(ions)>0:
|
||||
eps.randomsrc = 'ion0'
|
||||
#eps.randomsrc = 'ion0' # don't do randomsrc by default
|
||||
ips = particleset(
|
||||
name='ion0',
|
||||
)
|
||||
|
@ -2717,7 +2718,12 @@ class QmcpackInput(SimulationInput,Names):
|
|||
|
||||
if abs(net_spin) > 1e-1:
|
||||
if ddet!=None:
|
||||
ddet.occupation.spindataset = 1 #jtk mark check
|
||||
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
|
||||
|
@ -3175,7 +3181,8 @@ def generate_particlesets(electrons = 'e',
|
|||
ions = 'ion0',
|
||||
up = 'u',
|
||||
down = 'd',
|
||||
system = None
|
||||
system = None,
|
||||
randomsrc = False
|
||||
):
|
||||
if system is None:
|
||||
QmcpackInput.class_error('generate_particlesets argument system must not be None')
|
||||
|
@ -3219,7 +3226,9 @@ def generate_particlesets(electrons = 'e',
|
|||
)
|
||||
particlesets.append(eps)
|
||||
if len(ions)>0:
|
||||
eps.randomsrc = iname
|
||||
if randomsrc:
|
||||
eps.randomsrc = iname
|
||||
#end if
|
||||
ips = particleset(name=iname)
|
||||
groups = []
|
||||
for ion in ions:
|
||||
|
@ -4203,6 +4212,7 @@ def generate_basic_input(id = 'qmc',
|
|||
buffer = None,
|
||||
lr_dim_cutoff = 15,
|
||||
remove_cell = False,
|
||||
randomsrc = False,
|
||||
meshfactor = 1.0,
|
||||
precision = 'float',
|
||||
twistnum = None,
|
||||
|
@ -4258,7 +4268,8 @@ def generate_basic_input(id = 'qmc',
|
|||
|
||||
if system!=None:
|
||||
particlesets = generate_particlesets(
|
||||
system = system
|
||||
system = system,
|
||||
randomsrc = randomsrc or tuple(bconds)!=('p','p','p')
|
||||
)
|
||||
#end if
|
||||
|
||||
|
|
|
@ -346,9 +346,11 @@ class EnergyDensityAnalyzer(HDFAnalyzer):
|
|||
self.error('attempted load without data')
|
||||
#end if
|
||||
name = self.info.name
|
||||
self.data = QAdata()
|
||||
self.data = QAHDFdata()
|
||||
if name in data:
|
||||
self.data = data[name]
|
||||
hdfg = data[name]
|
||||
hdfg._remove_hidden(deep=False)
|
||||
self.data.transfer_from(hdfg)
|
||||
del data[name]
|
||||
else:
|
||||
self.info.should_remove = True
|
||||
|
@ -364,7 +366,7 @@ class EnergyDensityAnalyzer(HDFAnalyzer):
|
|||
data = self.data
|
||||
|
||||
#why is this called 3 times?
|
||||
print nbe
|
||||
#print nbe
|
||||
|
||||
#transfer hdf data
|
||||
sg_pattern = re.compile(r'spacegrid\d*')
|
||||
|
@ -449,7 +451,7 @@ class EnergyDensityAnalyzer(HDFAnalyzer):
|
|||
input = self.run_info.input
|
||||
xml = self.run_info.ordered_input
|
||||
ps = input.get('particlesets')
|
||||
if 'ion0' in ps and len(ps.ion0.groups)>1:
|
||||
if 'ion0' in ps and len(ps.ion0.groups)>1 and 'size' in ps.ion0:
|
||||
qsx = xml.simulation.qmcsystem
|
||||
if len(ps)==1:
|
||||
psx = qsx.particleset
|
||||
|
@ -488,11 +490,6 @@ class EnergyDensityAnalyzer(HDFAnalyzer):
|
|||
|
||||
#reorder the atomic data
|
||||
for sg in self.spacegrids:
|
||||
|
||||
#print 'qqa'
|
||||
#import code
|
||||
#code.interact(local=locals())
|
||||
|
||||
sg.reorder_atomic_data(imap)
|
||||
#end for
|
||||
#end if
|
||||
|
@ -1095,6 +1092,9 @@ class TracesAnalyzer(QAanalyzer):
|
|||
|
||||
|
||||
class DensityMatricesAnalyzer(HDFAnalyzer):
|
||||
|
||||
allowed_settings = ['save_data','jackknife','diagonal','occ_tol','coup_tol','stat_tol']
|
||||
|
||||
def __init__(self,name,nindent=0):
|
||||
HDFAnalyzer.__init__(self)
|
||||
self.info.name = name
|
||||
|
@ -1134,19 +1134,18 @@ class DensityMatricesAnalyzer(HDFAnalyzer):
|
|||
|
||||
|
||||
def analyze_local(self):
|
||||
|
||||
save_data = False
|
||||
jackknife = False
|
||||
diagonal = True
|
||||
|
||||
# 1) exclude states that do not contribute to the number trace
|
||||
# 2) exclude elements that are not statistically significant (1 sigma?)
|
||||
# 3) use remaining states to form filtered number and energy matrices
|
||||
# 4) perform jackknife sampling to get eigenvalue error bars
|
||||
# 5) consider using cross-correlations w/ excluded elements to reduce variance
|
||||
occ_tol = 1e-3
|
||||
coup_tol = 1e-4
|
||||
stat_tol = 2.0
|
||||
save_data = False
|
||||
jackknife = False
|
||||
diagonal = True
|
||||
occ_tol = 1e-3
|
||||
coup_tol = 1e-4
|
||||
stat_tol = 2.0
|
||||
|
||||
|
||||
nbe = QAanalyzer.method_info.nblocks_exclude
|
||||
self.info.nblocks_exclude = nbe
|
||||
|
|
|
@ -310,8 +310,10 @@ class Structure(Sobj):
|
|||
def group_atoms(self,folded=True):
|
||||
if len(self.elem)>0:
|
||||
order = self.elem.argsort()
|
||||
self.elem = self.elem[order]
|
||||
self.pos = self.pos[order]
|
||||
if (self.elem!=self.elem[order]).any():
|
||||
self.elem = self.elem[order]
|
||||
self.pos = self.pos[order]
|
||||
#end if
|
||||
#end if
|
||||
if self.folded_structure!=None and folded:
|
||||
self.folded_structure.group_atoms(folded)
|
||||
|
@ -2375,6 +2377,16 @@ class Crystal(Structure):
|
|||
basis = [[0,0,0],[.5,0,0]],
|
||||
basis_vectors = 'conventional'
|
||||
),
|
||||
('rocksalt','prim'):obj(
|
||||
lattice = 'cubic',
|
||||
cell = 'primitive',
|
||||
centering = 'F',
|
||||
constants = 5.64,
|
||||
units = 'A',
|
||||
atoms = ('Na','Cl'),
|
||||
basis = [[0,0,0],[.5,0,0]],
|
||||
basis_vectors = 'conventional'
|
||||
),
|
||||
('copper','prim'):obj(
|
||||
lattice = 'cubic',
|
||||
cell = 'primitive',
|
||||
|
|
Loading…
Reference in New Issue