Skip to content

Commit

Permalink
Refactoring GA
Browse files Browse the repository at this point in the history
  • Loading branch information
@ifilot
ifilot committed Feb 2, 2025
2 parents f0f2cc9 + 33d8cd8 commit b4b12c6
Show file tree
Hide file tree
Showing 28 changed files with 1,478 additions and 34 deletions.
1 change: 1 addition & 0 deletions MANIFEST.in
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
include sphecerix/charactertables/*.json
1 change: 1 addition & 0 deletions examples/.gitignore
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
figures/*
116 changes: 116 additions & 0 deletions examples/benzene.py
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,116 @@
# -*- coding: utf-8 -*-

import sys
import os
import numpy as np

# add a reference to load the Sphecerix library
sys.path.append(os.path.join(os.path.dirname(__file__), '..'))

from sphecerix import Molecule, BasisFunction, SymmetryOperations,\
visualize_matrices, CharacterTable, ProjectionOperator

def main():
mol = Molecule()
mol.add_atom('C', 0.0000000015, -1.3868467444, 0.0000000000, unit='angstrom')
mol.add_atom('C', 1.2010445126, -0.6934233709, 0.0000000000, unit='angstrom')
mol.add_atom('C', 1.2010445111, 0.6934233735, 0.0000000000, unit='angstrom')
mol.add_atom('C', -0.0000000015, 1.3868467444, 0.0000000000, unit='angstrom')
mol.add_atom('C', -1.2010445126, 0.6934233709, 0.0000000000, unit='angstrom')
mol.add_atom('C', -1.2010445111, -0.6934233735, 0.0000000000, unit='angstrom')
mol.add_atom('H', 0.0000000027, -2.4694205285, 0.0000000000, unit='angstrom')
mol.add_atom('H', 2.1385809117, -1.2347102619, 0.0000000000, unit='angstrom')
mol.add_atom('H', 2.1385809090, 1.2347102666, 0.0000000000, unit='angstrom')
mol.add_atom('H', -0.0000000027, 2.4694205285, 0.0000000000, unit='angstrom')
mol.add_atom('H', -2.1385809117, 1.2347102619, 0.0000000000, unit='angstrom')
mol.add_atom('H', -2.1385809090, -1.2347102666, 0.0000000000, unit='angstrom')

molset = {
'C': [BasisFunction(2,0,0),
# BasisFunction(2,1,1), # x
# BasisFunction(2,1,-1), # y
BasisFunction(2,1,0)], # z
'H': [BasisFunction(1,0,0)]
}
mol.build_basis(molset)

symops = SymmetryOperations(mol)

# E
symops.add('identity')

# 2C6
symops.add('rotation', '6+', np.array([0,0,1]), 2.0 * np.pi / 6)
symops.add('rotation', '6-', np.array([0,0,1]), -2.0 * np.pi / 6)

# 2C3
symops.add('rotation', '3+', np.array([0,0,1]), 2.0 * np.pi / 3)
symops.add('rotation', '3-', np.array([0,0,1]), -2.0 * np.pi / 3)

# C2
symops.add('rotation', '2', np.array([0,0,1]), np.pi)

# 3C2'
for i in range(0,3):
symops.add('rotation', '2,%i' % i, np.array([np.sin(2.0 * np.pi * i/6.),
-np.cos(2.0 * np.pi * i/6.),
0.0]), np.pi)

# 3C2''
for i in range(0,3):
symops.add('rotation', '2,%i' % i, np.array([np.sin(2.0 * np.pi * (i/6. + 1./12)),
-np.cos(2.0 * np.pi * (i/6. + 1./12)),
0.0]), np.pi)


# inversion
symops.add('inversion')

# 2S3
symops.add('improper', '3+', np.array([0,0,1]), 2.0 * np.pi / 3)
symops.add('improper', '3-', np.array([0,0,1]), -2.0 * np.pi / 3)

# 2S6
symops.add('improper', '6+', np.array([0,0,1]), 2.0 * np.pi / 6)
symops.add('improper', '6-', np.array([0,0,1]), -2.0 * np.pi / 6)

# sigma_h
symops.add('mirror', 'h(xy)', np.array([0,0,1]))

# sigma_d
for i in range(0,3):
symops.add('mirror', 'd,%i' % i, np.array([np.cos(2.0 * np.pi * (i/6. + 1./12)),
np.sin(2.0 * np.pi * (i/6. + 1./12)),
0.0]))

# sigma_v
for i in range(0,3):
symops.add('mirror', 'v,%i' % i, np.array([np.cos(2.0 * np.pi * i/6),
np.sin(2.0 * np.pi * i/6),
0.0]))

symops.run()

visualize_matrices(symops.operation_matrices,
[op.name for op in symops.operations],
[bf.name for bf in symops.mol.basis],
xlabelrot=90, figsize=(24,32), numcols=4)

# print result LOT
ct = CharacterTable('d6h')
print(ct.lot(np.trace(symops.operation_matrices, axis1=1, axis2=2)))

# apply projection operator
po = ProjectionOperator(ct, symops)

mos = po.build_mos(verbose=True)
newmats = [mos @ m @ mos.transpose() for m in symops.operation_matrices]

visualize_matrices(newmats,
[op.name for op in symops.operations],
['$\phi_{%i}$' % (i+1) for i in range(len(symops.mol.basis))],
figsize=(24,32), numcols=4,
highlight_groups=po.get_blocks())

if __name__ == '__main__':
main()
97 changes: 97 additions & 0 deletions examples/ccl4.py
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,97 @@
# -*- coding: utf-8 -*-

import sys
import os
import numpy as np
import matplotlib.pyplot as plt

# add a reference to load the Sphecerix library
sys.path.append(os.path.join(os.path.dirname(__file__), '..'))

from sphecerix import Molecule, BasisFunction, SymmetryOperations,\
visualize_matrices, CharacterTable, ProjectionOperator, \
plot_matrix

def main():
mol = Molecule()
mol.add_atom('C', 0.0, 0.0, 0.0, unit='angstrom')
mol.add_atom('Cl', 1,1,1, unit='angstrom')
mol.add_atom('Cl', 1,-1,-1, unit='angstrom')
mol.add_atom('Cl', -1,1,-1, unit='angstrom')
mol.add_atom('Cl', -1,-1,1, unit='angstrom')

molset = {
'C': [BasisFunction(1,0,0),
BasisFunction(2,0,0),
BasisFunction(2,1,1),
BasisFunction(2,1,-1),
BasisFunction(2,1,0)],
'Cl': [BasisFunction(1,0,0),
BasisFunction(2,1,1),
BasisFunction(2,1,-1),
BasisFunction(2,1,0)]
}
mol.build_basis(molset)

symops = SymmetryOperations(mol)
symops.add('identity')

# add C3 rotations
for i in range(0,4):
axis = mol.atoms[i+1][1]
axis /= np.linalg.norm(axis) # normalize axis
for j in range(0,2):
symops.add('rotation', '3,%i' % (i*2+j+1), axis, (-1)**j * 2.0 * np.pi / 3)

# C2 rotations
for i in range(0,3):
axis = np.zeros(3)
axis[i] = 1.0
symops.add('rotation', '2,%i' % (i+1), axis, np.pi)

# S4 rotations
for i in range(0,3):
axis = np.zeros(3)
axis[i] = 1.0
for j in range(0,2):
symops.add('improper', '4,%i' % (i+1), axis, (-1)**j * np.pi/2)

# sigma_d mirror planes
ctr = 0
for i in range(0,4):
axis1 = mol.atoms[i+1][1]
for j in range(i+1,4):
axis2 = mol.atoms[j+1][1]
normal = np.cross(axis1, axis2)
normal /= np.linalg.norm(normal)
ctr += 1
symops.add('mirror', ',d(%i)' % (ctr), normal)

symops.run()

visualize_matrices(symops.operation_matrices,
[op.name for op in symops.operations],
[bf.name for bf in symops.mol.basis],
xlabelrot=90, figsize=(36,24), numcols=6)

# print result LOT
ct = CharacterTable('td')
print(ct.lot(np.trace(symops.operation_matrices, axis1=1, axis2=2)))

# # apply projection operator
po = ProjectionOperator(ct, symops)
mos = po.build_mos(verbose=True)

fig, ax = plt.subplots(1,1,dpi=144,figsize=(20,20))
plot_matrix(ax, mos,[bf.name for bf in symops.mol.basis],None,0)

newmats = [mos @ m @ mos.transpose() for m in symops.operation_matrices]

visualize_matrices(newmats,
[op.name for op in symops.operations],
['$\phi_{%i}$' % (i+1) for i in range(len(symops.mol.basis))],
xlabelrot=90, figsize=(36,24), numcols=6,
highlight_groups=po.get_blocks())

if __name__ == '__main__':
main()
88 changes: 88 additions & 0 deletions examples/ch4.py
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,88 @@
# -*- coding: utf-8 -*-

import sys
import os
import numpy as np

# add a reference to load the Sphecerix library
sys.path.append(os.path.join(os.path.dirname(__file__), '..'))

from sphecerix import Molecule, BasisFunction, SymmetryOperations,\
visualize_matrices, CharacterTable, ProjectionOperator

def main():
mol = Molecule()
mol.add_atom('C', 0.0, 0.0, 0.0, unit='angstrom')
mol.add_atom('H', 1,1,1, unit='angstrom')
mol.add_atom('H', 1,-1,-1, unit='angstrom')
mol.add_atom('H', -1,1,-1, unit='angstrom')
mol.add_atom('H', -1,-1,1, unit='angstrom')

molset = {
'C': [BasisFunction(1,0,0),
BasisFunction(2,0,0),
BasisFunction(2,1,1),
BasisFunction(2,1,-1),
BasisFunction(2,1,0)],
'H': [BasisFunction(1,0,0)]
}
mol.build_basis(molset)

symops = SymmetryOperations(mol)
symops.add('identity')

# add C3 rotations
for i in range(0,4):
axis = mol.atoms[i+1][1]
axis /= np.linalg.norm(axis) # normalize axis
for j in range(0,2):
symops.add('rotation', '3,%i' % (i*2+j+1), axis, (-1)**j * 2.0 * np.pi / 3)

# C2 rotations
for i in range(0,3):
axis = np.zeros(3)
axis[i] = 1.0
symops.add('rotation', '2,%i' % (i+1), axis, np.pi)

# S4 rotations
for i in range(0,3):
axis = np.zeros(3)
axis[i] = 1.0
for j in range(0,2):
symops.add('improper', '4,%i' % (i+1), axis, (-1)**j * np.pi/2)

# sigma_d mirror planes
ctr = 0
for i in range(0,4):
axis1 = mol.atoms[i+1][1]
for j in range(i+1,4):
axis2 = mol.atoms[j+1][1]
normal = np.cross(axis1, axis2)
normal /= np.linalg.norm(normal)
ctr += 1
symops.add('mirror', ',d(%i)' % (ctr), normal)

symops.run()

visualize_matrices(symops.operation_matrices,
[op.name for op in symops.operations],
[bf.name for bf in symops.mol.basis],
xlabelrot=90, figsize=(18,12), numcols=6)

# print result LOT
ct = CharacterTable('td')
print(ct.lot(np.trace(symops.operation_matrices, axis1=1, axis2=2)))

# # apply projection operator
po = ProjectionOperator(ct, symops)
mos = po.build_mos(verbose=True)
newmats = [mos @ m @ mos.transpose() for m in symops.operation_matrices]

visualize_matrices(newmats,
[op.name for op in symops.operations],
['$\phi_{%i}$' % (i+1) for i in range(len(symops.mol.basis))],
xlabelrot=90, figsize=(18,12), numcols=6,
highlight_groups=po.get_blocks())

if __name__ == '__main__':
main()
Loading

0 comments on commit b4b12c6

Please sign in to comment.