diff --git a/pyqint/cgf.py b/pyqint/cgf.py index 4ec5c6c..3765c3c 100644 --- a/pyqint/cgf.py +++ b/pyqint/cgf.py @@ -1,5 +1,6 @@ from .gto import gto from .pyqint import PyCGF +from . import spherical_harmonics as sh class cgf: """ @@ -48,6 +49,17 @@ def add_gto(self, c, alpha, l, m, n): self.gtos.append(gto(c, self.p, alpha, l, m, n)) self.cgf.add_gto(c, alpha, l, m, n) + def add_spherical_gto(self, c, alpha, l, m): + """ + Add Spherical Gaussian Type Orbital to Contracted Gaussian Function. + l and m are the coefficients of the requested spherical harmonic function. + l must be <= 6 and -l <= m <= l. + """ + if not l <= 6 or not abs(m) <=l: + raise ValueError("l must be <= 6 and -l <= m <= l") + for gto in sh.spherical_harmonics[l][m]: + self.add_gto(gto[0] * c, alpha, gto[1][0], gto[1][1], gto[1][2]) + def get_amp_f(self, x, y, z): """ Get the amplitude of the wave function at position r diff --git a/pyqint/spherical_harmonics.py b/pyqint/spherical_harmonics.py new file mode 100644 index 0000000..e8e5a73 --- /dev/null +++ b/pyqint/spherical_harmonics.py @@ -0,0 +1,294 @@ +""" +coeffecients of spherical harmonics up to l=6 for cartesian gaussian orbitals. +The format is: +spherical_harmonics[l][m] = [[c, [l_x, l_y, l_z], ...] +where l and m are the spherical harmonic coefficents, c is the weight +and [l_x, l_y, l_z] are the cartesian exponents. +""" + +spherical_harmonics = { + 0: { + 0: [ + [1.0, [0, 0, 0]] + ] + }, + 1: { + -1: [ + [1.0, [0, 1, 0]] + ], + 0: [ + [1.0, [0, 0, 1]] + ], + 1: [ + [1.0, [1, 0, 0]] + ] + }, + 2: { + -2: [ + [1.0, [1, 1, 0]] + ], + -1: [ + [1.0, [0, 1, 1]] + ], + 0: [ + [-0.5, [2, 0, 0]], + [-0.5, [0, 2, 0]], + [1.0, [0, 0, 2]] + ], + 1: [ + [1.0, [1, 0, 1]] + ], + 2: [ + [0.8660254037844386, [2, 0, 0]], + [-0.8660254037844386, [0, 2, 0]] + ] + }, + 3: { + -3: [ + [1.0606601717798212, [2, 1, 0]], + [-0.7905694150420949, [0, 3, 0]] + ], + -2: [ + [1.0, [1, 1, 1]] + ], + -1: [ + [-0.27386127875258304, [2, 1, 0]], + [-0.6123724356957945, [0, 3, 0]], + [1.0954451150103321, [0, 1, 2]] + ], + 0: [ + [-0.6708203932499369, [2, 0, 1]], + [-0.6708203932499369, [0, 2, 1]], + [1.0, [0, 0, 3]] + ], + 1: [ + [-0.6123724356957945, [3, 0, 0]], + [-0.27386127875258304, [1, 2, 0]], + [1.0954451150103321, [1, 0, 2]] + ], + 2: [ + [0.8660254037844386, [2, 0, 1]], + [-0.8660254037844386, [0, 2, 1]] + ], + 3: [ + [0.7905694150420949, [3, 0, 0]], + [-1.0606601717798212, [1, 2, 0]] + ] + }, + 4: { + -4: [ + [1.118033988749895, [3, 1, 0]], + [-1.118033988749895, [1, 3, 0]] + ], + -3: [ + [1.0606601717798212, [2, 1, 1]], + [-0.7905694150420949, [0, 3, 1]] + ], + -2: [ + [-0.4225771273642583, [3, 1, 0]], + [-0.4225771273642583, [1, 3, 0]], + [1.1338934190276817, [1, 1, 2]] + ], + -1: [ + [-0.4008918628686366, [2, 1, 1]], + [-0.8964214570007952, [0, 3, 1]], + [1.1952286093343936, [0, 1, 3]] + ], + 0: [ + [0.375, [4, 0, 0]], + [0.21957751641341997, [2, 2, 0]], + [-0.8783100656536799, [2, 0, 2]], + [0.375, [0, 4, 0]], + [-0.8783100656536799, [0, 2, 2]], + [1.0, [0, 0, 4]] + ], + 1: [ + [-0.8964214570007952, [3, 0, 1]], + [-0.4008918628686366, [1, 2, 1]], + [1.1952286093343936, [1, 0, 3]] + ], + 2: [ + [-0.5590169943749475, [4, 0, 0]], + [0.9819805060619657, [2, 0, 2]], + [0.5590169943749475, [0, 4, 0]], + [-0.9819805060619657, [0, 2, 2]] + ], + 3: [ + [0.7905694150420949, [3, 0, 1]], + [-1.0606601717798212, [1, 2, 1]] + ], + 4: [ + [0.739509972887452, [4, 0, 0]], + [-1.299038105676658, [2, 2, 0]], + [0.739509972887452, [0, 4, 0]] + ] + }, + 5: { + -5: [ + [1.1692679333668567, [4, 1, 0]], + [-1.5309310892394863, [2, 3, 0]], + [0.701560760020114, [0, 5, 0]] + ], + -4: [ + [1.118033988749895, [3, 1, 1]], + [-1.118033988749895, [1, 3, 1]] + ], + -3: [ + [-0.5229125165837972, [4, 1, 0]], + [-0.22821773229381923, [2, 3, 0]], + [1.224744871391589, [2, 1, 2]], + [0.5229125165837972, [0, 5, 0]], + [-0.9128709291752769, [0, 3, 2]] + ], + -2: [ + [-0.6454972243679028, [3, 1, 1]], + [-0.6454972243679028, [1, 3, 1]], + [1.2909944487358056, [1, 1, 3]] + ], + -1: [ + [0.1613743060919757, [4, 1, 0]], + [0.21128856368212914, [2, 3, 0]], + [-0.5669467095138409, [2, 1, 2]], + [0.4841229182759271, [0, 5, 0]], + [-1.267731382092775, [0, 3, 2]], + [1.2909944487358056, [0, 1, 4]] + ], + 0: [ + [0.625, [4, 0, 1]], + [0.36596252735569995, [2, 2, 1]], + [-1.0910894511799618, [2, 0, 3]], + [0.625, [0, 4, 1]], + [-1.0910894511799618, [0, 2, 3]], + [1.0, [0, 0, 5]] + ], + 1: [ + [0.4841229182759271, [5, 0, 0]], + [0.21128856368212914, [3, 2, 0]], + [-1.267731382092775, [3, 0, 2]], + [0.1613743060919757, [1, 4, 0]], + [-0.5669467095138409, [1, 2, 2]], + [1.2909944487358056, [1, 0, 4]] + ], + 2: [ + [-0.8539125638299665, [4, 0, 1]], + [1.118033988749895, [2, 0, 3]], + [0.8539125638299665, [0, 4, 1]], + [-1.118033988749895, [0, 2, 3]] + ], + 3: [ + [-0.5229125165837972, [5, 0, 0]], + [0.22821773229381923, [3, 2, 0]], + [0.9128709291752769, [3, 0, 2]], + [0.5229125165837972, [1, 4, 0]], + [-1.224744871391589, [1, 2, 2]] + ], + 4: [ + [0.739509972887452, [4, 0, 1]], + [-1.299038105676658, [2, 2, 1]], + [0.739509972887452, [0, 4, 1]] + ], + 5: [ + [0.701560760020114, [5, 0, 0]], + [-1.5309310892394863, [3, 2, 0]], + [1.1692679333668567, [1, 4, 0]] + ] + }, + 6: { + -6: [ + [1.2151388809514738, [5, 1, 0]], + [-1.976423537605237, [3, 3, 0]], + [1.2151388809514738, [1, 5, 0]] + ], + -5: [ + [1.1692679333668567, [4, 1, 1]], + [-1.5309310892394863, [2, 3, 1]], + [0.701560760020114, [0, 5, 1]] + ], + -4: [ + [-0.5982930264130992, [5, 1, 0]], + [1.3055824196677337, [3, 1, 2]], + [0.5982930264130992, [1, 5, 0]], + [-1.3055824196677337, [1, 3, 2]] + ], + -3: [ + [-0.8192464664112216, [4, 1, 1]], + [-0.35754847096709713, [2, 3, 1]], + [1.4301938838683885, [2, 1, 3]], + [0.8192464664112216, [0, 5, 1]], + [-1.0660035817780522, [0, 3, 3]] + ], + -2: [ + [0.27308215547040715, [5, 1, 0]], + [0.26650089544451305, [3, 3, 0]], + [-0.9534625892455924, [3, 1, 2]], + [0.27308215547040715, [1, 5, 0]], + [-0.9534625892455924, [1, 3, 2]], + [1.4564381625088383, [1, 1, 4]] + ], + -1: [ + [0.2878538665448989, [4, 1, 1]], + [0.3768891807222045, [2, 3, 1]], + [-0.753778361444409, [2, 1, 3]], + [0.8635615996346968, [0, 5, 1]], + [-1.6854996561581053, [0, 3, 3]], + [1.381698559415515, [0, 1, 5]] + ], + 0: [ + [-0.3125, [6, 0, 0]], + [-0.1631978024584667, [4, 2, 0]], + [0.9791868147508004, [4, 0, 2]], + [-0.1631978024584667, [2, 4, 0]], + [0.5733530903673287, [2, 2, 2]], + [-1.3055824196677337, [2, 0, 4]], + [-0.3125, [0, 6, 0]], + [0.9791868147508004, [0, 4, 2]], + [-1.3055824196677337, [0, 2, 4]], + [1.0, [0, 0, 6]] + ], + 1: [ + [0.8635615996346968, [5, 0, 1]], + [0.3768891807222045, [3, 2, 1]], + [-1.6854996561581053, [3, 0, 3]], + [0.2878538665448989, [1, 4, 1]], + [-0.753778361444409, [1, 2, 3]], + [1.381698559415515, [1, 0, 5]] + ], + 2: [ + [0.45285552331841994, [6, 0, 0]], + [0.0788320279858614, [4, 2, 0]], + [-1.2613124477737825, [4, 0, 2]], + [-0.0788320279858614, [2, 4, 0]], + [1.2613124477737825, [2, 0, 4]], + [-0.45285552331841994, [0, 6, 0]], + [1.2613124477737825, [0, 4, 2]], + [-1.2613124477737825, [0, 2, 4]] + ], + 3: [ + [-0.8192464664112216, [5, 0, 1]], + [0.35754847096709713, [3, 2, 1]], + [1.0660035817780522, [3, 0, 3]], + [0.8192464664112216, [1, 4, 1]], + [-1.4301938838683885, [1, 2, 3]] + ], + 4: [ + [-0.49607837082461076, [6, 0, 0]], + [0.4317807998173484, [4, 2, 0]], + [0.8635615996346968, [4, 0, 2]], + [0.4317807998173484, [2, 4, 0]], + [-1.5169496905422948, [2, 2, 2]], + [-0.49607837082461076, [0, 6, 0]], + [0.8635615996346968, [0, 4, 2]] + ], + 5: [ + [0.701560760020114, [5, 0, 1]], + [-1.5309310892394863, [3, 2, 1]], + [1.1692679333668567, [1, 4, 1]] + ], + 6: [ + [0.6716932893813962, [6, 0, 0]], + [-1.7539019000502851, [4, 2, 0]], + [1.7539019000502851, [2, 4, 0]], + [-0.6716932893813962, [0, 6, 0]] + ] + } +} \ No newline at end of file diff --git a/tests/test_cgf.py b/tests/test_cgf.py index df47099..a379e18 100644 --- a/tests/test_cgf.py +++ b/tests/test_cgf.py @@ -1,5 +1,5 @@ import unittest -from pyqint import PyQInt, Molecule +from pyqint import PyQInt, Molecule, cgf from copy import deepcopy import numpy as np import os @@ -67,5 +67,24 @@ def testPlotGrid(self): ans = np.load(os.path.join(os.path.dirname(__file__), 'results', 'h2o_orb_1b2.npy')) np.testing.assert_almost_equal(res, ans, 6) + def testSphericalHarmonicsOnSite(self): + """ + Check if the overlap matrix of all spherical harmonics up to l=6 is the identity matrix + """ + basis_functions = [] + p0 = [0.0, 0.0, 0.0] + for l in range(7): + for m in range(-l, l+1): + orb = cgf(p0) + orb.add_spherical_gto(1.0, 1.0, l, m) + basis_functions.append(orb) + # build overlap matrix + om = np.zeros((len(basis_functions), len(basis_functions))) + integrator = PyQInt() + for i,orb1 in enumerate(basis_functions): + for j,orb2 in enumerate(basis_functions): + om[i,j] = integrator.overlap(orb1, orb2) + np.testing.assert_almost_equal(om, np.eye(om.shape[0]), 6) + if __name__ == '__main__': unittest.main()