<html>
<head>
<meta http-equiv="Content-Type" content="text/html; charset=iso-8859-1">
<style type="text/css" style="display:none"><!--P{margin-top:0;margin-bottom:0;} p
{margin-top:0;
margin-bottom:0}--></style>
</head>
<body dir="ltr" style="font-size:12pt;color:#000000;background-color:#FFFFFF;font-family:Calibri,Arial,Helvetica,sans-serif;">
<p>We are attempting to reproduce Figures 2 and 3 of <a href="http://dx.doi.org/10.1063/1.479214" target="_blank">
http://dx.doi.org/10.1063/1.479214</a> showing the CASSCF and MRCI potential energy surfaces (PESs) of H2S. Our input files are closely based on the h2s_diab1 example. The multi and adiabatic MRCI PES match the literature well. However, the diabatic MRCI
calculations do not match the paper and produce inconsistent results. With the script provided below, we obtain a conical intersection at 2.7 a0 with the MRCI calculation of the diabatized PES. Figure 2 of the included paper shows a lack of intersection
between the diabatic surfaces. Further, by reversing the order in which we take the non-fixed H-S bond length values we obtain a conic intersection at ~3.5 a0.
</p>
<p> </p>
<p>Similar issues arise with recreating Figure 3 of PES versus bond angle. Starting with small bond angles produces a conical intersection at ~110 degrees; starting with large bond angles has a conical intersection at ~95 degrees, relatively close to the literature.
<br>
</p>
<p><br>
</p>
<p>How should we fix this issue and consistently produce reliable diabatized MRCI surfaces?<br>
</p>
<p><br>
</p>
<p><br>
</p>
***,h2s Diabatization<br>
memory,3,m<br>
<br>
gprint,orbitals,civector<br>
<br>
symmetry,x<br>
orient,noorient !noorient should always be used for diabatization<br>
geometry={<br>
s;<br>
h1,s,r1;<br>
h2,s,r2,h1,theta}<br>
<br>
basis=aug-cc-pVQZ !This basis is too small for real application<br>
<br>
r1=2.7 !Reference geometry<br>
theta=92<br>
<br>
r=[ 3.6,3.55,3.5,3.45,3.4,3.35,3.3,3.25,3.2,3.15,3.1,3.05,3.0,2.95,2.9,2.85,2.8,2.75,2.7,2.65,2.6,2.55,2.5,2.45,2.4,2.35,2.3,2.25,2.2,2.15,2.1,2.05,2.0]<br>
<br>
reforb=2140.2 !Orbital dumprecord at reference geometry<br>
refci=6000.2 !MRCI record at reference geometry<br>
savci=6100.2 !MRCI record at displaced geometries<br>
<br>
text,compute wavefunction at reference geometry (C2v)<br>
r2=3.6<br>
<br>
{hf;occ,9,2;wf,18,2,4;<br>
orbital,2100.2}<br>
<br>
{multi;occ,9,2;closed,4,1;<br>
wf,18,2;state,2; !1B1 and 1A2 states<br>
natorb,reforb !Save reference orbitals on reforb<br>
noextra} !Dont use extra symmetries<br>
<br>
{ci;occ,9,2;closed,4,1; !MRCI at reference geometry<br>
wf,18,2,0;state,2; !1B1 and 1A2 states<br>
orbital,reforb !Use orbitals from previous CASSCF<br>
save,refci} !Save MRCI wavefunction<br>
<br>
Text,Displaced geometries<br>
<br>
do i=1,#r !Loop over different r values<br>
data,truncate,savci+1 !truncate dumpfile after reference<br>
r2=r(i) !Set current r2<br>
<br>
{multi;occ,9,2;closed,4,1;<br>
wf,18,2,0;state,2; !Wavefunction definition<br>
start,reforb !Starting orbitals<br>
orbital,3140.2; !Dump record for orbitals<br>
diab,reforb !Generate diabatic orbitals relative to reference geometry<br>
noextra} !Dont use extra symmetries<br>
<br>
{ci;occ,9,2;closed,4,1;<br>
wf,18,2,0;state,2; !1B1 and 1A2 states<br>
orbital,diabatic !Use diabatic orbitals<br>
save,savci} !Save MRCI for displaced geometries<br>
<br>
e1(i)=energy(1) !Save adiabatic energies<br>
e2(i)=energy(2)<br>
<br>
{ci;trans,savci,savci !Compute transition densities at R2<br>
dm,7000.2} !Save transition densities on this record<br>
{ci;trans,savci,refci; !Compute transition densities between R2 and R1<br>
dm,7100.2} !Save transition densities on this record<br>
<br>
{ddr<br>
density,7000.2,7100.2 !Densities for <R2||R2> and <R2||R1><br>
orbital,3140.2,2140.2 !Orbitals for <R2||R2> and <R2||R1><br>
energy,e1(i),e2(i) !Adiabatic energies<br>
mixing,1.2,2.2} !Compute mixing angle and diabatic energies<br>
<br>
mixci(i)=mixangci(1) !Mixing angle obtained from ci vectors only<br>
h11ci(i)=hdiaci(1) !Diabatic energies obtained from ci vectors only<br>
h21ci(i)=hdiaci(2)<br>
h22ci(i)=hdiaci(3)<br>
<br>
mixtot(i)=mixang(1) !Mixing angle from total overlap (including first-order correction)<br>
h11(i)=hdia(1) !Diabatic energies obtained from total overlap<br>
h21(i)=hdia(2)<br>
h22(i)=hdia(3)<br>
<br>
<br>
{table,r,e1,e2,h11ci,h22ci,h21ci,mixci<br>
title,Diabatic energies for H2S, obtained from CI-vectors<br>
format,'(f10.2,5f14.8,f12.2)'<br>
sort,1}<br>
<br>
{table,r,e1,e2,h11,h22,h21,mixtot<br>
title,Diabatic energies for H2S, obtained from CI-vectors and orbital correction<br>
format,'(f10.2,5f14.8,f12.2)'<br>
sort,1}<br>
<br>
enddo !end loop over i
<p><br>
</p>
</body>
</html>