9. fci — Full configuration interaction

Different FCI solvers are implemented to support different type of symmetry.
Symmetry

File Point group Spin singlet Real hermitian* Alpha/beta degeneracy direct_spin0_symm Yes Yes Yes Yes direct_spin1_symm Yes No Yes Yes direct_spin0 No Yes Yes Yes direct_spin1 No No Yes Yes direct_uhf No No Yes No direct_nosym No No No** Yes

  • Real hermitian Hamiltonian implies (ij|kl) = (ji|kl) = (ij|lk) = (ji|lk)

** Hamiltonian is real but not hermitian, (ij|kl) != (ji|kl) ...

9.1. direct CI

Full CI solver for spin-free Hamiltonian. This solver can be used to compute doublet, triplet,...

The CI wfn are stored as a 2D array [alpha,beta], where each row corresponds to an alpha string. For each row (alpha string), there are total-num-beta-strings of columns. Each column corresponds to a beta string.

Different FCI solvers are implemented to support different type of symmetry.
Symmetry

File Point group Spin singlet Real hermitian* Alpha/beta degeneracy direct_spin0_symm Yes Yes Yes Yes direct_spin1_symm Yes No Yes Yes direct_spin0 No Yes Yes Yes direct_spin1 No No Yes Yes direct_uhf No No Yes No direct_nosym No No No** Yes

  • Real hermitian Hamiltonian implies (ij|kl) = (ji|kl) = (ij|lk) = (ji|lk)

** Hamiltonian is real but not hermitian, (ij|kl) != (ji|kl) ...

pyscf.fci.direct_spin1.FCI

alias of FCISolver

class pyscf.fci.direct_spin1.FCISolver(mol=None)[source]

Full CI solver

Attributes:
verbose
: int
Print level. Default value equals to Mole.verbose.
max_cycle
: int
Total number of iterations. Default is 100
max_space
: tuple of int
Davidson iteration space size. Default is 14.
conv_tol
: float
Energy convergence tolerance. Default is 1e-10.
level_shift
: float
Level shift applied in the preconditioner to avoid singularity. Default is 1e-3
davidson_only
: bool
By default, the entire Hamiltonian matrix will be constructed and diagonalized if the system is small (see attribute pspace_size). Setting this parameter to True will enforce the eigenvalue problems being solved by Davidson subspace algorithm. This flag should be enabled when initial guess is given or particular spin symmetry or point-group symmetry is required because the initial guess or symmetry are completely ignored in the direct diagonlization.
pspace_size
: int
The dimension of Hamiltonian matrix over which Davidson iteration algorithm will be used for the eigenvalue problem. Default is 400. This is roughly corresponding to a (6e,6o) system.
nroots
: int
Number of states to be solved. Default is 1, the ground state.
spin
: int or None
Spin (2S = nalpha-nbeta) of the system. If this attribute is None, spin will be determined by the argument nelec (number of electrons) of the kernel function.
wfnsym
: str or int
Symmetry of wavefunction. It is used only in direct_spin1_symm and direct_spin0_symm solver.

Saved results

eci
: float or a list of float
FCI energy(ies)
ci
: nparray
FCI wfn vector(s)
converged
: bool (or a list of bool for multiple roots)
Whether davidson iteration is converged

Examples:

>>> from pyscf import gto, scf, ao2mo, fci
>>> mol = gto.M(atom='Li 0 0 0; Li 0 0 1', basis='sto-3g')
>>> mf = scf.RHF(mol).run()
>>> h1 = mf.mo_coeff.T.dot(mf.get_hcore()).dot(mf.mo_coeff)
>>> eri = ao2mo.kernel(mol, mf.mo_coeff)
>>> cisolver = fci.direct_spin1.FCI(mol)
>>> e, ci = cisolver.kernel(h1, eri, h1.shape[1], mol.nelec, ecore=mol.energy_nuc())
>>> print(e)
-14.4197890826
absorb_h1e(h1e, eri, norb, nelec, fac=1)[source]

Modify 2e Hamiltonian to include 1e Hamiltonian contribution.

contract_1e(f1e, fcivec, norb, nelec, link_index=None, **kwargs)[source]

Contract the 1-electron Hamiltonian with a FCI vector to get a new FCI vector.

contract_2e(eri, fcivec, norb, nelec, link_index=None, **kwargs)[source]

Contract the 2-electron Hamiltonian with a FCI vector to get a new FCI vector.

Note the input arg eri is NOT the 2e hamiltonian matrix, the 2e hamiltonian is

\[\begin{split}h2e &= eri_{pq,rs} p^+ q r^+ s \\ &= (pq|rs) p^+ r^+ s q - (pq|rs) \delta_{qr} p^+ s\end{split}\]

So eri is defined as

\[eri_{pq,rs} = (pq|rs) - (1/Nelec) \sum_q (pq|qs)\]

to restore the symmetry between pq and rs,

\[eri_{pq,rs} = (pq|rs) - (.5/Nelec) [\sum_q (pq|qs) + \sum_p (pq|rp)]\]

See also direct_spin1.absorb_h1e()

energy(h1e, eri, fcivec, norb, nelec, link_index=None)[source]

Compute the FCI electronic energy for given Hamiltonian and FCI vector.

get_init_guess(norb, nelec, nroots, hdiag)[source]

Initial guess is the single Slater determinant

make_hdiag(h1e, eri, norb, nelec)[source]

Diagonal Hamiltonian for Davidson preconditioner

make_rdm1(fcivec, norb, nelec, link_index=None)[source]

Spin-traced one-particle density matrix

dm1[p,q] = <q_alpha^dagger p_alpha> + <q_beta^dagger p_beta>

The convention is based on McWeeney’s book, Eq (5.4.20) The contraction between 1-particle Hamiltonian and rdm1 is E = einsum(‘pq,qp’, h1, rdm1)

make_rdm12(fcivec, norb, nelec, link_index=None, reorder=True)[source]

Spin traced 1- and 2-particle density matrices.

1pdm[p,q] = :math:`langle q_alpha^dagger p_alpha rangle +
langle q_beta^dagger p_beta rangle`;
2pdm[p,q,r,s] = :math:`langle p_alpha^dagger r_alpha^dagger s_alpha q_alpharangle +
langle p_beta^dagger r_alpha^dagger s_alpha q_betarangle + langle p_alpha^dagger r_beta^dagger s_beta q_alpharangle + langle p_beta^dagger r_beta^dagger s_beta q_betarangle`.

Energy should be computed as E = einsum(‘pq,qp’, h1, 1pdm) + 1/2 * einsum(‘pqrs,pqrs’, eri, 2pdm) where h1[p,q] = <p|h|q> and eri[p,q,r,s] = (pq|rs)

make_rdm12s(fcivec, norb, nelec, link_index=None, reorder=True)[source]

Spin separated 1- and 2-particle density matrices. The return values include two lists, a list of 1-particle density matrices and a list of 2-particle density matrices. The density matrices are: (alpha,alpha), (beta,beta) for 1-particle density matrices; (alpha,alpha,alpha,alpha), (alpha,alpha,beta,beta), (beta,beta,beta,beta) for 2-particle density matrices.

1pdm[p,q] = \(\langle q^\dagger p\rangle\); 2pdm[p,q,r,s] = \(\langle p^\dagger r^\dagger s q\rangle\).

Energy should be computed as E = einsum(‘pq,qp’, h1, 1pdm) + 1/2 * einsum(‘pqrs,pqrs’, eri, 2pdm) where h1[p,q] = <p|h|q> and eri[p,q,r,s] = (pq|rs)

make_rdm1s(fcivec, norb, nelec, link_index=None)[source]

Spin separated 1-particle density matrices. The return values include two density matrices: (alpha,alpha), (beta,beta)

dm1[p,q] = <q^dagger p>

The convention is based on McWeeney’s book, Eq (5.4.20). The contraction between 1-particle Hamiltonian and rdm1 is E = einsum(‘pq,qp’, h1, rdm1)

make_rdm2(fcivec, norb, nelec, link_index=None, reorder=True)[source]

Spin traced 2-particle density matrice

NOTE the 2pdm is \(\langle p^\dagger q^\dagger s r\rangle\) but stored as [p,r,q,s]

pspace(h1e, eri, norb, nelec, hdiag=None, np=400)[source]

pspace Hamiltonian to improve Davidson preconditioner. See, CPL, 169, 463

spin_square(fcivec, norb, nelec)[source]

Spin square for RHF-FCI CI wfn only (obtained from spin-degenerated Hamiltonian)

trans_rdm1(cibra, ciket, norb, nelec, link_index=None)[source]

Spin traced transition 1-particle transition density matrices.

1pdm[p,q] = :math:`langle q_alpha^dagger p_alpha rangle
  • langle q_beta^dagger p_beta rangle`
trans_rdm12(cibra, ciket, norb, nelec, link_index=None, reorder=True)[source]

Spin traced transition 1- and 2-particle transition density matrices.

1pdm[p,q] = \(\langle q^\dagger p\rangle\); 2pdm[p,q,r,s] = \(\langle p^\dagger r^\dagger s q\rangle\).

trans_rdm12s(cibra, ciket, norb, nelec, link_index=None, reorder=True)[source]

Spin separated 1- and 2-particle transition density matrices. The return values include two lists, a list of 1-particle transition density matrices and a list of 2-particle transition density matrices. The density matrices are: (alpha,alpha), (beta,beta) for 1-particle transition density matrices; (alpha,alpha,alpha,alpha), (alpha,alpha,beta,beta), (beta,beta,alpha,alpha), (beta,beta,beta,beta) for 2-particle transition density matrices.

1pdm[p,q] = \(\langle q^\dagger p\rangle\); 2pdm[p,q,r,s] = \(\langle p^\dagger r^\dagger s q\rangle\).

trans_rdm1s(cibra, ciket, norb, nelec, link_index=None)[source]

Spin separated transition 1-particle density matrices. The return values include two density matrices: (alpha,alpha), (beta,beta). See also function make_rdm1s()

1pdm[p,q] = \(\langle q^\dagger p \rangle\)

pyscf.fci.direct_spin1.absorb_h1e(h1e, eri, norb, nelec, fac=1)[source]

Modify 2e Hamiltonian to include 1e Hamiltonian contribution.

pyscf.fci.direct_spin1.contract_1e(f1e, fcivec, norb, nelec, link_index=None)[source]

Contract the 1-electron Hamiltonian with a FCI vector to get a new FCI vector.

pyscf.fci.direct_spin1.contract_2e(eri, fcivec, norb, nelec, link_index=None)[source]

Contract the 2-electron Hamiltonian with a FCI vector to get a new FCI vector.

Note the input arg eri is NOT the 2e hamiltonian matrix, the 2e hamiltonian is

\[\begin{split}h2e &= eri_{pq,rs} p^+ q r^+ s \\ &= (pq|rs) p^+ r^+ s q - (pq|rs) \delta_{qr} p^+ s\end{split}\]

So eri is defined as

\[eri_{pq,rs} = (pq|rs) - (1/Nelec) \sum_q (pq|qs)\]

to restore the symmetry between pq and rs,

\[eri_{pq,rs} = (pq|rs) - (.5/Nelec) [\sum_q (pq|qs) + \sum_p (pq|rp)]\]

See also direct_spin1.absorb_h1e()

pyscf.fci.direct_spin1.energy(h1e, eri, fcivec, norb, nelec, link_index=None)[source]

Compute the FCI electronic energy for given Hamiltonian and FCI vector.

pyscf.fci.direct_spin1.get_init_guess(norb, nelec, nroots, hdiag)[source]

Initial guess is the single Slater determinant

pyscf.fci.direct_spin1.make_hdiag(h1e, eri, norb, nelec)[source]

Diagonal Hamiltonian for Davidson preconditioner

pyscf.fci.direct_spin1.make_rdm1(fcivec, norb, nelec, link_index=None)[source]

Spin-traced one-particle density matrix

dm1[p,q] = <q_alpha^dagger p_alpha> + <q_beta^dagger p_beta>

The convention is based on McWeeney’s book, Eq (5.4.20) The contraction between 1-particle Hamiltonian and rdm1 is E = einsum(‘pq,qp’, h1, rdm1)

pyscf.fci.direct_spin1.make_rdm12(fcivec, norb, nelec, link_index=None, reorder=True)[source]

Spin traced 1- and 2-particle density matrices.

1pdm[p,q] = :math:`langle q_alpha^dagger p_alpha rangle +
langle q_beta^dagger p_beta rangle`;
2pdm[p,q,r,s] = :math:`langle p_alpha^dagger r_alpha^dagger s_alpha q_alpharangle +
langle p_beta^dagger r_alpha^dagger s_alpha q_betarangle + langle p_alpha^dagger r_beta^dagger s_beta q_alpharangle + langle p_beta^dagger r_beta^dagger s_beta q_betarangle`.

Energy should be computed as E = einsum(‘pq,qp’, h1, 1pdm) + 1/2 * einsum(‘pqrs,pqrs’, eri, 2pdm) where h1[p,q] = <p|h|q> and eri[p,q,r,s] = (pq|rs)

pyscf.fci.direct_spin1.make_rdm12s(fcivec, norb, nelec, link_index=None, reorder=True)[source]

Spin separated 1- and 2-particle density matrices. The return values include two lists, a list of 1-particle density matrices and a list of 2-particle density matrices. The density matrices are: (alpha,alpha), (beta,beta) for 1-particle density matrices; (alpha,alpha,alpha,alpha), (alpha,alpha,beta,beta), (beta,beta,beta,beta) for 2-particle density matrices.

1pdm[p,q] = \(\langle q^\dagger p\rangle\); 2pdm[p,q,r,s] = \(\langle p^\dagger r^\dagger s q\rangle\).

Energy should be computed as E = einsum(‘pq,qp’, h1, 1pdm) + 1/2 * einsum(‘pqrs,pqrs’, eri, 2pdm) where h1[p,q] = <p|h|q> and eri[p,q,r,s] = (pq|rs)

pyscf.fci.direct_spin1.make_rdm1s(fcivec, norb, nelec, link_index=None)[source]

Spin separated 1-particle density matrices. The return values include two density matrices: (alpha,alpha), (beta,beta)

dm1[p,q] = <q^dagger p>

The convention is based on McWeeney’s book, Eq (5.4.20). The contraction between 1-particle Hamiltonian and rdm1 is E = einsum(‘pq,qp’, h1, rdm1)

pyscf.fci.direct_spin1.pspace(h1e, eri, norb, nelec, hdiag=None, np=400)[source]

pspace Hamiltonian to improve Davidson preconditioner. See, CPL, 169, 463

pyscf.fci.direct_spin1.trans_rdm1(cibra, ciket, norb, nelec, link_index=None)[source]

Spin traced transition 1-particle transition density matrices.

1pdm[p,q] = :math:`langle q_alpha^dagger p_alpha rangle
  • langle q_beta^dagger p_beta rangle`
pyscf.fci.direct_spin1.trans_rdm12(cibra, ciket, norb, nelec, link_index=None, reorder=True)[source]

Spin traced transition 1- and 2-particle transition density matrices.

1pdm[p,q] = \(\langle q^\dagger p\rangle\); 2pdm[p,q,r,s] = \(\langle p^\dagger r^\dagger s q\rangle\).

pyscf.fci.direct_spin1.trans_rdm12s(cibra, ciket, norb, nelec, link_index=None, reorder=True)[source]

Spin separated 1- and 2-particle transition density matrices. The return values include two lists, a list of 1-particle transition density matrices and a list of 2-particle transition density matrices. The density matrices are: (alpha,alpha), (beta,beta) for 1-particle transition density matrices; (alpha,alpha,alpha,alpha), (alpha,alpha,beta,beta), (beta,beta,alpha,alpha), (beta,beta,beta,beta) for 2-particle transition density matrices.

1pdm[p,q] = \(\langle q^\dagger p\rangle\); 2pdm[p,q,r,s] = \(\langle p^\dagger r^\dagger s q\rangle\).

pyscf.fci.direct_spin1.trans_rdm1s(cibra, ciket, norb, nelec, link_index=None)[source]

Spin separated transition 1-particle density matrices. The return values include two density matrices: (alpha,alpha), (beta,beta). See also function make_rdm1s()

1pdm[p,q] = \(\langle q^\dagger p \rangle\)

Different FCI solvers are implemented to support different type of symmetry.
Symmetry

File Point group Spin singlet Real hermitian* Alpha/beta degeneracy direct_spin0_symm Yes Yes Yes Yes direct_spin1_symm Yes No Yes Yes direct_spin0 No Yes Yes Yes direct_spin1 No No Yes Yes direct_uhf No No Yes No direct_nosym No No No** Yes

  • Real hermitian Hamiltonian implies (ij|kl) = (ji|kl) = (ij|lk) = (ji|lk)

** Hamiltonian is real but not hermitian, (ij|kl) != (ji|kl) ...

FCI solver for Singlet state

Different FCI solvers are implemented to support different type of symmetry.
Symmetry

File Point group Spin singlet Real hermitian* Alpha/beta degeneracy direct_spin0_symm Yes Yes Yes Yes direct_spin1_symm Yes No Yes Yes direct_spin0 No Yes Yes Yes direct_spin1 No No Yes Yes direct_uhf No No Yes No direct_nosym No No No** Yes

  • Real hermitian Hamiltonian implies (ij|kl) = (ji|kl) = (ij|lk) = (ji|lk)

** Hamiltonian is real but not hermitian, (ij|kl) != (ji|kl) ...

pyscf.fci.direct_spin0.contract_1e(f1e, fcivec, norb, nelec, link_index=None)[source]

Contract the 1-electron Hamiltonian with a FCI vector to get a new FCI vector.

pyscf.fci.direct_spin0.contract_2e(eri, fcivec, norb, nelec, link_index=None)[source]

Contract the 2-electron Hamiltonian with a FCI vector to get a new FCI vector.

Note the input arg eri is NOT the 2e hamiltonian matrix, the 2e hamiltonian is

\[\begin{split}h2e &= eri_{pq,rs} p^+ q r^+ s \\ &= (pq|rs) p^+ r^+ s q - (pq|rs) \delta_{qr} p^+ s\end{split}\]

So eri is defined as

\[eri_{pq,rs} = (pq|rs) - (1/Nelec) \sum_q (pq|qs)\]

to restore the symmetry between pq and rs,

\[eri_{pq,rs} = (pq|rs) - (.5/Nelec) [\sum_q (pq|qs) + \sum_p (pq|rp)]\]

See also direct_spin1.absorb_h1e()

pyscf.fci.direct_spin0.make_hdiag(h1e, eri, norb, nelec)[source]

Diagonal Hamiltonian for Davidson preconditioner

pyscf.fci.direct_spin0.make_rdm1(fcivec, norb, nelec, link_index=None)[source]

Spin-traced one-particle density matrix

dm1[p,q] = <q_alpha^dagger p_alpha> + <q_beta^dagger p_beta>

The convention is based on McWeeney’s book, Eq (5.4.20) The contraction between 1-particle Hamiltonian and rdm1 is E = einsum(‘pq,qp’, h1, rdm1)

pyscf.fci.direct_spin0.make_rdm12(fcivec, norb, nelec, link_index=None, reorder=True)[source]

Spin traced 1- and 2-particle density matrices.

1pdm[p,q] = :math:`langle q_alpha^dagger p_alpha rangle +
langle q_beta^dagger p_beta rangle`;
2pdm[p,q,r,s] = :math:`langle p_alpha^dagger r_alpha^dagger s_alpha q_alpharangle +
langle p_beta^dagger r_alpha^dagger s_alpha q_betarangle + langle p_alpha^dagger r_beta^dagger s_beta q_alpharangle + langle p_beta^dagger r_beta^dagger s_beta q_betarangle`.

Energy should be computed as E = einsum(‘pq,qp’, h1, 1pdm) + 1/2 * einsum(‘pqrs,pqrs’, eri, 2pdm) where h1[p,q] = <p|h|q> and eri[p,q,r,s] = (pq|rs)

pyscf.fci.direct_spin0.make_rdm1s(fcivec, norb, nelec, link_index=None)[source]

Spin separated 1-particle density matrices. The return values include two density matrices: (alpha,alpha), (beta,beta)

dm1[p,q] = <q^dagger p>

The convention is based on McWeeney’s book, Eq (5.4.20). The contraction between 1-particle Hamiltonian and rdm1 is E = einsum(‘pq,qp’, h1, rdm1)

pyscf.fci.direct_spin0.trans_rdm1(cibra, ciket, norb, nelec, link_index=None)[source]

Spin traced transition 1-particle transition density matrices.

1pdm[p,q] = :math:`langle q_alpha^dagger p_alpha rangle
  • langle q_beta^dagger p_beta rangle`
pyscf.fci.direct_spin0.trans_rdm12(cibra, ciket, norb, nelec, link_index=None, reorder=True)[source]

Spin traced transition 1- and 2-particle transition density matrices.

1pdm[p,q] = \(\langle q^\dagger p\rangle\); 2pdm[p,q,r,s] = \(\langle p^\dagger r^\dagger s q\rangle\).

pyscf.fci.direct_spin0.trans_rdm1s(cibra, ciket, norb, nelec, link_index=None)[source]

Spin separated transition 1-particle density matrices. The return values include two density matrices: (alpha,alpha), (beta,beta). See also function make_rdm1s()

1pdm[p,q] = \(\langle q^\dagger p \rangle\)

Different FCI solvers are implemented to support different type of symmetry.
Symmetry

File Point group Spin singlet Real hermitian* Alpha/beta degeneracy direct_spin0_symm Yes Yes Yes Yes direct_spin1_symm Yes No Yes Yes direct_spin0 No Yes Yes Yes direct_spin1 No No Yes Yes direct_uhf No No Yes No direct_nosym No No No** Yes

  • Real hermitian Hamiltonian implies (ij|kl) = (ji|kl) = (ij|lk) = (ji|lk)

** Hamiltonian is real but not hermitian, (ij|kl) != (ji|kl) ...

Different FCI solvers are implemented to support different type of symmetry.
Symmetry

File Point group Spin singlet Real hermitian* Alpha/beta degeneracy direct_spin0_symm Yes Yes Yes Yes direct_spin1_symm Yes No Yes Yes direct_spin0 No Yes Yes Yes direct_spin1 No No Yes Yes direct_uhf No No Yes No direct_nosym No No No** Yes

  • Real hermitian Hamiltonian implies (ij|kl) = (ji|kl) = (ij|lk) = (ji|lk)

** Hamiltonian is real but not hermitian, (ij|kl) != (ji|kl) ...

9.2. cistring

pyscf.fci.cistring.addr2str(norb, nelec, addr)[source]

Convert CI determinant address to string

pyscf.fci.cistring.addrs2str(norb, nelec, addrs)[source]

Convert a list of CI determinant address to string

pyscf.fci.cistring.gen_cre_str_index(orb_list, nelec)[source]

linkstr_index to map between N electron string to N+1 electron string. It maps the given string to the address of the string which is generated by the creation operator.

For given string str0, index[str0] is nvir x 4 array. Each entry [i(cre),–,str1,sign] means starting from str0, creating i, to get str1.

pyscf.fci.cistring.gen_cre_str_index_o0(orb_list, nelec)[source]

Slow version of gen_cre_str_index function

pyscf.fci.cistring.gen_cre_str_index_o1(orb_list, nelec)[source]

C implementation of gen_cre_str_index function

pyscf.fci.cistring.gen_des_str_index(orb_list, nelec)[source]

linkstr_index to map between N electron string to N-1 electron string. It maps the given string to the address of the string which is generated by the annihilation operator.

For given string str0, index[str0] is nvir x 4 array. Each entry [–,i(des),str1,sign] means starting from str0, annihilating i, to get str1.

pyscf.fci.cistring.gen_des_str_index_o0(orb_list, nelec)[source]

Slow version of gen_des_str_index function

pyscf.fci.cistring.gen_des_str_index_o1(orb_list, nelec)[source]

C implementation of gen_des_str_index function

pyscf.fci.cistring.gen_linkstr_index(orb_list, nocc, strs=None, tril=False)[source]

Look up table, for the strings relationship in terms of a creation-annihilating operator pair.

For given string str0, index[str0] is (nocc+nocc*nvir) x 4 array. The first nocc rows [i(:occ),i(:occ),str0,sign] are occupied-occupied excitations, which do not change the string. The next nocc*nvir rows [a(:vir),i(:occ),str1,sign] are occupied-virtual exciations, starting from str0, annihilating i, creating a, to get str1.

pyscf.fci.cistring.gen_linkstr_index_trilidx(orb_list, nocc, strs=None)[source]

Generate linkstr_index with the assumption that \(p^+ q|0\rangle\) where \(p > q\). So the resultant link_index has the structure [pq, *, str1, sign]. It is identical to a call to reform_linkstr_index(gen_linkstr_index(...)).

pyscf.fci.cistring.gen_strings4orblist(orb_list, nelec)[source]

Generate string from the given orbital list.

Returns:
list of int64. One int64 element represents one string in binary format. The binary format takes the convention that the one bit stands for one orbital, bit-1 means occupied and bit-0 means unoccupied. The lowest (right-most) bit corresponds to the lowest orbital in the orb_list.

Exampels:

>>> [bin(x) for x in gen_strings4orblist((0,1,2,3),2)]
[0b11, 0b101, 0b110, 0b1001, 0b1010, 0b1100]
>>> [bin(x) for x in gen_strings4orblist((3,1,0,2),2)]
[0b1010, 0b1001, 0b11, 0b1100, 0b110, 0b101]
pyscf.fci.cistring.reform_linkstr_index(link_index)[source]

Compress the (a, i) pair index in linkstr_index to a lower triangular index, to match the 4-fold symmetry of integrals.

pyscf.fci.cistring.str2addr(norb, nelec, string)[source]

Convert string to CI determinant address

pyscf.fci.cistring.strs2addr(norb, nelec, strings)[source]

Convert a list of string to CI determinant address

pyscf.fci.cistring.sub_addrs(norb, nelec, orbital_indices, sub_nelec=0)[source]

The addresses of the determinants which include the specified orbital indices. The size of the returned addresses is equal to the number of determinants of (norb, nelec) system.

pyscf.fci.cistring.tn_strs(norb, nelec, n)[source]

Generate strings for Tn amplitudes. Eg n=1 (T1) has nvir*nocc strings, n=2 (T2) has nvir*(nvir-1)/2 * nocc*(nocc-1)/2 strings.

9.3. spin operator

pyscf.fci.spin_op.contract_ss(fcivec, norb, nelec)[source]

Contract spin square operator with FCI wavefunction \(S^2 |CI>\)

pyscf.fci.spin_op.local_spin(fcivec, norb, nelec, mo_coeff=None, ovlp=1, aolst=[])[source]

Local spin expectation value, which is defined as

<CI|(local S^2)|CI>

The local S^2 operator only couples the orbitals specified in aolst. The cross term which involves the interaction between the local part (in aolst) and non-local part (not in aolst) is not included. As a result, the value of local_spin is not additive. In other words, if local_spin is computed twice with the complementary aolst in the two runs, the summation does not equal to the S^2 of the entire system.

For a complete list of AOs, the value of local_spin is equivalent to <CI|S^2|CI>

pyscf.fci.spin_op.spin_square(fcivec, norb, nelec, mo_coeff=None, ovlp=1)[source]

General spin square operator.

... math:

<CI|S_+*S_-|CI> &= n_\alpha + \delta_{ik}\delta_{jl}Gamma_{i\alpha k\beta ,j\beta l\alpha } \\
<CI|S_-*S_+|CI> &= n_\beta + \delta_{ik}\delta_{jl}Gamma_{i\beta k\alpha ,j\alpha l\beta } \\
<CI|S_z*S_z|CI> &= \delta_{ik}\delta_{jl}(Gamma_{i\alpha k\alpha ,j\alpha l\alpha }
                 - Gamma_{i\alpha k\alpha ,j\beta l\beta }
                 - Gamma_{i\beta k\beta ,j\alpha l\alpha}
                 + Gamma_{i\beta k\beta ,j\beta l\beta})
                 + (n_\alpha+n_\beta)/4

Given the overlap betwen non-degenerate alpha and beta orbitals, this function can compute the expectation value spin square operator for UHF-FCI wavefunction

pyscf.fci.spin_op.spin_square0(fcivec, norb, nelec)[source]

Spin square for RHF-FCI CI wfn only (obtained from spin-degenerated Hamiltonian)

9.4. rdm

FCI 1, 2, 3, 4-particle density matrices.

Note the 1-particle density matrix has the same convention as the mean-field 1-particle density matrix (see McWeeney’s book Eq 5.4.20), which is

dm[p,q] = < q^+ p >
The contraction between 1-particle Hamiltonian and 1-pdm is
E = einsum(‘pq,qp’, h1, 1pdm)
Different conventions are used in the high order density matrices:
dm[p,q,r,s,...] = < p^+ r^+ ... s q >
pyscf.fci.rdm.make_dm123(fname, cibra, ciket, norb, nelec)[source]

Spin traced 1, 2 and 3-particle density matrices.

Note

In this function, 2pdm[p,q,r,s] is \(\langle p^\dagger q r^\dagger s\rangle\); 3pdm[p,q,r,s,t,u] is \(\langle p^\dagger q r^\dagger s t^\dagger u\rangle\).

After calling reorder_dm123, the 2pdm and 3pdm are transformed to the normal density matrices: 2pdm[p,r,q,s] = \(\langle p^\dagger q^\dagger s r\rangle\) 3pdm[p,s,q,t,r,u] = \(\langle p^\dagger q^\dagger r^\dagger u t s\rangle\).

pyscf.fci.rdm.make_dm1234(fname, cibra, ciket, norb, nelec)[source]

Spin traced 1, 2, 3 and 4-particle density matrices.

Note

In this function, 2pdm[p,q,r,s] is \(\langle p^\dagger q r^\dagger s\rangle\); 3pdm[p,q,r,s,t,u] is \(\langle p^\dagger q r^\dagger s t^\dagger u\rangle\); 4pdm[p,q,r,s,t,u,v,w] is \(\langle p^\dagger q r^\dagger s t^\dagger u v^\dagger w\rangle\).

After calling reorder_dm123, the 2pdm and 3pdm are transformed to the normal density matrices: 2pdm[p,r,q,s] = \(\langle p^\dagger q^\dagger s r\rangle\) 3pdm[p,s,q,t,r,u] = \(\langle p^\dagger q^\dagger r^\dagger u t s\rangle\). 4pdm[p,t,q,u,r,v,s,w] = \(\langle p^\dagger q^\dagger r^\dagger s^dagger w v u t\rangle\).

9.5. addons

pyscf.fci.addons.cre_a(ci0, norb, neleca_nelecb, ap_id)[source]

Construct (N+1)-electron wavefunction by adding an alpha electron in the N-electron wavefunction.

... math:

|N+1\rangle = \hat{a}^+_p |N\rangle
Args:
ci0
: 2D array
CI coefficients, row for alpha strings and column for beta strings.
norb
: int
Number of orbitals.
(neleca,nelecb)
: (int,int)
Number of (alpha, beta) electrons of the input CI function
ap_id
: int
Orbital index (0-based), for the creation operator
Returns:
2D array, row for alpha strings and column for beta strings. Note it has different number of rows to the input CI coefficients.
pyscf.fci.addons.cre_b(ci0, norb, neleca_nelecb, ap_id)[source]

Construct (N+1)-electron wavefunction by adding a beta electron in the N-electron wavefunction.

Args:
ci0
: 2D array
CI coefficients, row for alpha strings and column for beta strings.
norb
: int
Number of orbitals.
(neleca,nelecb)
: (int,int)
Number of (alpha, beta) electrons of the input CI function
ap_id
: int
Orbital index (0-based), for the creation operator
Returns:
2D array, row for alpha strings and column for beta strings. Note it has different number of columns to the input CI coefficients.
pyscf.fci.addons.cylindrical_init_guess(mol, norb, nelec, orbsym, wfnsym=0, singlet=True, nroots=1)[source]

FCI initial guess for system of cylindrical symmetry. (In testing)

Examples:

>>> mol = gto.M(atom='O; O 1 1.2', spin=2, symmetry=True)
>>> orbsym = [6,7,2,3]
>>> ci0 = fci.addons.cylindrical_init_guess(mol, 4, (3,3), orbsym, wfnsym=10)[0]
>>> print(ci0.reshape(4,4))
>>> ci0 = fci.addons.cylindrical_init_guess(mol, 4, (3,3), orbsym, wfnsym=10, singlet=False)[0]
>>> print(ci0.reshape(4,4))
pyscf.fci.addons.des_a(ci0, norb, neleca_nelecb, ap_id)[source]

Construct (N-1)-electron wavefunction by removing an alpha electron from the N-electron wavefunction.

... math:

|N-1\rangle = \hat{a}_p |N\rangle
Args:
ci0
: 2D array
CI coefficients, row for alpha strings and column for beta strings.
norb
: int
Number of orbitals.
(neleca,nelecb)
: (int,int)
Number of (alpha, beta) electrons of the input CI function
ap_id
: int
Orbital index (0-based), for the annihilation operator
Returns:
2D array, row for alpha strings and column for beta strings. Note it has different number of rows to the input CI coefficients
pyscf.fci.addons.des_b(ci0, norb, neleca_nelecb, ap_id)[source]

Construct (N-1)-electron wavefunction by removing a beta electron from N-electron wavefunction.

Args:
ci0
: 2D array
CI coefficients, row for alpha strings and column for beta strings.
norb
: int
Number of orbitals.
(neleca,nelecb)
: (int,int)
Number of (alpha, beta) electrons of the input CI function
ap_id
: int
Orbital index (0-based), for the annihilation operator
Returns:
2D array, row for alpha strings and column for beta strings. Note it has different number of columns to the input CI coefficients.
pyscf.fci.addons.det_overlap(string1, string2, norb, s=None)[source]

Determinants overlap on non-orthogonal one-particle basis

pyscf.fci.addons.energy(h1e, eri, fcivec, norb, nelec, link_index=None)[source]

Compute the FCI electronic energy for given Hamiltonian and FCI vector.

pyscf.fci.addons.fix_spin_(fciobj, shift=0.2, ss=None, **kwargs)[source]

If FCI solver cannot stay on spin eigenfunction, this function can add a shift to the states which have wrong spin.

\[(H + shift*S^2) |\Psi\rangle = E |\Psi\rangle\]
Args:
fciobj : An instance of FCISolver
Kwargs:
shift
: float
Level shift for states which have different spin
ss
: number
S^2 expection value == s*(s+1)
Returns
A modified FCI object based on fciobj.
pyscf.fci.addons.guess_wfnsym(ci, norb, nelec, orbsym)[source]

Guess the wavefunction symmetry based on the non-zero elements in the given CI coefficients.

Args:
ci
: 2D array
CI coefficients, row for alpha strings and column for beta strings.
norb
: int
Number of orbitals.
nelec
: int or 2-item list
Number of electrons, or 2-item list for (alpha, beta) electrons
orbsym
: list of int
The irrep ID for each orbital.
Returns:
Irrep ID
pyscf.fci.addons.initguess_triplet(norb, nelec, binstring)[source]

Generate a triplet initial guess for FCI solver

pyscf.fci.addons.large_ci(ci, norb, nelec, tol=0.1, return_strs=True)[source]

Search for the largest CI coefficients

pyscf.fci.addons.overlap(bra, ket, norb, nelec, s=None)[source]

Overlap between two CI wavefunctions

Args:
s
: 2D array or a list of 2D array
The overlap matrix of non-orthogonal one-particle basis
pyscf.fci.addons.reorder(ci, nelec, orbidxa, orbidxb=None)[source]

Reorder the CI coefficients, to adapt the reordered orbitals (The relation of the reordered orbitals and original orbitals is new = old[idx]). Eg.

The orbital ordering indices orbidx = [2,0,1] indicates the map old orbital a b c -> new orbital c a b. The strings are reordered as old-strings 0b011, 0b101, 0b110 == (1,2), (1,3), (2,3) <= apply orbidx to get orb-strings orb-strings (3,1), (3,2), (1,2) == 0B101, 0B110, 0B011 <= by gen_strings4orblist then argsort to translate the string representation to the address [2(=0B011), 0(=0B101), 1(=0B110)]

pyscf.fci.addons.symm_initguess(norb, nelec, orbsym, wfnsym=0, irrep_nelec=None)[source]

Generate CI wavefunction initial guess which has the given symmetry.

Args:
norb
: int
Number of orbitals.
nelec
: int or 2-item list
Number of electrons, or 2-item list for (alpha, beta) electrons
orbsym
: list of int
The irrep ID for each orbital.
Kwags:
wfnsym
: int
The irrep ID of target symmetry
irrep_nelec
: dict
Freeze occupancy for certain irreps
Returns:
CI coefficients 2D array which has the target symmetry.
pyscf.fci.addons.symmetrize_wfn(ci, norb, nelec, orbsym, wfnsym=0)[source]

Symmetrize the CI wavefunction by zeroing out the determinants which do not have the right symmetry.

Args:
ci
: 2D array
CI coefficients, row for alpha strings and column for beta strings.
norb
: int
Number of orbitals.
nelec
: int or 2-item list
Number of electrons, or 2-item list for (alpha, beta) electrons
orbsym
: list of int
The irrep ID for each orbital.
Kwags:
wfnsym
: int
The irrep ID of target symmetry
Returns:
2D array which is the symmetrized CI coefficients
pyscf.fci.addons.transform_ci_for_orbital_rotation(ci, norb, nelec, u)[source]

Transform CI coefficients to the representation in new one-particle basis. Solving CI problem for Hamiltonian h1, h2 defined in old basis, CI_old = fci.kernel(h1, h2, ...) Given orbital rotation u, the CI problem can be either solved by transforming the Hamiltonian, or transforming the coefficients. CI_new = fci.kernel(u^T*h1*u, ...) = transform_ci_for_orbital_rotation(CI_old, u)

Args:
u
: 2D array or a list of 2D array
the orbital rotation to transform the old one-particle basis to new one-particle basis