10. symm – Point group symmetry and spin symmetry¶
This module offers the functions to detect point group symmetry, basis symmetriziation, Clebsch-Gordon coefficients. This module works as a plugin of PySCF package. Symmetry is not hard coded in each method.
PySCF supports D2h symmetry and linear molecule symmetry (Dooh and Coov). For D2h, the direct production of representations are
D2h | A1g | B1g | B2g | B3g | A1u | B1u | B2u | B3u |
---|---|---|---|---|---|---|---|---|
A1g | A1g | |||||||
B1g | B1g | A1g | ||||||
B2g | B2g | B3g | A1g | |||||
B3g | B3g | B2g | B1g | A1g | ||||
A1u | A1u | B1u | B2u | B3u | A1g | |||
B1u | B1u | A1u | B3u | B2u | B1g | A1g | ||
B2u | B2u | B3u | A1u | B1u | B2g | B3g | A1g | |
B3u | B3u | B2u | B1u | A1u | B3g | B2g | B1g | A1g |
The multiplication table for XOR operator is
XOR | 000 | 001 | 010 | 011 | 100 | 101 | 110 | 111 |
---|---|---|---|---|---|---|---|---|
000 | 000 | |||||||
001 | 001 | 000 | ||||||
010 | 010 | 011 | 000 | |||||
011 | 011 | 010 | 001 | 000 | ||||
100 | 100 | 101 | 110 | 111 | 000 | |||
101 | 101 | 100 | 111 | 110 | 001 | 000 | ||
110 | 110 | 111 | 100 | 101 | 010 | 011 | 000 | |
111 | 111 | 110 | 101 | 100 | 011 | 010 | 001 | 000 |
Comparing the two table, we notice that the two tables can be changed to each other with the mapping
D2h | XOR | ID |
---|---|---|
A1g | 000 | 0 |
B1g | 001 | 1 |
B2g | 010 | 2 |
B3g | 011 | 3 |
A1u | 100 | 4 |
B1u | 101 | 5 |
B2u | 110 | 6 |
B3u | 111 | 7 |
The XOR operator and the D2h subgroups have the similar relationships.
We therefore use the XOR operator ID to assign the irreps (see
pyscf/symm/param.py
).
C2h | XOR | ID | C2v | XOR | ID | D2 | XOR | ID |
---|---|---|---|---|---|---|---|---|
Ag | 00 | 0 | A1 | 00 | 0 | A1 | 00 | 0 |
Bg | 01 | 1 | A2 | 01 | 1 | B1 | 01 | 1 |
Au | 10 | 2 | B1 | 10 | 2 | B2 | 10 | 2 |
Bu | 11 | 3 | B2 | 11 | 3 | B3 | 11 | 3 |
Cs | XOR | ID | Ci | XOR | ID | C2 | XOR | ID |
---|---|---|---|---|---|---|---|---|
A’ | 0 | 0 | Ag | 0 | 0 | A | 0 | 0 |
B” | 1 | 1 | Au | 1 | 1 | B | 1 | 1 |
To easily get the relationship between the linear molecule symmetry
and D2h/C2v, the ID for irreducible representations of linear molecule
symmetry are chosen as (see pyscf/symm/basis.py
)
\(D_{\infty h}\) | ID | \(D_{2h}\) | ID |
A1g | 0 | Ag | 0 |
A2g | 1 | B1g | 1 |
A1u | 5 | B1u | 5 |
A2u | 4 | Au | 4 |
E1gx | 2 | B2g | 2 |
E1gy | 3 | B3g | 3 |
E1uy | 6 | B2u | 6 |
E1ux | 7 | B3u | 7 |
E2gx | 10 | Ag | 0 |
E2gy | 11 | B1g | 1 |
E2ux | 15 | B1u | 5 |
E2uy | 14 | Au | 4 |
E3gx | 12 | B2g | 2 |
E3gy | 13 | B3g | 3 |
E3uy | 16 | B2u | 6 |
E3ux | 17 | B3u | 7 |
E4gx | 20 | Ag | 0 |
E4gy | 21 | B1g | 1 |
E4ux | 25 | B1u | 5 |
E4uy | 24 | Au | 4 |
E5gx | 22 | B2g | 2 |
E5gy | 23 | B3g | 3 |
E5uy | 26 | B2u | 6 |
E5ux | 27 | B3u | 7 |
and
\(C_{\infty v}\) | ID | \(C_{2v}\) | ID |
A1 | 0 | A1 | 0 |
A2 | 1 | A2 | 1 |
E1x | 2 | B1 | 2 |
E1y | 3 | B2 | 3 |
E2x | 10 | A1 | 0 |
E2y | 11 | A2 | 1 |
E3x | 12 | B1 | 2 |
E3y | 13 | B2 | 3 |
E4x | 20 | A1 | 0 |
E4y | 21 | A2 | 1 |
E5x | 22 | B1 | 2 |
E5y | 23 | B2 | 3 |
So that, the subduction from linear molecule symmetry to D2h/C2v can be
achieved by the modular operation %10
.
In many output messages, the irreducible representations are labeld with
the IDs instead of the irreps’ symbols. We can use
symm.addons.irrep_id2name()
to convert the ID to irreps’ symbol,
e.g.:
>>> from pyscf import symm
>>> [symm.irrep_id2name('Dooh', x) for x in [7, 6, 0, 10, 11, 0, 5, 3, 2, 5, 15, 14]]
['E1ux', 'E1uy', 'A1g', 'E2gx', 'E2gy', 'A1g', 'A1u', 'E1gy', 'E1gx', 'A1u', 'E2ux', 'E2uy']
10.1. Enabling symmetry in other module¶
SCF
To control the HF determinant symmetry, one can assign occupancy for particular irreps, e.g.
FCI
FCI wavefunction symmetry can be controlled by initial guess. Function
fci.addons.symm_initguess()
can generate the FCI initial guess with the right symmetry.MCSCF
The symmetry of active space in the CASCI/CASSCF calculations can controlled
MP2 and CCSD
Point group symmetry are not supported in CCSD, MP2.
10.1.1. Examples¶
Relevant examples
examples/symm/30-guess_symmetry.py
examples/symm/31-symmetrize_space.py
examples/symm/32-symmetrize_natural_orbital.py
examples/symm/33-lz_adaption.py
10.1.2. Program reference¶
10.2. geom¶
-
pyscf.symm.geom.
alias_axes
(axes, ref)[source]¶ Rename axes, make it as close as possible to the ref axes
-
pyscf.symm.geom.
detect_symm
(atoms, basis=None, verbose=2)[source]¶ Detect the point group symmetry for given molecule.
Return group name, charge center, and nex_axis (three rows for x,y,z)
10.3. Symmtry adapted basis¶
Generate symmetry adapted basis
10.4. addons¶
-
pyscf.symm.addons.
eigh
(h, orbsym)[source]¶ Solve eigenvalue problem based on the symmetry information for basis. See also pyscf/lib/linalg_helper.py
eigh_by_blocks()
Examples:
>>> from pyscf import gto, symm >>> mol = gto.M(atom='H 0 0 0; H 0 0 1', basis='ccpvdz', symmetry=True) >>> c = numpy.hstack(mol.symm_orb) >>> vnuc_so = reduce(numpy.dot, (c.T, mol.intor('int1e_nuc_sph'), c)) >>> orbsym = symm.label_orb_symm(mol, mol.irrep_name, mol.symm_orb, c) >>> symm.eigh(vnuc_so, orbsym) (array([-4.50766885, -1.80666351, -1.7808565 , -1.7808565 , -1.74189134, -0.98998583, -0.98998583, -0.40322226, -0.30242374, -0.07608981]), ...)
-
pyscf.symm.addons.
irrep_id2name
(gpname, irrep_id)[source]¶ Convert the internal irrep ID to irrep symbol
- Args:
- gpname : str
- The point group symbol
- irrep_id : int
- See IRREP_ID_TABLE in pyscf/symm/param.py
- Returns:
- Irrep sybmol, str
-
pyscf.symm.addons.
irrep_name2id
(gpname, symb)[source]¶ Convert the irrep symbol to internal irrep ID
- Args:
- gpname : str
- The point group symbol
- symb : str
- Irrep symbol
- Returns:
- Irrep ID, int
-
pyscf.symm.addons.
label_orb_symm
(mol, irrep_name, symm_orb, mo, s=None, check=True, tol=1e-09)[source]¶ Label the symmetry of given orbitals
irrep_name can be either the symbol or the ID of the irreducible representation. If the ID is provided, it returns the numeric code associated with XOR operator, see
symm.param.IRREP_ID_TABLE()
- Args:
mol : an instance of
Mole
- irrep_name : list of str or int
- A list of irrep ID or name, it can be either mol.irrep_id or mol.irrep_name. It can affect the return “label”.
- symm_orb : list of 2d array
- the symmetry adapted basis
- mo : 2d array
- the orbitals to label
- Returns:
- list of symbols or integers to represent the irreps for the given orbitals
Examples:
>>> from pyscf import gto, scf, symm >>> mol = gto.M(atom='H 0 0 0; H 0 0 1', basis='ccpvdz',verbose=0, symmetry=1) >>> mf = scf.RHF(mol) >>> mf.kernel() >>> symm.label_orb_symm(mol, mol.irrep_name, mol.symm_orb, mf.mo_coeff) ['Ag', 'B1u', 'Ag', 'B1u', 'B2u', 'B3u', 'Ag', 'B2g', 'B3g', 'B1u'] >>> symm.label_orb_symm(mol, mol.irrep_id, mol.symm_orb, mf.mo_coeff) [0, 5, 0, 5, 6, 7, 0, 2, 3, 5]
-
pyscf.symm.addons.
route
(target, nelec, orbsym)[source]¶ Pick orbitals to form a determinant which has the right symmetry. If solution is not found, return []
-
pyscf.symm.addons.
std_symb
(gpname)[source]¶ std_symb(‘d2h’) returns D2h; std_symb(‘D2H’) returns D2h
-
pyscf.symm.addons.
symmetrize_orb
(mol, mo, orbsym=None, s=None, check=False)[source]¶ Symmetrize the given orbitals.
This function is different to the
symmetrize_space()
: In this function, each orbital is symmetrized by removing non-symmetric components.symmetrize_space()
symmetrizes the entire space by mixing different orbitals.Note this function might return non-orthorgonal orbitals. Call
symmetrize_space()
to find the symmetrized orbitals that are close to the given orbitals.- Args:
- mo : 2D float array
- The orbital space to symmetrize
- Kwargs:
- orbsym : integer list
- Irrep id for each orbital. If not given, the irreps are guessed
by calling
label_orb_symm()
. - s : 2D float array
- Overlap matrix. If given, use this overlap than the the overlap of the input mol.
- Returns:
- 2D orbital coefficients
Examples:
>>> from pyscf import gto, symm, scf >>> mol = gto.M(atom = 'C 0 0 0; H 1 1 1; H -1 -1 1; H 1 -1 -1; H -1 1 -1', ... basis = 'sto3g') >>> mf = scf.RHF(mol).run() >>> mol.build(0, 0, symmetry='D2') >>> mo = symm.symmetrize_orb(mol, mf.mo_coeff) >>> print(symm.label_orb_symm(mol, mol.irrep_name, mol.symm_orb, mo)) ['A', 'A', 'B1', 'B2', 'B3', 'A', 'B1', 'B2', 'B3']
-
pyscf.symm.addons.
symmetrize_space
(mol, mo, s=None, check=True, tol=1e-07)[source]¶ Symmetrize the given orbital space.
This function is different to the
symmetrize_orb()
: In this function, the given orbitals are mixed to reveal the symmtery;symmetrize_orb()
projects out non-symmetric components for each orbital.- Args:
- mo : 2D float array
- The orbital space to symmetrize
- Kwargs:
- s : 2D float array
- Overlap matrix. If not given, overlap is computed with the input mol.
- Returns:
- 2D orbital coefficients
Examples:
>>> from pyscf import gto, symm, scf >>> mol = gto.M(atom = 'C 0 0 0; H 1 1 1; H -1 -1 1; H 1 -1 -1; H -1 1 -1', ... basis = 'sto3g') >>> mf = scf.RHF(mol).run() >>> mol.build(0, 0, symmetry='D2') >>> mo = symm.symmetrize_space(mol, mf.mo_coeff) >>> print(symm.label_orb_symm(mol, mol.irrep_name, mol.symm_orb, mo)) ['A', 'A', 'A', 'B1', 'B1', 'B2', 'B2', 'B3', 'B3']
10.5. Clebsch Gordon coefficients¶
10.6. Parameters¶
10.7. Spherical harmonics¶
Spherical harmonics
-
pyscf.symm.sph.
multipoles
(r, lmax, reorder_dipole=True)[source]¶ Compute all multipoles upto lmax
rad = numpy.linalg.norm(r, axis=1) ylms = real_ylm(r/rad.reshape(-1,1), lmax) pol = [rad**l*y for l, y in enumerate(ylms)]
- Kwargs:
- reorder_p : bool
- sort dipole to the order (x,y,z)
-
pyscf.symm.sph.
real2spinor_whole
(mol)¶ Transformation matrix that transforms real-spherical GTOs to spinor GTOs for all basis functions
Examples:
>>> from pyscf import gto >>> from pyscf.symm import sph >>> mol = gto.M(atom='H 0 0 0; F 0 0 1', basis='ccpvtz') >>> ca, cb = sph.sph2spinor_coeff(mol) >>> s0 = mol.intor('int1e_ovlp_spinor') >>> s1 = ca.conj().T.dot(mol.intor('int1e_ovlp_sph')).dot(ca) >>> s1+= cb.conj().T.dot(mol.intor('int1e_ovlp_sph')).dot(cb) >>> print(abs(s1-s0).max()) >>> 6.66133814775e-16
-
pyscf.symm.sph.
real_sph_vec
(r, lmax, reorder_p=False)[source]¶ Computes (all) real spherical harmonics up to the angular momentum lmax
-
pyscf.symm.sph.
sph2spinor_coeff
(mol)[source]¶ Transformation matrix that transforms real-spherical GTOs to spinor GTOs for all basis functions
Examples:
>>> from pyscf import gto >>> from pyscf.symm import sph >>> mol = gto.M(atom='H 0 0 0; F 0 0 1', basis='ccpvtz') >>> ca, cb = sph.sph2spinor_coeff(mol) >>> s0 = mol.intor('int1e_ovlp_spinor') >>> s1 = ca.conj().T.dot(mol.intor('int1e_ovlp_sph')).dot(ca) >>> s1+= cb.conj().T.dot(mol.intor('int1e_ovlp_sph')).dot(cb) >>> print(abs(s1-s0).max()) >>> 6.66133814775e-16
-
pyscf.symm.sph.
sph_pure2real
(l, reorder_p=True)[source]¶ Transformation matrix: from the pure spherical harmonic functions Y_m to the real spherical harmonic functions O_m.
O_m = sum Y_m’ * U(m’,m)Y(-1) = 1/sqrt(2){-iO(-1) + O(1)}; Y(1) = 1/sqrt(2){-iO(-1) - O(1)} Y(-2) = 1/sqrt(2){-iO(-2) + O(2)}; Y(2) = 1/sqrt(2){iO(-2) + O(2)} O(-1) = i/sqrt(2){Y(-1) + Y(1)}; O(1) = 1/sqrt(2){Y(-1) - Y(1)} O(-2) = i/sqrt(2){Y(-2) - Y(2)}; O(2) = 1/sqrt(2){Y(-2) + Y(2)}
- Kwargs:
- reorder_p (bool): Whether the p functions are in the (x,y,z) order.
- Returns:
- 2D array U_{complex,real}
10.8. Wigner rotation Dmatrix¶
Wigner rotation D-matrix for real spherical harmonics
-
pyscf.symm.Dmatrix.
Dmatrix
(l, alpha, beta, gamma, reorder_p=False)[source]¶ Wigner rotation D-matrix
D_{mm’} = <lm|R(alpha,beta,gamma)|lm’> alpha, beta, gamma are Euler angles (in z-y-z convention)
- Kwargs:
- reorder_p (bool): Whether to put the p functions in the (x,y,z) order.
-
pyscf.symm.Dmatrix.
dmatrix
(l, beta, reorder_p=False)[source]¶ Wigner small-d matrix (in z-y-z convention)
-
pyscf.symm.Dmatrix.
get_euler_angles
(c1, c2)[source]¶ Find the three Euler angles (alpha, beta, gamma in z-y-z convention) that rotates coordinates c1 to coordinates c2.
yp = numpy.einsum(‘j,kj->k’, c1[1], geom.rotation_mat(c1[2], beta)) tmp = numpy.einsum(‘ij,kj->ik’, c1 , geom.rotation_mat(c1[2], alpha)) tmp = numpy.einsum(‘ij,kj->ik’, tmp, geom.rotation_mat(yp , beta )) c2 = numpy.einsum(‘ij,kj->ik’, tmp, geom.rotation_mat(c2[2], gamma))
(For backward compatibility) if c1 and c2 are two points in the real space, the Euler angles define the rotation transforms the old coordinates to the new coordinates (new_x, new_y, new_z) in which c1 is identical to c2.
tmp = numpy.einsum(‘j,kj->k’, c1 , geom.rotation_mat((0,0,1), gamma)) tmp = numpy.einsum(‘j,kj->k’, tmp, geom.rotation_mat((0,1,0), beta) ) c2 = numpy.einsum(‘j,kj->k’, tmp, geom.rotation_mat((0,0,1), alpha))