Overview of adcc
Note
Work in progress. Should be expanded.
Note
For each logical section of adcc, some nice overview should be given in a separate pagce.
For documentation how to connect host programs to adcc, see Connecting host programs to adcc.
The adcc.adcN family of methods
- adcc.run_adc(data_or_matrix, n_states=None, kind='any', conv_tol=None, eigensolver=None, guesses=None, n_guesses=None, n_guesses_doubles=None, output=<_io.TextIOWrapper name='<stdout>' mode='w' encoding='utf-8'>, core_orbitals=None, frozen_core=None, frozen_virtual=None, method=None, n_singlets=None, n_triplets=None, n_spin_flip=None, environment=None, **solverargs)
Run an ADC calculation.
Main entry point to run an ADC calculation. The reference to build the ADC calculation upon is supplied using the data_or_matrix argument. adcc is pretty flexible here. Possible options include:
Hartree-Fock data from a host program, e.g. a molsturm SCF state, a pyscf SCF object or any class implementing the
adcc.HartreeFockProvider
interface. From this data all objects mentioned in (b) to (d) will be implicitly created and will become available in the returned state.A
adcc.ReferenceState
objectA
adcc.LazyMp
objectA
adcc.AdcMatrix
object
- Parameters
data_or_matrix – Data containing the SCF reference
n_states (int, optional) –
kind (str, optional) –
n_singlets (int, optional) –
n_triplets (int, optional) –
n_spin_flip (int, optional) – Specify the number and kind of states to be computed. Possible values for kind are “singlet”, “triplet”, “spin_flip” and “any”, which is the default. For unrestricted references clamping spin-pure singlets/triplets is currently not possible and kind has to remain as “any”. For restricted references kind=”singlets” or kind=”triplets” may be employed to enforce a particular excited states manifold. Specifying n_singlets is equivalent to setting kind=”singlet” and n_states=5. Similarly for n_triplets and n_spin_flip. n_spin_flip is only valid for unrestricted references.
conv_tol (float, optional) – Convergence tolerance to employ in the iterative solver for obtaining the ADC vectors (default: 1e-6 or 10 * SCF tolerance, whatever is larger)
eigensolver (str, optional) – The eigensolver algorithm to use.
n_guesses (int, optional) – Total number of guesses to compute. By default only guesses derived from the singles block of the ADC matrix are employed. See n_guesses_doubles for alternatives. If no number is given here n_guesses = min(4, 2 * number of excited states to compute) or a smaller number if the number of excitation is estimated to be less than the outcome of above formula.
n_guesses_doubles (int, optional) – Number of guesses to derive from the doubles block. By default none unless n_guesses as explicitly given or automatically determined is larger than the number of singles guesses, which can be possibly found.
guesses (list, optional) – Provide the guess vectors to be employed for the ADC run. Takes preference over n_guesses and n_guesses_doubles, such that these parameters are ignored.
output (stream, optional) – Python stream to which output will be written. If None all output is disabled.
core_orbitals (int or list or tuple, optional) – The orbitals to be put into the core-occupied space. For ways to define the core orbitals see the description in
adcc.ReferenceState
. Required if core-valence separation is applied and the input data is given as data from the host program (i.e. option (a) discussed above)frozen_core (int or list or tuple, optional) – The orbitals to select as frozen core orbitals (i.e. inactive occupied orbitals for both the MP and ADC methods performed). For ways to define these see the description in
adcc.ReferenceState
.frozen_virtual (int or list or tuple, optional) – The orbitals to select as frozen virtual orbitals (i.e. inactive virtuals for both the MP and ADC methods performed). For ways to define these see the description in
adcc.ReferenceState
.environment (bool or list or dict, optional) – The keywords to specify how coupling to an environment model, e.g. PE, is treated. For details see environment.
max_subspace (int, optional) – Maximal subspace size
max_iter (int, optional) – Maximal number of iterations
- Returns
An
adcc.ExcitedStates
object containing theadcc.AdcMatrix
, theadcc.LazyMp
ground state and theadcc.ReferenceState
as well as computed eigenpairs.- Return type
Examples
Run an ADC(2) calculation on top of a pyscf RHF reference of hydrogen flouride.
>>> from pyscf import gto, scf ... mol = gto.mole.M(atom="H 0 0 0; F 0 0 1.1", basis="sto-3g") ... mf = scf.RHF(mol) ... mf.conv_tol_grad = 1e-8 ... mf.kernel() ... ... state = adcc.run_adc(mf, method="adc2", n_singlets=3)
The same thing can also be achieved using the adcc.adcN family of short-hands (see e.g.
adcc.adc2()
,adcc.cvs_adc2x()
):>>> state = adcc.adc2(mf, n_singlets=3)
Run a CVS-ADC(3) calculation of O2 with one core-occupied orbital
>>> from pyscf import gto, scf ... mol = gto.mole.M(atom="O 0 0 0; O 0 0 1.2", basis="sto-3g") ... mf = scf.RHF(mol) ... mf.conv_tol_grad = 1e-8 ... mf.kernel() ... ... state = adcc.cvs_adc3(mf, core_orbitals=1, n_singlets=3)
- adcc.adc0(data_or_matrix, **kwargs)
Run an ADC(0) calculation. For more details see
adcc.run_adc()
.
- adcc.adc1(data_or_matrix, **kwargs)
Run an ADC(1) calculation. For more details see
adcc.run_adc()
.
- adcc.adc2(data_or_matrix, **kwargs)
Run an ADC(2) calculation. For more details see
adcc.run_adc()
.
- adcc.adc2x(data_or_matrix, **kwargs)
Run an ADC(2)-x calculation. For more details see
adcc.run_adc()
.
- adcc.adc3(data_or_matrix, **kwargs)
Run an ADC(3) calculation. For more details see
adcc.run_adc()
.
- adcc.cvs_adc0(data_or_matrix, **kwargs)
Run an CVS-ADC(0) calculation. For more details see
adcc.run_adc()
.
- adcc.cvs_adc1(data_or_matrix, **kwargs)
Run an CVS-ADC(1) calculation. For more details see
adcc.run_adc()
.
- adcc.cvs_adc2(data_or_matrix, **kwargs)
Run an CVS-ADC(2) calculation. For more details see
adcc.run_adc()
.
- adcc.cvs_adc2x(data_or_matrix, **kwargs)
Run an CVS-ADC(2)-x calculation. For more details see
adcc.run_adc()
.
- adcc.cvs_adc3(data_or_matrix, **kwargs)
Run an CVS-ADC(3) calculation. For more details see
adcc.run_adc()
.
- class adcc.ExcitedStates(data, method=None, property_method=None)
Construct an ExcitedStates class from some data obtained from an interative solver or another
ExcitedStates
object.The class provides access to the results from an ADC calculation as well as derived properties. Properties are computed lazily on the fly as requested by the user.
By default the ADC method is extracted from the data object and the property method in property_method is set equal to this method, except ADC(3) where property_method==”adc2”. This can be overwritten using the parameters.
- Parameters
data – Any kind of iterative solver state. Typically derived off a
solver.EigenSolverStateBase
. Can also be anExcitedStates
object.method (str, optional) – Provide an explicit method parameter if data contains none.
property_method (str, optional) – Provide an explicit method for property calculations to override the automatic selection.
- __init__(data, method=None, property_method=None)
Construct an ExcitedStates class from some data obtained from an interative solver or another
ExcitedStates
object.The class provides access to the results from an ADC calculation as well as derived properties. Properties are computed lazily on the fly as requested by the user.
By default the ADC method is extracted from the data object and the property method in property_method is set equal to this method, except ADC(3) where property_method==”adc2”. This can be overwritten using the parameters.
- Parameters
data – Any kind of iterative solver state. Typically derived off a
solver.EigenSolverStateBase
. Can also be anExcitedStates
object.method (str, optional) – Provide an explicit method parameter if data contains none.
property_method (str, optional) – Provide an explicit method for property calculations to override the automatic selection.
- describe(oscillator_strengths=True, rotatory_strengths=False, state_dipole_moments=False, transition_dipole_moments=False, block_norms=True)
Return a string providing a human-readable description of the class
- Parameters
oscillator_strengths (bool optional) – Show oscillator strengths, by default
True
.rotatory_strengths (bool optional) – Show rotatory strengths, by default
False
.state_dipole_moments (bool, optional) – Show state dipole moments, by default
False
.transition_dipole_moments (bool, optional) – Show state dipole moments, by default
False
.block_norms (bool, optional) – Show the norms of the (1p1h, 2p2h, …) blocks of the excited states, by default
True
.
- describe_amplitudes(tolerance=0.01, index_format=None)
Return a string describing the dominant amplitudes of each excitation vector in human-readable form. The
kwargs
are forFormatExcitationVector
.- Parameters
tolerance (float, optional) – Minimal absolute value of the excitation amplitudes considered.
index_format (NoneType or str or FormatIndexBase, optional) – Formatter to use for displaying tensor indices. Valid are
"adcc"
to keep the adcc-internal indexing,"hf"
to select the HFProvider indexing,"homolumo"
to index relative on the HOMO / LUMO / HOCO orbitals. IfNone
an automatic selection will be made.
- property excitation_property_keys
Extracts the property keys which are marked as excitation property with
mark_excitation_property()
.
- property excitations
Provides a list of Excitations, i.e., a view to all individual excited states and their properties. Still under heavy development.
- property state_diffdm
List of difference density matrices of all computed states
- property state_dipole_moment
List of state dipole moments
- property state_dm
List of state density matrices of all computed states
- to_dataframe(**kwargs)
Exports the ExcitedStates object as
pandas.DataFrame
. Atomic units are used for all values.
- to_qcvars(properties=False, recurse=False)
Return a dictionary with property keys compatible to a Psi4 wavefunction or a QCEngine Atomicresults object.
- property transition_dm
List of transition density matrices of all computed states
Visualisation
- class adcc.visualisation.ExcitationSpectrum(energy, intensity)
Construct an ExcitationSpectrum object
- Parameters
energy – Energies for plotting the spectrum
intensity – Intensities for plotting the spectrum
- __init__(energy, intensity)
Construct an ExcitationSpectrum object
- Parameters
energy – Energies for plotting the spectrum
intensity – Intensities for plotting the spectrum
- broaden_lines(width=None, shape='lorentzian', xmin=None, xmax=None)
Apply broadening to the current spectral data and return the broadened spectrum.
- Parameters
shape (str or callable, optional) – The shape of the broadening to use (lorentzian or gaussian), by default lorentzian broadening is used. This can be a callable to directly specify the function with which each line of the spectrum is convoluted.
width (float, optional) – The width to use for the broadening (stddev for the gaussian, gamma parameter for the lorentzian). Optional if shape is a callable.
xmin (float, optional) – Explicitly set the minimum value of the x-axis for broadening
xmax (float, optional) – Explicitly set the maximum value of the x-axis for broadening
- copy()
Return a consistent copy of the object.
- plot(**kwargs)
Plot the Spectrum represented by this class.
Parameters not listed below are passed to the matplotlib plot function.
- Parameters
style (str, optional) – Use some default setup of matplotlib parameters for certain types of spectra commonly plotted. Valid are “discrete” and “continuous”. By default no special style is chosen.
Adc Middle layer
- class adcc.AdcMatrix(method, hf_or_mp, block_orders=None, intermediates=None, diagonal_precomputed=None)
Initialise an ADC matrix.
- Parameters
method (str or AdcMethod) – Method to use.
hf_or_mp (adcc.ReferenceState or adcc.LazyMp) – HF reference or MP ground state
block_orders (optional) – The order of perturbation theory to employ for each matrix block. If not set, defaults according to the selected ADC method are chosen.
intermediates (adcc.Intermediates or NoneType) – Allows to pass intermediates to re-use to this class.
diagonal_precomputed (adcc.AmplitudeVector) – Allows to pass a pre-computed diagonal, for internal use only.
- __init__(method, hf_or_mp, block_orders=None, intermediates=None, diagonal_precomputed=None)
Initialise an ADC matrix.
- Parameters
method (str or AdcMethod) – Method to use.
hf_or_mp (adcc.ReferenceState or adcc.LazyMp) – HF reference or MP ground state
block_orders (optional) – The order of perturbation theory to employ for each matrix block. If not set, defaults according to the selected ADC method are chosen.
intermediates (adcc.Intermediates or NoneType) – Allows to pass intermediates to re-use to this class.
diagonal_precomputed (adcc.AmplitudeVector) – Allows to pass a pre-computed diagonal, for internal use only.
- property axis_blocks
Return the blocks used along one of the axes of the ADC matrix (e.g. [‘ph’, ‘pphh’]).
- block_apply(block, tensor)
Compute the application of a block of the ADC matrix with another AmplitudeVector or Tensor. Non-matching blocks in the AmplitudeVector will be ignored.
- block_spaces(block)
- block_view(block)
Return a view into the AdcMatrix that represents a single block of the matrix. Currently only diagonal blocks are supported.
- property blocks
- compute_apply(block, tensor)
- compute_matvec(ampl)
Compute the matrix-vector product of the ADC matrix with an excitation amplitude and return the result.
- construct_symmetrisation_for_blocks()
Construct the symmetrisation functions, which need to be applied to relevant blocks of an AmplitudeVector in order to symmetrise it to the right symmetry in order to be used with the various matrix-vector-products of this function.
Most importantly the returned functions antisymmetrise the occupied and virtual parts of the doubles parts if this is sensible for the method behind this adcmatrix.
Returns a dictionary block identifier -> function
- default_block_orders = {'adc0': {'ph_ph': 0, 'ph_pphh': None, 'pphh_ph': None, 'pphh_pphh': None}, 'adc1': {'ph_ph': 1, 'ph_pphh': None, 'pphh_ph': None, 'pphh_pphh': None}, 'adc2': {'ph_ph': 2, 'ph_pphh': 1, 'pphh_ph': 1, 'pphh_pphh': 0}, 'adc2x': {'ph_ph': 2, 'ph_pphh': 1, 'pphh_ph': 1, 'pphh_pphh': 1}, 'adc3': {'ph_ph': 3, 'ph_pphh': 2, 'pphh_ph': 2, 'pphh_pphh': 1}}
- dense_basis(axis_blocks=None, ordering='adcc')
Return the list of indices and their values of the dense basis representation
ordering: adcc, spin, spatial
- diagonal(block=None)
Return the diagonal of the ADC matrix
- has_block(block)
- matvec(*args, **kwargs)
Compute the matrix-vector product of the ADC matrix with an excitation amplitude and return the result.
- rmatvec(v)
- to_ndarray(out=None)
Return the ADC matrix object as a dense numpy array. Converts the sparse internal representation of the ADC matrix to a dense matrix and return as a numpy array.
Notes
This method is only intended to be used for debugging and visualisation purposes as it involves computing a large amount of matrix-vector products and the returned array consumes a considerable amount of memory.
The resulting matrix has no spin symmetry imposed, which means that its eigenspectrum may contain non-physical excitations (e.g. with linear combinations of α->β and α->α components in the excitation vector).
This function has not been sufficiently tested to be considered stable.
- class adcc.AdcMethod(method)
- at_level(newlevel)
Return an equivalent method, where only the level is changed (e.g. calling this on a CVS method returns a CVS method)
- available_methods = ['adc0', 'adc1', 'adc2', 'adc2x', 'adc3', 'cvs-adc0', 'cvs-adc1', 'cvs-adc2', 'cvs-adc2x', 'cvs-adc3']
- property base_method
The base (full) method, i.e. with all approximations such as CVS stripped off.
- property name
- property property_method
The name of the canonical method to use for computing properties for this ADC method. This only differs from the name property for the ADC(2)-x family of methods.
- class adcc.ReferenceState(hfdata, core_orbitals=None, frozen_core=None, frozen_virtual=None, symmetry_check_on_import=False, import_all_below_n_orbs=10)
Construct a ReferenceState holding information about the employed SCF reference.
The constructed object is lazy and will at construction only setup orbital energies and coefficients. Fock matrix blocks and electron-repulsion integral blocks are imported as needed.
Orbital subspace selection: In order to specify frozen_core, core_orbitals and frozen_virtual, adcc allows a range of specifications including
A number: Just put this number of alpha orbitals and this number of beta orbitals into the respective space. For frozen core and core orbitals these are counted from below, for frozen virtual orbitals, these are counted from above. If both frozen core and core orbitals are specified like this, the lowest-energy, occupied orbitals will be put into frozen core.
A range: The orbital indices given by this range will be put into the orbital subspace.
An explicit list of orbital indices to be placed into the subspace.
A pair of (a) to (c): If the orbital selection for alpha and beta orbitals should differ, a pair of ranges, or a pair of index lists or a pair of numbers can be specified.
- Parameters
hfdata – Object with Hartree-Fock data (e.g. a molsturm scf state, a pyscf SCF object or any class implementing the
adcc.HartreeFockProvider
interface or in fact any python object representing a pointer to a C++ object derived off theadcc::HartreeFockSolution_i
.core_orbitals (int or list or tuple, optional) – The orbitals to be put into the core-occupied space. For ways to define the core orbitals see the description above.
frozen_core (int or list or tuple, optional) – The orbitals to be put into the frozen core space. For ways to define the core orbitals see the description above. For an automatic selection of the frozen core space one may also specify
frozen_core=True
.frozen_virtuals (int or list or tuple, optional) – The orbitals to be put into the frozen virtual space. For ways to define the core orbitals see the description above.
symmetry_check_on_import (bool, optional) – Should symmetry of the imported objects be checked explicitly during the import process. This massively slows down the import and has a dramatic impact on memory usage. Thus one should enable this only for debugging (e.g. for testing import routines from the host programs). Do not enable this unless you know what you are doing.
import_all_below_n_orbs (int, optional) – For small problem sizes lazy make less sense, since the memory requirement for storing the ERI tensor is neglibile and thus the flexiblity gained by having the full tensor in memory is advantageous. Below the number of orbitals specified by this parameter, the class will thus automatically import all ERI tensor and Fock matrix blocks.
Examples
To start a calculation with the 2 lowest alpha and beta orbitals in the core occupied space, construct the class as
>>> ReferenceState(hfdata, core_orbitals=2)
or
>>> ReferenceState(hfdata, core_orbitals=range(2))
or
>>> ReferenceState(hfdata, core_orbitals=[0, 1])
or
>>> ReferenceState(hfdata, core_orbitals=([0, 1], [0, 1]))
There is no restriction to choose the core occupied orbitals from the bottom end of the occupied orbitals. For example to select the 2nd and 3rd orbital setup the class as
>>> ReferenceState(hfdata, core_orbitals=range(1, 3))
or
>>> ReferenceState(hfdata, core_orbitals=[1, 2])
If different orbitals should be placed in the alpha and beta orbitals, this can be achievd like so
>>> ReferenceState(hfdata, core_orbitals=([1, 2], [0, 1]))
which would place the 2nd and 3rd alpha and the 1st and second beta orbital into the core space.
- __init__(hfdata, core_orbitals=None, frozen_core=None, frozen_virtual=None, symmetry_check_on_import=False, import_all_below_n_orbs=10)
Construct a ReferenceState holding information about the employed SCF reference.
The constructed object is lazy and will at construction only setup orbital energies and coefficients. Fock matrix blocks and electron-repulsion integral blocks are imported as needed.
Orbital subspace selection: In order to specify frozen_core, core_orbitals and frozen_virtual, adcc allows a range of specifications including
A number: Just put this number of alpha orbitals and this number of beta orbitals into the respective space. For frozen core and core orbitals these are counted from below, for frozen virtual orbitals, these are counted from above. If both frozen core and core orbitals are specified like this, the lowest-energy, occupied orbitals will be put into frozen core.
A range: The orbital indices given by this range will be put into the orbital subspace.
An explicit list of orbital indices to be placed into the subspace.
A pair of (a) to (c): If the orbital selection for alpha and beta orbitals should differ, a pair of ranges, or a pair of index lists or a pair of numbers can be specified.
- Parameters
hfdata – Object with Hartree-Fock data (e.g. a molsturm scf state, a pyscf SCF object or any class implementing the
adcc.HartreeFockProvider
interface or in fact any python object representing a pointer to a C++ object derived off theadcc::HartreeFockSolution_i
.core_orbitals (int or list or tuple, optional) – The orbitals to be put into the core-occupied space. For ways to define the core orbitals see the description above.
frozen_core (int or list or tuple, optional) – The orbitals to be put into the frozen core space. For ways to define the core orbitals see the description above. For an automatic selection of the frozen core space one may also specify
frozen_core=True
.frozen_virtuals (int or list or tuple, optional) – The orbitals to be put into the frozen virtual space. For ways to define the core orbitals see the description above.
symmetry_check_on_import (bool, optional) – Should symmetry of the imported objects be checked explicitly during the import process. This massively slows down the import and has a dramatic impact on memory usage. Thus one should enable this only for debugging (e.g. for testing import routines from the host programs). Do not enable this unless you know what you are doing.
import_all_below_n_orbs (int, optional) – For small problem sizes lazy make less sense, since the memory requirement for storing the ERI tensor is neglibile and thus the flexiblity gained by having the full tensor in memory is advantageous. Below the number of orbitals specified by this parameter, the class will thus automatically import all ERI tensor and Fock matrix blocks.
Examples
To start a calculation with the 2 lowest alpha and beta orbitals in the core occupied space, construct the class as
>>> ReferenceState(hfdata, core_orbitals=2)
or
>>> ReferenceState(hfdata, core_orbitals=range(2))
or
>>> ReferenceState(hfdata, core_orbitals=[0, 1])
or
>>> ReferenceState(hfdata, core_orbitals=([0, 1], [0, 1]))
There is no restriction to choose the core occupied orbitals from the bottom end of the occupied orbitals. For example to select the 2nd and 3rd orbital setup the class as
>>> ReferenceState(hfdata, core_orbitals=range(1, 3))
or
>>> ReferenceState(hfdata, core_orbitals=[1, 2])
If different orbitals should be placed in the alpha and beta orbitals, this can be achievd like so
>>> ReferenceState(hfdata, core_orbitals=([1, 2], [0, 1]))
which would place the 2nd and 3rd alpha and the 1st and second beta orbital into the core space.
- property backend
The identifier of the back end used for the SCF calculation.
- property cached_eri_blocks
Get or set the list of momentarily cached ERI tensor blocks
Setting this property allows to drop ERI tensor blocks if they are no longer needed to save memory.
- property cached_fock_blocks
Get or set the list of momentarily cached Fock matrix blocks
Setting this property allows to drop fock matrix blocks if they are no longer needed to save memory.
- property conv_tol
SCF convergence tolererance
- property density
Return the Hartree-Fock density in the MO basis
- property dipole_moment
Return the HF dipole moment of the reference state (that is the sum of the electronic and the nuclear contribution.)
- property energy_scf
Final total SCF energy
- eri(self: libadcc.ReferenceState, arg0: str) libadcc::Tensor
Return the ERI (electron-repulsion integrals) tensor block corresponding to the provided space.
- flush_hf_cache(self: libadcc.ReferenceState) None
Tell the contained HartreeFockSolution_i object (which was passed upon construction), that a larger amount of import operations is done and that the next request for further imports will most likely take some time, such that intermediate caches can now be flushed to save some memory or other resources.
- fock(self: libadcc.ReferenceState, arg0: str) libadcc::Tensor
Return the Fock matrix block corresponding to the provided space.
- property has_core_occupied_space
Is a core occupied space setup, such that a core-valence separation can be applied.
- import_all(self: libadcc.ReferenceState) None
Normally the class only imports the Fock matrix blocks and electron-repulsion integrals of a particular space combination when this is requested by a call to above fock() or eri() functions. This function call, however, instructs the class to immediately import all such blocks. Typically you do not want to do this.
- property irreducible_representation
Reference state irreducible representation
- property is_aufbau_occupation
Returns whether the molecular orbital occupation in this reference is according to the Aufbau principle (lowest-energy orbitals are occupied)
- property mospaces
The MoSpaces object supplied on initialisation
- property n_alpha
Number of alpha electrons
- property n_beta
Number of beta electrons
- property n_orbs
Number of molecular orbitals
- property n_orbs_alpha
Number of alpha orbitals
- property n_orbs_beta
Number of beta orbitals
- property nuclear_dipole
- property nuclear_total_charge
- orbital_coefficients(self: libadcc.ReferenceState, arg0: str) libadcc::Tensor
Return the molecular orbital coefficients corresponding to the provided space (alpha and beta coefficients are returned)
- orbital_coefficients_alpha(self: libadcc.ReferenceState, arg0: str) libadcc::Tensor
Return the alpha molecular orbital coefficients corresponding to the provided space
- orbital_coefficients_beta(self: libadcc.ReferenceState, arg0: str) libadcc::Tensor
Return the beta molecular orbital coefficients corresponding to the provided space
- orbital_energies(self: libadcc.ReferenceState, arg0: str) libadcc::Tensor
Return the orbital energies corresponding to the provided space
- property restricted
Return whether the reference is restricted or not.
- property spin_multiplicity
Return the spin multiplicity of the reference state. 0 indicates that the spin cannot be determined or is not integer (e.g. UHF)
- property timer
Obtain the timer object of this class.
- to_qcvars(properties=False, recurse=False)
Return a dictionary with property keys compatible to a Psi4 wavefunction or a QCEngine Atomicresults object.
- class adcc.LazyMp(hf)
Initialise the class dealing with the Møller-Plesset ground state.
- __init__(hf)
Initialise the class dealing with the Møller-Plesset ground state.
- density(level=2)
Return the MP density in the MO basis with all corrections up to the specified order of perturbation theory
- df(space)
Delta Fock matrix
- dipole_moment(level=2)
Return the MP dipole moment at the specified level of perturbation theory.
- energy(level=2)
Obtain the total energy (SCF energy plus all corrections) at a particular level of perturbation theory.
- energy_correction(level=2)
Obtain the MP energy correction at a particular level
- property mp2_density
- property mp2_diffdm
Return the MP2 differensce density in the MO basis.
- property mp2_dipole_moment
- t2(space)
T2 amplitudes
- t2eri(space, contraction)
Return the T2 tensor with ERI tensor contraction intermediates. These are called pi1 to pi7 in libadc.
- td2(space)
Return the T^D_2 term
- to_qcvars(properties=False, recurse=False, maxlevel=2)
Return a dictionary with property keys compatible to a Psi4 wavefunction or a QCEngine Atomicresults object.
Tensor and symmetry interface
- class adcc.Tensor(sym_or_mo, space=None, permutations=None, spin_block_maps=None, spin_blocks_forbidden=None)
Construct an uninitialised Tensor from an
MoSpaces
or aSymmetry
object.More information about the last four, symmetry-related parameters see the documentation of the
Symmetry
object.- Parameters
sym_or_mo – Symmetry or MoSpaces object
spaces (str, optional) – Space of the tensor, can be None if the first argument is a
Symmetry
object.permutations (list, optional) – List of permutational symmetries of the Tensor.
spin_block_maps (list, optional) – List of mappings between spin blocks
spin_blocks_forbidden (list, optional) – List of forbidden (i.e. forced-to-zero) spin blocks.
Notes
An
MoSpaces
object is contained in many datastructures of adcc, including theAdcMatrix
, theLazyMp
, theReferenceState
and any solver or ADC results state.Examples
Construct a symmetric tensor in the “o1o1” (occupied-occupied) spaces:
>>> Tensor(mospaces, "o1o1", permutations=["ij", "ji"])
Construct an anti-symmetric tensor in the “v1v1” spaces:
>>> Tensor(mospaces, "v1v1", permutations=["ab", "-ba"])
Construct a tensor in “o1v1”, which maps the alpha-alpha block anti-symmetrically to the beta-beta block and which has the other spin blocks set to zero:
>>> Tensor(mospaces, "o1v1", spin_block_maps=[("aa", "bb", -1)], ... spin_blocks_forbidden=["ab", "ba"])
- property T
- __init__(sym_or_mo, space=None, permutations=None, spin_block_maps=None, spin_blocks_forbidden=None)
Construct an uninitialised Tensor from an
MoSpaces
or aSymmetry
object.More information about the last four, symmetry-related parameters see the documentation of the
Symmetry
object.- Parameters
sym_or_mo – Symmetry or MoSpaces object
spaces (str, optional) – Space of the tensor, can be None if the first argument is a
Symmetry
object.permutations (list, optional) – List of permutational symmetries of the Tensor.
spin_block_maps (list, optional) – List of mappings between spin blocks
spin_blocks_forbidden (list, optional) – List of forbidden (i.e. forced-to-zero) spin blocks.
Notes
An
MoSpaces
object is contained in many datastructures of adcc, including theAdcMatrix
, theLazyMp
, theReferenceState
and any solver or ADC results state.Examples
Construct a symmetric tensor in the “o1o1” (occupied-occupied) spaces:
>>> Tensor(mospaces, "o1o1", permutations=["ij", "ji"])
Construct an anti-symmetric tensor in the “v1v1” spaces:
>>> Tensor(mospaces, "v1v1", permutations=["ab", "-ba"])
Construct a tensor in “o1v1”, which maps the alpha-alpha block anti-symmetrically to the beta-beta block and which has the other spin blocks set to zero:
>>> Tensor(mospaces, "o1v1", spin_block_maps=[("aa", "bb", -1)], ... spin_blocks_forbidden=["ab", "ba"])
- antisymmetrise(*args, **kwargs)
Overloaded function.
antisymmetrise(self: libadcc.Tensor, arg0: list) -> libadcc.Tensor
antisymmetrise(self: libadcc.Tensor, *args) -> libadcc.Tensor
- copy(self: libadcc.Tensor) libadcc.Tensor
Returns a deep copy of the tensor.
- describe_expression(*args, **kwargs)
Overloaded function.
describe_expression(self: libadcc.Tensor, arg0: str) -> str
Return a string providing a hopefully descriptive representation of the tensor expression stored inside the object.
describe_expression(self: libadcc.Tensor) -> str
- describe_symmetry(self: libadcc.Tensor) str
Return a string providing a hopefully descriptive representation of the symmetry information stored inside the tensor.
- diagonal(self: libadcc.Tensor, *args) libadcc.Tensor
- dot(*args, **kwargs)
Overloaded function.
dot(self: libadcc.Tensor, arg0: libadcc.Tensor) -> float
dot(self: libadcc.Tensor, arg0: list) -> numpy.ndarray[numpy.float64]
- empty_like(self: libadcc.Tensor) libadcc.Tensor
- evaluate(self: libadcc.Tensor) libadcc.Tensor
Ensure the tensor to be fully evaluated and resilient in memory. Usually happens automatically when needed. Might be useful for fine-tuning, however.
- property flags
- is_allowed(self: libadcc.Tensor, arg0: tuple) bool
Is a particular index allowed by symmetry
- property mutable
- property ndim
- property needs_evaluation
Does the tensor need evaluation or is it fully evaluated and resilient in memory.
- nosym_like(self: libadcc.Tensor) libadcc.Tensor
- ones_like(self: libadcc.Tensor) libadcc.Tensor
- select_below_absmax(tolerance)
Select the absolute maximal values in the tensor, which are below the given tolerance.
- select_n_absmax(self: libadcc.Tensor, arg0: int) list
Select the n absolute maximal elements.
- select_n_absmin(self: libadcc.Tensor, arg0: int) list
Select the n absolute minimal elements.
- select_n_max(self: libadcc.Tensor, arg0: int) list
Select the n maximal elements.
- select_n_min(self: libadcc.Tensor, arg0: int) list
Select the n minimal elements.
- set_from_ndarray(*args, **kwargs)
Overloaded function.
set_from_ndarray(self: libadcc.Tensor, arg0: numpy.ndarray) -> libadcc.Tensor
Set all tensor elements from a standard np::ndarray by making a copy. Provide an optional tolerance argument to increase the tolerance for the check for symmetry consistency.
set_from_ndarray(self: libadcc.Tensor, arg0: numpy.ndarray[numpy.float64], arg1: float) -> libadcc.Tensor
Set all tensor elements from a standard np::ndarray by making a copy. Provide an optional tolerance argument to increase the tolerance for the check for symmetry consistency.
- set_immutable(self: libadcc.Tensor) None
Set the tensor as immutable, allowing some optimisations to be performed.
- set_mask(self: libadcc.Tensor, arg0: str, arg1: float) None
Set all elements corresponding to an index mask, which is given by a string eg. ‘iijkli’ sets elements T_{iijkli}
- set_random(self: libadcc.Tensor) libadcc.Tensor
Set all tensor elements to random data, adhering to the internal symmetry.
- property shape
- property size
- property space
- property subspaces
- symmetrise(*args, **kwargs)
Overloaded function.
symmetrise(self: libadcc.Tensor, arg0: list) -> libadcc.Tensor
symmetrise(self: libadcc.Tensor, *args) -> libadcc.Tensor
- to_ndarray(self: libadcc.Tensor) numpy.ndarray[numpy.float64]
Export the tensor data to a standard np::ndarray by making a copy.
- transpose(*args, **kwargs)
Overloaded function.
transpose(self: libadcc.Tensor) -> libadcc.Tensor
transpose(self: libadcc.Tensor, arg0: tuple) -> libadcc.Tensor
- zeros_like(self: libadcc.Tensor) libadcc.Tensor
- class adcc.Symmetry(*args, **kwargs)
Overloaded function.
__init__(self: libadcc.Symmetry, arg0: libadcc.MoSpaces, arg1: str) -> None
Construct a Symmetry class from an MoSpaces object and the identifier for the space (e.g. o1o1, v1o1, o3v2o1v1, …). Python binding to
libadcc::Symmetry
.__init__(self: libadcc.Symmetry, arg0: libadcc.MoSpaces, arg1: str, arg2: Dict[str, Tuple[int, int]]) -> None
Construct a Symmetry class from an MoSpaces object, a space string and a map to supply the number of orbitals for some additional axes. For the additional axis the pair contains either two numbers (for the number of alpha and beta orbitals in this axis) or only one number and a zero (for an axis, which as only one spin kind, alpha or beta).
This is an advanced constructor. Use only if you know what you do.
- clear(self: libadcc.Symmetry) None
Clear the symmetry.
- describe(self: libadcc.Symmetry) str
Return a descriptive string.
- property empty
Is the symmetry empty (i.e. noy symmetry setup)
- property irreps_allowed
The list of irreducible representations, for which the tensor shall be non-zero. If this is not set, i.e. an empty list, all irreps will be allowed.
- property mospaces
Return the MoSpaces object supplied on initialisation
- property ndim
Return the number of dimensions.
- property permutations
The list of index permutations, which do not change the tensor. A minus may be used to indicate anti-symmetric permutations with respect to the first (reference) permutation.
For example the list [“ij”, “ji”] defines a symmetric matrix and [“ijkl”, “-jikl”, “-ijlk”, “klij”] the symmetry of the ERI tensor. Not all permutations need to be given to fully describe the symmetry. Beware that the check for errors and conflicts is only rudimentary at the moment.
- property shape
Return the shape of tensors constructed from this symmetry.
- property space
Return the space supplied on initialisation.
- property spin_block_maps
A list of tuples of the form (“aaaa”, “bbbb”, -1.0), i.e. two spin blocks followed by a factor. This maps the second onto the first with a factor of -1.0 between them.
- property spin_blocks_forbidden
List of spin-blocks, which are marked forbidden (i.e. enforce them to stay zero). Blocks are given as a string in the letters ‘a’ and ‘b’, e.g. [“aaba”, “abba”]
- class adcc.AmplitudeVector(*args, **kwargs)
Construct an AmplitudeVector. Typical use cases are
AmplitudeVector(ph=tensor_singles, pphh=tensor_doubles)
.- __init__(*args, **kwargs)
Construct an AmplitudeVector. Typical use cases are
AmplitudeVector(ph=tensor_singles, pphh=tensor_doubles)
.
- property blocks_ph
Return the blocks which are used inside the vector. Note: This is a temporary name. The attribute will be removed in 0.16.0.
- clear() None. Remove all items from D.
- copy()
Return a copy of the AmplitudeVector
- dot(other)
Return the dot product with another AmplitudeVector or the dot products with a list of AmplitudeVectors. In the latter case a np.ndarray is returned.
- empty_like()
Return an empty AmplitudeVector of the same shape and symmetry
- fromkeys(value=None, /)
Create a new dictionary with keys from iterable and values set to value.
- get(key, default=None, /)
Return the value for key if key is in the dictionary, else default.
- items() a set-like object providing a view on D's items
- keys() a set-like object providing a view on D's keys
- nosym_like()
Return an empty AmplitudeVector of the same shape and symmetry
- ones_like()
Return an empty AmplitudeVector of the same shape and symmetry
- pop(k[, d]) v, remove specified key and return the corresponding value.
If key is not found, default is returned if given, otherwise KeyError is raised
- popitem()
Remove and return a (key, value) pair as a 2-tuple.
Pairs are returned in LIFO (last-in, first-out) order. Raises KeyError if the dict is empty.
- setdefault(key, default=None, /)
Insert key with a value of default if key is not in the dictionary.
Return the value for key if key is in the dictionary, else default.
- update([E, ]**F) None. Update D from dict/iterable E and F.
If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]
- values() an object providing a view on D's values
- zeros_like()
Return an AmplitudeVector of the same shape and symmetry with all elements set to zero
- adcc.copy(a)
Return a copy of the input tensor.
- adcc.dot(a, b)
Form the scalar product between two tensors.
- adcc.einsum(subscripts, *operands, optimise='auto')
Evaluate Einstein summation convention for the operands similar to numpy’s einsum function. Uses opt_einsum and libadcc to perform the contractions.
Using this function does not evaluate, but returns a contraction expression tree where contractions are queued in optimal order.
- Parameters
subscripts (str) – Specifies the subscripts for summation.
*operands (list of tensors) – These are the arrays for the operation.
optimise (str, list or bool, optional (default:
auto
)) – Choose the type of the path optimisation, see opt_einsum.contract for details.
- adcc.lincomb(coefficients, tensors, evaluate=False)
Form a linear combination from a list of tensors.
If coefficients is a 1D array, just form a single linear combination, else return a list of vectors representing the linear combination by reading the coefficients row-by-row.
- Parameters
coefficients (list) – Coefficients for the linear combination
tensors (list) – Tensors for the linear combination
evaluate (bool) – Should the linear combination be evaluated (True) or should just a lazy expression be formed (False). Notice that lincomb(…, evaluate=True) is identical to lincomb(…, evaluate=False).evaluate(), but the former is generally faster.
- adcc.empty_like(a)
Return an empty tensor of the same shape and symmetry as the input tensor.
- adcc.nosym_like(a)
Return tensor of the same shape, but without the symmetry setup of the input tensor.
- adcc.ones_like(a)
Return tensor of the same shape and symmetry as the input tensor, but initialised to 1, that is the canonical blocks are 1 and the other ones are symmetry-equivalent (-1 or 0)
- adcc.zeros_like(a)
Return a zero tensor of the same shape and symmetry as the input tensor.
- adcc.evaluate(a)
Force full evaluation of a tensor expression
- adcc.direct_sum(subscripts, *operands)
Solvers
- adcc.solver.davidson.davidson_iterations(matrix, state, max_subspace, max_iter, n_ep, is_converged, which, callback=None, preconditioner=None, preconditioning_method='Davidson', debug_checks=False, residual_min_norm=None, explicit_symmetrisation=None)
Drive the davidson iterations
- Parameters
matrix – Matrix to diagonalise
state – DavidsonState containing the eigenvector guess
max_subspace (int or NoneType, optional) – Maximal subspace size
max_iter (int, optional) – Maximal number of iterations
n_ep (int or NoneType, optional) – Number of eigenpairs to be computed
is_converged – Function to test for convergence
callback (callable, optional) – Callback to run after each iteration
which (str, optional) – Which eigenvectors to converge to. Needs to be chosen such that it agrees with the selected preconditioner.
preconditioner – Preconditioner (type or instance)
preconditioning_method (str, optional) – Precondititoning method. Valid values are “Davidson” or “Sleijpen-van-der-Vorst”
debug_checks (bool, optional) – Enable some potentially costly debug checks (Loss of orthogonality etc.)
residual_min_norm (float or NoneType, optional) – Minimal norm a residual needs to have in order to be accepted as a new subspace vector (defaults to 2 * len(matrix) * machine_expsilon)
explicit_symmetrisation – Explicit symmetrisation to apply to new subspace vectors before adding them to the subspace. Allows to correct for loss of index or spin symmetries (type or instance)
- adcc.solver.davidson.default_print(state, identifier, file=<_io.TextIOWrapper name='<stdout>' mode='w' encoding='utf-8'>)
A default print function for the davidson callback
- adcc.solver.davidson.eigsh(matrix, guesses, n_ep=None, max_subspace=None, conv_tol=1e-09, which='SA', max_iter=70, callback=None, preconditioner=None, preconditioning_method='Davidson', debug_checks=False, residual_min_norm=None, explicit_symmetrisation=<class 'adcc.solver.explicit_symmetrisation.IndexSymmetrisation'>)
Davidson eigensolver for ADC problems
- Parameters
matrix – ADC matrix instance
guesses (list) – Guess vectors (fixes also the Davidson block size)
n_ep (int or NoneType, optional) – Number of eigenpairs to be computed
max_subspace (int or NoneType, optional) – Maximal subspace size
conv_tol (float, optional) – Convergence tolerance on the l2 norm squared of residuals to consider them converged
which (str, optional) – Which eigenvectors to converge to (e.g. LM, LA, SM, SA)
max_iter (int, optional) – Maximal number of iterations
callback (callable, optional) – Callback to run after each iteration
preconditioner – Preconditioner (type or instance)
preconditioning_method (str, optional) – Precondititoning method. Valid values are “Davidson” or “Sleijpen-van-der-Vorst”
explicit_symmetrisation – Explicit symmetrisation to apply to new subspace vectors before adding them to the subspace. Allows to correct for loss of index or spin symmetries (type or instance)
debug_checks (bool, optional) – Enable some potentially costly debug checks (Loss of orthogonality etc.)
residual_min_norm (float or NoneType, optional) – Minimal norm a residual needs to have in order to be accepted as a new subspace vector (defaults to 2 * len(matrix) * machine_expsilon)
- adcc.solver.lanczos.amend_true_residuals(state, subspace, rvals, rvecs, epair_mask)
Compute the true residuals and residual norms (and not the ones estimated from the Lanczos subspace) and amend the state accordingly.
- adcc.solver.lanczos.compute_true_residuals(subspace, rvals, rvecs, epair_mask)
Compute the true residuals and residual norms (and not the ones estimated from the Lanczos subspace).
- adcc.solver.lanczos.default_print(state, identifier, file=<_io.TextIOWrapper name='<stdout>' mode='w' encoding='utf-8'>)
A default print function for the lanczos callback
- adcc.solver.lanczos.lanczos(matrix, guesses, n_ep, max_subspace=None, conv_tol=1e-09, which='LM', max_iter=100, callback=None, debug_checks=False, explicit_symmetrisation=<class 'adcc.solver.explicit_symmetrisation.IndexSymmetrisation'>, min_subspace=None)
Lanczos eigensolver for ADC problems
- Parameters
matrix – ADC matrix instance
guesses (list) – Guess vectors (fixes also the Lanczos block size)
n_ep (int) – Number of eigenpairs to be computed
max_subspace (int or NoneType, optional) – Maximal subspace size
conv_tol (float, optional) – Convergence tolerance on the l2 norm squared of residuals to consider them converged
which (str, optional) – Which eigenvectors to converge to (e.g. LM, LA, SM, SA)
max_iter (int, optional) – Maximal number of iterations
callback (callable, optional) – Callback to run after each iteration
debug_checks (bool, optional) – Enable some potentially costly debug checks (Loss of orthogonality etc.)
explicit_symmetrisation (optional) – Explicit symmetrisation to use after orthogonalising the subspace vectors. Allows to correct for loss of index or spin symmetries during orthogonalisation (type or instance).
min_subspace (int or NoneType, optional) – Subspace size to collapse to when performing a thick restart.
- adcc.solver.lanczos.lanczos_iterations(iterator, n_ep, min_subspace, max_subspace, conv_tol=1e-09, which='LA', max_iter=100, callback=None, debug_checks=False, state=None)
Drive the Lanczos iterations
- Parameters
iterator (LanczosIterator) – Iterator generating the Lanczos subspace (contains matrix, guess, residual, Ritz pairs from restart, symmetrisation and orthogonalisation)
n_ep (int) – Number of eigenpairs to be computed
min_subspace (int) – Subspace size to collapse to when performing a thick restart.
max_subspace (int) – Maximal subspace size
conv_tol (float, optional) – Convergence tolerance on the l2 norm squared of residuals to consider them converged
which (str, optional) – Which eigenvectors to converge to (e.g. LM, LA, SM, SA)
max_iter (int, optional) – Maximal number of iterations
callback (callable, optional) – Callback to run after each iteration
debug_checks (bool, optional) – Enable some potentially costly debug checks (Loss of orthogonality etc.)
- adcc.solver.power_method.default_print(state, identifier, file=<_io.TextIOWrapper name='<stdout>' mode='w' encoding='utf-8'>)
A default print function for the power_method callback
- adcc.solver.power_method.power_method(A, guess, conv_tol=1e-09, max_iter=70, callback=None, explicit_symmetrisation=<class 'adcc.solver.explicit_symmetrisation.IndexSymmetrisation'>)
Use the power iteration to solve for the largest eigenpair of A.
The power method is a very simple diagonalisation method, which solves for the (by magnitude) largest eigenvalue of the matrix A.
- Parameters
A – Matrix object. Only the @ operator needs to be implemented.
guess – Matrix used as a guess
conv_tol (float) – Convergence tolerance on the l2 norm of residuals to consider them converged.
max_iter (int) – Maximal numer of iterations
callback – Callback function called after each iteration
explicit_symmetrisation – Explicit symmetrisation to perform during iteration to ensure obtaining an eigenvector with matching symmetry criteria.
Adc Equations
- adcc.adc_pp.modified_transition_moments(method, ground_state, dipole_operator=None, intermediates=None)
Compute the modified transition moments (MTM) for the provided ADC method with reference to the passed ground state.
- Parameters
method (adc.Method) – Provide a method at which to compute the MTMs
ground_state (adcc.LazyMp) – The MP ground state
dipole_operator (adcc.OneParticleOperator or list, optional) – Only required if different dipole operators than the standard dipole operators in the MO basis should be used.
intermediates (adcc.Intermediates) – Intermediates from the ADC calculation to reuse
- Returns
- Return type
adcc.AmplitudeVector or list of adcc.AmplitudeVector
- adcc.adc_pp.state2state_transition_dm(method, ground_state, amplitude_from, amplitude_to, intermediates=None)
Compute the state to state transition density matrix state in the MO basis using the intermediate-states representation.
- Parameters
method (str, AdcMethod) – The method to use for the computation (e.g. “adc2”)
ground_state (LazyMp) – The ground state upon which the excitation was based
amplitude_from (AmplitudeVector) – The amplitude vector of the state to start from
amplitude_to (AmplitudeVector) – The amplitude vector of the state to excite to
intermediates (adcc.Intermediates) – Intermediates from the ADC calculation to reuse
- adcc.adc_pp.state_diffdm(method, ground_state, amplitude, intermediates=None)
Compute the one-particle difference density matrix of an excited state in the MO basis.
- Parameters
method (str, AdcMethod) – The method to use for the computation (e.g. “adc2”)
ground_state (LazyMp) – The ground state upon which the excitation was based
amplitude (AmplitudeVector) – The amplitude vector
intermediates (adcc.Intermediates) – Intermediates from the ADC calculation to reuse
- adcc.adc_pp.transition_dm(method, ground_state, amplitude, intermediates=None)
Compute the one-particle transition density matrix from ground to excited state in the MO basis.
- Parameters
method (str, AdcMethod) – The method to use for the computation (e.g. “adc2”)
ground_state (LazyMp) – The ground state upon which the excitation was based
amplitude (AmplitudeVector) – The amplitude vector
intermediates (adcc.Intermediates) – Intermediates from the ADC calculation to reuse
- adcc.adc_pp.matrix.block(ground_state, spaces, order, variant=None, intermediates=None)
Gets ground state, potentially intermediates, spaces (ph, pphh and so on) and the perturbation theory order for the block, variant is “cvs” or sth like that.
It is assumed largely, that CVS is equivalent to mp.has_core_occupied_space, while one would probably want in the long run that one can have an “o2” space, but not do CVS.
Utilities
- adcc.banner(colour=True)
Return a nice banner describing adcc and its components
The returned string contains version information, maintainer emails and references.
- Parameters
colour (bool) – Should colour be used in the print out