IP/EA-ADC up to 3rd order

[fedy9]

• added IP/EA ADC folders like adc_pp
• added type attribute to ADC matrix to distinguish between IP/EA/PP (e.g. matrix.type="ip"). Accessible wherever it's necessary
• IP/EA ADC matrices up to 3rd order yielding excitation energies
• added pole strengths attribute to ExcitedStates class in analogy to oscillator strengths
• IP/EA pole strengths up to 3rd order computing f amplitudes as intermediate
• added state_diffdm's and s2s transition dm's up to 2nd order for IP/EA
• new guess setup for IP/EA: appropriate estimate of guesses, symmetry setup, guess construction on C++ side
• generalisation of multiple PP specific routines, e.g. the singles block is selected dynamically (p/h/ph) and no longer automatically "ph"
• New and adapted error handling
• modified 'validate_state_parameters()' function in workflow to check for valid IP/EA inputs and returns a new Boolean 'is_alpha' which is None for PP and specifies whether an alpha or beta electron is attached/detached for EA/IP. Default value is "True". Added keyword 'n_doublets' specifying the doublet states for IP/EA
• spin_change calculation was updated for IP/EA case in workflow
• New enforce_spin_kind routines for IP/EA
• adapted output: during calc. (solver printing) and adaptations for the ExcitedStates methods 'describe()' and 'describe_amplitudes()'
• added total dipole moment property in ExcitedStates for all calculations (PP/IP/EA)
• Warning for state dipole moments of IP/EA states since dipole moments of charged molecules are gauge dependent

New methods on user side: ip_adc{0..3}() ea_adc{0..3}() New keywords on user side: n_doublets, is_alpha e.g. ip_adc3(scf_result, n_doublets=5, is_alpha=True)

Remarks: Energies and properties were verified against Q-Chem calculations.

Pole strengths of 3rd order were implemented for completeness. They are currently not used in the calculation since all properties are computed at most at 2nd order even for ADC(3) calculations. Hence they were not verified yet.

No simultaneous calculation of alpha and beta states to simplify the code and error handling and to keep consistency with PP calculations. If alpha and beta would be computed in the same adcc instance, intermediates could be reused. However, these intermediates are not the bottleneck of the calculation and thus recomputed for the other spin.

Quartets not yet implemented since they are pure doubles states which are inaccurate up to 3rd order.

To Do Testing procedures for IP/EA ADC Verify 3rd order pole strengths

Observed issue True for PP/IP/EA If 'kind="any"' and the user requests a significant amount of states getting somewhere near the number of total accessible states, the Davidson will quickly converge to 0-vector. This is not true if kind is specified. Example input

from pyscf import gto, scf
import adcc

# Run SCF in pyscf
mol = gto.M(
    atom='O 0.0 0.0 0.0;'
         'H 0.0 0.0 0.957;'
         'H 0.927 0.0 -0.240',
    basis='sto-3g',
    spin=0
)
scfres = scf.RHF(mol)
scfres.conv_tol = 1e-13
scfres.kernel()

state = adcc.adc2(scfres, n_states=65)

[fedy9]

Added ip/ea adc2x for consistency with pp basemethods

Added tests for IP/EA-ADC:

• workflow
• state_densities
• consistent properties
• guesses

Modified adcc testdata generation: the adcc hdf5 reference data have to be regenerated so the pp tests pass, since I moved "available_kinds" to the "adc" level

Modified read-in of reference data and creating MockAdcStates

Generated ip/ea reference data for h2o_sto3g and cn_sto3g Testing alpha and beta attachment/detachment for unrestricted cn_sto3g Only alpha for restricted h2o_sto3g

Removed some intermediates in ip and ea matrix.py to avoid conflicts and reuse code in pp matrix.py

Introduced a blacklist for excitation_properties that depend on the type (pp vs ip/ea)

All former tests are running (except: pcm-lr/ptlr using pyscf that fails with pyscf version 2.3.0, scf energy fails, pcm with ddx solver since not installed)

[fedy9]

Looks nice I think! If you can come up with some minimal test for each method, that'd be great!

And thanks for reviewing :)

Since I already added some tests, I'm not sure what you have in mind here. The obvious missing tests are testing against Q-Chem. However I struggled with that on the one hand with the interface and on the other hand due to differing different SCF solutions (psi4/pyscf vs qchem energies) that occured for some of the small set of open-shell test molecules leading to error propagation using more or less default inputs. But that should be doable.

[maxscheurer]

Since I already added some tests, I'm not sure what you have in mind here. The obvious missing tests are testing against Q-Chem. However I struggled with that on the one hand with the interface and on the other hand due to differing different SCF solutions (psi4/pyscf vs qchem energies) that occured for some of the small set of open-shell test molecules leading to error propagation using more or less default inputs. But that should be doable.

I was talking about testing against Q-Chem. That would be great, because then one can be sure nothing gets broken with the IP/EA feature.

The adcc-testdata setup is somewhat complicated, right. For this reason, for the PE-ADC checks against Q-Chem, I wrote this file which essentially just invokes Q-Chem, runs PE-ADC, and dumps the results to a YAML file. Maybe you could do the same for IP/EA-ADC, which would be totally fine for me. Then we'd have some minimal functionality tests that run on CI.

[fedy9]

I reuse some of the code that appears in all three matrix files by adding a new function covering that which can be called in the other matrix files

Further I added minimal testing against Q-Chem reference data, i.e. h2o-sto3g excitation energies and pole strengths for the first 2-3 states of IP/EA-ADC(2). I dumped the reference values in the qchem yml file

[jonasleitner]

I don't understand the current separation of methods between ElectronicTransition and ExcitedStates. Some methods on both classes are PP-ADC specific while others are general for PP/IP/EA-ADC. For instance, describe and the amplitude printout on ExcitedStates or transition_dipole_moment and oscillator_strength on ElectronicTransition are PP specific. I think it would make more sense to put all methods that are general for PP/IP/EA-ADC on the ElectronicTransition class and then have a PP-ADC specific ExcitedStates child class and another IonizedStates child class for IP/EA. And similarly, FormatExcitationVector should also be refactored using a generic base class and specific child classes for PP/IP/EA. What do you think @apapapostolou @maxscheurer @fedy9?

[apapapostolou]

I don't understand the current separation of methods between ElectronicTransition and ExcitedStates. Some methods on both classes are PP-ADC specific while others are general for PP/IP/EA-ADC. For instance, describe and the amplitude printout on ExcitedStates or transition_dipole_moment and oscillator_strength on ElectronicTransition are PP specific. I think it would make more sense to put all methods that are general for PP/IP/EA-ADC on the ElectronicTransition class and then have a PP-ADC specific ExcitedStates child class and another IonizedStates child class for IP/EA. And similarly, FormatExcitationVector should also be refactored using a generic base class and specific child classes for PP/IP/EA. What do you think @apapapostolou @maxscheurer @fedy9?

I agree with @jonasleitner. I also don't understand the current separation between ElectronicTransition and ExcitedStates. How was it decided which methods are implemented where @maxscheurer? I think it would make sense to have a parent class for methods that are needed for all variants (PP/IP/EA) and then more specific child classes.

[maxscheurer]

From the top of my head, it can be rationalized as follows: Electronic excitation refers to a single excitation and the corresponding properties, whereas ExcitedStates is kind of a collection of ElectronicExcitations (list of structs vs struct of lists). That way, specific excitations can be addressed individually, for example when computing the gradient or so. Let's put this on the agenda for the next meeting, and I will look thru the code again. Please try to play around with the features of both classes also. I think there is a simple function on ExcitedStates that converts it to a list of ElectronicExcitations.

[jonasleitner]

From the top of my head, it can be rationalized as follows: Electronic excitation refers to a single excitation and the corresponding properties, whereas ExcitedStates is kind of a collection of ElectronicExcitations (list of structs vs struct of lists). That way, specific excitations can be addressed individually, for example when computing the gradient or so. Let's put this on the agenda for the next meeting, and I will look thru the code again. Please try to play around with the features of both classes also. I think there is a simple function on ExcitedStates that converts it to a list of ElectronicExcitations.

I think you are confusing ElectronicTransition with Excitation. But Excitation will certainly also be affected when we restructure ElectronicTransition and ExcitedStates. So good point. Depending on the use case it might be a good idea to refactor Excitation anyway, because currently it computes e.g., the diff_dm of all states and then extracts the result of the state of interest, which is very convenient. But we don't gain any performance if we are only interested in the properties of a single state, which was kind of my initial expectation.

[maxscheurer]

Ups sorry 😂 I was away from my computer, so didn't check again... I need to think again why it was split up this way, it's been a while.