E_tot_final=`grep "SCF Done" final.out | awk -F = '{print $2}'| awk '{print $1}'`
E_T_final=`grep "ET=" final.out | awk -F = '{print $2}' | awk '{print $1}'`
E_x_final=`grep "ENTVJ=" final.out | awk -F = '{print $3}' | awk '{print $1}'`
E_c_final=`grep "ENTVJ=" final.out | awk -F = '{print $4}' | awk '{print $1}'`
...
E_disp_final=`grep "Dispersion energy" final.out | awk -F = '{print $2}'| awk '{print $1}'`
...
E_tot_frz=`grep " E=" final.out | head -1 | awk '{print $2}'`Before running these and subsequent lines in the sobEDA.sh, you only need to replace the values after "SCF Done", "ET=", "ENTVJ=" in final.out with those in final2.out.
]]>I see at the bottom of the sobEDA.sh script that the only energy values that depend on the final energy are the total (dE_tot), dispersion (dE_disp), Pauli repulsion energy (dE_rep), and orbital energy (dE_orb). For reference, these are the lines of code I am looking at in the sobEDA.sh script.
##### Show summary of interaction #####
dE_tot=$E_tot_final
dE_els=$E_els_promol
dE_x=$E_x_promol
dE_c=$E_c_promol
dE_disp=$E_disp_final
for ((i=1;i<=$nfrag;i=i+1))
do
dE_tot=`echo "$dE_tot-(${E_tot[$i]})" | bc -l`
dE_els=`echo "$dE_els-(${E_els[$i]})" | bc -l`
dE_x=`echo "$dE_x-(${E_x[$i]})" | bc -l`
dE_c=`echo "$dE_c-(${E_c[$i]})" | bc -l`
dE_disp=`echo "$dE_disp-(${E_disp[$i]})" | bc -l`
done
dE_orb=`echo "$E_tot_final-($E_tot_frz)" | bc -l`
dE_rep=`echo "$E_tot_frz-($E_tot_promol)" | bc -l`
The final total energy (E_tot_final) looks for the energy in the line that says "SCF Done" in final.out, dispersion energy (E_disp_final) looks for the energy in the line that says "Dispersion energy" in final.out, frozen energy (E_tot_frz) looks for the energy shown in the first cycle in final.out (which final.gjf reads from the promol.chk file), and E_tot_promol looks for the energy in the line that says "SCF Done" in promol.out.
After running a single point calculation (let's call this new single point calculation "final2") using the same route section as final.gjf, but reading the guess from the stable wavefunction, I receive an output that finishes after only 1 cycle (I guess because it is reading from an already optimized, stable wavefunction). With the knowledge of the way the parameters of EDA are calculated by the program, I compare between final and final2 the dispersion energy (it does not change, and therefore dE_disp would not change) and the final energy printed after "SCF Done" (it becomes more negative because of the stable wavefunction).
If I try to perform the EDA on final2.out, obviously I will get erroneous results because the frozen energy is the same as the final energy (it finishes in one cycle). So, my next thought is to modify the line in the original final.out that says "SCF Done" and replace it with the energy from that same line in final2.out. (i.e. change "SCF Done: E(UPBE1PBE) = XXX.XXX" to " SCF Done: E(UPBE1PBE) = YYY.YYY") That way the final energy represents the energy from the stable wavefunction. Although, this is assuming the frozen energy would remain the same.
Please let me know if my logic is correct and if it is reasonable to replace the energy value shown after "SCF Done" in final.out with the energy from final2.out.
]]>