Patch 0.17.0 (#36)

* Adding construction of .abo files for managlyph

* Removing dependency

* Removing MacOS support

.github/workflows/build_conda.yml

@@ -150,61 +150,3 @@ jobs:
         run: |
           export ANACONDA_API_TOKEN=$INPUT_ANACONDATOKEN
           anaconda upload packages/*.tar.bz2
-
-#-------------------------------------------------------------------------------
-# Anaconda / MacOSx
-#-------------------------------------------------------------------------------
-  build-anaconda-macos:
-    needs: check-version-strings
-    runs-on: macos-latest
-
-    steps:
-      - name: Checkout repo
-        uses: actions/checkout@v3
-      - name: Set-up miniconda
-        uses: conda-incubator/setup-miniconda@v2
-        with:
-          activate-environment: test
-          environment-file: environment.yml
-          python-version: 3.8
-          auto-activate-base: false
-      - name: Build
-        run: |
-          conda install conda-build
-          conda build . --no-include-recipe
-      - name: Archive packages
-        uses: actions/upload-artifact@v3
-        with:
-          name: anaconda-macos-packages
-          path: /usr/local/miniconda/conda-bld/osx-64/pyqint-*.tar.bz2
-
-  publish-anaconda-macos:
-    name: Publish Anaconda / MacOS
-    if: startsWith(github.ref, 'refs/tags/v')
-    needs: build-anaconda-macos
-    runs-on: ubuntu-latest
-    environment:
-      name: anaconda
-      url: https://anaconda.org/ifilot/pyqint
-    steps:
-      - name: Checkout repo
-        uses: actions/checkout@v3
-      - name: Set-up miniconda
-        uses: conda-incubator/setup-miniconda@v2
-        with:
-          activate-environment: test
-          environment-file: environment.yml
-          python-version: 3.8
-          auto-activate-base: false
-      - name: Retrieve packages
-        uses: actions/download-artifact@v3
-        with:
-          name: anaconda-macos-packages
-          path: packages
-      - name: publish-to-conda
-        shell: bash -l {0}
-        env:
-          INPUT_ANACONDATOKEN: ${{ secrets.ANACONDA_TOKEN }}
-        run: |
-          export ANACONDA_API_TOKEN=$INPUT_ANACONDATOKEN
-          anaconda upload packages/*.tar.bz2

.gitignore

@@ -149,3 +149,6 @@ test.py
 nose_debug
 anaconda-upload/*
 *.vscode
+
+# managlyph .abo files
+*.abo

COMPILATION.md

@@ -11,7 +11,7 @@ pipx run cibuildwheel --only cp310-manylinux_x86_64
 To install the `whl` file
 
 ```bash
-VERSION=0.16.2 pip3 install wheelhouse/pyqint-$VERSION-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl
+pip3 install wheelhouse/pyqint-0.17.0-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl
 ```
 
 and to locally test
@@ -20,6 +20,12 @@ and to locally test
 pytest-3 tests/*.py
 ```
 
+For skipping unit testing when compiling with `cibuildwheel`, run
+
+```bash
+CIBW_TEST_SKIP="cp310-manylinux_x86_64" pipx run cibuildwheel --only cp310-manylinux_x86_64
+```
+
 ### Uploading to PyPi
 
 This will place wheels in the `dist` folder. To upload these wheels

examples/abobuilder.py

@@ -0,0 +1,34 @@
+from pyqint import MoleculeBuilder, GeometryOptimization, AboBuilder, FosterBoys
+import numpy as np
+
+#
+# Plot the isosurfaces for the CO molecule
+#
+
+def main():
+    # calculate sto-3g coefficients for co
+    print('Calculating CH4')
+    res = optimize_ch4()
+
+    print('Building ABO')
+    abo = AboBuilder(res)
+    abo.build_abo('ch4_canonical.abo', isovalue=0.01, alpha=0.9)
+
+    fb = FosterBoys(res)
+    res = fb.run()
+
+    print('Building ABO')
+    abo = AboBuilder(res)
+    abo.build_abo('ch4_fb.abo', isovalue=0.01, alpha=0.9)
+
+def optimize_ch4():
+    """
+    Optimization function for scipy.optimize.minimize
+    """
+    mol = MoleculeBuilder().from_name('CH4')
+    res = GeometryOptimization().run(mol, 'sto3g')
+
+    return res['data']
+
+if __name__ == '__main__':
+    main()

meta.yaml

@@ -1,6 +1,6 @@
 package:
   name: "pyqint"
-  version: "0.16.2"
+  version: "0.17.0"
 
 source:
   path: .

pyqint/__init__.py

@@ -8,5 +8,6 @@
 from .molecule_builder import MoleculeBuilder
 from .geometry_optimization import GeometryOptimization
 from .blenderrender import BlenderRender
+from .abobuilder import AboBuilder
 
 from ._version import __version__

pyqint/_version.py

@@ -1,2 +1,2 @@
-__version__ = "0.16.2"
+__version__ = "0.17.0"
 

pyqint/abobuilder.py

@@ -0,0 +1,157 @@
+import os
+import numpy as np
+from pyqint import PyQInt
+
+# try to import PyTessel but do not throw an error if it cannot be loaded
+try:
+    from pytessel import PyTessel
+except ModuleNotFoundError:
+    print('Cannot find module PyTessel')
+
+try:
+    from mendeleev import element
+except ModuleNotFoundError:
+    print('Cannot find module mendeleev')
+
+class AboBuilder:
+    def __init__(self, res):
+        # copy objects from Hartree-Fock or FosterBoys result dictionary
+        self.orbc_canonical = res['orbc']
+        self.orbe_canonical = res['orbe']
+        self.mol = res['mol']
+        self.nelec = res['nelec']
+        self.cgfs = res['cgfs']
+
+    def build_abo(self, outfile, isovalue=0.01, alpha=0.8, sz=5, npts=100):
+        """
+        Generate an ABO file
+
+        arguments:
+            outfile: path to write .abo file to
+            isovalue: isovalue used for the generation of the isosurface
+            alpha: transparency setting for the isosurface
+            sz: size of the box; note that the box dimensions in one direction are twice the size
+            npts: number of sampling points inside the box
+        """
+        frame_contents = []
+        atoms = self.mol.get_atoms()
+
+        # specify colors for occupied and virtual orbitals
+        colors = [
+            np.array([0.592, 0.796, 0.369, alpha], dtype=np.float32),
+            np.array([0.831, 0.322, 0.604, alpha], dtype=np.float32),
+            np.array([1.000, 0.612, 0.000, alpha], dtype=np.float32),
+            np.array([0.400, 0.831, 0.706, alpha], dtype=np.float32)
+        ]
+
+        # build output file
+        f = open(outfile, 'wb')
+
+        # write number of frames (number of atomic orbitals + 1 for the base geometry)
+        f.write(int(len(self.cgfs) + 1).to_bytes(2, byteorder='little'))
+        print('Generating %i frames for ABO file' % int(len(self.cgfs) + 1))
+
+        #
+        # First write the bare geometry of the molecule
+        #
+
+        # write frame_idx
+        f.write(int(1).to_bytes(2, byteorder='little'))
+
+        descriptor = 'Geometry'
+
+        f.write(len(descriptor).to_bytes(2, byteorder='little'))
+        f.write(bytearray(descriptor, encoding='utf8'))
+
+        # write nr_atoms
+        f.write(len(self.mol).to_bytes(2, byteorder='little'))
+        for atom in atoms:
+            f.write(element(atom[0]).atomic_number.to_bytes(1, byteorder='little'))
+            f.write(np.array(atom[1] * 0.529177, dtype=np.float32).tobytes())
+
+        # write number of models (always zero for pure geometry)
+        f.write(int(0).to_bytes(1, byteorder='little'))
+        f.write(int(0).to_bytes(1, byteorder='little'))
+
+        #
+        # Build the molecular orbitals including the geometry
+        #
+        for i in range(0, len(self.cgfs)):
+            # write frame_idx
+            f.write((i+1).to_bytes(2, byteorder='little'))
+
+            descriptor = "MO %02i (E = %.1f Ht)" % (i+1, self.orbe_canonical[i])
+
+            f.write(len(descriptor).to_bytes(2, byteorder='little'))
+            f.write(bytearray(descriptor, encoding='utf8'))
+
+            # write nr_atoms
+            f.write(len(self.mol).to_bytes(2, byteorder='little'))
+            for atom in atoms:
+                f.write(element(atom[0]).atomic_number.to_bytes(1, byteorder='little'))
+                f.write(np.array(atom[1] * 0.529177, dtype=np.float32).tobytes())
+
+            print('Writing MO #%02i' % (i+1))
+
+            # write number of models (always two, positive and negative lobe)
+            f.write(int(2).to_bytes(2, byteorder='little'))
+
+            # generate isosurfaces
+            posorb, negorb = self.__build_isosurface(self.cgfs, self.orbc_canonical[:,i], 
+                                                        isovalue, sz=5, npts=100)
+
+            # loop over lobes
+            for k in range(0, 2):
+                vertices = posorb[0] if k == 0 else negorb[0]
+                normals = posorb[1] if k == 0 else negorb[1]
+                indices = posorb[2] if k == 0 else negorb[2]
+                vertices_normals = np.hstack([vertices * 0.529177, normals])
+
+                # write model idx
+                f.write(int(k).to_bytes(2, byteorder='little'))
+
+                # write model color
+                color = np.array(colors[k]) if i < (self.nelec//2) else np.array(colors[k+2])
+                f.write(color.tobytes())
+
+                # write number of vertices
+                f.write(vertices_normals.shape[0].to_bytes(4, byteorder='little'))
+
+                # write vertices
+                f.write(vertices_normals.tobytes())
+
+                # write number of indices
+                f.write(int(len(indices)/3).to_bytes(4, byteorder='little'))
+
+                # write indices
+                f.write(indices.tobytes())
+
+                if k == 0:
+                    print('    Writing positive lobe: %i vertices and %i facets' % (vertices_normals.shape[0], indices.shape[0] / 3))
+                else:
+                    print('    Writing negative lobe: %i vertices and %i facets' % (vertices_normals.shape[0], indices.shape[0] / 3))
+
+        f.close()
+
+        # report filesize
+        print("Creating file: %s" % outfile)
+        print("Size: %f MB" % (os.stat(outfile).st_size / (1024*1024)))
+
+    def __build_isosurface(self, cgfs, coeff, isovalue, sz=5, npts=100):
+        """
+        Construct isosurfaces from PyQInt output
+        """
+        # generate some data
+        isovalue = np.abs(isovalue)
+        integrator = PyQInt()
+        grid = integrator.build_rectgrid3d(-sz, sz, npts)
+        scalarfield = np.reshape(integrator.plot_wavefunction(grid, coeff, cgfs), (npts, npts, npts))
+        unitcell = np.diag(np.ones(3) * 2 * sz)
+
+        # write out
+        # (vertices, normals, indices)
+        pytessel = PyTessel()
+        posorb = pytessel.marching_cubes(scalarfield.flatten(), scalarfield.shape, unitcell.flatten(), isovalue)
+        negorb = pytessel.marching_cubes(scalarfield.flatten(), scalarfield.shape, unitcell.flatten(), -isovalue)
+
+        return (posorb, negorb)

pyqint/foster_boys.py

@@ -72,7 +72,9 @@ def __single_runner(self):
             'nriter': nriter,
             'mol': self.mol,
             'r2start': self.__calculate_r2(self.orbc_canonical),
-            'r2final': self.__calculate_r2(orbc)
+            'r2final': self.__calculate_r2(orbc),
+            'nelec': self.nelec,
+            'cgfs': self.cgfs,
         }
 
         return result

pyqint/molecule.py

@@ -23,6 +23,9 @@ def __str__(self):
             res += " %s (%f,%f,%f)\n" % (atom[0], atom[1][0], atom[1][1], atom[1][2])
 
         return res
+
+    def __len__(self):
+        return len(self.__atoms)
 
     def get_nelec(self):
         """