pyscf.mp package

Submodules

pyscf.mp.dfgmp2 module

density fitting GMP2, 3-center integrals incore.

class pyscf.mp.dfgmp2.DFGMP2(mf, frozen=None, mo_coeff=None, mo_occ=None)[source]

Bases: DFMP2

ao2mo(mo_coeff=None)[source]
init_amps(mo_energy=None, mo_coeff=None, eris=None, with_t2=True)[source]
loop_ao2mo(mo_coeff, nocc, orbspin)[source]
make_rdm1(t2=None, ao_repr=False)[source]

Spin-traced one-particle density matrix. The occupied-virtual orbital response is not included.

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

The convention of 1-pdm 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)

Kwargs:
ao_reprboolean

Whether to transfrom 1-particle density matrix to AO representation.

make_rdm2(t2=None, ao_repr=False)[source]

Spin-traced two-particle density matrix in MO basis

dm2[p,q,r,s] = sum_{sigma,tau} <p_sigma^dagger r_tau^dagger s_tau q_sigma>

Note the contraction between ERIs (in Chemist’s notation) and rdm2 is E = einsum(‘pqrs,pqrs’, eri, rdm2)

pyscf.mp.dfgmp2.kernel(mp, mo_energy=None, mo_coeff=None, eris=None, with_t2=True, verbose=3)[source]

pyscf.mp.dfmp2 module

density fitting MP2, 3-center integrals incore.

class pyscf.mp.dfmp2.DFMP2(mf, frozen=None, mo_coeff=None, mo_occ=None)[source]

Bases: MP2

ao2mo(mo_coeff=None)[source]
init_amps(mo_energy=None, mo_coeff=None, eris=None, with_t2=True)[source]
loop_ao2mo(mo_coeff, nocc)[source]
make_rdm1(t2=None, ao_repr=False)[source]

Spin-traced one-particle density matrix. The occupied-virtual orbital response is not included.

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

The convention of 1-pdm 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)

Kwargs:
ao_reprboolean

Whether to transfrom 1-particle density matrix to AO representation.

make_rdm2(t2=None, ao_repr=False)[source]

Spin-traced two-particle density matrix in MO basis

dm2[p,q,r,s] = sum_{sigma,tau} <p_sigma^dagger r_tau^dagger s_tau q_sigma>

Note the contraction between ERIs (in Chemist’s notation) and rdm2 is E = einsum(‘pqrs,pqrs’, eri, rdm2)

nuc_grad_method()[source]
reset(mol=None)[source]
update_amps(t2, eris)[source]

Update non-canonical MP2 amplitudes

pyscf.mp.dfmp2.MP2

alias of DFMP2

pyscf.mp.dfmp2.kernel(mp, mo_energy=None, mo_coeff=None, eris=None, with_t2=True, verbose=3)[source]

pyscf.mp.dfmp2_native module

native implementation of DF-MP2/RI-MP2 with an RHF reference

exception pyscf.mp.dfmp2_native.BatchSizeError[source]

Bases: Exception

pyscf.mp.dfmp2_native.DFMP2

alias of DFRMP2

class pyscf.mp.dfmp2_native.DFRMP2(mf, frozen=None, auxbasis=None)[source]

Bases: StreamObject

native implementation of DF-MP2/RI-MP2 with an RHF reference

calculate_energy()[source]

Calculates the MP2 correlation energy.

calculate_integrals_()[source]

Calculates the three center integrals for MP2.

delete()[source]

Delete the temporary file(s).

dump_flags(logger=None)[source]

Prints selected information.

Args:

logger : Logger object

property e_tot

total energy (SCF + MP2)

property has_ints
kernel()[source]

Alias for the MP2 energy calculation. Does not need to be called to calculate the 1-RDM only.

make_natorbs(rdm1_mo=None, relaxed=False)[source]

Calculate natural orbitals. Note: the most occupied orbitals come first (left)

and the least occupied orbitals last (right).

Args:
rdm1_mo1-RDM in MO basis

the function calculates a density matrix if none is provided

relaxed : calculated relaxed or unrelaxed density matrix

Returns:

natural occupation numbers, natural orbitals

make_rdm1(relaxed=False, ao_repr=False)[source]

Calculates the MP2 1-RDM. - The relaxed density matrix can be used to calculate properties of systems

for which MP2 is well-behaved.

  • The unrelaxed density is less suited to calculate properties accurately, but it can be used to calculate CASSCF starting orbitals.

Args:

relaxed : relaxed density if True, unrelaxed density if False ao_repr : density in AO or in MO basis

Returns:

the 1-RDM

make_rdm1_relaxed(ao_repr=False)[source]
make_rdm1_unrelaxed(ao_repr=False)[source]
nuc_grad_method()[source]
pyscf.mp.dfmp2_native.MP2

alias of DFRMP2

pyscf.mp.dfmp2_native.RMP2

alias of DFRMP2

pyscf.mp.dfmp2_native.SCSDFMP2

alias of SCSDFRMP2

class pyscf.mp.dfmp2_native.SCSDFRMP2(mf, ps=1.2, pt=0.3333333333333333, *args, **kwargs)[source]

Bases: DFRMP2

RHF-DF-MP2 with spin-component scaling S. Grimme, J. Chem. Phys. 118 (2003), 9095 https://doi.org/10.1063/1.1569242

dump_flags(logger=None)[source]

Prints selected information.

Args:

logger : Logger object

pyscf.mp.dfmp2_native.SCSMP2

alias of SCSDFRMP2

pyscf.mp.dfmp2_native.SCSRMP2

alias of SCSDFRMP2

pyscf.mp.dfmp2_native.emp2_rhf(intsfile, mo_energy, frozen_mask, logger, ps=1.0, pt=1.0)[source]

Calculates the DF-MP2 energy with an RHF reference.

Args:

intsfile : contains the three center integrals in MO basis mo_energy : energies of the molecular orbitals frozen_mask : boolean mask for frozen orbitals logger : Logger instance ps : SCS factor for opposite-spin contributions pt : SCS factor for same-spin contributions

Returns:

the MP2 correlation energy

pyscf.mp.dfmp2_native.fock_response_rhf(mf, dm, full=True)[source]

Calculate the Fock response function for a given density matrix: sum_pq [ 4 (ai|pq) - (ap|iq) - (aq|ip) ] dm[p, q]

Args:

mf : RHF instance dm : density matrix in MO basis full : full MO density matrix if True, virtual x occupied if False

Returns:

Fock response in MO basis. Shape: virtual x occupied.

pyscf.mp.dfmp2_native.ints3c_cholesky(mol, auxmol, mo_coeff1, mo_coeff2, max_memory, logger)[source]

Calculate the three center electron repulsion integrals in MO basis multiplied with the Cholesky factor of the inverse Coulomb metric matrix. Only integrals in MO basis are stored on disk; integral-direct with regard to AO integrals.

Args:

mol : Mole instance auxmol : Mole instance with auxiliary basis mo_coeff1 : MO coefficient matrix for the leading MO index, typically occupied mo_coeff2 : MO coefficient matrix for the secondary MO index, typically virtual max_memory : memory threshold in MB logger : Logger instance

Returns:

A HDF5 temporary file containing the integrals in the dataset “ints_cholesky”. Indexing order: [mo1, aux, mo2]

pyscf.mp.dfmp2_native.make_rdm1(mp2, relaxed, logger=None)[source]

Calculates the unrelaxed or relaxed MP2 density matrix.

Args:

mp2 : DFRMP2 instance relaxed : relaxed density if True, unrelaxed density if False logger : Logger instance

Returns:

the 1-RDM in MO basis

pyscf.mp.dfmp2_native.orbgrad_from_Gamma(mol, auxmol, Gamma, mo_coeff, frozen_mask, max_memory, logger)[source]

Calculates the orbital gradient of the two-electron term in the Hylleraas functional.

Args:

mol : Mole object auxmol : Mole object for the auxiliary functions Gamma : h5py dataset with the 3c2e density, order: [occ. orbs., aux. fcns., virt. orbs.] mo_coeff : molecular orbital coefficients frozen_mask : boolean mask for frozen orbitals max_memory : memory limit in MB logger : Logger object

Returns:

orbital gradient in shape: virt. orbitals x occ. orbitals, orbital gradient in shape: froz. orbitals x occ. orbitals

pyscf.mp.dfmp2_native.rmp2_densities_contribs(intsfile, mo_energy, frozen_mask, max_memory, logger, calcGamma=False, auxmol=None, ps=1.0, pt=1.0)[source]

Calculates the unrelaxed DF-MP2 density matrix contribution with an RHF reference. Note: this is the difference density, i.e. without HF contribution. Also calculates the three-center two-particle density if requested.

Args:

intsfile : contains the three center integrals mo_energy : molecular orbital energies frozen_mask : boolean mask for frozen orbitals max_memory : memory threshold in MB logger : Logger instance calcGamma : if True, calculate 3c2e density auxmol : required if relaxed is True ps : SCS factor for opposite-spin contributions pt : SCS factor for same-spin contributions

Returns:

matrix containing the 1-RDM contribution, file with 3c2e density if requested

pyscf.mp.dfmp2_native.shellBatchGenerator(mol, nao_max)[source]

Generates sets of shells with a limited number of functions.

Args:

mol : the molecule object nao_max : maximum number of AOs in each set

Returns:

generator yields ((first shell, last shell+1), (first AO, last AO+1))

pyscf.mp.dfmp2_native.solve_cphf_rhf(mf, Lvo, max_cycle, tol, logger)[source]

Solve the CPHF equations. (e[i] - e[a]) zvo[a, i] - sum_bj [ 4 (ai|bj) - (ab|ij) - (aj|ib) ] zvo[b, j] = Lvo[a, i]

Args:

mf : an RHF object Lvo : right-hand side the the response equation max_cycle : number of iterations for the CPHF solver tol : convergence tolerance for the CPHF solver logger : Logger object

pyscf.mp.dfump2_native module

native implementation of DF-MP2/RI-MP2 with a UHF reference

pyscf.mp.dfump2_native.DFMP2

alias of DFUMP2

class pyscf.mp.dfump2_native.DFUMP2(mf, frozen=None, auxbasis=None)[source]

Bases: DFRMP2

native implementation of DF-MP2/RI-MP2 with a UHF reference

calculate_energy()[source]

Calculates the MP2 correlation energy.

calculate_integrals_()[source]

Calculates the three center integrals for MP2.

delete()[source]

Delete the temporary file(s).

dump_flags(logger=None)[source]

Prints selected information.

Args:

logger : Logger object

make_natorbs(rdm1_mo=None, relaxed=False)[source]

Calculate natural orbitals. Note: the most occupied orbitals come first (left)

and the least occupied orbitals last (right).

Args:
rdm1_mo1-RDM in MO basis

the function calculates a density matrix if none is provided

relaxed : calculated relaxed or unrelaxed density matrix

Returns:

natural occupation numbers, natural orbitals

make_rdm1(relaxed=False, ao_repr=False)[source]

Calculates the MP2 1-RDM. - The relaxed density matrix can be used to calculate properties of systems

for which MP2 is well-behaved.

  • The unrelaxed density is less suited to calculate properties accurately, but it can be used to calculate CASSCF starting orbitals.

Args:

relaxed : relaxed density if True, unrelaxed density if False ao_repr : density in AO or in MO basis

Returns:

the 1-RDM

nuc_grad_method()[source]
pyscf.mp.dfump2_native.MP2

alias of DFUMP2

pyscf.mp.dfump2_native.SCSDFMP2

alias of SCSDFUMP2

class pyscf.mp.dfump2_native.SCSDFUMP2(mf, ps=1.2, pt=0.3333333333333333, *args, **kwargs)[source]

Bases: DFUMP2

UHF-DF-MP2 with spin-component scaling S. Grimme, J. Chem. Phys. 118 (2003), 9095 https://doi.org/10.1063/1.1569242

dump_flags(logger=None)[source]

Prints selected information.

Args:

logger : Logger object

pyscf.mp.dfump2_native.SCSMP2

alias of SCSDFUMP2

pyscf.mp.dfump2_native.SCSUMP2

alias of SCSDFUMP2

pyscf.mp.dfump2_native.UMP2

alias of DFUMP2

pyscf.mp.dfump2_native.emp2_uhf(intsfiles, mo_energy, frozen_mask, logger, ps=1.0, pt=1.0)[source]

Calculates the DF-MP2 energy with an UHF reference.

Args:

intsfiles : contains the three center integrals in MO basis mo_energy : energies of the molecular orbitals frozen_mask : boolean mask for frozen orbitals logger : Logger instance ps : SCS factor for opposite-spin contributions pt : SCS factor for same-spin contributions

Returns:

the MP2 correlation energy

pyscf.mp.dfump2_native.fock_response_uhf(mf, dm, full=True)[source]

Calculate the unrestricted Fock response function for a given density matrix.

Args:

mf : UHF instance dm : density matrix in MO basis full : full MO density matrix if True, [virt. x occ., virt. x occ.] if False

Returns:

Fock response in MO basis. Shape: [virt. x occ., virt. x occ.]

pyscf.mp.dfump2_native.make_rdm1(mp2, relaxed, logger=None)[source]

Calculates the unrelaxed or relaxed MP2 density matrix.

Args:

mp2 : DFUMP2 instance relaxed : relaxed density if True, unrelaxed density if False logger : Logger instance

Returns:

the 1-RDM in MO basis

pyscf.mp.dfump2_native.solve_cphf_uhf(mf, Lvo, max_cycle, tol, logger)[source]

Solve the CPHF equations.

Args:

mf : a UHF object Lvo : right-hand side the the response equation max_cycle : number of iterations for the CPHF solver tol : convergence tolerance for the CPHF solver logger : Logger object

pyscf.mp.dfump2_native.ump2_densities_contribs(intsfiles, mo_energy, frozen_mask, max_memory, logger, calcGamma=False, auxmol=None, ps=1.0, pt=1.0)[source]

Calculates the unrelaxed DF-MP2 density matrix contribution with a UHF reference. Note: this is the difference density, i.e. without HF contribution.A lso calculates the three-center two-particle density if requested.

Args:

intsfile : contains the three center integrals mo_energy : molecular orbital energies frozen_mask : boolean mask for frozen orbitals max_memory : memory threshold in MB logger : Logger instance calcGamma : if True, calculate 3c2e density auxmol : required if relaxed is True ps : SCS factor for opposite-spin contributions pt : SCS factor for same-spin contributions

Returns:

matrix containing the 1-RDM contribution, file with 3c2e density if requested

pyscf.mp.gmp2 module

GMP2 in spin-orbital form E(MP2) = 1/4 <ij||ab><ab||ij>/(ei+ej-ea-eb)

class pyscf.mp.gmp2.GMP2(mf, frozen=None, mo_coeff=None, mo_occ=None)[source]

Bases: MP2

ao2mo(mo_coeff=None)[source]
density_fit(auxbasis=None, with_df=None)[source]
energy(t2, eris)

MP2 energy

init_amps(mo_energy=None, mo_coeff=None, eris=None, with_t2=True)[source]
make_rdm1(t2=None, ao_repr=False)

One-particle density matrix in the molecular spin-orbital representation (the occupied-virtual blocks from the orbital response contribution are not included).

dm1[p,q] = <q^dagger p> (p,q are spin-orbitals)

The convention of 1-pdm 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(t2=None, ao_repr=False)

Two-particle density matrix in the molecular spin-orbital representation

dm2[p,q,r,s] = <p^dagger r^dagger s q>

where p,q,r,s are spin-orbitals. p,q correspond to one particle and r,s correspond to another particle. The contraction between ERIs (in Chemist’s notation) and rdm2 is E = einsum(‘pqrs,pqrs’, eri, rdm2)

nuc_grad_method()[source]
update_amps(t2, eris)

Update non-canonical MP2 amplitudes

pyscf.mp.gmp2.MP2

alias of GMP2

pyscf.mp.gmp2.energy(mp, t2, eris)[source]

MP2 energy

pyscf.mp.gmp2.kernel(mp, mo_energy=None, mo_coeff=None, eris=None, with_t2=True, verbose=None)[source]
pyscf.mp.gmp2.make_rdm1(mp, t2=None, ao_repr=False)[source]

One-particle density matrix in the molecular spin-orbital representation (the occupied-virtual blocks from the orbital response contribution are not included).

dm1[p,q] = <q^dagger p> (p,q are spin-orbitals)

The convention of 1-pdm 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.mp.gmp2.make_rdm2(mp, t2=None, ao_repr=False)[source]

Two-particle density matrix in the molecular spin-orbital representation

dm2[p,q,r,s] = <p^dagger r^dagger s q>

where p,q,r,s are spin-orbitals. p,q correspond to one particle and r,s correspond to another particle. The contraction between ERIs (in Chemist’s notation) and rdm2 is E = einsum(‘pqrs,pqrs’, eri, rdm2)

pyscf.mp.gmp2.update_amps(mp, t2, eris)[source]

Update non-canonical MP2 amplitudes

pyscf.mp.mp2 module

RMP2

class pyscf.mp.mp2.MP2(mf, frozen=None, mo_coeff=None, mo_occ=None)[source]

Bases: StreamObject

restricted MP2 with canonical HF and non-canonical HF reference

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default value equals to Mole.max_memory

conv_tolfloat

For non-canonical MP2, converge threshold for MP2 correlation energy. Default value is 1e-7.

conv_tol_normtfloat

For non-canonical MP2, converge threshold for norm(t2). Default value is 1e-5.

max_cycleint

For non-canonical MP2, max number of MP2 iterations. Default value is 50.

diis_spaceint

For non-canonical MP2, DIIS space size in MP2 iterations. Default is 6.

level_shiftfloat

A shift on virtual orbital energies to stablize the MP2 iterations.

frozenint or list

If integer is given, the inner-most orbitals are excluded from MP2 amplitudes. Given the orbital indices (0-based) in a list, both occupied and virtual orbitals can be frozen in MP2 calculation.

>>> mol = gto.M(atom = 'H 0 0 0; F 0 0 1.1', basis = 'ccpvdz')
>>> mf = scf.RHF(mol).run()
>>> # freeze 2 core orbitals
>>> pt = mp.MP2(mf).set(frozen = 2).run()
>>> # auto-generate the number of core orbitals to be frozen (1 in this case)
>>> pt = mp.MP2(mf).set_frozen().run()
>>> # freeze 2 core orbitals and 3 high lying unoccupied orbitals
>>> pt.set(frozen = [0,1,16,17,18]).run()

Saved results

e_corrfloat

MP2 correlation correction

e_corr_ss/osfloat

Same-spin and opposite-spin component of the MP2 correlation energy

e_totfloat

Total MP2 energy (HF + correlation)

t2 :

T amplitudes t2[i,j,a,b] (i,j in occ, a,b in virt)

Gradients(*args, **kwargs)
ao2mo(mo_coeff=None)[source]
as_scanner()

Generating a scanner/solver for MP2 PES.

The returned solver is a function. This function requires one argument “mol” as input and returns total MP2 energy.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the MP2 and the underlying SCF objects (conv_tol, max_memory etc) are automatically applied in the solver.

Note scanner has side effects. It may change many underlying objects (_scf, with_df, with_x2c, …) during calculation.

Examples:

>>> from pyscf import gto, scf, mp
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> mp2_scanner = mp.MP2(scf.RHF(mol)).as_scanner()
>>> e_tot = mp2_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot = mp2_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
conv_tol = 1e-07
conv_tol_normt = 1e-05
density_fit(auxbasis=None, with_df=None)[source]
dump_flags(verbose=None)[source]
property e_tot
property e_tot_scs
property emp2
property emp2_scs
energy(t2, eris)

MP2 energy

get_e_hf(mo_coeff=None)
get_frozen_mask()

Get boolean mask for the restricted reference orbitals.

In the returned boolean (mask) array of frozen orbital indices, the element is False if it corresonds to the frozen orbital.

get_nmo()
get_nocc()
init_amps(mo_energy=None, mo_coeff=None, eris=None, with_t2=True)[source]
kernel(mo_energy=None, mo_coeff=None, eris=None, with_t2=True)[source]
Args:
with_t2bool

Whether to generate and hold t2 amplitudes in memory.

make_fno(thresh=1e-06, pct_occ=None, nvir_act=None, t2=None)

Frozen natural orbitals

Attributes:
threshfloat

Threshold on NO occupation numbers. Default is 1e-6 (very conservative).

pct_occfloat

Percentage of total occupation number. Default is None. If present, overrides thresh.

nvir_actint

Number of virtual NOs to keep. Default is None. If present, overrides thresh and pct_occ.

Returns:
frozenlist or ndarray

List of orbitals to freeze

no_coeffndarray

Semicanonical NO coefficients in the AO basis

make_rdm1(t2=None, eris=None, ao_repr=False)

Spin-traced one-particle density matrix. The occupied-virtual orbital response is not included.

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

The convention of 1-pdm 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)

Kwargs:
ao_reprboolean

Whether to transfrom 1-particle density matrix to AO representation.

make_rdm2(t2=None, eris=None, ao_repr=False)

Spin-traced two-particle density matrix in MO basis

dm2[p,q,r,s] = sum_{sigma,tau} <p_sigma^dagger r_tau^dagger s_tau q_sigma>

Note the contraction between ERIs (in Chemist’s notation) and rdm2 is E = einsum(‘pqrs,pqrs’, eri, rdm2)

max_cycle = 50
property nmo
property nocc
nuc_grad_method()[source]
reset(mol=None)[source]
set_frozen(method='auto', window=(-1000.0, 1000.0))[source]
update_amps(t2, eris)

Update non-canonical MP2 amplitudes

class pyscf.mp.mp2.MP2_Scanner(mp)[source]

Bases: SinglePointScanner

pyscf.mp.mp2.RMP2

alias of MP2

pyscf.mp.mp2.as_scanner(mp)[source]

Generating a scanner/solver for MP2 PES.

The returned solver is a function. This function requires one argument “mol” as input and returns total MP2 energy.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the MP2 and the underlying SCF objects (conv_tol, max_memory etc) are automatically applied in the solver.

Note scanner has side effects. It may change many underlying objects (_scf, with_df, with_x2c, …) during calculation.

Examples:

>>> from pyscf import gto, scf, mp
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> mp2_scanner = mp.MP2(scf.RHF(mol)).as_scanner()
>>> e_tot = mp2_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot = mp2_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
pyscf.mp.mp2.energy(mp, t2, eris)[source]

MP2 energy

pyscf.mp.mp2.get_e_hf(mp, mo_coeff=None)[source]
pyscf.mp.mp2.get_frozen_mask(mp)[source]

Get boolean mask for the restricted reference orbitals.

In the returned boolean (mask) array of frozen orbital indices, the element is False if it corresonds to the frozen orbital.

pyscf.mp.mp2.get_nmo(mp)[source]
pyscf.mp.mp2.get_nocc(mp)[source]
pyscf.mp.mp2.kernel(mp, mo_energy=None, mo_coeff=None, eris=None, with_t2=True, verbose=None)[source]
pyscf.mp.mp2.make_fno(mp, thresh=1e-06, pct_occ=None, nvir_act=None, t2=None)[source]

Frozen natural orbitals

Attributes:
threshfloat

Threshold on NO occupation numbers. Default is 1e-6 (very conservative).

pct_occfloat

Percentage of total occupation number. Default is None. If present, overrides thresh.

nvir_actint

Number of virtual NOs to keep. Default is None. If present, overrides thresh and pct_occ.

Returns:
frozenlist or ndarray

List of orbitals to freeze

no_coeffndarray

Semicanonical NO coefficients in the AO basis

pyscf.mp.mp2.make_rdm1(mp, t2=None, eris=None, ao_repr=False)[source]

Spin-traced one-particle density matrix. The occupied-virtual orbital response is not included.

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

The convention of 1-pdm 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)

Kwargs:
ao_reprboolean

Whether to transfrom 1-particle density matrix to AO representation.

pyscf.mp.mp2.make_rdm2(mp, t2=None, eris=None, ao_repr=False)[source]

Spin-traced two-particle density matrix in MO basis

dm2[p,q,r,s] = sum_{sigma,tau} <p_sigma^dagger r_tau^dagger s_tau q_sigma>

Note the contraction between ERIs (in Chemist’s notation) and rdm2 is E = einsum(‘pqrs,pqrs’, eri, rdm2)

pyscf.mp.mp2.update_amps(mp, t2, eris)[source]

Update non-canonical MP2 amplitudes

pyscf.mp.mp2f12_slow module

MP2-F12 (In testing)

Refs: * JCC 32, 2492 (2011); DOI:10.1002/jcc.21825 * JCP 139, 084112 (2013); DOI:10.1063/1.4818753

With strong orthogonalization ansatz 2

pyscf.mp.mp2f12_slow.energy_f12(mf, auxmol, zeta)[source]
pyscf.mp.mp2f12_slow.find_cabs(mol, auxmol, lindep=1e-08)[source]
pyscf.mp.mp2f12_slow.trans(eri, mos)[source]

pyscf.mp.ump2 module

UMP2 with spatial integals

pyscf.mp.ump2.MP2

alias of UMP2

class pyscf.mp.ump2.UMP2(mf, frozen=None, mo_coeff=None, mo_occ=None)[source]

Bases: MP2

Gradients(*args, **kwargs)
ao2mo(mo_coeff=None)[source]
energy(t2, eris)

MP2 energy

get_frozen_mask()

Get boolean mask for the unrestricted reference orbitals.

In the returned boolean (mask) array of frozen orbital indices, the element is False if it corresonds to the frozen orbital.

get_nmo()
get_nocc()
init_amps(mo_energy=None, mo_coeff=None, eris=None, with_t2=True)[source]
make_fno(thresh=1e-06, pct_occ=None, t2=None, eris=None)

Frozen natural orbitals

Returns:
frozenlist or ndarray

Length-2 list of orbitals to freeze

no_coeffndarray

Length-2 list of semicanonical NO coefficients in the AO basis

make_rdm1(t2=None, ao_repr=False)

One-particle spin density matrices dm1a, dm1b in MO basis (the occupied-virtual blocks due to the orbital response contribution are not included).

dm1a[p,q] = <q_alpha^dagger p_alpha> dm1b[p,q] = <q_beta^dagger p_beta>

The convention of 1-pdm is based on McWeeney’s book, Eq (5.4.20).

make_rdm2(t2=None, ao_repr=False)

Two-particle spin density matrices dm2aa, dm2ab, dm2bb in MO basis

dm2aa[p,q,r,s] = <q_alpha^dagger s_alpha^dagger r_alpha p_alpha> dm2ab[p,q,r,s] = <q_alpha^dagger s_beta^dagger r_beta p_alpha> dm2bb[p,q,r,s] = <q_beta^dagger s_beta^dagger r_beta p_beta>

(p,q correspond to one particle and r,s correspond to another particle) Two-particle density matrix should be contracted to integrals with the pattern below to compute energy

E = numpy.einsum(‘pqrs,pqrs’, eri_aa, dm2_aa) E+= numpy.einsum(‘pqrs,pqrs’, eri_ab, dm2_ab) E+= numpy.einsum(‘pqrs,rspq’, eri_ba, dm2_ab) E+= numpy.einsum(‘pqrs,pqrs’, eri_bb, dm2_bb)

where eri_aa[p,q,r,s] = (p_alpha q_alpha | r_alpha s_alpha ) eri_ab[p,q,r,s] = ( p_alpha q_alpha | r_beta s_beta ) eri_ba[p,q,r,s] = ( p_beta q_beta | r_alpha s_alpha ) eri_bb[p,q,r,s] = ( p_beta q_beta | r_beta s_beta )

nuc_grad_method()[source]
update_amps(t2, eris)

Update non-canonical MP2 amplitudes

pyscf.mp.ump2.energy(mp, t2, eris)[source]

MP2 energy

pyscf.mp.ump2.get_frozen_mask(mp)[source]

Get boolean mask for the unrestricted reference orbitals.

In the returned boolean (mask) array of frozen orbital indices, the element is False if it corresonds to the frozen orbital.

pyscf.mp.ump2.get_nmo(mp)[source]
pyscf.mp.ump2.get_nocc(mp)[source]
pyscf.mp.ump2.kernel(mp, mo_energy=None, mo_coeff=None, eris=None, with_t2=True, verbose=None)[source]
pyscf.mp.ump2.make_fno(mp, thresh=1e-06, pct_occ=None, t2=None, eris=None)[source]

Frozen natural orbitals

Returns:
frozenlist or ndarray

Length-2 list of orbitals to freeze

no_coeffndarray

Length-2 list of semicanonical NO coefficients in the AO basis

pyscf.mp.ump2.make_rdm1(mp, t2=None, ao_repr=False)[source]

One-particle spin density matrices dm1a, dm1b in MO basis (the occupied-virtual blocks due to the orbital response contribution are not included).

dm1a[p,q] = <q_alpha^dagger p_alpha> dm1b[p,q] = <q_beta^dagger p_beta>

The convention of 1-pdm is based on McWeeney’s book, Eq (5.4.20).

pyscf.mp.ump2.make_rdm2(mp, t2=None, ao_repr=False)[source]

Two-particle spin density matrices dm2aa, dm2ab, dm2bb in MO basis

dm2aa[p,q,r,s] = <q_alpha^dagger s_alpha^dagger r_alpha p_alpha> dm2ab[p,q,r,s] = <q_alpha^dagger s_beta^dagger r_beta p_alpha> dm2bb[p,q,r,s] = <q_beta^dagger s_beta^dagger r_beta p_beta>

(p,q correspond to one particle and r,s correspond to another particle) Two-particle density matrix should be contracted to integrals with the pattern below to compute energy

E = numpy.einsum(‘pqrs,pqrs’, eri_aa, dm2_aa) E+= numpy.einsum(‘pqrs,pqrs’, eri_ab, dm2_ab) E+= numpy.einsum(‘pqrs,rspq’, eri_ba, dm2_ab) E+= numpy.einsum(‘pqrs,pqrs’, eri_bb, dm2_bb)

where eri_aa[p,q,r,s] = (p_alpha q_alpha | r_alpha s_alpha ) eri_ab[p,q,r,s] = ( p_alpha q_alpha | r_beta s_beta ) eri_ba[p,q,r,s] = ( p_beta q_beta | r_alpha s_alpha ) eri_bb[p,q,r,s] = ( p_beta q_beta | r_beta s_beta )

pyscf.mp.ump2.update_amps(mp, t2, eris)[source]

Update non-canonical MP2 amplitudes

Module contents

Moller-Plesset perturbation theory

pyscf.mp.GMP2(mf, frozen=None, mo_coeff=None, mo_occ=None)[source]
pyscf.mp.MP2(mf, frozen=None, mo_coeff=None, mo_occ=None)[source]

restricted MP2 with canonical HF and non-canonical HF reference

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default value equals to Mole.max_memory

conv_tolfloat

For non-canonical MP2, converge threshold for MP2 correlation energy. Default value is 1e-7.

conv_tol_normtfloat

For non-canonical MP2, converge threshold for norm(t2). Default value is 1e-5.

max_cycleint

For non-canonical MP2, max number of MP2 iterations. Default value is 50.

diis_spaceint

For non-canonical MP2, DIIS space size in MP2 iterations. Default is 6.

level_shiftfloat

A shift on virtual orbital energies to stablize the MP2 iterations.

frozenint or list

If integer is given, the inner-most orbitals are excluded from MP2 amplitudes. Given the orbital indices (0-based) in a list, both occupied and virtual orbitals can be frozen in MP2 calculation.

>>> mol = gto.M(atom = 'H 0 0 0; F 0 0 1.1', basis = 'ccpvdz')
>>> mf = scf.RHF(mol).run()
>>> # freeze 2 core orbitals
>>> pt = mp.MP2(mf).set(frozen = 2).run()
>>> # auto-generate the number of core orbitals to be frozen (1 in this case)
>>> pt = mp.MP2(mf).set_frozen().run()
>>> # freeze 2 core orbitals and 3 high lying unoccupied orbitals
>>> pt.set(frozen = [0,1,16,17,18]).run()

Saved results

e_corrfloat

MP2 correlation correction

e_corr_ss/osfloat

Same-spin and opposite-spin component of the MP2 correlation energy

e_totfloat

Total MP2 energy (HF + correlation)

t2 :

T amplitudes t2[i,j,a,b] (i,j in occ, a,b in virt)

pyscf.mp.RMP2(mf, frozen=None, mo_coeff=None, mo_occ=None)[source]

restricted MP2 with canonical HF and non-canonical HF reference

Attributes:
verboseint

Print level. Default value equals to Mole.verbose

max_memoryfloat or int

Allowed memory in MB. Default value equals to Mole.max_memory

conv_tolfloat

For non-canonical MP2, converge threshold for MP2 correlation energy. Default value is 1e-7.

conv_tol_normtfloat

For non-canonical MP2, converge threshold for norm(t2). Default value is 1e-5.

max_cycleint

For non-canonical MP2, max number of MP2 iterations. Default value is 50.

diis_spaceint

For non-canonical MP2, DIIS space size in MP2 iterations. Default is 6.

level_shiftfloat

A shift on virtual orbital energies to stablize the MP2 iterations.

frozenint or list

If integer is given, the inner-most orbitals are excluded from MP2 amplitudes. Given the orbital indices (0-based) in a list, both occupied and virtual orbitals can be frozen in MP2 calculation.

>>> mol = gto.M(atom = 'H 0 0 0; F 0 0 1.1', basis = 'ccpvdz')
>>> mf = scf.RHF(mol).run()
>>> # freeze 2 core orbitals
>>> pt = mp.MP2(mf).set(frozen = 2).run()
>>> # auto-generate the number of core orbitals to be frozen (1 in this case)
>>> pt = mp.MP2(mf).set_frozen().run()
>>> # freeze 2 core orbitals and 3 high lying unoccupied orbitals
>>> pt.set(frozen = [0,1,16,17,18]).run()

Saved results

e_corrfloat

MP2 correlation correction

e_corr_ss/osfloat

Same-spin and opposite-spin component of the MP2 correlation energy

e_totfloat

Total MP2 energy (HF + correlation)

t2 :

T amplitudes t2[i,j,a,b] (i,j in occ, a,b in virt)

pyscf.mp.UMP2(mf, frozen=None, mo_coeff=None, mo_occ=None)[source]