Skip to content

Commit

Permalink
Feature/spherical harmonics (#39)
Browse files Browse the repository at this point in the history 
* add coefficients for spherical harmonics

* add spherica gto

* also add test

* more docs

* test for valid input

---------

Co-authored-by: Moritz Gubler <moritz.gubler@e.email>
  • Loading branch information
@moritzgubler
moritzgubler and Moritz Gubler authored Sep 11, 2024
1 parent e46527f commit 4a6f3f6
Show file tree
Hide file tree
Showing 3 changed files with 326 additions and 1 deletion.
12 changes: 12 additions & 0 deletions pyqint/cgf.py
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
from .gto import gto
from .pyqint import PyCGF
from . import spherical_harmonics as sh

class cgf:
"""
Expand Down Expand Up @@ -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
Expand Down
294 changes: 294 additions & 0 deletions pyqint/spherical_harmonics.py
View file
Original file line number Diff line number Diff line change
@@ -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]]
]
}
}
21 changes: 20 additions & 1 deletion tests/test_cgf.py
View file
Original file line number Diff line number Diff line change
@@ -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
Expand Down Expand Up @@ -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()

0 comments on commit 4a6f3f6

Please sign in to comment.