[molpro-user] Curves not smooth

[Bhargava Anusuri]

Dear molpro users,
                      I've been computing two-state potential energy curves
of same symmetry for Li+ ion and NO system. For some r(NO) values the the
second state is not smooth with sharp inflections. I circumvented this
problem in some cases by increasing n value in 'state,n;'. But in some
cases it still persists. I tried incresing 'nstati' in CI but then the 3rd
root is getting involved and CI is not converging. Any help is greatly
appreciated. Given below is my input:
***,LiNO+  Diabatization and NACME calculation, noorient and x is removed,
rno fixed at 2.4 as in paper
memory,300,m
gthresh,energy=0.32d-6
gprint,orbitals


rli=16.0d0
conv=0.52917725d0
rliang=rli*conv
rnang=1.28d0*conv
roang=-1.12d0*conv

geomtyp=xyz

geometry
   3
LiNO+ 15 degrees
Li      rliang*sin(15)  0.00000000       rliang*cos(15)
N       0.00000000      0.00000000       rnang
O       0.00000000      0.00000000       roang
end

basis=cc-pVTZ

r=2.8d0
dr=[0,0.0002,-0.0002]                                              !Small
displacements for finite difference NACME calculation

reforb1=2140.2                                                     !Orbital
dumprecord at reference geometry
refci=6000.2                                                       !MRCI
record at reference geometry
savci=6100.2                                                       !MRCI
record at displaced geometries

text,compute wavefunction at reference geometry (C2v)

hf;wf,17,1,1;orbital,2100.2

multi;
maxiter,40;
pspace,0.05;
wf,17,1,1;state,2;                                    !X 1A1 and A 1A1
states,optimize two states of symmetry 1
natorb,reforb1                                                     !Save
reference orbitals on reforb1
noextra                                                            !Dont
use extra symmetries

ci;
option,maxit=40                                                   !MRCI at
reference geometry
option,maxiti=3000;
!option,nstati=5;
closed,5;
occ,10,3;
wf,17,1,1;state,2;
    !X 1A1 and A 1A1 states
orbital,reforb1
   !Use orbitals from previous CASSCF
save,refci;
   !Save MRCI wavefunction

Text,Displaced geometries

do i=1,42
   !Loop over different r values
data,truncate,savci+1
   !truncate dumpfile after reference
reforb=reforb1

do j=1,3
    !Loop over small displacements for NACME
rli=r+dr(j)
    !Set current rh
!--------------------------------------------------------------------
conv=0.52917725d0
rliang=rli*conv
rnang=1.28d0*conv
roang=-1.12d0*conv

geomtyp=xyz


geometry
   3
LiNO+ 15 degrees
Li      rliang*sin(15)  0.00000000       rliang*cos(15)
N       0.00000000      0.00000000       rnang
O       0.00000000      0.00000000       roang
end
!--------------------------------------------------------------------
multi;
maxiter,40;
pspace,0.7;
wf,17,1,1;state,2;                                         !Wavefunction
definition
start,reforb
 !Starting orbitals
orbital,3140.2+j;
!Dumprecord for orbitals
diab,reforb
!Generate diabatic orbitals relative to reference geometry
noextra
!Dont use extra symmetries
reforb=3141.2
!Use orbitals for j=1 as reference for j=2,3

ci;
option,maxit=40
option,maxiti=6000;
!option,nstati=5;
closed,5;
occ,10,3;
wf,17,1,1;state,2;
orbital,diabatic
!Use diabatic orbitals
save,savci+j;
 !Save MRCI for displaced geometries

eadia=energy
!Save adiabatic energies for use in ddr
if(j.eq.1) then
e1(i)=energy(1)
 !Save adiabatic energies for table printing
e2(i)=energy(2)
end if

ci;trans,savci+j,savci+j;
 !Compute transition densities at R2+DR(j)
dm,7000.2+j;
!Save transition densities on this record
ci;trans,savci+j,refci;
 !Compute transition densities between R2+DR(j) and R1
dm,7100.2+j;
!Save transition densities on this record
ci;trans,savci+j,savci+1;
 !Compute transition densities between R and R2+DR(j)
dm,7200.2+j;
!Save transition densities on this record

ddr
density,7000.2+j,7100.2+j
 !Densities for <R2+DR||R2+DR> and <R2+DR||R1>
orbital,3140.2+j,2140.2
 !Orbitals for <R2+DR||R2+DR> and <R2+DR||R1>
energy,eadia(1),eadia(2)
!Adiabatic energies
mixing,1.1,2.1
!Compute mixing angle and diabatic energies

if(j.eq.1) then
 !Store diabatic energies for R2  (DR(1)=0)
  mixci(i)=mixangci(1)
!Mixing angle obtained from ci vectors only
  h11ci(i)=hdiaci(1)
!Diabatic energies obtained from ci vectors only
  h21ci(i)=hdiaci(2)
!HDIA contains the lower triangle of the diabatic hamiltonian
  h22ci(i)=hdiaci(3)
  mixtot(i)=mixang(1)                                        !Mixing angle
from total overlap (including first-order correction)
  h11(i)=hdia(1)
!Diabatic energies obtained from total overlap
  h21(i)=hdia(2)
  h22(i)=hdia(3)
end if

mix(j)=mixang(1)
!Store mixing angles for R2+DR(j)

enddo
 !End loop over j

dchi(i)=(mix(3)-mix(2))/(dr(3)-dr(2))*pi/180
!Finite difference derivative of mixing angle

ddr
density,7201.2,7202.2,7203.2
 !Compute NACME using 3-point formula
orbital,3141.2,3142.2,3143.2
states,2.1,1.1
nacmeci(i)=nacme
capr(i)=r

c=2.4
p(i)=c*capr(i)/capr(i)

if(i.lt.7) then
r=r+0.2d0
elseif(i.ge.7).and.(i.lt.19) then
r=r+0.1d0
elseif(i.ge.19).and.(i.lt.33) then
r=r+0.2d0
elseif(i.ge.33).and.(i.lt.37) then
r=r+0.5d0
else
r=r+1.0d0
endif


enddo

{table,p,capr,mixci,mixtot,dchi,nacmeci
Title,Mixing angles and non-adiabatic coupling matrix elements for LiNO+
save,2.4n.tab}

table,p,capr,e1,e2,h11ci,h22ci,h21ci
Title,Diabatic energies for LiNO+, obtained from CI-vectors
!format,'(f10.2,6f14.8)'
!sort,1

{table,p,capr,e1,e2,h11,h22,h21
title,Diabatic energies for LiNO+, obtained from CI-vectors and orbital
correction
save,2.4e.tab}

Regards,
Bhargava Anusuri,
C/O Prof. Sanjay Kumar,
Chemistry Dept.,
IIT Madras.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.molpro.net/pipermail/molpro-user/attachments/20130914/3007b7ca/attachment.html>