Full adcc reference

Note

Work in progress. Many function do not yet follow the numpy standard in their documentation!

This page contains a structured overview of the python API of adcc. See also the Index. 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, solver_method=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, n_core_orbitals=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:

  1. 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.
  2. A adcc.ReferenceState object
  3. A adcc.LazyMp object
  4. A 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 SCF tolerance / 100, whatever is larger)
  • solver_method (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.
Other Parameters:
 
  • max_subspace (int, optional) – Maximal subspace size
  • max_iter (int, optional) – Maximal number of iterations
Returns:

An adcc.ExcitedStates object containing the adcc.AdcMatrix, the adcc.LazyMp ground state and the adcc.ReferenceState as well as computed eigenpairs.

Return type:

ExcitedStates

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)
__init__(data, method=None, property_method=None)

Construct an ExcitedStates class from some data obtained from an interative solver.

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.
  • 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()

Return a string providing a human-readable description of the class

oscillator_strengths

List of oscillator strengths of all computed states

plot_spectrum(broadening='lorentzian', xaxis='eV', yaxis='cross_section', width=0.01, **kwargs)

One-shot plotting function for the spectrum generated by all states known to this class.

Makes use of the adcc.visualisation.ExcitationSpectrum class in order to generate and format the spectrum to be plotted, using many sensible defaults.

Parameters:
  • broadening (str or None or callable, optional) – The broadening type to used for the computed excitations. A value of None disables broadening any other value is passed straight to adcc.visualisation.ExcitationSpectrum.broaden_lines().
  • xaxis (str) – Energy unit to be used on the x-Axis. Options: [“eV”, “au”, “nm”, “cm-1”]
  • yaxis (str) – Quantity to plot on the y-Axis. Options are “cross_section”, “osc_strength”, “dipole” (plots norm of transition dipole).
  • width (float, optional) – Gaussian broadening standard deviation or Lorentzian broadening gamma parameter. The value should be given in atomic units and will be converted to the unit of the energy axis.
property_method

The method used to evaluate ADC properties

state_diffdms

List of difference density matrices of all computed states

state_dipole_moments

List of state dipole moments

state_dms

List of state density matrices of all computed states

timer

Return a cumulative timer collecting timings from the calculation

transition_dipole_moments

List of transition dipole moments of all computed states

transition_dms

List of transition density matrices of all computed states

Visualisation

class adcc.visualisation.ExcitationSpectrum(energy, intensity)
__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(*args, style=None, **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, mp_results)
__init__(method, mp_results)

Initialise an ADC matrix from a method, the reference_state and appropriate MP results.

block_spaces(self: libadcc.AdcMatrix, arg0: str) → List[str]
blocks
compute_apply(self: libadcc.AdcMatrix, arg0: str, arg1: libadcc.Tensor, arg2: libadcc.Tensor) → None
compute_matvec(in_ampl, out_ampl=None)

Compute the matrix-vector product of the ADC matrix with an excitation amplitude and return the result in the out_ampl if it is given, else the result will be returned.

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

dense_basis(blocks=None)

Return the list of indices and their values of the dense basis representation

diagonal(self: libadcc.AdcMatrix, arg0: str) → libadcc.Tensor
ground_state
has_block(self: libadcc.AdcMatrix, arg0: str) → bool
intermediates
is_core_valence_separated
matvec(v)
mospaces
ndim
reference_state
rmatvec(v)
shape
timer

Obtain the timer object of this class.

to_dense_matrix(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']
name
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)
__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

  1. 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.
  2. A range: The orbital indices given by this range will be put into the orbital subspace.
  3. An explicit list of orbital indices to be placed into the subspace.
  4. 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 the adcc::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.

backend

The identifier of the back end used for the SCF calculation.

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.

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.

conv_tol

SCF convergence tolererance

density

Return the Hartree-Fock density in the MO basis

dipole_moment

Return the HF dipole moment of the reference state (that is the sum of the electronic and the nuclear contribution.)

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.

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.

irreducible_representation

Reference state irreducible representation

mospaces

The MoSpaces object supplied on initialisation

n_alpha

Number of alpha electrons

n_beta

Number of beta electrons

n_orbs

Number of molecular orbitals

n_orbs_alpha

Number of alpha orbitals

n_orbs_beta

Number of beta orbitals

nuclear_dipole
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

restricted

Return whether the reference is restricted or not.

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)

timer

Obtain the timer object of this class.

class adcc.LazyMp(hf, caching_policy=<adcc.caching_policy.DefaultCachingPolicy object>)
__init__(hf, caching_policy=<adcc.caching_policy.DefaultCachingPolicy object>)

Initialise the class dealing with the M/oller-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(self: libadcc.LazyMp, arg0: str) → libadcc.Tensor

Obtain the Delta Fock matrix.

dipole_moment(level=2)

Return the MP dipole moment at the specified level of perturbation theory.

energy_correction(self: libadcc.LazyMp, arg0: int) → float

Obtain the appropriate MP energy correction.

has_core_occupied_space
mospaces
mp2_density
mp2_diffdm

Return the MP2 differensce density in the MO basis.

mp2_dipole_moment
reference_state
set_caching_policy(self: libadcc.LazyMp, arg0: CachingPolicy_i) → None
set_t2(self: libadcc.LazyMp, arg0: str, arg1: libadcc.Tensor) → None

Set the T2 amplitudes (invalidates dependent data automatically.

t2(self: libadcc.LazyMp, arg0: str) → libadcc.Tensor

Obtain the T2 amplitudes.

t2eri(self: libadcc.LazyMp, arg0: str, arg1: str) → libadcc.Tensor

Obtain a cached T2 tensor with ERI tensor contraction.

td2(self: libadcc.LazyMp, arg0: str) → libadcc.Tensor

Obtain the T^D_2 term.

timer

Obtain the timer object of this class.

Tensor and symmetry interface

class adcc.Tensor(sym_or_mo, space=None, irreps_allowed=None, permutations=None, spin_block_maps=None, spin_blocks_forbidden=None)
__init__(sym_or_mo, space=None, irreps_allowed=None, permutations=None, spin_block_maps=None, spin_blocks_forbidden=None)

Construct an uninitialised Tensor from an MoSpaces or a Symmetry 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.
  • irreps_allowed (list, optional) – List of allowed irreducible representations.
  • 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 the AdcMatrix, the LazyMp, the ReferenceState 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 transforms like the irrep “A”, 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", irreps_allowed=["A"],
...        spin_block_maps=[("aa", "bb", -1)],
...        spin_blocks_forbidden=["ab", "ba"])
add_linear_combination(self: libadcc.Tensor, arg0: numpy.ndarray[float64], arg1: list) → libadcc.Tensor

Add a linear combination of tensors to this tensor

antisymmetrise_to(self: libadcc.Tensor, arg0: libadcc.Tensor, arg1: list) → None
copy(self: libadcc.Tensor) → libadcc.Tensor

Returns a deep copy of the tensor.

copy_to(self: libadcc.Tensor, arg0: libadcc.Tensor) → None

Writes a deep copy of the tensor to another tensor

describe_symmetry(self: libadcc.Tensor) → str

Return a string providing a hopefully discriptive rerpesentation of the symmetry information stored inside the tensor.

dot(*args, **kwargs)

Overloaded function.

  1. dot(self: libadcc.Tensor, arg0: libadcc.Tensor) -> float
  2. dot(self: libadcc.Tensor, arg0: list) -> array
empty_like(self: libadcc.Tensor) → libadcc.Tensor
is_allowed(self: libadcc.Tensor, arg0: tuple) → bool

Is a particular index allowed by symmetry

mutable
ndim
nosym_like(self: libadcc.Tensor) → libadcc.Tensor
ones_like(self: libadcc.Tensor) → libadcc.Tensor
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.

  1. set_from_ndarray(self: libadcc.Tensor, arg0: array) -> None

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.

  1. set_from_ndarray(self: libadcc.Tensor, arg0: numpy.ndarray[float64], arg1: float) -> None

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) → None

Set all tensor elements to random data, adhering to the internal symmetry.

shape
size
symmetrise_to(self: libadcc.Tensor, arg0: libadcc.Tensor, arg1: list) → None
to_ndarray(self: libadcc.Tensor) → numpy.ndarray[float64]

Export the tensor data to a standard np::ndarray by making a copy.

transpose(*args, **kwargs)

Overloaded function.

  1. transpose(self: libadcc.Tensor) -> libadcc.Tensor
  2. transpose(self: libadcc.Tensor, arg0: tuple) -> libadcc.Tensor
zeros_like(self: libadcc.Tensor) → libadcc.Tensor
class adcc.Symmetry(mospaces, space, irreps_allowed=None, permutations=None, spin_block_maps=None, spin_blocks_forbidden=None)
clear(self: libadcc.Symmetry) → None

Clear the symmetry.

describe(self: libadcc.Symmetry) → str

Return a descriptive string.

empty

Is the symmetry empty (i.e. noy symmetry setup)

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.

mospaces

Return the MoSpaces object supplied on initialisation

ndim

Return the number of dimensions.

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.

shape

Return the shape of tensors constructed from this symmetry.

space

Return the space supplied on initialisation.

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.

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(*tensors)
__init__(*tensors)

Initialise an AmplitudeVector from some blocks

add_linear_combination(scalars, others)

Return an AmplitudeVector of the same shape and symmetry with all elements set to zero

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

nosym_like()

Return an empty AmplitudeVector of the same shape and symmetry

ones_like()

Return an empty AmplitudeVector of the same shape and symmetry

to_cpp()

Return the C++ equivalent of this object. This is needed at the interface to the C++ code.

zeros_like()

Return an AmplitudeVector of the same shape and symmetry with all elements set to zero

adcc.add(a, b, out=None)

Return the elementwise sum of two objects If out is given the result will be written to the latter tensor.

adcc.contract(contraction, a, b, out=None)

Form a single, einsum-like contraction, that is contract tensor a and be to form out via a contraction defined by the first argument string, e.g. “ab,bc->ac” or “abc,bcd->ad”.

Note: The contract function is experimental. Its interface can change
and the function may disappear in the future.
adcc.copy(a)

Return a copy of the input tensor.

adcc.divide(a, b, out=None)

Return the elementwise division of two objects If out is given the result will be written to the latter tensor.

Note: If out is not given, the symmetry of the contained objects will be destroyed!

adcc.dot(a, b)

Form the scalar product between two tensors.

adcc.empty_like(a)

Return an empty tensor of the same shape and symmetry as the input tensor.

adcc.linear_combination(coefficients, tensors)

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.

adcc.multiply(a, b, out=None)

Return the elementwise product of two objects If out is given the result will be written to the latter tensor.

Note: If out is not given, the symmetry of the contained objects will be destroyed!

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.subtract(a, b, out=None)

Return the elementwise difference of two objects If out is given the result will be written to the latter tensor.

adcc.transpose(a, axes=None)

Return the transpose of a tensor as a copy. If axes is not given all axes are reversed. Else the axes are expect as a tuple of indices, e.g. (1,0,2,3) will permute first two axes in the returned tensor.

adcc.zeros_like(a)

Return a zero tensor of the same shape and symmetry as the input tensor.

Solvers

adcc.solver.conjugate_gradient.conjugate_gradient(matrix, rhs, x0=None, conv_tol=1e-09, max_iter=100, callback=None, Pinv=None, cg_type='polak_ribiere', explicit_symmetrisation=<class 'adcc.solver.explicit_symmetrisation.IndexSymmetrisation'>)

An implementation of the conjugate gradient algorithm.

This algorithm implements the “flexible” conjugate gradient using the Polak-Ribière formula, but allows to employ the “traditional” Fletcher-Reeves formula as well. It solves matrix @ x = rhs for x by minimising the residual matrix @ x - rhs.

Parameters:
  • matrix – Matrix object. Should be an ADC matrix.
  • rhs – Right-hand side, source.
  • x0 – Initial guess
  • conv_tol (float) – Convergence tolerance on the l2 norm of residuals to consider them converged.
  • max_iter (int) – Maximum number of iterations
  • callback – Callback to call after each iteration
  • Pinv – Preconditioner to A, typically an estimate for A^{-1}
  • cg_type (string) – Identifier to select between polak_ribiere and fletcher_reeves
  • explicit_symmetrisation – Explicit symmetrisation to perform during iteration to ensure obtaining an eigenvector with matching symmetry criteria.
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)

@param matrix Matrix to diagonalise @param state DavidsonState containing the eigenvector guess

to propagate

@param max_subspace Maximal subspace size @param max_iter Maximal numer of iterations @param n_ep Number of eigenpairs to be computed @param is_converged Function to test for convergence @param callback Callback to run after each iteration @param which Which eigenvectors to converge to.

Needs to be compatible with the selected preconditioner.

@param preconditioner Preconditioner (type or instance) @param preconditioning_method Precondititoning method. Valid values are

“Davidson” or “Sleijpen-van-der-Vorst”
@param debug_checks Enable some potentially costly debug checks
(loss of orthogonality in subspace etc)
@param residual_min_norm Minimal norm a residual needs to have in order
to be accepted as a new subspace vector (defaults to 2 * len(matrix) * machine_expsilon)
@param explicit_symmetrisation Explicit symmetrisation to perform
on new subspace vectors before adding them to the subspace.
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

@param matrix ADC matrix instance @param guesses Guess vectors (fixes the block size) @param n_ep Number of eigenpairs to be computed @param max_subspace Maximal subspace size @param conv_tol Convergence tolerance on the l2 norm of residuals

to consider them converged
@param which Which eigenvectors to converge to.
Needs to be chosen such that it agrees with the selected preconditioner.

@param max_iter Maximal numer of iterations @param callback Callback to run after each iteration @param preconditioner Preconditioner (type or instance) @param preconditioning_method Precondititoning method. Valid values are

“Davidson” or “Sleijpen-van-der-Vorst”
@param 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)
@param debug_checks Enable some potentially costly debug checks
(Loss of orthogonality etc.)
@param residual_min_norm 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.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.

Properties

adcc.modified_transition_moments.compute_modified_transition_moments(gs_or_matrix, dipole_operator, method=None)

Compute the modified transition moments (MTM) for the provided ADC method with reference to the passed ground state and the appropriate dipole integrals in the MO basis.

Parameters:
  • gs_or_matrix – The MP ground state or ADC matrix for which level of theory the modified transition moments are to be computed
  • dipole_operator (OneParticleOperator) – Electric dipole operator
  • method (optional) – Provide an explicit method to override the automatic selection. Only required if LazyMp is provided as gs_or_matrix
Returns:

Return type:

AmplitudeVector

State analysis

TODO

Other stuff and utilities

adcc.banner(colour=True, show_doi=True, show_website=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
  • show_doi (bool) – Should DOI and publication information be printed.
  • show_website (bool) – Should a website for each project be printed.