10.24. rt — Real-time time-dependent density functional theory

The rt module implements the real-time time-dependent mean-field methods.

10.24.1. Program reference

class pyscf.rt.tdscf.RTTDSCF(ks, prm=None, output='log.dat', auxbas='weigend')[source]

RT-TDSCF base object. Other types of propagations may inherit from this. Calling this class starts the propagation Attributes:

verbose: int

Print level. Default value equals to ks.verbose

conv_tol: float

converge threshold. Default value equals to ks.conv_tol

auxbas: str

auxilliary basis for 2c/3c eri. Default is weigend

prm: str

string object with |variable value| on each line

Saved results

output: str

name of the file with result of propagation

auxmol_set(mol, auxbas='weigend')[source]

Generate 2c/3c electron integral (eri2c,eri3c) Generate ovlp matrix (S), and AO to Lowdin AO matrix transformation matrix (X)

Args:
mol: Mole class

Default is ks.mol

Kwargs:
auxbas: str

auxilliary basis for 2c/3c eri. Default is weigend

Returns:
eri3c: float

3 center eri. shape: (AO,AO,AUX)

eri2c: float

2 center eri. shape: (AUX,AUX)

dipole(rho, c_am)[source]
Args:
c_am: float or complex

Transformation Matrix |AO><MO|

rho: complex

MO density matrix.

Returns:
dipole: float

xyz component of dipole of a molecule. [x y z]

energy(dm_lao, fmat, jmat, kmat)[source]
Args:
dm_lao: complex

Density in LAO basis.

fmat: complex

Fock matrix in Lowdin AO basis

jmat: complex

Coulomb matrix in AO basis

kmat: complex

Exact Exchange in AO basis

Returns:
e_tot: float

Total Energy of a system

fockbuild(dm_lao, it=-1)[source]

Updates Fock matrix

Args:
dm_lao: float or complex

Lowdin AO density matrix.

it: int

iterator for SCF DIIS

Returns:
fmat: float or complex

Fock matrix in Lowdin AO basis

jmat: float or complex

Coulomb matrix in AO basis

kmat: float or complex

Exact Exchange in AO basis

get_j(dm)[source]

Update Coulomb Matrix

Args:
dm: float or complex

AO density matrix.

Returns:
jmat: float or complex

Coulomb matrix in AO basis

get_jk(dm)[source]

Update Coulomb and Exact Exchange Matrix

Args:
dm: float or complex

AO density matrix.

Returns:
jmat: float or complex

Coulomb matrix in AO basis

kmat: float or complex

Exact Exchange in AO basis

get_k(dm)[source]

Update Exact Exchange Matrix

Args:
dm: float or complex

AO density matrix.

Returns:
kmat: float or complex

Exact Exchange in AO basis

get_vxc(dm)[source]

Update exchange matrices and energy Args:

dm: float or complex

AO density matrix.

Returns:
vxc: float or complex

exchange-correlation matrix in AO basis

excsum: float

exchange-correlation energy

kmat: float or complex

Exact Exchange in AO basis

initfockbuild()[source]

Using Roothan’s equation to build a Initial Fock matrix and Transformation Matrices

Returns:
fmat: float or complex

Fock matrix in Lowdin AO basis

c_am: float

Transformation Matrix |AO><MO|

v_lm: float

Transformation Matrix |LAO><MO|

initialcondition(prm)[source]

Prepare the variables/Matrices needed for propagation The SCF is done here to make matrices that are not accessable from pyscf.scf Args:

prm: str

string object with |variable value| on each line

Returns:
fmat: float or complex

Fock matrix in Lowdin AO basis

c_am: float

Transformation Matrix |AO><MO|

v_lm: float

Transformation Matrix |LAO><MO|

rho: float or complex

Initial MO density matrix.

loginstant(rho, c_am, v_lm, fmat, jmat, kmat, tnow, it)[source]

time is logged in atomic units. Args:

rho: complex

MO density matrix.

c_am: complex

Transformation Matrix |AO><MO|

v_lm: complex

Transformation Matrix |LAO><MO|

fmat: complex

Fock matrix in Lowdin AO basis

jmat: complex

Coulomb matrix in AO basis

kmat: complex

Exact Exchange in AO basis

tnow: float

Current time in propagation in A.U.

it: int

Number of iteration of propagation

Returns:
tore: str

|t, dipole(x,y,z), energy|

prop(fmat, c_am, v_lm, rho, output)[source]

The main tdscf propagation loop. Args:

fmat: complex

Fock matrix in Lowdin AO basis

c_am: complex

Transformation Matrix |AO><MO|

v_lm: complex

Transformation Matrix |LAO><MO|

rho: complex

MO density matrix.

output: str

name of the file with result of propagation

Saved results:
f: file

output file with |t, dipole(x,y,z), energy|

readparams(prm)[source]

Set Defaults, Read the file and fill the params dictionary

Args:
prm: str

string object with |variable value| on each line

tddftstep(fmat, c_am, v_lm, rho, rhom12, tnow)[source]

Take dt step in propagation updates matrices and rho to next timestep Args:

fmat: float or complex

Fock matrix in Lowdin AO basis

c_am: float or complex

Transformation Matrix |AO><MO|

v_lm: float or complex

Transformation Matrix |LAO><MO|

rho: complex

MO density matrix.

rhom12: complex tnow: float

current time in A.U.

Returns:
n_rho: complex

MO density matrix.

n_rhom12: complex n_c_am: complex

Transformation Matrix |AO><MO|

n_v_lm: complex

Transformation Matrix |LAO><MO|

n_fmat: complex

Fock matrix in Lowdin AO basis

n_jmat: complex

Coulomb matrix in AO basis

n_kmat: complex

Exact Exchange in AO basis

pyscf.rt.tdscf.matrixpower(A, p, PrintCondition=False)[source]

Raise a Hermitian Matrix to a possibly fractional power.

class pyscf.rt.tdfields.FIELDS(rt, prm, field_tol=1e-09)[source]

A class which manages field perturbations.

Attributes:
dip_ints: float

Dipole matrix of a molecule

pert_xyz: float

Strength of perturbation in |x,y,z| direction

dip0: float

Initial Dipole moment of a molecule in |x y z| direction

field_tol: float

Field amplitude value at which it is considered no field

applyfield(a_mat, c_am, tnow)[source]
Args:
a_mat: float or complex

an MO matrix

c_am: float or complex

Transformation Matrix |AO><MO|

tnow: float

current time.

Returns:
a_mat_field: float or complex

an MO matrix with the field added

ison: bool

On whether field is on or off

expectation(rho, c_am)[source]
Args:
rho: float or complex

current MO density.

c_am: float or complex

Transformation Matrix |AO><MO|

Returns:
mu: float or complex

dipole moment in |x y z| direction

impulseamp(tnow)[source]

Apply impulsive wave to the system Args:

tnow: float

Current time in A.U.

Returns:
amp: float

Amplitude of field at time

ison: bool

On whether field is on or off

initializeexpectation(rho0, c_am)[source]

Calculate the initial dipole moment