nexus: improve handling of twist input

git-svn-id: https://subversion.assembla.com/svn/qmcdev/trunk@7058 e5b18d87-469d-4833-9cc0-8cdfa06e9491

nexus/library/qmcpack.py

@@ -60,19 +60,11 @@ class Qmcpack(Simulation):
 
 
     def post_init(self):
-        #jtk mark
-        #  may need to put this back
-        #  removed because particleset is not required by qmcpack
-        #   and thus the system cannot be determined without access to the h5file
-        #if self.system is None:
-        #    self.error('system must be specified to determine type of run')
-        ##end if
-
         generic_input = self.has_generic_input()
 
         if self.system is None:
             if not generic_input:
-                self.warn('system must be specified to determine whether to twist average\n  proceeding under the assumption of no twist averaging')
+                self.warn('system must be specified to determine whether to twist average\nproceeding under the assumption of no twist averaging')
             #end if
             self.should_twist_average = False
         else:
@@ -166,7 +158,7 @@ class Qmcpack(Simulation):
                     orb_elem.tilematrix = array(system.structure.tmatrix)
                 #end if
                 defs = obj(
-                    twistnum   = 0,
+                    #twistnum   = 0,
                     meshfactor = 1.0
                     )
                 for var,val in defs.iteritems():
@@ -174,6 +166,11 @@ class Qmcpack(Simulation):
                         orb_elem[var] = val
                     #end if
                 #end for
+                has_twist    = 'twist' in orb_elem
+                has_twistnum = 'twistnum' in orb_elem
+                if  not has_twist and not has_twistnum:
+                    orb_elem.twistnum = 0
+                #end if
 
                 system = self.system
                 structure = system.structure
@@ -189,7 +186,7 @@ class Qmcpack(Simulation):
                 twistnums = range(len(structure.kpoints))
                 if self.should_twist_average:
                     self.twist_average(twistnums)
-                elif orb_elem.twistnum is None:
+                elif not has_twist and orb_elem.twistnum is None:
                     orb_elem.twistnum = twistnums[0]
                 #end if
 

nexus/library/qmcpack_input.py

@@ -4163,7 +4163,7 @@ def generate_bspline_builder(type           = 'bspline',
         bsb.buffer = buffer
     #end if
     if twist!=None:
-        bsb.twist = tuple(twist)
+        bsb.twistnum = system.structure.select_twist(twist)
     elif twistnum!=None:
         bsb.twistnum = twistnum
     elif len(system.structure.kpoints)==1:
@@ -4287,11 +4287,13 @@ def generate_determinantset_old(type           = 'bspline',
             )
         )
     if twist!=None:
-        dset.twist = tuple(twist)
+        dset.twistnum = system.structure.select_twist(twist)
     elif twistnum!=None:
         dset.twistnum = twistnum
-    #else:
-    #    dset.twistnum = 0
+    elif len(system.structure.kpoints)==1:
+        dset.twistnum = 0
+    else:
+        dset.twistnum = None
     #end if
     return dset
 #end def generate_determinantset_old

nexus/library/structure.py

@@ -2870,6 +2870,38 @@ class Structure(Sobj):
     #end def kmap
 
 
+    def select_twist(self,selector='smallest',tol=1e-6):
+        index = None
+        invalid_selector = False
+        if isinstance(selector,str):
+            if selector=='smallest':
+                index = (self.kpoints**2).sum(1).argmin()
+            else:
+                invalid_selector = True
+            #end if
+        elif isinstance(selector,(tuple,list,ndarray)):
+            ku_sel = array(selector,dtype=float)
+            n = 0
+            for ku in self.kpoints_unit():
+                if norm(ku-ku_sel)<tol:
+                    index = n
+                    break
+                #end if
+                n+=1
+            #end for
+            if index is None:
+                self.error('cannot identify twist number\ntwist requested: {0}\ntwists present: {1}'.format(ku_sel,sorted([tuple(k) for k in self.kpoints_unit()])))
+            #end if
+        else:
+            invalid_selector = True
+        #end if
+        if invalid_selector:
+            self.error('cannot identify twist number\ninvalid selector provided: {0}\nvalid string inputs for selector: smallest\nselector can also be a length 3 tuple, list or array (a twist vector)'.format(selector))
+        #end if
+        return index
+    #end def select_twist
+
+
     def pos_unit(self,pos=None):
         if pos is None:
             pos = self.pos