[molpro-user] icpks convergence problem

[Andreas Hesselmann]

Dear Glen,

the default convergence tolerance for CPKS
is set to 1d-6 for the density matrix.
I have submit a patch (sapt_icpks) with which
you can adjust this value (sapt_cpksthr) as
well as the maximum number of iterations
in the CPKS (which by default is set to 50)
in the input.

I would be sceptical about the result of a 
CPKS calculation if the convergence criteria
of 1d-6 is not met within 50 iterations.
You might check your result for the induction energy
by comparing it with the finite field result.

Here is an example for a finite field calculation
of the induction energy for the water dimer:


gthresh,energy=1d-11,orbital=1d-12
geometry={nosym;noorient;angstroms;
O1,, 0.0000000000, 0.0000000000, 0.0000000000
H1,, 0.0000000000, 0.0000000000, 0.9621100000
H2,, 0.0000000000, 0.9368337299,-0.2465760535
O2,, 0.0000000000, 2.7788154372,-0.9708520107
H3,,-0.7632231375, 2.9252149352,-1.5405904257
H4,, 0.7632231375, 2.9252149352,-1.5405904257}
basis=avdz
int

field=[0,0.001,-0.001]

!---Monomer A-----------------
dummy,o2,h3,h4

{hf; maxit,20; save,2151.2}

{matrop
load,d,den,2151.2
load,epot,epot
load,h0,h0
coul,J,d
add,v,1,epot,2,J
save,v,1211.2,square}

!---Monomer B-----------------
int
dummy,o1,h1,h2

{hf; start,atden; maxit,20; save,2152.2}
e(1)=energy

{matrop
load,d,den,2152.2
load,epot,epot
load,h0,h0
coul,J,d
add,v,1,epot,2,J
save,v,1212.2,square}


delete,1
!--->Induktion A->B
{matrop
load,v,square,1211.2
load,h0,h0
add,h01,h0,field(2),v
save,h01,,h0}

{hf; maxit,20; save,2153.2}
e(2)=energy

{matrop
load,d,den,2153.2
load,v,square,1211.2
trace,v1,d,v
set,v1,0.5*v1}

delete,1
{matrop
load,v,square,1211.2
load,h0,h0
add,h01,h0,field(3),v
save,h01,,h0}

{hf; maxit,20; save,2154.2}
e(3)=energy

{matrop
load,d,den,2154.2
load,v,square,1211.2
trace,v2,d,v
set,v2,0.5*v2}

show,v*
table,e
digits,12
!values in mH
e10_elst_A=-(e(2)-e(3))/(field(2)-field(3))*1000
e20_ind_A= -0.5*(e(2)+e(3)-2*e(1))/((field(2)-field(1))*(field(3)-field(1)))*1000
e20_ind_A=(v1-v2)/(field(2)-field(3))*1000

text,Energies in mH
table,e10_elst_A,e20_ind_A
digits,12

!--->Induktion B->A
int
dummy,o2,h3,h4

{hf; start,atden; maxit,20; save,2155.2}
e(1)=energy

delete,1
{matrop
load,v,square,1212.2
load,h0,h0
add,h01,h0,field(2),v
save,h01,,h0}

{hf; maxit,20; save,2156.2}
e(2)=energy

{matrop
load,d,den,2156.2
load,v,square,1212.2
trace,v1,d,v
set,v1,0.5*v1}

delete,1
{matrop
load,v,square,1212.2
load,h0,h0
add,h01,h0,field(3),v
save,h01,,h0}

{hf; maxit,20; save,2157.2}
e(3)=energy

{matrop
load,d,den,2157.2
load,v,square,1212.2
trace,v2,d,v
set,v2,0.5*v2}


table,e
digits,12
!values in mH
e10_elst_B=-(e(2)-e(3))/(field(2)-field(3))*1000
e20_ind_B= -0.5*(e(2)+e(3)-2*e(1))/((field(2)-field(1))*(field(3)-field(1)))*1000
e20_ind_B=(v1-v2)/(field(2)-field(3))*1000
e20_ind=e20_ind_A+e20_ind_B

text,Energies in mH
table,e20_ind_A,e20_ind_B,e20_ind
digits,8,8,8


Of course the convergence thresholds for the SCF
should be even higher in the finite-field case.

Best regards,
Andreas








--------------------------------------------------
Andreas Hesselmann
Institut für Physikalische und Theoretische Chemie
Universität Erlangen
Egerlandstraße 3
91058 Erlangen / Germany
Phone:  +49 9131/85-25021
E-Mail: andreas.hesselmann at chemie.uni-erlangen.de
-------------------------------------------------