Plotting Spectrum from ADCC/Psi4 Calculation

[sbembenek]

I’m running the example shown here: https://github.com/psi4/psi4/blob/master/samples/adcc/h2o-cvs-adc2/input.dat

It runs with no issues.

But I am also interested in plotting the spectrum (as described, from the adcc side, here: https://github.com/psi4/psi4/blob/master/samples/adcc/h2o-cvs-adc2/input.dat) which is said to be accessible in psi4 via the object ‘adcc_state’ as described here: ADC: Ab Initio Polarization Propagator

Pulling it altogether, I added to the example this bit:

# Plot spectrum
import matplotlib.pyplot as plt
import adcc
plt.figure(figsize=(4, 5))
adcc_state.plot_spectrum()
plt.tight_layout()
plt.savefig('spectrum.png')

Giving the overall code:

#! CVS-ADC(2)/cc-pvdz calculation of 10 water singlet excited states
#! tackling the Oxygen 1s edge core exitations

molecule h2o {
    O 0 0 0
    H 0 0 1.795239827225189
    H 1.693194615993441 0 -0.599043184453037
    symmetry c1
    units au
}

set {
    reference rhf
    basis cc-pvdz
    guess core
    num_core_orbitals 1
    roots_per_irrep [10]
}

set_num_threads(62)

# Run normal calculation of excitation energies (no properties)
energy_adc, wfn = energy('cvs-adc(2)', return_wfn=True)

# Alternative: Run computation of properties as well
properties('cvs-adc(2)', properties=["oscillator_strength", "dipole"])


# Plot spectrum
import matplotlib.pyplot as plt
import adcc
plt.figure(figsize=(4, 5))
adcc_state.plot_spectrum()
plt.tight_layout()
plt.savefig('spectrum.png')

But I get the error:

Printing out the relevant lines from the Psithon --> Python processed input file:
    energy_adc, wfn = energy('cvs-adc(2)', return_wfn=True)
    properties('cvs-adc(2)', properties=["oscillator_strength", "dipole"])
    import matplotlib.pyplot as plt
    import adcc
    plt.figure(figsize=(4, 5))
--> adcc_state.plot_spectrum()
    plt.tight_layout()
    plt.savefig('spectrum.png')

!-----------------------------------!
!                                   !
!  name 'adcc_state' is not defined !
!                                   !
!-----------------------------------!

Which I don’t understand since:

–>The following attribute is set on returned wavefunctions:

• adcc_state: The adcc.ExcitedStates object used by adcc to store the ADC(n) excitation energies and all precomputed data in the format used by adcc. Provides direct access to analysis and plotting capabilities from adcc. For example adcc_state.plot_spectrum() plots a broadened excited states spectrum in matplotlib. See the adcc calculations documentation for details.

What am I missing?

[jmisiewicz]

What you’re missing is the first line of the documentation snippet you posted. adcc_state is an attribute on a returned wavefunction. When the documentation says adcc_state.plot_spsectrum(), it’s shorthand for whatever-you-named-your-wfn-object.adcc_state.plot_spectrum(). It doesn’t create some free-floating variable named adcc_state.

This is standard shorthand for an attribute of an object that documentation hasn’t named.

[sbembenek]

The input file is:

import adcc
import matplotlib.pyplot as plt

molecule h2o {
O     0.000000     0.000000     0.065441
H    -0.763619    -0.000000    -0.519300
H     0.763619     0.000000    -0.519300
    symmetry c1
    units au
}

set {
    guess core
    num_core_orbitals 1
    roots_per_irrep [10]
}

set scf_type df
set basis aug-cc-pVTZ
set reference rhf
set s_tolerance 1e-9
set_num_threads(62)

# Run normal calculation of excitation energies
scf_e, wfn = psi4.energy('scf', return_wfn=True)

# Alternative: Run computation of properties
properties('cvs-adc(2)', properties=["oscillator_strength", "dipole"])

plt.figure(figsize=(4, 5))
wfn.adcc_state.plot_spectrum()
plt.tight_layout()
plt.savefig('spectrum.png')

So, the wavefunction object is ‘wfn’ – yes? So, then →

wfn.adcc_state.plot_spectrum()

I must still be missing something … didn’t work for me.

What did work was:

state = adcc.cvs_adc3(wfn, n_singlets=10, core_orbitals=1)
state.plot_spectrum()

Thanks for the help!

[jmisiewicz]

wfn in your example comes from the SCF computaiton. It of course knows nothing about ADC.

[sbembenek]

My bad – I was bouncing back and forth between the Python API.

Yes, that works, here’s the whole thing:

#! tackling the Oxygen 1s edge core excitations

import matplotlib.pyplot as plt

molecule h2o {
O     0.000000     0.000000     0.065441
H    -0.763619    -0.000000    -0.519300
H     0.763619     0.000000    -0.519300
    symmetry c1
    units au
}

set {
    guess core
    num_core_orbitals 1
    roots_per_irrep [10]
}

set scf_type df
set basis aug-cc-pVTZ
set reference rhf
set s_tolerance 1e-9
set_num_threads(64)

# Run normal calculation of excitation energies
scf_e, wfn = energy('cvs-adc(2)', return_wfn=True)

# Plot spectrum
plt.figure(figsize=(4, 5))
wfn.adcc_state.plot_spectrum()
plt.tight_layout()
plt.savefig('spectrum.png')

And then leveraging the Python API →

import adcc
import matplotlib.pyplot as plt

molecule h2o {
O     0.000000     0.000000     0.065441
H    -0.763619    -0.000000    -0.519300
H     0.763619     0.000000    -0.519300
    symmetry c1
    units au
}

set scf_type df
set basis aug-cc-pVTZ
set reference rhf
set s_tolerance 1e-9
set_num_threads(64)

print(adcc.banner())

# Run normal SCF calculation for ADCC
scf_e, wfn = psi4.energy('scf', return_wfn=True)

# Run state calculation
state = adcc.cvs_adc2(wfn, n_singlets=10, core_orbitals=1)

# Plot spectrum
plt.figure(figsize=(4, 5))
state.plot_spectrum()
plt.tight_layout()
plt.savefig('spectrum.png')

Thanks!