How to compute Two-photon absorption cross-sections with ADC(2)/ADC(3) methods?

Hi all,

I am currently using Qchem 6.0 and according to the manual it is possible to compute two-photon absorption cross-sections (TPAs) with ADC methods through the matrix inversion technique using the keyword ADC_PROP_TPA, however in the output I don’t find this information. If I include the keyword ADC_PROP_ES2ES I am able to recover TPAs but computed with the Sum over states expressions which is not the quantity I am interested in.

This is the input that I am using:

$molecule
0 1
O 0.000000 0.000000 -0.069903
H 0.000000 0.757532 0.518435
H 0.000000 -0.757532 0.518435
$end
$rem
METHOD adc(2)
BASIS aug-cc-pvdz
AUX_BASIS rimp2-aug-cc-pvdz
ADC_PROP_ES True
ADC_PROP_TPA True
EE_STATES [1,1,0,1]
$end
and this is the part of the output related to one of excited states:
Excited state 1 (singlet, B2) [converged]

Term symbol: 1 (1) B2 R^2 = 7.51604e-07

Total energy: -76.0072675789 a.u.
Excitation energy: 6.967051 eV

Osc. strength: 0.054586
Trans. dip. moment [a.u.]: [ -0.000000, -0.565506, -0.000000]
<i|r^2|0> [a.u.]: [ 0.000000, -0.000000, -0.000000]

Dip. moment [a.u.]: [ 0.000000, 0.000000, 0.690673]
Total dipole [Debye]: 1.755515
<r^2> [a.u.]: [ 15.147806, 10.478238, 13.978771]
Total <r^2> [a.u.]: 39.604816

V1^2 = 0.9292, V2^2 = 0.0708

Important amplitudes:
occ i occ j vir a vir b v
---------------------------------------------------------
1 (B2) A 4 (A1) A 0.6139
1 (B2) A 6 (A1) A 0.2067
---------------------------------------------------------

As you can see there is no trace of TPAs, do you know if I am doing anything wrong?

Best Regards,
Carmelo

I am not familiar with this methodology but I found the following sample input, whose output prints the two-photon abs. matrix without the sum-over-states message that appears in other jobs:

$comment
TPA cross-sections of H2O (direct)
$end

$molecule
0 1
O         0.000000    0.000000    0.146192
H        -0.748943    0.000000   -0.461081
H         0.748943    0.000000   -0.461081
$end

$rem
jobtype = sp
method = adc(2)
n_frozen_core = 0
basis = cc-pvdz
ee_singlets = [2,0,0,0]
adc_prop_tpa = true
$end

Output looks like:

  Excited state 1 (singlet, A1)                                    [converged]
  ----------------------------------------------------------------------------
  Term symbol:  2 (1) A1                                    R^2 =  1.96980e-08

  Total energy:                                            -75.8347960040 a.u.
  Excitation energy:                                              10.780969 eV

  Osc. strength:                                                      0.101837
  Trans. dip. moment [a.u.]:           [  -0.000000,    0.000000,    0.620933]
  <i|r^2|0> [a.u.]:                    [  -0.534828,   -0.104117,   -0.396852]
  Two-photon absorption cross-section [a.u.]:                       594.495795
  Two-photon absorption matrix [a.u.]:
                                       |   -4.896661   -0.000000   -0.000000 |
                                       |   -0.000000   -0.786115    0.000000 |
                                       |   -0.000000    0.000000   -6.795354 |

Hi jhelbert,

thank you very much for your useful answer! Indeed in this way it works, it seems that the problem is the employment of the resolution of identity which prevents the printing of these informations. Do you think it is a bug or maybe they haven’t implemented it yet?

Morevoer the same keywords don’t work for adc(3), which in principle should be possible. Did you find an example also for adc(3)?

Best
Carmelo

That is the only example that I found, in the samples directory that should ship with Q-Chem (sp_adc_tpa.in). Could be this is not yet implemented with RI or for ADC(3), please check the manual and if it’s not clear from there you could contact Q-Chem support.

Hi jhelbert, I am checking again this problem and I realized that even for adc2 there is something a bit weird for the two-photon absorption, in fact even in the output that you sent me the TPA is very large and not coherent with other values in the literature (see the paper of the implementation of this method on qchem Calculations of nonlinear response properties using the intermediate state representation and the algebraic-diagrammatic construction polarization propagator approach: Two-photon absorption spectra | The Journal of Chemical Physics | AIP Publishing). I am wondering if it is a conversion issue or a bug, do you know to whom I should talk about it?

Best,
Carmelo

Please email Q-Chem customer support.

My understanding form previously raising such queries is that TPA is only implemented (currently) for ADC(2), and (as mentioned above) not for Resolution of Identity calculations

Doing TPA with resolution of the identity does not work (you can only use the approximate sum-over-state version via ADC_PROP_ES2ES).

TPA cross-sections of organic molecules are on the order of several hundred [a.u.]. This is also seen on Table V of the paper you referred to. So, the job the John shared looks alright.

To convert the number in a.u. to a cross-sections in GM you have to do

GM = 0.00313184 * TPA / 30 * dE

Here TPA is the cross section printed (in this case 594.496) and dE is the excitation energy in eV.

Hi Felix, thanks for your comments! However, this value doesn’t look correct to me. In fact if you check in this same paper in Table III for water you can see that the values of delta^{TPA} are always lower than 100 for state 2A1 which is different from the output of John. I repeated the calculations with the same structure and aug-cc-pvdz to directly compare with the values in the Table in the paper to show to you:

Excited state 3 (singlet, A1) [converged]

Term symbol: 2 (1) A1 R^2 = 2.69827e-08

Total energy: -75.9196501516 a.u.
Excitation energy: 9.351234 eV

Osc. strength: 0.098405
Trans. dip. moment [a.u.]: [ -0.000000, 0.000000, 0.655382]
<i|r^2|0> [a.u.]: [ 0.388866, -0.248583, 0.458807]
Two-photon absorption cross-section [a.u.]: 413.857055

So as you can see the energy is well reproduced but the TPA is not, and it is the same for the other states. So I believe there is some kind of issue with the units or some error in the paper… Could you please let me know what do you think about that?

Be aware that ADC(2) two-photon cross section (a.u.) reported by the program does not include a necessary division by the factor 30, which you will need if comparing with the literature and other methods (e.g. EOM-CCSD in QChem)

It may help to compare against the McClain formalism for converting the TPA tensor into a rotationally averaged cross section (J.Chem. Phys 53 (1970) 29) which is used in these calculations. Its not too difficult to convert the TPA tensor S_ij directly to the cross-section (a.u) using the formulae there (which explicitly show the division by 30). This will use quantities Df - sum{S_ii_S_jj}/30 and Dg=sum{S_ijS_ij}/30.
Then for linear polarization sigma=2Df+4Dg (a.u)

If you compare with e.g. an EOM-CCSD two photon calculation in QChem you will see that there the linear polarization cross-section is computed correctly with the factor 1/30 included (although then there was a trivial error in the formula used for circular polarization cross-section in version <5.4 which was quoted as sigma=-Df+6Dg but according to the McClain reference should be sigma=-2Df+6*Dg ).

But at least you will need to divide the printed cross-section from ADC(2) calculation by 30 to compare with the literature values.