- Basic trexio support in Fragment, read and write. Molecular info on…

…ly for now

- nucchargestoelems function

ash/functions/functions_elstructure.py

@@ -2333,6 +2333,7 @@ def read_trexio_file(filename="trexio", back_end_type="text"):
 
     return charges_r, coords_r, labels_r
 
+#Very basic trexio file writer, only molecule info for the moment
 def write_trexio_file(fragment, filename="trexio", back_end_type="text"):
     try:
         import trexio
@@ -2364,6 +2365,10 @@ def write_trexio_file(fragment, filename="trexio", back_end_type="text"):
 
     open_file.close()
 
+    print("Wrote trexio file:", filename)
+
+
+
 def yoshimine_sort(a,b,c,d):
     if a > b:
         ab = a*(a+1)/2 + b

ash/modules/module_coords.py

@@ -83,7 +83,7 @@ def calculate_reaction_energy(self):
 # ASH Fragment class
 class Fragment:
     def __init__(self, coordsstring=None, fragfile=None, databasefile=None, xyzfile=None, pdbfile=None, grofile=None,
-                 amber_inpcrdfile=None, amber_prmtopfile=None, smiles=None,
+                 amber_inpcrdfile=None, amber_prmtopfile=None, trexiofile=None, smiles=None,
                  chemshellfile=None, coords=None, elems=None, connectivity=None, atom=None, diatomic=None, diatomic_bondlength=None,
                  bondlength=None,
                  atomcharges=None, atomtypes=None, conncalc=False, scale=None, tol=None, printlevel=2, charge=None,
@@ -215,6 +215,26 @@ def __init__(self, coordsstring=None, fragfile=None, databasefile=None, xyzfile=
         elif fragfile is not None:
             self.label = fragfile.split('/')[-1].split('.')[0]
             self.read_fragment_from_file(fragfile)
+        # Trexio
+        elif trexiofile is not None:
+            from ash.functions.functions_elstructure import read_trexio_file
+            print("Reading TREXIO file:", trexiofile)
+            if 'h5' in trexiofile or 'hdf5' in trexiofile:
+                print("Assuming TREXIO HDF5 format based on file suffix")
+                charges, coords, labels = read_trexio_file(filename=trexiofile, back_end_type="hdf5")
+            elif '.text' in trexiofile:
+                print("Assuming TREXIO TEXT format based on file suffix")
+                charges, coords, labels = read_trexio_file(filename=trexiofile, back_end_type="text")
+            else:
+                print("No recognized suffix in TREXIO file found. Exiting")
+                ashexit()
+
+            self.coords = np.array(coords)
+            elems = [reformat_element(el) for el in labels]
+            # Converting charges to elements
+            self.elems = nucchargestoelems(charges)
+            # NOTE: Labels not currently used
+
         # Reading an XYZ-file from the ASH database
         elif databasefile is not None:
             databasepath=ashpath+"/databases/fragments/"
@@ -231,7 +251,7 @@ def __init__(self, coordsstring=None, fragfile=None, databasefile=None, xyzfile=
                 exit()
         # If all else fails, exit
         else:
-            ashexit(errormessage="Fragment requires some kind of valid coordinates input!")
+            ashexit(errormessage="Fragment requires some kind of valid coordinate input!")
         # Label for fragment (string). Useful for distinguishing different fragments
         # This overrides label-definitions above (self.label=xyzfile etc)
         if label is not None:
@@ -243,7 +263,6 @@ def __init__(self, coordsstring=None, fragfile=None, databasefile=None, xyzfile=
         if mult != None:
             self.mult = mult
 
-
         # Now update attributes after defining coordinates, getting charge, mult
         self.update_attributes()
         if conncalc is True:
@@ -844,6 +863,10 @@ def write_XYZ_for_atoms(self,xyzfilename="Fragment-subset.xyz", atoms=None):
                 line = "{:4} {:>12.6f} {:>12.6f} {:>12.6f}".format(el, c[0], c[1], c[2])
                 ofile.write(line + '\n')
 
+    def write_trexio(self,filename=None,format="text"):
+        from ash.functions.functions_elstructure import write_trexio_file
+        write_trexio_file(self, filename=filename, back_end_type=format)
+
     #Function to get subset-coordinates with linkatoms
     def get_subset_coords_with_linkatoms(self,qmatoms):
         conn_scale = ash.settings_ash.settings_dict["scale"]
@@ -2318,6 +2341,14 @@ def elemstonuccharges(ellist):
         nuccharges.append(atcharge)
     return nuccharges
 
+def nucchargestoelems(chargelist):
+    elems = []
+    for c in chargelist:
+        for key, value in elematomnumbers.items():
+            if value == c:
+                elems.append(key.upper())
+    return elems
+
 
 # Calculate molecular mass from list of atoms
 def totmasslist(ellist):