Fixed scoring for multi-elements

GRSlib/GRS.py

@@ -149,14 +149,12 @@ def get_score(self,data):
         between file types (xyz=lammps-data, ase.Atoms, etc)
         """
         #Pass data to, and do something with the functs of scoring
-        #if self.config.sections['TARGET'].target_fname is None:
-            #print("Provided target descriptors superceed target data file")
         try:
             self.descriptors['target'] = np.load(self.config.sections['TARGET'].target_fdesc)    
+            print("Provided target descriptors superceed target data file")
         except:
             self.descriptors['target'] = self.convert_to_desc(self.config.sections['TARGET'].target_fname)
         self.descriptors['current'] = self.convert_to_desc(data)
-
         try: # If prior is empty when getting a score, set to starting structure. 
             if self.descriptors.get('prior',None)==None: 
                 self.set_prior([self.config.sections['TARGET'].start_fname])
@@ -243,7 +241,7 @@ def genetic_move(self,data):
         self.score = Scoring(self.pt, self.config, self.loss_func, self.descriptors)  # Set scoring class to assign scores to moves
         self.genmove = Optimize(self.pt, self.config, self.score, self.convert) #Set desired motion class with scoring attached
 
-        gen_scores = self.genmove.unique_selection(data)
+        gen_scores = self.genmove.advance_generations(data)
         #self.write_output()
         best_candidate = sorted(gen_scores,key=lambda x: x[3])[0][2]
         return gen_scores, best_candidate

GRSlib/basis/construct.py

@@ -24,7 +24,7 @@ def make_ACE_functions(elements,ranks,lmax,lmin,nmax):
 
         # radial function hyperparameters are defined per bond type:
         bonds = [p for p in itertools.product(elements,elements)]
-        print('bonds',bonds)
+        #print('bonds',bonds)
         rc_range,rc_default,lmb_default,rcin_default = get_default_settings(elements,nshell=2,return_range=True,apply_shift=False)
         rcutfac = [float(k) for k in rc_default.split()[2:]]
         lmbda = [float(k) for k in lmb_default.split()[2:]]
@@ -66,6 +66,6 @@ def make_ACE_functions(elements,ranks,lmax,lmin,nmax):
                 store_generalized(ccs, coupling_type='wig',L_R=L_R)
         Apot = AcePot(elements,reference_ens,ranks,nmax,lmax,nradbase,rcutfac,lmbda,rcinner,drcinner,lmin=lmin, **{'ccs':ccs[M_R]})
         Apot.write_pot('coupling_coefficients')
-        print(Apot.nus)
+        #print(Apot.nus)
     except ModuleNotFoundError:
         raise ModuleNotFoundError("Need to have fitsnap3 module in your python path to automatically generate ACE bases")

GRSlib/converters/convert.py

@@ -35,7 +35,7 @@ def lammps_to_ase(self,data):
             ase_data = read(data,format='lammps-data',Z_of_type=Z_of_type)
         except:
             ase_data = read(data+".lammps-data",format='lammps-data',Z_of_type=Z_of_type)
-        print('in GRSlib/converters/convert.py lammps_to_ase ',ase_data,types)
+        #print('in GRSlib/converters/convert.py lammps_to_ase ',ase_data,types)
         return ase_data
 
     def lammps_ace(self,data):

GRSlib/converters/sections/lammps_ace.py

@@ -40,7 +40,10 @@ def run_lammps_single(self,data):
 
     def _collect_lammps_single(self):
         num_atoms = self._lmp.extract_global("natoms")
-        num_types = self.config.sections['BASIS'].numtypes
+        try:
+            types = self._lmp.numpy.extract_atom(name="type", nelem=num_atoms).ravel()
+        except:
+            types = self._lmp.numpy.extract_atom_iarray(name="type", nelem=num_atoms).ravel()         
         with open('coupling_coefficients.yace','r') as readcoeff:
             lines = readcoeff.readlines()
             elemline = [line for line in lines if 'elements' in line][0]
@@ -51,18 +54,27 @@ def _collect_lammps_single(self):
             elems = elemstr4.split()
             nelements = len(elems)
             desclines = [line for line in lines if 'mu0' in line]
+
+#        if nelements != int(len(set(types))):
+#            print("Target structure doesnt have the full number of desired element types")
 
-        ncols_pace = int(len(desclines)/nelements) +1
-        #ncols_pace = int(len(desclines)/nelements) + nelements #with bzero?
+#        ncols_pace = int(len(desclines)/nelements) +1
+        ncols_pace = int(len(desclines)) + 1 #need to bring in full columns and reassign based on element type
         nrows_pace = num_atoms
         lmp_pace = _extract_compute_np(self._lmp, "pace", 0, 2, (nrows_pace, ncols_pace))
-
         if (np.isinf(lmp_pace)).any() or (np.isnan(lmp_pace)).any():
             self.pt.single_print('WARNING! Applying np.nan_to_num()')
             lmp_pace = np.nan_to_num(lmp_pace)
         if (np.isinf(lmp_pace)).any() or (np.isnan(lmp_pace)).any():
             raise ValueError('Nan in computed data of file')
-#        print("Got descriptors, returning")
-#        print(np.shape(lmp_pace))
-        return lmp_pace
+        ndesc = int(len(desclines)/nelements)
+        if nelements > 1:                
+            for atoms in range(num_atoms):
+                try:
+                    perelem_desc= np.r_[[lmp_pace[atoms,(types[atoms]-1)*ndesc:(types[atoms])*ndesc]], perelem_desc]
+                except:
+                    perelem_desc = [lmp_pace[atoms,(types[atoms]-1)*ndesc:(types[atoms])*ndesc]]
+            return perelem_desc
+        else:
+            return lmp_pace
 

GRSlib/converters/sections/lammps_base.py

@@ -114,6 +114,7 @@ def _extract_compute_np(lmp, name, compute_style, result_type, array_shape):
 
     if array_shape is None:
         array_np = lmp.numpy.extract_compute(name,compute_style, result_type)
+        return array_np
     else:
         ptr = lmp.extract_compute(name, compute_style, result_type)
         if result_type == 0:
@@ -130,7 +131,7 @@ def _extract_compute_np(lmp, name, compute_style, result_type, array_shape):
         array_np.shape = array_shape
         reshaped_array_np = np.delete(array_np, np.s_[-1:], axis=1)
         #LAMMPS is returning an extra descriptor at the end of the list, deleting by reshaping
-    return reshaped_array_np
+        return reshaped_array_np
 
 def _extract_commands(string):
     return [x for x in string.splitlines() if x.strip() != '']

GRSlib/io/sections/basis.py

@@ -34,8 +34,6 @@ def __init__(self, name, config, pt, infile, args):
             self.rcinner = self.get_value("BASIS","rcinner",'0.0').split()
             self.drcinner = self.get_value("BASIS","drcinner",'0.01').split()
             self.elements = self.get_value("BASIS", "elements", "H").split()
-            print('elements',self.elements)
-            print('rcinner',self.rcinner)
             self.mumax = len(self.elements)
             self.numtypes = len(self.elements)
             #self.erefs = self.get_value("ACE", "erefs", "0.0").split() 
@@ -115,10 +113,10 @@ def _generate_b_list(self):
             nus.sort(key = lambda x : mu0s[nus_unsort.index(x)],reverse = False)
             byattyp = srt_by_attyp(nus)
             #config.nus = [item for sublist in list(byattyp.values()) for item in sublist]
-            print('byatttyp',byattyp)
-            print('num type',self.numtypes)
+#            print('byatttyp',byattyp)
+#            print('num type',self.numtypes)
             for atype in range(self.numtypes):
-                print('atype',atype)
+#                print('atype',atype)
                 nus = byattyp[str(atype)]
                 for nu in nus:
                     i += 1
@@ -146,7 +144,7 @@ def decorated_write_couple():
                 bondinds=range(len(self.elements))
                 bonds = [b for b in itertools.product(bondinds,bondinds)]
                 bondstrs = ['[%d, %d]' % b for b in bonds]
-                print('bond strs basis.py',bondstrs)
+#                print('bond strs basis.py',bondstrs)
                 assert len(self.lmbda) == len(bondstrs), "must provide rc, lambda, for each BOND type" 
                 assert len(self.rcutfac) == len(bondstrs), "must provide rc, lambda, for each BOND type" 
                 if len(self.lmbda) == 1:
@@ -186,7 +184,7 @@ def decorated_write_couple():
                         #store them for later so they don't need to be recalculated
                         store_generalized(ccs, coupling_type='wig',L_R=L_R)
 
-                print('rcinner vals basis.py', rcinnervals)
+#                print('rcinner vals basis.py', rcinnervals)
 
                 apot = AcePot(self.elements, reference_ens, [int(k) for k in self.ranks], [int(k) for k in self.nmax],  [int(k) for k in self.lmax], self.nmaxbase, rcvals, lmbdavals, rcinnervals, drcinnervals, [int(k) for k in self.lmin], self.b_basis, **{'ccs':ccs[M_R]})
                 apot.write_pot('coupling_coefficients')

GRSlib/io/sections/genetic.py

@@ -8,19 +8,23 @@ def __init__(self, name, config, pt,infile, args):
         super().__init__(name, config, pt, infile,args)
         self.allowedkeys = ['mutation_rate', 'mutation_types', 'population_size', 'ngenerations', 'start_type',
                             'max_atoms', 'min_atoms', 'max_length_aspect', 'max_angle_aspect', 'density_ratio',
-                            'composition_constraint', 'composition', 'start_type', 'lattice_type', 'structure_template', 'frac_pop_per_ref', 'dev_per_ref','reference_phases', 'change_ele_tol']
+                            'composition_constraint', 'composition', 'start_type', 'lattice_type', 'structure_template', 
+                            'frac_pop_per_ref', 'dev_per_ref','reference_phases', 'change_ele_tol', 'tournament_size',
+                            'replacements', 'multi_mutate','randomseed']
         self._check_section()
         self.mutation_rate = self.get_value("GENETIC", "mutation_rate", 0.5, interpreter="float")
-        self.frac_pop_per_ref = [float(k) for k in self.get_value("GENETIC", "frac_pop_per_ref", "1.0").split()]
-        self.dev_per_ref =[float(kk) for kk in self.get_value("GENETIC", "dev_per_ref", "0.2").split()]
-        self.reference_phases = self.get_value("GENETIC", "reference_phases", "bcc").split()
+        #self.reference_phases = self.get_value("GENETIC", "reference_phases", "bcc").split()
         self.mutation_types = self.get_value("GENETIC", "mutation_types", {"perturb": 0.5, "change_ele": 0.0, "atom_count" : 0.1, "volume" : 0.2, "minimize" : 0.2}, interpreter="dict")
         self.population_size = self.get_value("GENETIC", "population_size", 20, interpreter="int")
         self.ngenerations = self.get_value("GENETIC", "ngenerations", 10, interpreter="int")
+        self.tournament_size = self.get_value("GENETIC", "tournament_size", 10, interpreter="int")
+        self.replacements = self.get_value("GENETIC", "replacements", False)
+        self.multi_mutate = self.get_value("GENETIC", "multi_mutate", True)
+        self.randomseed = self.get_value("GENETIC", "randomseed", 12345, interpreter="int")
         self.max_atoms = self.get_value("GENETIC", "max_atoms", 100, interpreter="int")
         self.change_ele_tol = self.get_value("GENETIC", "change_ele_tol", 0.1, interpreter="float")
         self.min_atoms = self.get_value("GENETIC", "min_atoms", 10, interpreter="int")
-        self.reference_phases = self.get_value("GENETIC","reference_phases", "bcc").split()
+        #self.reference_phases = self.get_value("GENETIC","reference_phases", "bcc").split()
         self.max_length_aspect = self.get_value("GENETIC", "max_length_aspec", 3.0, interpreter="float")
         self.max_angle_aspect = self.get_value("GENETIC", "max_angle_aspec", 3.0, interpreter="float")
         self.density_ratio = self.get_value("GENETIC", "density_ratio", 1.3, interpreter="float") #This will allow for 30% changes in either direction of density
@@ -35,7 +39,9 @@ def __init__(self, name, config, pt,infile, args):
             self.lattice_type = self.get_value("GENETIC", "lattice_type", None)
         elif self.start_type == "template":
             self.structure_template = self.get_value("GENETIC", "structure_template", None)
-#        elif self.start_type == "phases":
- #           self.reference_phases = self.get_value("GENETIC", "reference_phases", None, interpreter="list")
+        elif self.start_type == "phases":
+            self.reference_phases = self.get_value("GENETIC","reference_phases", "bcc").split()
+            self.frac_pop_per_ref = [float(k) for k in self.get_value("GENETIC", "frac_pop_per_ref", "1.0").split()]
+            self.dev_per_ref =[float(kk) for kk in self.get_value("GENETIC", "dev_per_ref", "0.2").split()]
 
         self.delete()

GRSlib/io/sections/scoring.py

@@ -24,6 +24,9 @@ def __init__(self, name, config, pt,infile, args):
             #Score reduction if exact moment value is matched
             self.smartmask = self.get_value("SCORING", "smartmask", 0, interpreter="int")
             #Use a subset of descriptors of given count
+        elif self.score_type == "wasserstein":
+            self.smartmask = self.get_value("SCORING", "smartmask", 0, interpreter="int")
+            #Use a subset of descriptors of given count
 
         elif self.score_type == "entropy":
             self.internal_entropy = self.get_value("SCORING", "internal_entropy", 1.0)

GRSlib/motion/create.py

@@ -94,9 +94,9 @@ def from_template(self,*args):
         random_shuffle = True
         #More of a super function that will call a bunch of the ones below
 #        print("Starting population using provided template")
-        population = []
         duplicate = self.convert.lammps_to_ase(args[0][0])
-        for candidate in range(self.config.sections["GENETIC"].population_size):
+        population = [duplicate]
+        for candidate in range(self.config.sections["GENETIC"].population_size-1):
             tmp_atoms = duplicate.copy()
             if random_shuffle:
                 duplicate_syms = [atom.symbol for atom in tmp_atoms]
@@ -113,8 +113,7 @@ def from_lattice(self,*args):
 
     def from_phases(self,*args): # similar to template but start from a known lattice
         verbose = False
-        #value to compress/expand cells by (see TODO below)
-        compex = 0.0
+        compex = 0.0 #value to compress/expand cells by (see TODO below)
         population = []
         references_to_try = self.config.sections["GENETIC"].reference_phases #= ['hpc','fcc','bcc']
         pop_size = self.config.sections["GENETIC"].population_size
@@ -174,7 +173,7 @@ def from_phases(self,*args): # similar to template but start from a known lattic
 
     def from_random(self,*args):
         #More of a super function that will call a bunch of the ones below
-        #print("Starting population of random low energy structures of provided elements")
+#        print("Starting population of random low energy structures of provided elements")
         population = []
         # From types, find cell
         num_ele = len(self.config.sections["BASIS"].elements)