diff --git a/bin/lowdin b/bin/lowdin index 3a709f9c..da82a9d3 100755 --- a/bin/lowdin +++ b/bin/lowdin @@ -204,6 +204,9 @@ if [ $extFile="lowdin" ]; then else if(toupper($i)==toupper("m")){ printf("\tInputParticle_mass = %s\n",toupper($(i+2)) ) } + else if(toupper($i)==toupper("eta")){ + printf("\tInputParticle_eta = %s\n",toupper($(i+2)) ) + } else if(toupper($i)==toupper("omega")){ printf("\tInputParticle_omega = %s\n",toupper($(i+2)) ) } @@ -381,10 +384,10 @@ if [ $extFile="lowdin" ]; then ' | $SED "s/,/./g" >> $nameFile.aux ########################################### - # Check custom basis sets in the input + # Check custom basis sets/potentials in the input ########################################### - BASIS_NAMES=(`gawk '($1~/BASIS/){print toupper($2)}' $nameFile`) + BASIS_NAMES=(`gawk '($1~/^BASIS$/){print toupper($2)}' $nameFile`) if [ ${#BASIS_NAMES[@]} -gt "0" ] then for BASIS_NAME in ${BASIS_NAMES[@]} @@ -401,6 +404,23 @@ if [ $extFile="lowdin" ]; then done fi + POTENTIALS_NAMES=(`gawk '($1~/^POTENTIAL$/){print toupper($2)}' $nameFile`) + if [ ${#POTENTIALS_NAMES[@]} -gt "0" ] + then + for POTENTIALS_NAME in ${POTENTIALS_NAMES[@]} + do + if [ -e $LOWDIN_DATA/basis/$POTENTIALS_NAME ] + then + echo "## ERROR: ## The custom potential file already exists in " $LOWDIN_DATA/potentials/$POTENTIALS_NAME + echo "Modify the POTENTIALS block in your input and select a different name" + exit 1 + fi + gawk '($1~/^POTENTIAL$/ && toupper($2)~/^'$POTENTIALS_NAME'$/){flag=1; next} + ($0~/END/){flag=0}; + (flag==1){print toupper($0)}' $nameFile > $POTENTIALS_NAME.$PID + done + fi + ########################################### # Exec lowdin.x @@ -421,6 +441,7 @@ if [ $extFile="lowdin" ]; then cp $nameFile*.vec $LOWDIN_SCRATCH/$nameFile &> /dev/null cp $nameFile*.plainvec $LOWDIN_SCRATCH/$nameFile &> /dev/null + cp $nameFile*.fchk $LOWDIN_SCRATCH/$nameFile &> /dev/null cp $nameFile*.val $LOWDIN_SCRATCH/$nameFile &> /dev/null cp $nameFile*.dens $LOWDIN_SCRATCH/$nameFile &> /dev/null cp $nameFile*.sup $LOWDIN_SCRATCH/$nameFile &> /dev/null @@ -432,7 +453,9 @@ if [ $extFile="lowdin" ]; then mv $nameFile*.over $LOWDIN_SCRATCH/$nameFile &> /dev/null mv $nameFile*.kin $LOWDIN_SCRATCH/$nameFile &> /dev/null mv $nameFile*.coeff $LOWDIN_SCRATCH/$nameFile &> /dev/null - #PID to avoid basis duplicates in simultaneous calculations + cp $nameFile*.gms.bs $LOWDIN_SCRATCH/$nameFile &> /dev/null + + #PID to avoid basis/potentials duplicates in simultaneous calculations if [ ${#BASIS_NAMES[@]} -gt "0" ] then for BASIS_NAME in ${BASIS_NAMES[@]} @@ -440,12 +463,14 @@ if [ $extFile="lowdin" ]; then mv $BASIS_NAME.$PID $LOWDIN_SCRATCH/$nameFile/$BASIS_NAME &> /dev/null done fi - - if [ -e $nameFile.gms.bs ] + if [ ${#POTENTIALS_NAMES[@]} -gt "0" ] then - cp $nameFile.gms.bs $LOWDIN_SCRATCH/$nameFile + for POTENTIALS_NAME in ${POTENTIALS_NAMES[@]} + do + mv $POTENTIALS_NAME.$PID $LOWDIN_SCRATCH/$nameFile/$POTENTIALS_NAME &> /dev/null + done fi - + # setting default number of cores for OpenMP if [ -z "$OMP_NUM_THREADS" ]; then @@ -515,6 +540,7 @@ if [ $extFile="lowdin" ]; then mv $LOWDIN_SCRATCH/$nameFile/$nameFile*.47 $currentPath &> 2 mv $LOWDIN_SCRATCH/$nameFile/*.vec $currentPath &> 2 mv $LOWDIN_SCRATCH/$nameFile/*.plainvec $currentPath &> 2 + mv $LOWDIN_SCRATCH/$nameFile/*.fchk $currentPath &> 2 # mv $LOWDIN_SCRATCH/$nameFile/*.val $currentPath &> 2 mv $LOWDIN_SCRATCH/$nameFile/*.NOCI.coords $currentPath &> 2 mv $LOWDIN_SCRATCH/$nameFile/*.NOCI.s* $currentPath &> 2 diff --git a/lib/basis/H2O-1S1P1D b/lib/basis/H2O-1S1P1D deleted file mode 100644 index cc7027ba..00000000 --- a/lib/basis/H2O-1S1P1D +++ /dev/null @@ -1,44 +0,0 @@ -O-HYDROGEN H_1 (1S) BASIS TYPE: 2 -3 -1 0 1 -14.509888498676842 1.0 -2 1 1 -6.885507269761004 1.0 -3 2 1 -9.023681376783887 1.0 - -O-HYDROGEN H-A_1 (1S) BASIS TYPE: 2 -3 -1 0 1 -14.509888498676842 1.0 -2 1 1 -6.885507269761004 1.0 -3 2 1 -9.023681376783887 1.0 - -O-HYDROGEN H-B_1 (1S) BASIS TYPE: 2 -3 -1 0 1 -14.509888498676842 1.0 -2 1 1 -6.885507269761004 1.0 -3 2 1 -9.023681376783887 1.0 - -O-H-TIP X0.5+ (1S) BASIS TYPE: 2 -3 -1 0 1 -14.509888498676842 1.0 -2 1 1 -6.885507269761004 1.0 -3 2 1 -9.023681376783887 1.0 - -O-H-TIP Y0.5+ (1S) BASIS TYPE: 2 -3 -1 0 1 -14.509888498676842 1.0 -2 1 1 -6.885507269761004 1.0 -3 2 1 -9.023681376783887 1.0 diff --git a/lib/basis/PSX-DZ b/lib/basis/PSX-DZ index 35b403de..0d5600b6 100644 --- a/lib/basis/PSX-DZ +++ b/lib/basis/PSX-DZ @@ -1,5 +1,53 @@ -O-POSITRON E+ (5S) BASIS TYPE: 3 -# (5S)-[5S] +O-POSITRON E+ (5S3P2D) BASIS TYPE: 3 +# (5S3P2D)-[5S3P2D] +10 +1 0 1 +.0189693659 1.0 +1 0 1 +.05186863351164733038 1.0 +1 0 1 +.14182630861506996771 1.0 +1 0 1 +.38780088183468771642 1.0 +1 0 1 +1.06037818667291718709 1.0 +1 1 1 +.0590955656 1.0 +1 1 1 +.16127967470846848928 1.0 +1 1 1 +.44015372744092004227 1.0 +1 2 1 +.11654481880000000000 1.0 +1 2 1 +.31287113568838127412 1.0 + +O-POSITRON E+A (5S3P2D) BASIS TYPE: 3 +# (5S3P2D)-[5S3P2D] +10 +1 0 1 +.0189693659 1.0 +1 0 1 +.05186863351164733038 1.0 +1 0 1 +.14182630861506996771 1.0 +1 0 1 +.38780088183468771642 1.0 +1 0 1 +1.06037818667291718709 1.0 +1 1 1 +.0590955656 1.0 +1 1 1 +.16127967470846848928 1.0 +1 1 1 +.44015372744092004227 1.0 +1 2 1 +.11654481880000000000 1.0 +1 2 1 +.31287113568838127412 1.0 + +O-POSITRON E+B (5S3P2D) BASIS TYPE: 3 +# (5S3P2D)-[5S3P2D] 10 1 0 1 .0189693659 1.0 diff --git a/lib/dataBases/constantsOfCoupling.lib b/lib/dataBases/constantsOfCoupling.lib index bda5b291..1f940075 100644 --- a/lib/dataBases/constantsOfCoupling.lib +++ b/lib/dataBases/constantsOfCoupling.lib @@ -4525,19 +4525,3 @@ LAMBDA = 2.0 PARTICLESFRACTION = 0.5 / -&SPECIE - NAME = "HA-TIP" - SYMBOL = "X0.5+" - KAPPA = -1.0 - ETA = 1.0 - LAMBDA = 1.0 - PARTICLESFRACTION = 1 -/ -&SPECIE - NAME = "HB-TIP" - SYMBOL = "Y0.5+" - KAPPA = -1.0 - ETA = 1.0 - LAMBDA = 1.0 - PARTICLESFRACTION = 1 -/ diff --git a/lib/dataBases/elementalParticles.lib b/lib/dataBases/elementalParticles.lib index d6630a9b..54376343 100644 --- a/lib/dataBases/elementalParticles.lib +++ b/lib/dataBases/elementalParticles.lib @@ -171,7 +171,7 @@ SYMBOL = "HEA3" CATEGORY = "LEPTON" CHARGE = 1 - MASS = 5494.892576965 + MASS = 5495.8851 SPIN = 0.5 / &PARTICLE @@ -179,7 +179,7 @@ SYMBOL = "HEB3" CATEGORY = "LEPTON" CHARGE = 1 - MASS = 5494.892576965 + MASS = 5495.8851 SPIN = 0.5 / &PARTICLE @@ -187,7 +187,7 @@ SYMBOL = "HES3" CATEGORY = "LEPTON" CHARGE = 1 - MASS = 5494.892576965 + MASS = 5495.8851 SPIN = 0.5 / &PARTICLE @@ -195,7 +195,7 @@ SYMBOL = "HEA4" CATEGORY = "LEPTON" CHARGE = 1 - MASS = 7292.327967297 + MASS = 7294.2994 SPIN = 0.5 / &PARTICLE @@ -203,7 +203,7 @@ SYMBOL = "HEB4" CATEGORY = "LEPTON" CHARGE = 1 - MASS = 7292.327967297 + MASS = 7294.2994 SPIN = 0.5 / &PARTICLE @@ -211,31 +211,7 @@ SYMBOL = "HES4" CATEGORY = "LEPTON" CHARGE = 1 - MASS = 7292.327967297 - SPIN = 0.5 -/ -&PARTICLE - NAME = "HA-TIP" - SYMBOL = "X0.5+" - CATEGORY = "LEPTON" - CHARGE = 0.5564 - MASS = 1836.15267247 - SPIN = 0.5 -/ -&PARTICLE - NAME = "HB-TIP" - SYMBOL = "Y0.5+" - CATEGORY = "LEPTON" - CHARGE = 0.5564 - MASS = 1836.15267247 - SPIN = -0.5 -/ -&PARTICLE - NAME = "M-TIP" - SYMBOL = "X1.1-" - CATEGORY = "LEPTON" - CHARGE = -1.1128 - MASS = 1836.15267247 + MASS = 7294.2994 SPIN = 0.5 / &PARTICLE diff --git a/lib/potentials/VHH-CCSDT b/lib/potentials/VHH-CCSDT deleted file mode 100644 index 138f0597..00000000 --- a/lib/potentials/VHH-CCSDT +++ /dev/null @@ -1,567 +0,0 @@ -#Fitted from CCSD(T)/def2-TZVPPD with constant ROH -O-H_1H_1 -26 -1 0 -10.000000000 1.313644267 -0.0 0.0 0.0 -2 0 -7.498942093 3.312208762 -0.0 0.0 0.0 -3 0 -5.623413252 -4.290599715 -0.0 0.0 0.0 -4 0 -4.216965034 -2.400490962 -0.0 0.0 0.0 -5 0 -3.162277660 6.160643793 -0.0 0.0 0.0 -6 0 -2.371373706 2.406172066 -0.0 0.0 0.0 -7 0 -1.778279410 -6.852617537 -0.0 0.0 0.0 -8 0 -1.333521432 -1.043562213 -0.0 0.0 0.0 -9 0 -1.000000000 7.556790649 -0.0 0.0 0.0 -10 0 -0.749894209 0.777037995 -0.0 0.0 0.0 -11 0 -0.562341325 -8.976316536 -0.0 0.0 0.0 -12 0 -0.421696503 0.575296040 -0.0 0.0 0.0 -13 0 -0.316227766 8.672468936 -0.0 0.0 0.0 -14 0 -0.237137371 -0.855388714 -0.0 0.0 0.0 -15 0 -0.177827941 -6.271012179 -0.0 0.0 0.0 -16 0 -0.133352143 1.152328020 -0.0 0.0 0.0 -17 0 -0.100000000 3.403807143 -0.0 0.0 0.0 -18 0 -0.074989421 -1.745300340 -0.0 0.0 0.0 -19 0 -0.056234133 -1.349326021 -0.0 0.0 0.0 -20 0 -0.042169650 1.605561812 -0.0 0.0 0.0 -21 0 -0.031622777 -0.267553445 -0.0 0.0 0.0 -22 0 -0.023713737 0.034566428 -0.0 0.0 0.0 -23 0 -0.017782794 -0.185356609 -0.0 0.0 0.0 -24 0 -0.013335214 -0.103128899 -0.0 0.0 0.0 -25 0 -0.010000000 0.128087663 -0.0 0.0 0.0 -26 0 -0.000000000 0.085213897 -0.0 0.0 0.0 - -O-H-A_1H-A_1 -26 -1 0 -10.000000000 1.313644267 -0.0 0.0 0.0 -2 0 -7.498942093 3.312208762 -0.0 0.0 0.0 -3 0 -5.623413252 -4.290599715 -0.0 0.0 0.0 -4 0 -4.216965034 -2.400490962 -0.0 0.0 0.0 -5 0 -3.162277660 6.160643793 -0.0 0.0 0.0 -6 0 -2.371373706 2.406172066 -0.0 0.0 0.0 -7 0 -1.778279410 -6.852617537 -0.0 0.0 0.0 -8 0 -1.333521432 -1.043562213 -0.0 0.0 0.0 -9 0 -1.000000000 7.556790649 -0.0 0.0 0.0 -10 0 -0.749894209 0.777037995 -0.0 0.0 0.0 -11 0 -0.562341325 -8.976316536 -0.0 0.0 0.0 -12 0 -0.421696503 0.575296040 -0.0 0.0 0.0 -13 0 -0.316227766 8.672468936 -0.0 0.0 0.0 -14 0 -0.237137371 -0.855388714 -0.0 0.0 0.0 -15 0 -0.177827941 -6.271012179 -0.0 0.0 0.0 -16 0 -0.133352143 1.152328020 -0.0 0.0 0.0 -17 0 -0.100000000 3.403807143 -0.0 0.0 0.0 -18 0 -0.074989421 -1.745300340 -0.0 0.0 0.0 -19 0 -0.056234133 -1.349326021 -0.0 0.0 0.0 -20 0 -0.042169650 1.605561812 -0.0 0.0 0.0 -21 0 -0.031622777 -0.267553445 -0.0 0.0 0.0 -22 0 -0.023713737 0.034566428 -0.0 0.0 0.0 -23 0 -0.017782794 -0.185356609 -0.0 0.0 0.0 -24 0 -0.013335214 -0.103128899 -0.0 0.0 0.0 -25 0 -0.010000000 0.128087663 -0.0 0.0 0.0 -26 0 -0.000000000 0.085213897 -0.0 0.0 0.0 - -O-H-A_1H-B_1 -26 -1 0 -10.000000000 1.313644267 -0.0 0.0 0.0 -2 0 -7.498942093 3.312208762 -0.0 0.0 0.0 -3 0 -5.623413252 -4.290599715 -0.0 0.0 0.0 -4 0 -4.216965034 -2.400490962 -0.0 0.0 0.0 -5 0 -3.162277660 6.160643793 -0.0 0.0 0.0 -6 0 -2.371373706 2.406172066 -0.0 0.0 0.0 -7 0 -1.778279410 -6.852617537 -0.0 0.0 0.0 -8 0 -1.333521432 -1.043562213 -0.0 0.0 0.0 -9 0 -1.000000000 7.556790649 -0.0 0.0 0.0 -10 0 -0.749894209 0.777037995 -0.0 0.0 0.0 -11 0 -0.562341325 -8.976316536 -0.0 0.0 0.0 -12 0 -0.421696503 0.575296040 -0.0 0.0 0.0 -13 0 -0.316227766 8.672468936 -0.0 0.0 0.0 -14 0 -0.237137371 -0.855388714 -0.0 0.0 0.0 -15 0 -0.177827941 -6.271012179 -0.0 0.0 0.0 -16 0 -0.133352143 1.152328020 -0.0 0.0 0.0 -17 0 -0.100000000 3.403807143 -0.0 0.0 0.0 -18 0 -0.074989421 -1.745300340 -0.0 0.0 0.0 -19 0 -0.056234133 -1.349326021 -0.0 0.0 0.0 -20 0 -0.042169650 1.605561812 -0.0 0.0 0.0 -21 0 -0.031622777 -0.267553445 -0.0 0.0 0.0 -22 0 -0.023713737 0.034566428 -0.0 0.0 0.0 -23 0 -0.017782794 -0.185356609 -0.0 0.0 0.0 -24 0 -0.013335214 -0.103128899 -0.0 0.0 0.0 -25 0 -0.010000000 0.128087663 -0.0 0.0 0.0 -26 0 -0.000000000 0.085213897 -0.0 0.0 0.0 - -O-H-B_1H-B_1 -26 -1 0 -10.000000000 1.313644267 -0.0 0.0 0.0 -2 0 -7.498942093 3.312208762 -0.0 0.0 0.0 -3 0 -5.623413252 -4.290599715 -0.0 0.0 0.0 -4 0 -4.216965034 -2.400490962 -0.0 0.0 0.0 -5 0 -3.162277660 6.160643793 -0.0 0.0 0.0 -6 0 -2.371373706 2.406172066 -0.0 0.0 0.0 -7 0 -1.778279410 -6.852617537 -0.0 0.0 0.0 -8 0 -1.333521432 -1.043562213 -0.0 0.0 0.0 -9 0 -1.000000000 7.556790649 -0.0 0.0 0.0 -10 0 -0.749894209 0.777037995 -0.0 0.0 0.0 -11 0 -0.562341325 -8.976316536 -0.0 0.0 0.0 -12 0 -0.421696503 0.575296040 -0.0 0.0 0.0 -13 0 -0.316227766 8.672468936 -0.0 0.0 0.0 -14 0 -0.237137371 -0.855388714 -0.0 0.0 0.0 -15 0 -0.177827941 -6.271012179 -0.0 0.0 0.0 -16 0 -0.133352143 1.152328020 -0.0 0.0 0.0 -17 0 -0.100000000 3.403807143 -0.0 0.0 0.0 -18 0 -0.074989421 -1.745300340 -0.0 0.0 0.0 -19 0 -0.056234133 -1.349326021 -0.0 0.0 0.0 -20 0 -0.042169650 1.605561812 -0.0 0.0 0.0 -21 0 -0.031622777 -0.267553445 -0.0 0.0 0.0 -22 0 -0.023713737 0.034566428 -0.0 0.0 0.0 -23 0 -0.017782794 -0.185356609 -0.0 0.0 0.0 -24 0 -0.013335214 -0.103128899 -0.0 0.0 0.0 -25 0 -0.010000000 0.128087663 -0.0 0.0 0.0 -26 0 -0.000000000 0.085213897 -0.0 0.0 0.0 - -O-X0.5+X0.5+ -26 -1 0 -10.000000000 1.313644267 -0.0 0.0 0.0 -2 0 -7.498942093 3.312208762 -0.0 0.0 0.0 -3 0 -5.623413252 -4.290599715 -0.0 0.0 0.0 -4 0 -4.216965034 -2.400490962 -0.0 0.0 0.0 -5 0 -3.162277660 6.160643793 -0.0 0.0 0.0 -6 0 -2.371373706 2.406172066 -0.0 0.0 0.0 -7 0 -1.778279410 -6.852617537 -0.0 0.0 0.0 -8 0 -1.333521432 -1.043562213 -0.0 0.0 0.0 -9 0 -1.000000000 7.556790649 -0.0 0.0 0.0 -10 0 -0.749894209 0.777037995 -0.0 0.0 0.0 -11 0 -0.562341325 -8.976316536 -0.0 0.0 0.0 -12 0 -0.421696503 0.575296040 -0.0 0.0 0.0 -13 0 -0.316227766 8.672468936 -0.0 0.0 0.0 -14 0 -0.237137371 -0.855388714 -0.0 0.0 0.0 -15 0 -0.177827941 -6.271012179 -0.0 0.0 0.0 -16 0 -0.133352143 1.152328020 -0.0 0.0 0.0 -17 0 -0.100000000 3.403807143 -0.0 0.0 0.0 -18 0 -0.074989421 -1.745300340 -0.0 0.0 0.0 -19 0 -0.056234133 -1.349326021 -0.0 0.0 0.0 -20 0 -0.042169650 1.605561812 -0.0 0.0 0.0 -21 0 -0.031622777 -0.267553445 -0.0 0.0 0.0 -22 0 -0.023713737 0.034566428 -0.0 0.0 0.0 -23 0 -0.017782794 -0.185356609 -0.0 0.0 0.0 -24 0 -0.013335214 -0.103128899 -0.0 0.0 0.0 -25 0 -0.010000000 0.128087663 -0.0 0.0 0.0 -26 0 -0.000000000 0.085213897 -0.0 0.0 0.0 - -O-X0.5+Y0.5+ -26 -1 0 -10.000000000 1.313644267 -0.0 0.0 0.0 -2 0 -7.498942093 3.312208762 -0.0 0.0 0.0 -3 0 -5.623413252 -4.290599715 -0.0 0.0 0.0 -4 0 -4.216965034 -2.400490962 -0.0 0.0 0.0 -5 0 -3.162277660 6.160643793 -0.0 0.0 0.0 -6 0 -2.371373706 2.406172066 -0.0 0.0 0.0 -7 0 -1.778279410 -6.852617537 -0.0 0.0 0.0 -8 0 -1.333521432 -1.043562213 -0.0 0.0 0.0 -9 0 -1.000000000 7.556790649 -0.0 0.0 0.0 -10 0 -0.749894209 0.777037995 -0.0 0.0 0.0 -11 0 -0.562341325 -8.976316536 -0.0 0.0 0.0 -12 0 -0.421696503 0.575296040 -0.0 0.0 0.0 -13 0 -0.316227766 8.672468936 -0.0 0.0 0.0 -14 0 -0.237137371 -0.855388714 -0.0 0.0 0.0 -15 0 -0.177827941 -6.271012179 -0.0 0.0 0.0 -16 0 -0.133352143 1.152328020 -0.0 0.0 0.0 -17 0 -0.100000000 3.403807143 -0.0 0.0 0.0 -18 0 -0.074989421 -1.745300340 -0.0 0.0 0.0 -19 0 -0.056234133 -1.349326021 -0.0 0.0 0.0 -20 0 -0.042169650 1.605561812 -0.0 0.0 0.0 -21 0 -0.031622777 -0.267553445 -0.0 0.0 0.0 -22 0 -0.023713737 0.034566428 -0.0 0.0 0.0 -23 0 -0.017782794 -0.185356609 -0.0 0.0 0.0 -24 0 -0.013335214 -0.103128899 -0.0 0.0 0.0 -25 0 -0.010000000 0.128087663 -0.0 0.0 0.0 -26 0 -0.000000000 0.085213897 -0.0 0.0 0.0 - -O-Y0.5+Y0.5+ -26 -1 0 -10.000000000 1.313644267 -0.0 0.0 0.0 -2 0 -7.498942093 3.312208762 -0.0 0.0 0.0 -3 0 -5.623413252 -4.290599715 -0.0 0.0 0.0 -4 0 -4.216965034 -2.400490962 -0.0 0.0 0.0 -5 0 -3.162277660 6.160643793 -0.0 0.0 0.0 -6 0 -2.371373706 2.406172066 -0.0 0.0 0.0 -7 0 -1.778279410 -6.852617537 -0.0 0.0 0.0 -8 0 -1.333521432 -1.043562213 -0.0 0.0 0.0 -9 0 -1.000000000 7.556790649 -0.0 0.0 0.0 -10 0 -0.749894209 0.777037995 -0.0 0.0 0.0 -11 0 -0.562341325 -8.976316536 -0.0 0.0 0.0 -12 0 -0.421696503 0.575296040 -0.0 0.0 0.0 -13 0 -0.316227766 8.672468936 -0.0 0.0 0.0 -14 0 -0.237137371 -0.855388714 -0.0 0.0 0.0 -15 0 -0.177827941 -6.271012179 -0.0 0.0 0.0 -16 0 -0.133352143 1.152328020 -0.0 0.0 0.0 -17 0 -0.100000000 3.403807143 -0.0 0.0 0.0 -18 0 -0.074989421 -1.745300340 -0.0 0.0 0.0 -19 0 -0.056234133 -1.349326021 -0.0 0.0 0.0 -20 0 -0.042169650 1.605561812 -0.0 0.0 0.0 -21 0 -0.031622777 -0.267553445 -0.0 0.0 0.0 -22 0 -0.023713737 0.034566428 -0.0 0.0 0.0 -23 0 -0.017782794 -0.185356609 -0.0 0.0 0.0 -24 0 -0.013335214 -0.103128899 -0.0 0.0 0.0 -25 0 -0.010000000 0.128087663 -0.0 0.0 0.0 -26 0 -0.000000000 0.085213897 -0.0 0.0 0.0 diff --git a/lib/potentials/VOH-CCSDT b/lib/potentials/VOH-CCSDT deleted file mode 100644 index d3b76b4d..00000000 --- a/lib/potentials/VOH-CCSDT +++ /dev/null @@ -1,420 +0,0 @@ -#Fitted from CCSD(T)/def2-TZVPPD with constant RHH -O-H-A_1 -27 -1 0 -13.892336977 9.839539381 -0.0 0.0 0.0 -2 0 -10.000000000 15.916259785 -0.0 0.0 0.0 -3 0 -7.498942093 -18.199008237 -0.0 0.0 0.0 -4 0 -5.623413252 6.723316202 -0.0 0.0 0.0 -5 0 -4.216965034 4.671467422 -0.0 0.0 0.0 -6 0 -3.162277660 1.703771884 -0.0 0.0 0.0 -7 0 -2.371373706 0.297676262 -0.0 0.0 0.0 -8 0 -1.778279410 0.967600283 -0.0 0.0 0.0 -9 0 -1.333521432 1.034932150 -0.0 0.0 0.0 -10 0 -1.000000000 0.657621305 -0.0 0.0 0.0 -11 0 -0.749894209 -0.327381267 -0.0 0.0 0.0 -12 0 -0.562341325 0.221528399 -0.0 0.0 0.0 -13 0 -0.421696503 -0.255694044 -0.0 0.0 0.0 -14 0 -0.316227766 0.097853405 -0.0 0.0 0.0 -15 0 -0.237137371 0.744489008 -0.0 0.0 0.0 -16 0 -0.177827941 -0.750090202 -0.0 0.0 0.0 -17 0 -0.133352143 -0.309090719 -0.0 0.0 0.0 -18 0 -0.100000000 -0.998800172 -0.0 0.0 0.0 -19 0 -0.074989421 0.718795407 -0.0 0.0 0.0 -20 0 -0.056234133 0.228002288 -0.0 0.0 0.0 -21 0 -0.042169650 0.573744431 -0.0 0.0 0.0 -22 0 -0.031622777 -0.056969092 -0.0 0.0 0.0 -23 0 -0.023713737 -0.392241803 -0.0 0.0 0.0 -24 0 -0.017782794 -0.569609681 -0.0 0.0 0.0 -25 0 -0.013335214 0.759323484 -0.0 0.0 0.0 -26 0 -0.010000000 -0.194525832 -0.0 0.0 0.0 -27 0 -0.000000000 0.138912412 -0.0 0.0 0.0 - -O-H-B_1 -27 -1 0 -13.892336977 9.839539381 -0.0 0.0 0.0 -2 0 -10.000000000 15.916259785 -0.0 0.0 0.0 -3 0 -7.498942093 -18.199008237 -0.0 0.0 0.0 -4 0 -5.623413252 6.723316202 -0.0 0.0 0.0 -5 0 -4.216965034 4.671467422 -0.0 0.0 0.0 -6 0 -3.162277660 1.703771884 -0.0 0.0 0.0 -7 0 -2.371373706 0.297676262 -0.0 0.0 0.0 -8 0 -1.778279410 0.967600283 -0.0 0.0 0.0 -9 0 -1.333521432 1.034932150 -0.0 0.0 0.0 -10 0 -1.000000000 0.657621305 -0.0 0.0 0.0 -11 0 -0.749894209 -0.327381267 -0.0 0.0 0.0 -12 0 -0.562341325 0.221528399 -0.0 0.0 0.0 -13 0 -0.421696503 -0.255694044 -0.0 0.0 0.0 -14 0 -0.316227766 0.097853405 -0.0 0.0 0.0 -15 0 -0.237137371 0.744489008 -0.0 0.0 0.0 -16 0 -0.177827941 -0.750090202 -0.0 0.0 0.0 -17 0 -0.133352143 -0.309090719 -0.0 0.0 0.0 -18 0 -0.100000000 -0.998800172 -0.0 0.0 0.0 -19 0 -0.074989421 0.718795407 -0.0 0.0 0.0 -20 0 -0.056234133 0.228002288 -0.0 0.0 0.0 -21 0 -0.042169650 0.573744431 -0.0 0.0 0.0 -22 0 -0.031622777 -0.056969092 -0.0 0.0 0.0 -23 0 -0.023713737 -0.392241803 -0.0 0.0 0.0 -24 0 -0.017782794 -0.569609681 -0.0 0.0 0.0 -25 0 -0.013335214 0.759323484 -0.0 0.0 0.0 -26 0 -0.010000000 -0.194525832 -0.0 0.0 0.0 -27 0 -0.000000000 0.138912412 -0.0 0.0 0.0 - -O-H_1 -27 -1 0 -13.892336977 9.839539381 -0.0 0.0 0.0 -2 0 -10.000000000 15.916259785 -0.0 0.0 0.0 -3 0 -7.498942093 -18.199008237 -0.0 0.0 0.0 -4 0 -5.623413252 6.723316202 -0.0 0.0 0.0 -5 0 -4.216965034 4.671467422 -0.0 0.0 0.0 -6 0 -3.162277660 1.703771884 -0.0 0.0 0.0 -7 0 -2.371373706 0.297676262 -0.0 0.0 0.0 -8 0 -1.778279410 0.967600283 -0.0 0.0 0.0 -9 0 -1.333521432 1.034932150 -0.0 0.0 0.0 -10 0 -1.000000000 0.657621305 -0.0 0.0 0.0 -11 0 -0.749894209 -0.327381267 -0.0 0.0 0.0 -12 0 -0.562341325 0.221528399 -0.0 0.0 0.0 -13 0 -0.421696503 -0.255694044 -0.0 0.0 0.0 -14 0 -0.316227766 0.097853405 -0.0 0.0 0.0 -15 0 -0.237137371 0.744489008 -0.0 0.0 0.0 -16 0 -0.177827941 -0.750090202 -0.0 0.0 0.0 -17 0 -0.133352143 -0.309090719 -0.0 0.0 0.0 -18 0 -0.100000000 -0.998800172 -0.0 0.0 0.0 -19 0 -0.074989421 0.718795407 -0.0 0.0 0.0 -20 0 -0.056234133 0.228002288 -0.0 0.0 0.0 -21 0 -0.042169650 0.573744431 -0.0 0.0 0.0 -22 0 -0.031622777 -0.056969092 -0.0 0.0 0.0 -23 0 -0.023713737 -0.392241803 -0.0 0.0 0.0 -24 0 -0.017782794 -0.569609681 -0.0 0.0 0.0 -25 0 -0.013335214 0.759323484 -0.0 0.0 0.0 -26 0 -0.010000000 -0.194525832 -0.0 0.0 0.0 -27 0 -0.000000000 0.138912412 -0.0 0.0 0.0 - -O-X0.5+ -27 -1 0 -13.892336977 9.839539381 -0.0 0.0 0.0 -2 0 -10.000000000 15.916259785 -0.0 0.0 0.0 -3 0 -7.498942093 -18.199008237 -0.0 0.0 0.0 -4 0 -5.623413252 6.723316202 -0.0 0.0 0.0 -5 0 -4.216965034 4.671467422 -0.0 0.0 0.0 -6 0 -3.162277660 1.703771884 -0.0 0.0 0.0 -7 0 -2.371373706 0.297676262 -0.0 0.0 0.0 -8 0 -1.778279410 0.967600283 -0.0 0.0 0.0 -9 0 -1.333521432 1.034932150 -0.0 0.0 0.0 -10 0 -1.000000000 0.657621305 -0.0 0.0 0.0 -11 0 -0.749894209 -0.327381267 -0.0 0.0 0.0 -12 0 -0.562341325 0.221528399 -0.0 0.0 0.0 -13 0 -0.421696503 -0.255694044 -0.0 0.0 0.0 -14 0 -0.316227766 0.097853405 -0.0 0.0 0.0 -15 0 -0.237137371 0.744489008 -0.0 0.0 0.0 -16 0 -0.177827941 -0.750090202 -0.0 0.0 0.0 -17 0 -0.133352143 -0.309090719 -0.0 0.0 0.0 -18 0 -0.100000000 -0.998800172 -0.0 0.0 0.0 -19 0 -0.074989421 0.718795407 -0.0 0.0 0.0 -20 0 -0.056234133 0.228002288 -0.0 0.0 0.0 -21 0 -0.042169650 0.573744431 -0.0 0.0 0.0 -22 0 -0.031622777 -0.056969092 -0.0 0.0 0.0 -23 0 -0.023713737 -0.392241803 -0.0 0.0 0.0 -24 0 -0.017782794 -0.569609681 -0.0 0.0 0.0 -25 0 -0.013335214 0.759323484 -0.0 0.0 0.0 -26 0 -0.010000000 -0.194525832 -0.0 0.0 0.0 -27 0 -0.000000000 0.138912412 -0.0 0.0 0.0 - -O-Y0.5+ -27 -1 0 -13.892336977 9.839539381 -0.0 0.0 0.0 -2 0 -10.000000000 15.916259785 -0.0 0.0 0.0 -3 0 -7.498942093 -18.199008237 -0.0 0.0 0.0 -4 0 -5.623413252 6.723316202 -0.0 0.0 0.0 -5 0 -4.216965034 4.671467422 -0.0 0.0 0.0 -6 0 -3.162277660 1.703771884 -0.0 0.0 0.0 -7 0 -2.371373706 0.297676262 -0.0 0.0 0.0 -8 0 -1.778279410 0.967600283 -0.0 0.0 0.0 -9 0 -1.333521432 1.034932150 -0.0 0.0 0.0 -10 0 -1.000000000 0.657621305 -0.0 0.0 0.0 -11 0 -0.749894209 -0.327381267 -0.0 0.0 0.0 -12 0 -0.562341325 0.221528399 -0.0 0.0 0.0 -13 0 -0.421696503 -0.255694044 -0.0 0.0 0.0 -14 0 -0.316227766 0.097853405 -0.0 0.0 0.0 -15 0 -0.237137371 0.744489008 -0.0 0.0 0.0 -16 0 -0.177827941 -0.750090202 -0.0 0.0 0.0 -17 0 -0.133352143 -0.309090719 -0.0 0.0 0.0 -18 0 -0.100000000 -0.998800172 -0.0 0.0 0.0 -19 0 -0.074989421 0.718795407 -0.0 0.0 0.0 -20 0 -0.056234133 0.228002288 -0.0 0.0 0.0 -21 0 -0.042169650 0.573744431 -0.0 0.0 0.0 -22 0 -0.031622777 -0.056969092 -0.0 0.0 0.0 -23 0 -0.023713737 -0.392241803 -0.0 0.0 0.0 -24 0 -0.017782794 -0.569609681 -0.0 0.0 0.0 -25 0 -0.013335214 0.759323484 -0.0 0.0 0.0 -26 0 -0.010000000 -0.194525832 -0.0 0.0 0.0 -27 0 -0.000000000 0.138912412 -0.0 0.0 0.0 diff --git a/src/core/ConstantsOfCoupling.f90 b/src/core/ConstantsOfCoupling.f90 index 3e3bba2c..0c190802 100644 --- a/src/core/ConstantsOfCoupling.f90 +++ b/src/core/ConstantsOfCoupling.f90 @@ -71,59 +71,53 @@ subroutine ConstantsOfCoupling_load( this, symbolSelected ) !! Looking for library inquire(file=trim(CONTROL_instance%DATA_DIRECTORY)//"/dataBases/constantsOfCoupling.lib", exist=existFile) - if ( existFile ) then - - !! Open library - open(unit=10, file=trim(CONTROL_instance%DATA_DIRECTORY)//"/dataBases/constantsOfCoupling.lib", status="old", form="formatted" ) - - !! Read information - symbol = "NONE" - stat = 0 - - do while(trim(symbol) /= trim(symbolSelected)) - - !! Setting defaults - name = "NONE" - kappa = 0 - eta = 0 - particlesFraction = 1 - - if (stat == -1 ) then - - call ConstantsOfCoupling_exception( ERROR, "Elemental particle: "//trim(symbolSelected)//" NOT found!!", "In ConstantsOfCoupling at load function.") - this%isInstanced = .false. + if ( .not. existFile ) call ConstantsOfCoupling_exception( ERROR, "LOWDIN library not found!! please export lowdinvars.sh file.", "In ConstantsOfCoupling at load function.") + + !! Open library + open(unit=10, file=trim(CONTROL_instance%DATA_DIRECTORY)//"/dataBases/constantsOfCoupling.lib", status="old", form="formatted" ) - end if - - read(10,NML=specie, iostat=stat) - - if (stat > 0 ) then - - call ConstantsOfCoupling_exception( ERROR, "Failed reading ConstantsOfCouplings.lib file!! please check this file.", "In ConstantsOfCoupling at load function.") - - end if + !! Read information + symbol = "NONE" + stat = 0 - end do + do while(trim(symbol) /= trim(symbolSelected)) - !! Set object variables - this%name = name - this%symbol = symbol - this%kappa = kappa - this%eta = eta - this%lambda = lambda - this%particlesFraction = particlesFraction - - !! Debug information. - !! call ConstantsOfCoupling_show(this) - - close(10) + !! Setting defaults + name = "NONE" + kappa = 0 + eta = 0 + particlesFraction = 1 - else + if (stat == -1 ) then - call ConstantsOfCoupling_exception( ERROR, "LOWDIN library not found!! please export lowdinvars.sh file.", "In ConstantsOfCoupling at load function.") + ! call ConstantsOfCoupling_exception( WARNING, "Elemental particle: "//trim(symbolSelected)//" NOT found!!", "Setting default values") + this%isInstanced=.false. + exit + end if - end if + read(10,NML=specie, iostat=stat) + if (stat > 0 ) then + + call ConstantsOfCoupling_exception( ERROR, "Failed reading ConstantsOfCouplings.lib file!! please check this file.", "In ConstantsOfCoupling at load function.") + + end if + + end do + + !! Set object variables + this%name = name + this%symbol = symbol + this%kappa = kappa + this%eta = eta + this%lambda = lambda + this%particlesFraction = particlesFraction + + !! Debug information. + ! call ConstantsOfCoupling_show(this) + + close(10) + !! Done end subroutine ConstantsOfCoupling_load diff --git a/src/core/ElementalParticle.f90 b/src/core/ElementalParticle.f90 index ce1b69b9..5b961321 100644 --- a/src/core/ElementalParticle.f90 +++ b/src/core/ElementalParticle.f90 @@ -40,6 +40,7 @@ module ElementalParticle_ real(8) :: mass real(8) :: charge real(8) :: spin + logical :: custom end type ElementalParticle public :: & @@ -59,6 +60,7 @@ subroutine ElementalParticle_load( this, symbolSelected ) character(*) :: symbolSelected logical :: existFile + logical :: custom integer :: stat integer :: i @@ -81,59 +83,64 @@ subroutine ElementalParticle_load( this, symbolSelected ) !! Looking for library inquire(file=trim(CONTROL_instance%DATA_DIRECTORY)//trim(CONTROL_instance%ELEMENTAL_PARTICLES_DATABASE), exist=existFile) - if ( existFile ) then + if ( .not. existFile ) call ElementalParticle_exception( ERROR, "LOWDIN library not found!! please export lowdinvars.sh file.", "In ElementalParticle at load function.") - !! Open library - open(unit=10, file=trim(CONTROL_instance%DATA_DIRECTORY)//trim(CONTROL_instance%ELEMENTAL_PARTICLES_DATABASE), status="old", form="formatted" ) - - !! Read information - symbol = "NONE" - stat = 0 - - do while(trim(symbol) /= trim(symbolSelected)) + !! Open library + open(unit=10, file=trim(CONTROL_instance%DATA_DIRECTORY)//trim(CONTROL_instance%ELEMENTAL_PARTICLES_DATABASE), status="old", form="formatted" ) + + !! Read information + symbol = "NONE" + stat = 0 + + do while(trim(symbol) /= trim(symbolSelected)) + + !! Setting defaults + name = "NONE" + category = "NONE" + mass = -1 + charge = 0 + spin = 0 + custom = .false. - !! Setting defaults - name = "NONE" - category = "NONE" - mass = -1 - charge = 0 - spin = 0 - - if (stat == -1 ) then - - call ElementalParticle_exception( ERROR, "Elemental particle: "//trim(symbolSelected)//" NOT found!!", "In ElementalParticle at load function.") + if (stat == -1 ) then - end if + ! call ElementalParticle_exception( WARNING, "Elemental particle: "//trim(symbolSelected)//" NOT found in ElementalParticles.lib", "Setting default values") + name = trim(symbolSelected) + symbol = trim(symbolSelected) + category = "FERMION" + mass = 1.0 + charge = 1.0 + spin = 0.5 + custom = .true. + + exit - read(10,NML=particle, iostat=stat) - - if (stat > 0 ) then - - call ElementalParticle_exception( ERROR, "Failed reading ElementalParticles.lib file!! please check this file.", "In ElementalParticle at load function.") - - end if + end if - end do - - !! Set object variables - this%name = name - this%symbol = symbol - this%category = category - this%mass = mass - this%charge = charge - this%spin = spin - - !! Debug information. - !! call ElementalParticle_show(this) - - close(10) - - else + read(10,NML=particle, iostat=stat) + + if (stat > 0 ) then + + call ElementalParticle_exception( ERROR, "Failed reading ElementalParticles.lib file!! please check this file.", "In ElementalParticle at load function.") - call ElementalParticle_exception( ERROR, "LOWDIN library not found!! please export lowdinvars.sh file.", "In ElementalParticle at load function.") + end if - end if + end do + !! Set object variables + this%name = name + this%symbol = symbol + this%category = category + this%mass = mass + this%charge = charge + this%spin = spin + this%custom = custom + + !! Debug information. + ! call ElementalParticle_show(this) + + close(10) + !! Done end subroutine ElementalParticle_load diff --git a/src/core/ExternalPotential.f90 b/src/core/ExternalPotential.f90 deleted file mode 100644 index f6459bfa..00000000 --- a/src/core/ExternalPotential.f90 +++ /dev/null @@ -1,378 +0,0 @@ -!!****************************************************************************** -!! This code is part of LOWDIN Quantum chemistry package -!! -!! this program has been developed under direction of: -!! -!! Prof. A REYES' Lab. Universidad Nacional de Colombia -!! http://www.qcc.unal.edu.co -!! Prof. R. FLORES' Lab. Universidad de Guadajara -!! http://www.cucei.udg.mx/~robertof -!! -!! Todos los derechos reservados, 2013 -!! -!!****************************************************************************** - -!> @brief This module contains all the routines to handle external potentials -!! @author E. F. Posada (efposadac@unal.edu.co) -!! @version 2.0 -!! Creation data : 06-08-10 -!! -!! History change: -!! -!! - 06-08-10 : Sergio A. Gonzalez ( sergmonic@gmail.com ) -!! -# Creacioon del modulo y metodos basicos -!! - 2011-02-14 : Fernando Posada ( efposadac@unal.edu.co ) -!! -# Reescribe y adapta el módulo para su inclusion en Lowdin -module ExternalPotential_ - use ContractedGaussian_ - use String_ - use Matrix_ - use Units_ - use Exception_ - implicit none - - type, public :: ExternalPot - character(20) :: name - character(50) :: specie - character(50) :: ttype - character(50) :: units - integer :: numOfComponents - integer :: iter - type(ContractedGaussian), allocatable :: gaussianComponents(:) - end type - - type, public :: ExternalPotential - integer :: ssize - type(ExternalPot), allocatable :: potentials(:) - logical :: isInstanced - end type - - type(ExternalPotential), public, target :: ExternalPotential_instance - -contains - - - !> - !! @brief Constructor by default - !! @param this - subroutine ExternalPotential_constructor(numberOfPotentials) - implicit none - - integer :: numberOfPotentials - - ExternalPotential_instance%ssize = numberOfPotentials - allocate(ExternalPotential_instance%potentials(numberOfPotentials)) - ExternalPotential_instance%isInstanced = .true. - - end subroutine ExternalPotential_constructor - - !> - !! @brief Destroys the object - !! @param this - subroutine ExternalPotential_destructor() - implicit none - - integer :: i - - do i = 1, ExternalPotential_instance%ssize - if (allocated(ExternalPotential_instance%potentials(i)%gaussianComponents)) deallocate(ExternalPotential_instance%potentials(i)%gaussianComponents) - end do - - if (allocated(ExternalPotential_instance%potentials) ) deallocate(ExternalPotential_instance%potentials) - ExternalPotential_instance%isInstanced=.false. - - end subroutine ExternalPotential_destructor - - !> - !! @brief Shows information of the object - !! @param this - subroutine ExternalPotential_show() - implicit none - type(ExternalPot), pointer :: this - integer :: potId, i - - do potId = 1, ExternalPotential_instance%ssize - this => ExternalPotential_instance%potentials(potId) - - print *,"" - print *,"=======" - print *, "External Potential for ", trim(this%specie), " : ", trim(this%name) - print *, "Type : ", trim(this%ttype) - write(6,"(T10,A20,A10,A10,A10,A10,A20)") "Exponent", "l", "R_x", "R_y", "R_z", "Factor" - - do i=1,this%numOfComponents - write(6,"(T10,F20.10,I10,F10.5,F10.5,F10.5,F20.10)") & - this%gaussianComponents(i)%orbitalExponents(1), & - this%gaussianComponents(i)%angularMoment, this%gaussianComponents(i)%origin(:), & - this%gaussianComponents(i)%contractionCoefficients(1) - end do - end do - - end subroutine ExternalPotential_show - - !> - !! @brief loads information from the input file - !! @param this - !! @author E. F. Posada, 2013 - subroutine ExternalPotential_load(potId, name, species) - implicit none - integer :: potId - character(*) :: name - character(*) :: species - - type(ExternalPot), pointer :: this - integer :: status, i, j - character(150) :: fileName - character(20) :: token - character(10) :: symbol - logical :: existFile, found - - this => ExternalPotential_instance%potentials(potId) - - this%name= trim(name) - this%specie= trim(species) - this%ttype="" - this%units="bohr" - this%numOfComponents=0 - this%iter=1 - - fileName = trim( trim( CONTROL_instance%DATA_DIRECTORY ) // & - trim(CONTROL_instance%POTENTIALS_DATABASE)// String_getUppercase(trim(name))) - - inquire(file=trim(fileName), exist = existFile) - if(existFile) then - - !! Open File - open(unit=30, file=trim(fileName), status="old",form="formatted") - rewind(30) - - found = .false. - - !! Open element and Find the proper potential - do while(found .eqv. .false.) - read(30,*, iostat=status) token - symbol = token(3:) - - !! Some debug information in case of error! - if (status > 0 ) then - - call ExternalPotential_exception(ERROR, & - "ERROR reading ExternalPotential file: "//trim(this%name)//& - " Please check that file!","ExternalPotential module at Load function.") - - end if - - if (status == -1 ) then - - call ExternalPotential_exception(ERROR, & - "The ExternalPotential: "//trim(this%name)//& - " for: "//trim(species)//& - " was not found!","ExternalPotential module at Load function.") - - end if - - if(trim(token(1:2)) == "O-") then - if(trim(symbol) == trim(species)) then - found = .true. - - end if - - end if - - end do - - !! Neglect any comment - token = "#" - do while(trim(token(1:1)) == "#") - - read(30,*) token - - end do - - !! Start reading Potential - backspace(30) - - read(30,*, iostat=status) this%numOfComponents - - !! Some debug information in case of error! - if (status > 0 ) then - - call ExternalPotential_exception(ERROR, & - "ERROR reading ExternalPotential file: "//trim(this%name)//& - " Please check that file!","ExternalPotential module at Load function.") - - end if - - allocate(this%gaussianComponents(this%numOfComponents)) - - do i = 1, this%numOfComponents - - read(30,*,iostat=status) this%gaussianComponents(i)%id, & - this%gaussianComponents(i)%angularMoment - this%gaussianComponents(i)%length = 1 - - !! Some debug information in case of error! - if (status > 0 ) then - - call ExternalPotential_exception(ERROR, & - "ERROR reading ExternalPotential file: "//trim(this%name)//& - " Please check that file!","ExternalPotential module at Load function.") - - end if - - allocate(this%gaussianComponents(i)%orbitalExponents(this%gaussianComponents(i)%length)) - allocate(this%gaussianComponents(i)%contractionCoefficients(this%gaussianComponents(i)%length)) - - do j = 1, this%gaussianComponents(i)%length - - read(30,*,iostat=status) this%gaussianComponents(i)%orbitalExponents(j), & - this%gaussianComponents(i)%contractionCoefficients(j) - read(30,*,iostat=status) this%gaussianComponents(i)%origin - - !! Some debug information in case of error! - if (status > 0 ) then - - call ExternalPotential_exception(ERROR, & - "ERROR reading ExternalPotential file: "//trim(this%name)//& - " Please check that file!","ExternalPotential module at Load function.") - - end if - - end do - - - !! Calculates the number of Cartesian orbitals, by dimensionality - select case(CONTROL_instance%DIMENSIONALITY) - case(3) - this%gaussianComponents(i)%numCartesianOrbital = ( ( this%gaussianComponents(i)%angularMoment + 1_8 )*( this%gaussianComponents(i)%angularMoment + 2_8 ) ) / 2_8 - case(2) - this%gaussianComponents(i)%numCartesianOrbital = ( ( this%gaussianComponents(i)%angularMoment + 1_8 ) ) - case(1) - this%gaussianComponents(i)%numCartesianOrbital = 1 - case default - call ExternalPotential_exception( ERROR, & - "Class object ExternalPotential in load function",& - "This Dimensionality is not available") - end select - - !! Normalize - allocate(this%gaussianComponents(i)%contNormalization(this%gaussianComponents(i)%numCartesianOrbital)) - allocate(this%gaussianComponents(i)%primNormalization(this%gaussianComponents(i)%length, & - this%gaussianComponents(i)%length*this%gaussianComponents(i)%numCartesianOrbital)) - - this%gaussianComponents(i)%contNormalization = 1.0_8 - this%gaussianComponents(i)%primNormalization = 1.0_8 - - call ContractedGaussian_normalizePrimitive(this%gaussianComponents(i)) - call ContractedGaussian_normalizeContraction(this%gaussianComponents(i)) - - !! DEBUG - ! call ContractedGaussian_showInCompactForm(ExternalPotential_instance%potentials(potId)%gaussianComponents(i)) - - end do - - close(30) - - !!DONE - - else - - call ExternalPotential_exception(ERROR, & - "The ExternalPotential file: "//trim(name)//& - " was not found!","ExternalPotential module at Load function.") - - end if - - end subroutine ExternalPotential_load - -! !> -! !! @brief -! !! @param this -! function ExternalPotential_getPotential( this, coords ) result(output) -! implicit none -! type(ExternalPotential) :: this -! real(8) :: coords(3) -! real(8) :: output - -! ! integer :: i - -! ! output=0.0 - -! ! do i=1, this%gaussianComponents%length -! ! output = output+( this%gaussianComponents%contractionCoefficients(i)* & -! ! exp(-this%gaussianComponents%primitives(i)%orbitalExponent*( dot_product(coords,coords) ) ) ) -! ! end do - -! end function ExternalPotential_getPotential - - -! ! function ExternalPotential_getInteractionMtx(this, contractions) result(output) -! subroutine ExternalPotential_getInteractionMtx( this, contractions ) -! implicit none -! type(ExternalPotential) :: this -! type(ContractedGaussian) :: contractions(:) -! type(Matrix) :: output - -! ! integer :: i, j, k, l, m, a, b -! ! integer :: numContractions -! ! real(8), allocatable :: auxVal(:) -! ! type(ContractedGaussian) :: auxContract - -! ! do i = 1, size(contractions) -! ! numContractions = numContractions + contractions(i)%numCartesianOrbital -! ! end do - -! ! call Matrix_constructor(output,int(numContractions,8),int(numContractions,8)) - -! ! a = 1 -! ! b = 1 - -! ! do i=1, size(contractions) - -! ! call ContractedGaussian_product(contractions(i), & -! ! this%gaussianComponents, auxContract) - -! ! do j=1, size(contractions) - -! ! call ContractedGaussian_overlapIntegral( auxContract, contractions(j), auxVal) - -! ! m = 0 -! ! do k = a, auxContract%numCartesianOrbital - 1 -! ! do l = b, contractions(j)%numCartesianOrbital - 1 -! ! m = m + 1 -! ! output%values(k,l)= auxVal(m) -! ! end do -! ! end do -! ! b = b + contractions(j)%numCartesianOrbital -! ! end do -! ! a = a + auxContract%numCartesianOrbital -! ! call ContractedGaussian_destructor(auxContract) -! ! end do - -! ! call Matrix_show(output) - -! ! stop "ExternalPotential_getInteractionMtx" - -! ! ! end function ExternalPotential_getInteractionMtx -! end subroutine ExternalPotential_getInteractionMtx - - - !> - !! @brief Maneja excepciones de la clase - subroutine ExternalPotential_exception( typeMessage, description, debugDescription) - implicit none - integer :: typeMessage - character(*) :: description - character(*) :: debugDescription - - type(Exception) :: ex - - call Exception_constructor( ex , typeMessage ) - call Exception_setDebugDescription( ex, debugDescription ) - call Exception_setDescription( ex, description ) - call Exception_show( ex ) - call Exception_destructor( ex ) - - end subroutine ExternalPotential_exception - -end module ExternalPotential_ diff --git a/src/core/GTFPotential.f90 b/src/core/GTFPotential.f90 new file mode 100644 index 00000000..4ac4b93a --- /dev/null +++ b/src/core/GTFPotential.f90 @@ -0,0 +1,309 @@ +!!****************************************************************************** +!! This code is part of LOWDIN Quantum chemistry package +!! +!! this program has been developed under direction of: +!! +!! Prof. A REYES' Lab. Universidad Nacional de Colombia +!! http://www.qcc.unal.edu.co +!! Prof. R. FLORES' Lab. Universidad de Guadajara +!! http://www.cucei.udg.mx/~robertof +!! +!! Todos los derechos reservados, 2013 +!! +!!****************************************************************************** + +!> @brief This module contains all the routines to handle external and interal GTF potentials +!! @author E. F. Posada (efposadac@unal.edu.co) +!! @version 2.0 +!! Creation data : 06-08-10 +!! +!! History change: +!! +!! - 06-08-10 : Sergio A. Gonzalez ( sergmonic@gmail.com ) +!! -# Creacioon del modulo y metodos basicos +!! - 2011-02-14 : Fernando Posada ( efposadac@unal.edu.co ) +!! -# Reescribe y adapta el módulo para su inclusion en Lowdin +!! - 2024-11-26 : Felix +!! -# Merges ExternalPotential and InternalPotentials modules into a single file (GTFPotential) +module GTFPotential_ + use ContractedGaussian_ + use String_ + use Matrix_ + use Units_ + use Exception_ + implicit none + + type, public :: GaussPot + character(20) :: name + character(50) :: species + character(50) :: otherSpecies + character(50) :: ttype + character(50) :: units + integer :: numOfComponents + integer :: iter + type(ContractedGaussian), allocatable :: gaussianComponents(:) + end type GaussPot + + type, public :: GTFPotential + integer :: ssize + type(GaussPot), allocatable :: potentials(:) + character(50) :: type + logical :: isInstanced + end type GTFPotential + + type(GTFPotential), public, target :: ExternalPotential_instance, InterPotential_instance + +contains + + !> + !! @brief Initializes the class + !! @param this, n + !! @author E. F. Posada, 2013 + subroutine GTFPotential_constructor(this,numberOfPotentials,type) + implicit none + type(GTFPotential) :: this + integer :: numberOfPotentials + character(*) :: type + + + this%ssize = numberOfPotentials + allocate(this%potentials(numberOfPotentials)) + this%isInstanced = .true. + this%type = type + + end subroutine GTFPotential_constructor + + !> + !! @brief Destroys the object + !! @param this + subroutine GTFPotential_destructor(this) + implicit none + type(GTFPotential) :: this + + integer :: i + + do i = 1, this%ssize + if (allocated(this%potentials(i)%gaussianComponents)) deallocate(this%potentials(i)%gaussianComponents) + end do + + if (allocated(this%potentials) ) deallocate(this%potentials) + this%isInstanced=.false. + + end subroutine GTFPotential_destructor + + !> + !! @brief loads information from the input file + !! @param this + !! @author E. F. Posada, 2013 + subroutine GTFPotential_load(this, potId, name, species, otherSpecies) + implicit none + type(GTFPotential) :: this + integer :: potId + character(*) :: name + character(*) :: species + character(*), optional :: otherSpecies + + if(present(otherSpecies)) then + call GaussPot_load(this%potentials(potId), potId, name, species, otherSpecies) + else + call GaussPot_load(this%potentials(potId), potId, name, species) + end if + end subroutine GTFPotential_load + + !> + !! @brief Shows information of the object + !! @param this + subroutine GTFPotential_show(this) + implicit none + type(GTFPotential) :: this + integer :: i, j + + do i=1,this%ssize + if( this%potentials(i)%ttype .eq. "INTERNAL") then + print *,"" + print *,"=======" + write(*,"(A30,A)") "GTF Interparticle potential: ", trim(this%potentials(i)%name) + write(*,"(A4,A10,A5,A10)") "for ", trim(this%potentials(i)%species) ," and ", trim(this%potentials(i)%otherSpecies) + write(*,"(T10,A10,A10)") "Units:", trim(this%potentials(i)%units) + write(*,"(T10,A16,A16)") "Exponent", "Factor" + do j=1,this%potentials(i)%numOfComponents + write(*,"(T10,E16.8,E16.8)") this%potentials(i)%gaussianComponents(j)%orbitalExponents, & + this%potentials(i)%gaussianComponents(j)%contractionCoefficients(1) + end do + else if( this%potentials(i)%ttype .eq. "EXTERNAL") then + print *,"" + print *,"=======" + write(*,"(A25,A20,A5,A10)") "GTF External potential: ", trim(this%potentials(i)%name), " for ", trim(this%potentials(i)%species) + write(*,"(T10,A10,A10)") "Units:", trim(this%potentials(i)%units) + write(*,"(T10,A16,A10,A10,A10,A16)") "Exponent", "R_x", "R_y", "R_z", "Factor" + + do j=1,this%potentials(i)%numOfComponents + write(*,"(T10,E16.8,F10.5,F10.5,F10.5,E16.8)") & + this%potentials(i)%gaussianComponents(j)%orbitalExponents(1), & + this%potentials(i)%gaussianComponents(j)%origin(:), & + this%potentials(i)%gaussianComponents(j)%contractionCoefficients(1) + end do + end if + end do + + end subroutine GTFPotential_show + + !> + !! @brief loads information from the input file + !! @param this + !! @author Felix, 2024 + subroutine GaussPot_load(this, potId, name, species, otherSpecies) + type(GaussPot) :: this + integer :: potId + character(*) :: name + character(*) :: species + character(*), optional :: otherSpecies + + integer :: status, i, j + character(150) :: fileName + character(20) :: token + character(50) :: symbol + logical :: existFile, found + + this%name= trim(name) + this%species= trim(species) + this%units="BOHR" + this%numOfComponents=0 + this%iter=1 + + this%ttype="EXTERNAL" + this%otherSpecies="" + if(present(otherSpecies) ) then + this%ttype="INTERNAL" + this%otherSpecies=otherSpecies + end if + + fileName = trim( trim( CONTROL_instance%DATA_DIRECTORY ) // & + trim(CONTROL_instance%POTENTIALS_DATABASE)// String_getUppercase(trim(name))) + + !! Open Potential file from library + inquire(file=trim(fileName), exist = existFile) + if(existFile) then + open(unit=30, file=trim(fileName), status="old",form="formatted") + else + !! Open Potential file from directory + inquire(file=trim(this%name), exist = existFile) + if(existFile) then + open(unit=30, file=trim(this%name), status="old",form="formatted") + else + !! File not found + call Exception_stopError("The GTFPotential file: "//trim(this%name)//& + " was not found!","GTFPotential module at Load function.") + end if + end if + + !! Open File + rewind(30) + + found = .false. + + !! Open element and Find the proper potential + do while(found .eqv. .false.) + read(30,*, iostat=status) token + symbol = token(3:) + + !! Some debug information in case of error! + if (status > 0 ) call Exception_stopError("ERROR reading InterPotential file: "//trim(this%name)//& + " Please check that file!","GTFPotential module at Load function.") + + if (status == -1 ) call Exception_stopError("The InterPotential: "//trim(this%name)//& + " for: "//trim(species)//trim(otherSpecies)//& + " was not found!","GTFPotential module at Load function.") + + if(trim(token(1:2)) == "O-") then + if(this%ttype=="EXTERNAL" .and. trim(symbol) == trim(species)) found = .true. + if(this%ttype=="INTERNAL" .and. trim(symbol) == trim(species)//trim(otherSpecies)) found = .true. + end if + end do + + !! Neglect any comment + token = "#" + do while(trim(token(1:1)) == "#") + read(30,*) token + end do + + !! Start reading Potential + backspace(30) + read(30,*, iostat=status) this%numOfComponents + + !! Some debug information in case of error! + if (status > 0 ) call Exception_stopError("ERROR reading InternalPotential file: "//trim(this%name)//& + " Please check that file!","GTFPotential module at Load function.") + + allocate(this%gaussianComponents(this%numOfComponents)) + + do i = 1, this%numOfComponents + read(30,*,iostat=status) this%gaussianComponents(i)%id, & + this%gaussianComponents(i)%angularMoment + + if(this%gaussianComponents(i)%angularMoment .gt. 0) then + print *, "Warning! you provided a non-zero angular momentum for a GTFpotential ", this%name ,"this feature is not yet implemented, will be ignored and set to zeo" + this%gaussianComponents(i)%angularMoment=0 + end if + + this%gaussianComponents(i)%length = 1 + + !! Some debug information in case of error! + if (status > 0 ) call Exception_stopError("ERROR reading InternalPotential file: "//trim(this%name)//& + " Please check that file!","GTFPotential module at Load function.") + + allocate(this%gaussianComponents(i)%orbitalExponents(this%gaussianComponents(i)%length)) + allocate(this%gaussianComponents(i)%contractionCoefficients(this%gaussianComponents(i)%length)) + + do j = 1, this%gaussianComponents(i)%length + + read(30,*,iostat=status) this%gaussianComponents(i)%orbitalExponents(j), & + this%gaussianComponents(i)%contractionCoefficients(j) + read(30,*,iostat=status) this%gaussianComponents(i)%origin + + !! Some debug information in case of error! + if (status > 0 ) call Exception_stopError("ERROR reading InternalPotential file: "//trim(this%name)//& + " Please check that file!","GTFPotential module at Load function.") + + end do + + if(this%ttype=="INTERNAL" .and. sum(this%gaussianComponents(i)%origin(:)**2) .gt. CONTROL_instance%DOUBLE_ZERO_THRESHOLD) then + print *, "Warning! you provided a non-zero origin for interpotential ", this%name ,"this feature is not yet implemented, will be ignored and set to zeo" + this%gaussianComponents(i)%origin=0.0_8 + end if + + !! Calculates the number of Cartesian orbitals, by dimensionality + select case(CONTROL_instance%DIMENSIONALITY) + case(3) + this%gaussianComponents(i)%numCartesianOrbital = ( ( this%gaussianComponents(i)%angularMoment + 1_8 )*( this%gaussianComponents(i)%angularMoment + 2_8 ) ) / 2_8 + case(2) + this%gaussianComponents(i)%numCartesianOrbital = ( ( this%gaussianComponents(i)%angularMoment + 1_8 ) ) + case(1) + this%gaussianComponents(i)%numCartesianOrbital = 1 + case default + call Exception_stopError("Class object InternalPotential in load function",& + "This Dimensionality is not available") + end select + + !! Normalize + allocate(this%gaussianComponents(i)%contNormalization(this%gaussianComponents(i)%numCartesianOrbital)) + allocate(this%gaussianComponents(i)%primNormalization(this%gaussianComponents(i)%length, & + this%gaussianComponents(i)%length*this%gaussianComponents(i)%numCartesianOrbital)) + + this%gaussianComponents(i)%contNormalization = 1.0_8 + this%gaussianComponents(i)%primNormalization = 1.0_8 + + call ContractedGaussian_normalizePrimitive(this%gaussianComponents(i)) + call ContractedGaussian_normalizeContraction(this%gaussianComponents(i)) + + !! DEBUG + ! call ContractedGaussian_showInCompactForm(InterPotential_instance%potentials(potId)%gaussianComponents(i)) + + end do + + close(30) + + !!DONE + end subroutine GaussPot_load + +end module GTFPotential_ diff --git a/src/core/InputManager.f90 b/src/core/InputManager.f90 index e5a0abde..506a89fb 100644 --- a/src/core/InputManager.f90 +++ b/src/core/InputManager.f90 @@ -23,8 +23,7 @@ module InputManager_ use Exception_ use Particle_ use MolecularSystem_ - use InterPotential_ - use ExternalPotential_ + use GTFPotential_ implicit none @@ -344,6 +343,7 @@ subroutine InputManager_loadGeometry() real(8):: InputParticle_origin(3) real(8) :: InputParticle_charge real(8) :: InputParticle_mass + integer :: InputParticle_eta real(8) :: InputParticle_omega character(15):: InputParticle_qdoCenterOf character(3):: InputParticle_fixedCoordinates @@ -359,6 +359,7 @@ subroutine InputManager_loadGeometry() InputParticle_basisSetName, & InputParticle_charge, & InputParticle_mass, & + InputParticle_eta, & InputParticle_omega, & InputParticle_qdoCenterOf, & InputParticle_origin, & @@ -538,6 +539,7 @@ subroutine InputManager_loadGeometry() InputParticle_basisSetName = "NONE" InputParticle_charge=0.0_8 InputParticle_mass=0.0_8 + InputParticle_eta=0 InputParticle_omega=0.0_8 InputParticle_qdoCenterOf = "NONE" InputParticle_origin=0.0_8 @@ -620,7 +622,8 @@ subroutine InputManager_loadGeometry() spin="ALPHA", & id = particlesID(speciesID), & charge = InputParticle_charge, & - mass = InputParticle_mass, & + mass = InputParticle_mass, & + eta = InputParticle_eta, & omega = InputParticle_omega ) !!BETA SET @@ -644,6 +647,7 @@ subroutine InputManager_loadGeometry() id = particlesID(speciesID), & charge = InputParticle_charge, & mass = InputParticle_mass, & + eta = InputParticle_eta, & omega = InputParticle_omega ) else @@ -671,6 +675,7 @@ subroutine InputManager_loadGeometry() id = particlesID(speciesID), & charge = InputParticle_charge, & mass = InputParticle_mass, & + eta = InputParticle_eta, & omega = InputParticle_omega ) end if @@ -750,7 +755,7 @@ subroutine InputManager_loadPotentials() ! Load interpotentials if(CONTROL_instance%IS_THERE_INTERPARTICLE_POTENTIAL) then - call InterPotential_constructor(Input_instance%numberOfInterPots) + call GTFPotential_constructor(InterPotential_instance, Input_instance%numberOfInterPots,"INTERNAL") !! Reload input file rewind(4) @@ -760,10 +765,10 @@ subroutine InputManager_loadPotentials() read(4,NML=InterPot, iostat=stat) if( stat > 0 ) then - call InputManager_exception( ERROR, "check the TASKS block in your input file", "InputManager loadTask function" ) + call InputManager_exception( ERROR, "check the INTERPOTENTIAL block in your input file", "InputManager loadTask function" ) end if - call InterPotential_load(potId, trim(InterPot_name), trim(InterPot_specie), trim(InterPot_otherSpecie)) + call GTFPotential_load(InterPotential_instance, potId, trim(InterPot_name), trim(InterPot_specie), trim(InterPot_otherSpecie)) end do @@ -772,7 +777,7 @@ subroutine InputManager_loadPotentials() ! Load External Potentials if(CONTROL_instance%IS_THERE_EXTERNAL_POTENTIAL) then - call ExternalPotential_constructor(Input_instance%numberOfExternalPots) + call GTFPotential_constructor(ExternalPotential_instance, Input_instance%numberOfExternalPots,"EXTERNAL") !! Reload input file rewind(4) @@ -782,10 +787,10 @@ subroutine InputManager_loadPotentials() read(4,NML=ExternalPot, iostat=stat) if( stat > 0 ) then - call InputManager_exception( ERROR, "check the TASKS block in your input file", "InputManager loadTask function" ) + call InputManager_exception( ERROR, "check the EXTERPOTENTIAL block in your input file", "InputManager loadTask function" ) end if - call ExternalPotential_load(potId, trim(ExternalPot_name), trim(ExternalPot_specie)) + call GTFPotential_load(ExternalPotential_instance, potId, trim(ExternalPot_name), trim(ExternalPot_specie)) end do end if diff --git a/src/core/InterPotential.f90 b/src/core/InterPotential.f90 deleted file mode 100644 index f30008e7..00000000 --- a/src/core/InterPotential.f90 +++ /dev/null @@ -1,298 +0,0 @@ -!!****************************************************************************** -!! This code is part of LOWDIN Quantum chemistry package -!! -!! this program has been developed under direction of: -!! -!! Prof. A REYES' Lab. Universidad Nacional de Colombia -!! http://www.qcc.unal.edu.co -!! Prof. R. FLORES' Lab. Universidad de Guadajara -!! http://www.cucei.udg.mx/~robertof -!! -!! Todos los derechos reservados, 2013 -!! -!!****************************************************************************** - -!> @brief This module contains all the routines to handle inter particle potentials -!! @author E. F. Posada (efposadac@unal.edu.co) -!! @version 2.0 - -module InterPotential_ - use ContractedGaussian_ - use String_ - use Exception_ - implicit none - - type, public :: InterPot - character(20) :: name - character(50) :: specie - character(50) :: otherSpecie - character(50) :: ttype - character(50) :: units - integer :: numOfComponents - integer :: iter - type(ContractedGaussian), allocatable :: gaussianComponents(:) - end type - - type, public :: InterPotential - integer :: ssize - type(InterPot), allocatable :: Potentials(:) - logical :: isInstanced - end type - - type(InterPotential), public, target :: InterPotential_instance - -contains - - !> - !! @brief Initializes the class - !! @param this - !! @author E. F. Posada, 2013 - subroutine InterPotential_constructor(numberOfPotentials) - implicit none - integer :: numberOfPotentials - - InterPotential_instance%ssize = numberOfPotentials - allocate(InterPotential_instance%potentials(numberOfPotentials)) - InterPotential_instance%isInstanced = .true. - - end subroutine InterPotential_constructor - - !> - !! @brief destroy the class - !! @param this - !! @author E. F. Posada, 2013 - subroutine InterPotential_destructor() - implicit none - - integer :: i - - do i = 1, InterPotential_instance%ssize - if (allocated(InterPotential_instance%potentials(i)%gaussianComponents) ) deallocate(InterPotential_instance%potentials(i)%gaussianComponents) - end do - - if (allocated(InterPotential_instance%potentials) )deallocate(InterPotential_instance%potentials) - - end subroutine InterPotential_destructor - - !> - !! @brief Shows information of the object - !! @param this - subroutine InterPotential_show() - implicit none - integer :: i, j - type(InterPot), pointer :: this - - do i=1,InterPotential_instance%ssize - this => InterPotential_instance%potentials(i) - print *,"" - print *,"=======" - print *, "InterParticle Potential for ", trim(this%specie) ," and ", trim(this%otherSpecie), " : ", trim(this%name) - print *, "Type : ", trim(this%ttype) - write(6,"(T10,A20,A10,A10,A10,A20)") "Exponent", "l", "Factor" - do j=1,this%numOfComponents - write(6,"(T10,F20.15,I10,F20.15)") this%gaussianComponents(j)%orbitalExponents, & - this%gaussianComponents(j)%angularMoment, this%gaussianComponents(j)%contractionCoefficients(1) - end do - end do - - end subroutine InterPotential_show - - - !> - !! @brief loads information from the input file - !! @param this - !! @author E. F. Posada, 2015 - subroutine InterPotential_load(potId, name, species, otherSpecies) - implicit none - integer :: potId - character(*) :: name - character(*) :: species - character(*) :: otherSpecies - - type(InterPot), pointer :: this - integer :: status, i, j - character(150) :: fileName - character(20) :: token - character(20) :: symbol - logical :: existFile, found - - this => InterPotential_instance%potentials(potId) - - this%name= trim(name) - this%specie= trim(species) - this%otherSpecie= trim(otherSpecies) - this%ttype="" - this%units="bohr" - this%numOfComponents=0 - this%iter=1 - - fileName = trim( trim( CONTROL_instance%DATA_DIRECTORY ) // & - trim(CONTROL_instance%POTENTIALS_DATABASE)// String_getUppercase(trim(name))) - - - inquire(file=trim(fileName), exist = existFile) - if(existFile) then - - !! Open File - open(unit=30, file=trim(fileName), status="old",form="formatted") - rewind(30) - - found = .false. - - !! Open element and Find the proper potential - do while(found .eqv. .false.) - read(30,*, iostat=status) token - symbol = token(3:) - - !! Some debug information in case of error! - if (status > 0 ) then - - call InterPotential_exception(ERROR, & - "ERROR reading InterPotential file: "//trim(this%name)//& - " Please check that file!","InternalPotential module at Load function.") - - end if - - if (status == -1 ) then - - call InterPotential_exception(ERROR, & - "The InterPotential: "//trim(this%name)//& - " for: "//trim(species)//trim(otherSpecies)//& - " was not found!","InternalPotential module at Load function.") - - end if - - if(trim(token(1:2)) == "O-") then - if(trim(symbol) == trim(species)//trim(otherSpecies)) then - found = .true. - - end if - - end if - - end do - - !! Neglect any comment - token = "#" - do while(trim(token(1:1)) == "#") - - read(30,*) token - - end do - - !! Start reading Potential - backspace(30) - - read(30,*, iostat=status) this%numOfComponents - - !! Some debug information in case of error! - if (status > 0 ) then - - call InterPotential_exception(ERROR, & - "ERROR reading InternalPotential file: "//trim(this%name)//& - " Please check that file!","InternalPotential module at Load function.") - - end if - - allocate(this%gaussianComponents(this%numOfComponents)) - - do i = 1, this%numOfComponents - - read(30,*,iostat=status) this%gaussianComponents(i)%id, & - this%gaussianComponents(i)%angularMoment - this%gaussianComponents(i)%length = 1 - - !! Some debug information in case of error! - if (status > 0 ) then - - call InterPotential_exception(ERROR, & - "ERROR reading InternalPotential file: "//trim(this%name)//& - " Please check that file!","InternalPotential module at Load function.") - - end if - - allocate(this%gaussianComponents(i)%orbitalExponents(this%gaussianComponents(i)%length)) - allocate(this%gaussianComponents(i)%contractionCoefficients(this%gaussianComponents(i)%length)) - - do j = 1, this%gaussianComponents(i)%length - - read(30,*,iostat=status) this%gaussianComponents(i)%orbitalExponents(j), & - this%gaussianComponents(i)%contractionCoefficients(j) - read(30,*,iostat=status) this%gaussianComponents(i)%origin - - !! Some debug information in case of error! - if (status > 0 ) then - - call InterPotential_exception(ERROR, & - "ERROR reading InternalPotential file: "//trim(this%name)//& - " Please check that file!","InternalPotential module at Load function.") - - end if - - end do - - - !! Calculates the number of Cartesian orbitals, by dimensionality - select case(CONTROL_instance%DIMENSIONALITY) - case(3) - this%gaussianComponents(i)%numCartesianOrbital = ( ( this%gaussianComponents(i)%angularMoment + 1_8 )*( this%gaussianComponents(i)%angularMoment + 2_8 ) ) / 2_8 - case(2) - this%gaussianComponents(i)%numCartesianOrbital = ( ( this%gaussianComponents(i)%angularMoment + 1_8 ) ) - case(1) - this%gaussianComponents(i)%numCartesianOrbital = 1 - case default - call InterPotential_exception( ERROR, & - "Class object InternalPotential in load function",& - "This Dimensionality is not available") - end select - - !! Normalize - allocate(this%gaussianComponents(i)%contNormalization(this%gaussianComponents(i)%numCartesianOrbital)) - allocate(this%gaussianComponents(i)%primNormalization(this%gaussianComponents(i)%length, & - this%gaussianComponents(i)%length*this%gaussianComponents(i)%numCartesianOrbital)) - - this%gaussianComponents(i)%contNormalization = 1.0_8 - this%gaussianComponents(i)%primNormalization = 1.0_8 - - call ContractedGaussian_normalizePrimitive(this%gaussianComponents(i)) - call ContractedGaussian_normalizeContraction(this%gaussianComponents(i)) - - !! DEBUG - ! call ContractedGaussian_showInCompactForm(InterPotential_instance%potentials(potId)%gaussianComponents(i)) - - end do - - close(30) - - !!DONE - - else - - call InterPotential_exception(ERROR, & - "The InternalPotential file: "//trim(name)//& - " was not found!","InternalPotential module at Load function.") - - end if - - end subroutine InterPotential_load - - - !> - !! @brief Handles class exceptions - !< - subroutine InterPotential_exception( typeMessage, description, debugDescription) - implicit none - integer :: typeMessage - character(*) :: description - character(*) :: debugDescription - - type(Exception) :: ex - - call Exception_constructor( ex , typeMessage ) - call Exception_setDebugDescription( ex, debugDescription ) - call Exception_setDescription( ex, description ) - call Exception_show( ex ) - call Exception_destructor( ex ) - - end subroutine InterPotential_exception -end module InterPotential_ diff --git a/src/core/MolecularSystem.f90 b/src/core/MolecularSystem.f90 index b5a766fe..32dced1b 100644 --- a/src/core/MolecularSystem.f90 +++ b/src/core/MolecularSystem.f90 @@ -38,8 +38,7 @@ module MolecularSystem_ use Matrix_ use Vector_ use InternalCoordinates_ - use ExternalPotential_ - use InterPotential_ + use GTFPotential_ implicit none type , public :: MolecularSystem @@ -343,9 +342,10 @@ subroutine MolecularSystem_showInformation(this) print *," MOLECULAR SYSTEM: ",trim(system%name) print *,"-----------------" print *,"" - write (6,"(T5,A16,A)") "DESCRIPTION : ", trim( system%description ) - write (6,"(T5,A16,I3)") "CHARGE : ",system%charge - write (6,"(T5,A16,A4)") "PUNTUAL GROUP : ", "NONE" + write (6,"(T5,A16,A)") "DESCRIPTION : ", trim( system%description ) + write (6,"(T5,A16,I3)") "CHARGE : ",system%charge + write (6,"(T5,A16,F12.4)") "MASS (m_e) : ", MolecularSystem_getTotalMass(system) + write (6,"(T5,A16,A4)") "PUNTUAL GROUP : ", "NONE" print *,"" @@ -382,9 +382,9 @@ subroutine MolecularSystem_showParticlesInformation(this) !! print *,"" print *," INFORMATION OF QUANTUM SPECIES " - write (6,"(T5,A70)") "---------------------------------------------------------------------------------------------" + write (6,"(T5,A85)") "------------------------------------------------------------------------------------------------------------" write (6,"(T10,A2,A4,A8,A12,A4,A5,A6,A5,A6,A5,A4,A5,A12)") "ID", " ","Symbol", " ","mass", " ","charge", " ","omega","","spin","","multiplicity" - write (6,"(T5,A70)") "---------------------------------------------------------------------------------------------" + write (6,"(T5,A85)") "------------------------------------------------------------------------------------------------------------" do i = 1, system%numberOfQuantumSpecies write (6,'(T8,I3.0,A5,A10,A5,F10.4,A5,F5.2,A5,F5.2,A5,F5.2,A5,F5.2)') & @@ -415,16 +415,16 @@ subroutine MolecularSystem_showParticlesInformation(this) print *,"" print *," BASIS SET FOR SPECIES " - write (6,"(T7,A60)") "------------------------------------------------------------" - write (6,"(T10,A8,A10,A8,A5,A12,A5,A9)") "Symbol", " ","N. Basis", " ","N. Particles"," ","Basis Set" - write (6,"(T7,A60)") "------------------------------------------------------------" + write (6,"(T7,A70)") "----------------------------------------------------------------------" + write (6,"(T10,A8,A10,A8,A5,A12,A5,A20)") "Symbol", " ","N. Basis", " ","N. Particles"," ","Basis Set" + write (6,"(T7,A70)") "----------------------------------------------------------------------" !! Only shows the basis-set name of the first particle by specie. do i = 1, system%numberOfQuantumSpecies if( system%species(i)%isElectron .and. CONTROL_instance%IS_OPEN_SHELL ) then - write (6,'(T10,A10,A5,I8,A5,I12,A5,A10)') & + write (6,'(T10,A10,A5,I8,A5,I12,A5,A20)') & trim(system%species(i)%symbol)," ",& !MolecularSystem_getTotalNumberOfContractions(i)," ",& system%species(i)%basisSetSize," ",& @@ -434,7 +434,7 @@ subroutine MolecularSystem_showParticlesInformation(this) else - write (6,'(T10,A10,A5,I8,A5,I12,A5,A10)') & + write (6,'(T10,A10,A5,I8,A5,I12,A5,A20)') & trim(system%species(i)%symbol)," ",& !MolecularSystem_getTotalNumberOfContractions(i)," ",& system%species(i)%basisSetSize," ",& @@ -501,7 +501,7 @@ subroutine MolecularSystem_showParticlesInformation(this) if(CONTROL_instance%IS_THERE_EXTERNAL_POTENTIAL) then print *,"" print *," INFORMATION OF EXTERNAL POTENTIALS " - call ExternalPotential_show() + call GTFPotential_show(ExternalPotential_instance) print *,"" print *," END INFORMATION OF EXTERNAL POTENTIALS" print *,"" @@ -510,7 +510,7 @@ subroutine MolecularSystem_showParticlesInformation(this) if(CONTROL_instance%IS_THERE_INTERPARTICLE_POTENTIAL) then print *,"" print *," INFORMATION OF INTER-PARTICLE POTENTIALS " - call InterPotential_show() + call GTFPotential_show(InterPotential_instance) print *,"" print *," END INFORMATION OF INTER-PARTICLE POTENTIALS" print *,"" @@ -671,7 +671,7 @@ subroutine MolecularSystem_saveToFile(targetFilePrefix) do i = 1, ExternalPotential_instance%ssize write(40,*) i write(40,*) ExternalPotential_instance%potentials(i)%name - write(40,*) ExternalPotential_instance%potentials(i)%specie + write(40,*) ExternalPotential_instance%potentials(i)%species end do end if @@ -681,8 +681,8 @@ subroutine MolecularSystem_saveToFile(targetFilePrefix) do i = 1, InterPotential_instance%ssize write(40,*) i write(40,*) InterPotential_instance%potentials(i)%name - write(40,*) InterPotential_instance%potentials(i)%specie - write(40,*) InterPotential_instance%potentials(i)%otherSpecie + write(40,*) InterPotential_instance%potentials(i)%species + write(40,*) InterPotential_instance%potentials(i)%otherSpecies end do end if @@ -950,16 +950,14 @@ subroutine MolecularSystem_loadFromFile( form, targetPrefix ) if(CONTROL_instance%IS_THERE_EXTERNAL_POTENTIAL) then read(40,*) auxValue - call ExternalPotential_constructor(auxValue) - - !! FELIX TODO: create function to get potential ID + call GTFPotential_constructor(ExternalPotential_instance,auxValue,"EXTERNAL") do j = 1, ExternalPotential_instance%ssize read(40,*) i read(40,*) name read(40,*) species - call ExternalPotential_load(i, name, species) + call GTFPotential_load(ExternalPotential_instance, i, name, species) end do @@ -968,7 +966,7 @@ subroutine MolecularSystem_loadFromFile( form, targetPrefix ) if(CONTROL_instance%IS_THERE_INTERPARTICLE_POTENTIAL) then read(40,*) auxValue - call InterPotential_constructor(auxValue) + call GTFPotential_constructor(InterPotential_instance,auxValue,"INTERNAL") do j = 1, InterPotential_instance%ssize read(40,*) i @@ -976,7 +974,7 @@ subroutine MolecularSystem_loadFromFile( form, targetPrefix ) read(40,*) species read(40,*) otherSpecies - call InterPotential_load(i, name, species, otherSpecies) + call GTFPotential_load(InterPotential_instance, i, name, species, otherSpecies) end do @@ -1979,13 +1977,14 @@ end subroutine MolecularSystem_changeOrbitalOrder !> !! @brief Lee la matriz de densidad y los orbitales de un archivo fchk tipo Gaussian - subroutine MolecularSystem_readFchk( fileName, coefficients, densityMatrix, nameOfSpecies ) + subroutine MolecularSystem_readFchk( fileName, coefficients, densityMatrix, nameOfSpecies, readSuccess ) implicit none character(*), intent(in) :: fileName type(Matrix), intent(inout) :: coefficients type(Matrix), intent(inout) :: densityMatrix character(*) :: nameOfSpecies + logical, optional :: readSuccess integer :: speciesID integer :: numberOfContractions @@ -2006,9 +2005,15 @@ subroutine MolecularSystem_readFchk( fileName, coefficients, densityMatrix, name speciesID=MolecularSystem_getSpecieID(nameOfSpecies) numberOfContractions=MolecularSystem_getTotalnumberOfContractions( speciesID ) inquire(FILE = trim(fileName), EXIST = existFchk ) - if ( .not. existFchk ) call MolecularSystem_exception( ERROR, "I did not find any .fchk coefficients file", "At MolecularSystem_readFchk") + if ( .not. existFchk .and. present(readSuccess)) then + readSuccess=.false. + call MolecularSystem_exception( WARNING, "I did not find the "//trim(filename)//" coefficients file for "//trim(nameOfSpecies), "At MolecularSystem_readFchk") + return + end if + if ( .not. existFchk) then + call MolecularSystem_exception( ERROR, "I did not find the "//trim(filename)//" coefficients file for "//trim(nameOfSpecies), "At MolecularSystem_readFchk") + end if - fchkUnit = 50 open(unit=fchkUnit, file=filename, status="old", form="formatted", access='sequential', action='read') @@ -2097,8 +2102,8 @@ subroutine MolecularSystem_readFchk( fileName, coefficients, densityMatrix, name ! print *, "density matrix from orbitals read" ! call Matrix_show(densityMatrix) - close(fchkUnit) + if(present(readSuccess)) readSuccess=.true. end subroutine MolecularSystem_readFchk @@ -2783,6 +2788,9 @@ function MolecularSystem_getTotalMass( this, unid ) result( output ) case("AMU") output = output * AMU + case("DALTON") + output = output * DALTON + case default end select diff --git a/src/core/Particle.f90 b/src/core/Particle.f90 index 7a845404..f927d006 100644 --- a/src/core/Particle.f90 +++ b/src/core/Particle.f90 @@ -32,6 +32,7 @@ module Particle_ use Exception_ use CONTROL_ + use Units_ use PhysicalConstants_ use AtomicElement_ use ElementalParticle_ @@ -49,6 +50,7 @@ module Particle_ real(8) :: origin(3) !< Posicion espacial real(8) :: charge !< Carga asociada a la particula. real(8) :: mass !< Masa asociada a la particula. + real(8) :: eta !< Particles per orbital real(8) :: omega !< harmonic oscillator frequency real(8) :: spin !< Especifica el espin de la particula real(8) :: totalCharge !< Carga total asociada a la particula. @@ -83,7 +85,7 @@ module Particle_ !! -Re-written and Verified, 2013. E. F. Posada !! @version 2.0 subroutine Particle_load( this, name, baseName, origin, fix, multiplicity, & - addParticles, subsystem, translationCenter, rotationPoint, rotateAround, spin, id, charge, mass, omega, qdoCenterOf ) + addParticles, subsystem, translationCenter, rotationPoint, rotateAround, spin, id, charge, mass, eta, omega, qdoCenterOf ) implicit none type(particle) :: this character(*), intent(in) :: name @@ -101,6 +103,7 @@ subroutine Particle_load( this, name, baseName, origin, fix, multiplicity, & integer, intent(in) :: id real(8), intent(in), optional :: charge real(8), intent(in), optional :: mass + integer, intent(in), optional :: eta real(8), intent(in), optional :: omega type(AtomicElement) :: element @@ -120,6 +123,7 @@ subroutine Particle_load( this, name, baseName, origin, fix, multiplicity, & real(8) :: auxOrigin(3) real(8) :: auxCharge real(8) :: auxMass + integer :: auxEta real(8) :: auxOmega real(8) :: auxMultiplicity logical :: isDummy @@ -137,6 +141,9 @@ subroutine Particle_load( this, name, baseName, origin, fix, multiplicity, & auxMass=0.0_8 if ( present(mass) ) auxMass= mass + auxEta=0 + if ( present(eta) ) auxEta= eta + auxOmega=0.0_8 if ( present(omega) ) auxOmega= omega @@ -219,6 +226,7 @@ subroutine Particle_load( this, name, baseName, origin, fix, multiplicity, & origin=auxOrigin, & charge=auxCharge, & mass=auxMass, & + eta=auxEta, & omega=auxOmega, & basisSetName=trim(baseName), & elementSymbol=trim(elementSymbol), & @@ -309,6 +317,7 @@ subroutine Particle_load( this, name, baseName, origin, fix, multiplicity, & origin = auxOrigin,& charge = auxCharge,& mass = auxMass,& + eta=auxEta, & omega=auxOmega, & basisSetName = trim(baseName), & elementSymbol = trim(elementSymbol), & @@ -330,10 +339,10 @@ subroutine Particle_load( this, name, baseName, origin, fix, multiplicity, & this%charge = element%atomicNumber end if - if ( auxMass == 0.0_8) then - this%mass = element%massicNumber * PhysicalConstants_NEUTRON_MASS & + if ( auxMass .eq. 0.0_8 .and. element%atomicWeight .eq. 0.0_8 ) & + this%mass = element%massicNumber * PhysicalConstants_NEUTRON_MASS & + element%atomicNumber * (PhysicalConstants_PROTON_MASS - PhysicalConstants_NEUTRON_MASS) - end if + if ( auxMass .eq. 0.0_8 ) this%mass = element%atomicWeight*DALTON - element%atomicNumber*PhysicalConstants_ELECTRON_MASS this%totalCharge = element%atomicNumber this%internalSize = 1 + auxAdditionOfParticles @@ -393,8 +402,14 @@ subroutine Particle_load( this, name, baseName, origin, fix, multiplicity, & this%charge = element%atomicNumber this%totalCharge = element%atomicNumber end if - this%mass = element%massicNumber * PhysicalConstants_NEUTRON_MASS & - + element%atomicNumber * (PhysicalConstants_PROTON_MASS - PhysicalConstants_NEUTRON_MASS) + + + if ( element%atomicWeight .eq. 0.0_8 ) then + this%mass = element%massicNumber * PhysicalConstants_NEUTRON_MASS & + + element%atomicNumber * (PhysicalConstants_PROTON_MASS - PhysicalConstants_NEUTRON_MASS) + else + this%mass = element%atomicWeight*DALTON - element%atomicNumber*PhysicalConstants_ELECTRON_MASS + end if this%internalSize = 1 + auxAdditionOfParticles this%id = id @@ -439,8 +454,14 @@ subroutine Particle_load( this, name, baseName, origin, fix, multiplicity, & if ( auxCharge == 0.0_8) then this%charge = element%atomicNumber end if - this%mass = element%massicNumber * PhysicalConstants_NEUTRON_MASS & - + element%atomicNumber * (PhysicalConstants_PROTON_MASS - PhysicalConstants_NEUTRON_MASS) + + if ( element%atomicWeight .eq. 0.0_8 ) then + this%mass = element%massicNumber * PhysicalConstants_NEUTRON_MASS & + + element%atomicNumber * (PhysicalConstants_PROTON_MASS - PhysicalConstants_NEUTRON_MASS) + else + this%mass = element%atomicWeight*DALTON - element%atomicNumber*PhysicalConstants_ELECTRON_MASS + end if + this%totalCharge = element%atomicNumber this%internalSize = 1 + auxAdditionOfParticles this%id = id @@ -468,12 +489,19 @@ subroutine Particle_load( this, name, baseName, origin, fix, multiplicity, & else if ( present(baseName) .and. trim(baseName) /= "DIRAC" .and. trim(baseName) /= "MM") then call ElementalParticle_load( eparticle, trim(name) ) + + if(eparticle%custom .and. auxMass .ne. 0.0_8 .and. auxCharge .ne. 0.0_8) & + print *, "I'm loading a custom particle from the input", trim(name), "with mass", auxMass, "and charge", auxCharge + if(eparticle%custom .and. auxMass .eq. 0.0_8 .and. auxCharge .eq. 0.0_8) & + call Particle_exception( ERROR, "Elemental particle: "//trim(name)//& + " NOT found in ElementalParticles.lib. If you want to use a custom particle, provide at least its mass and charge in the input", "In Particle_load routine") call Particle_build( this = this, & isQuantum = .true., & origin = auxOrigin, & charge = auxCharge, & mass = auxMass, & + eta=auxEta, & omega = auxOmega, & basisSetName = trim(baseName), & elementSymbol = trim(elementSymbol), & @@ -586,7 +614,7 @@ end subroutine Particle_load !! @param this Quantum particle or point charge !! @author S. A. Gonzalez (before known as Particle_constructor) subroutine Particle_build( this, name, symbol, basisSetName, elementSymbol, nickname, & - origin, mass, omega, qdoCenterOf, charge, totalCharge, spin, & + origin, mass, eta, omega, qdoCenterOf, charge, totalCharge, spin, & owner, subsystem, translationCenter, rotationPoint, rotateAround, massNumber, isQuantum, isDummy) implicit none @@ -600,6 +628,7 @@ subroutine Particle_build( this, name, symbol, basisSetName, elementSymbol, nick character(*), optional, intent(in) :: qdoCenterOf real(8), optional, intent(in) :: origin(3) real(8), optional, intent(in) :: mass + integer, optional, intent(in) :: eta real(8), optional, intent(in) :: omega real(8), optional, intent(in) :: charge real(8), optional, intent(in) :: totalCharge @@ -621,6 +650,7 @@ subroutine Particle_build( this, name, symbol, basisSetName, elementSymbol, nick this%isQuantum = .false. this%origin = 0.0_8 this%mass = PhysicalConstants_ELECTRON_MASS + this%eta = 0 this%omega = 0.0_8 this%qdoCenterOf = "NONE" this%charge = PhysicalConstants_ELECTRON_CHARGE @@ -646,6 +676,7 @@ subroutine Particle_build( this, name, symbol, basisSetName, elementSymbol, nick !! Loads optional information if ( present(origin) ) this%origin = origin if ( present( mass ) ) this%mass = mass + if ( present( eta ) ) this%eta = eta if ( present( omega ) ) this%omega = omega if ( present( qdoCenterOf ) ) this%qdoCenterOf = qdoCenterOf if ( present(charge) ) this%charge = charge @@ -722,6 +753,7 @@ subroutine Particle_destroy( this ) this%origin = 0.0_8 this%charge = 0 this%mass = 0.0_8 + this%eta = 0 this%omega = 0.0_8 this%spin = 0.0_8 this%totalCharge = 0 @@ -760,6 +792,7 @@ subroutine Particle_show( this ) write (6,"(T10,A16,I8)") "Owner : ",this%owner write (6,"(T10,A16,F8.2)") "Charge : ",this%charge write (6,"(T10,A16,F8.2)") "Mass : ",this%mass + write (6,"(T10,A16,I8)") "Eta : ",this%eta write (6,"(T10,A16,F8.2)") "Omega : ",this%omega write (6,"(T10,A16,F8.2)") "QDO center of : ",this%qdoCenterOf write (6,"(T10,A16,F8.2)") "Spin : ",this%spin @@ -828,6 +861,7 @@ subroutine Particle_saveToFile( this, unit ) write(unit,*) this%origin write(unit,*) this%charge write(unit,*) this%mass + write(unit,*) this%eta write(unit,*) this%omega write(unit,*) this%qdoCenterOf write(unit,*) this%spin @@ -888,6 +922,7 @@ subroutine Particle_loadFromFile( this, unit ) read(unit,*) this%origin read(unit,*) this%charge read(unit,*) this%mass + read(unit,*) this%eta read(unit,*) this%omega read(unit,*) this%qdoCenterOf read(unit,*) this%spin diff --git a/src/core/Species.f90 b/src/core/Species.f90 index 061113e8..e4d0de4c 100644 --- a/src/core/Species.f90 +++ b/src/core/Species.f90 @@ -93,9 +93,9 @@ subroutine Species_setSpecie( this, speciesID) if(trim(this%statistics) == "FERMION") then this%kappa = -1.0_8 - this%eta = 2.0_8 - this%lambda = 2.0_8 - this%particlesFraction = 0.5_8 + this%eta = 1.0_8 + this%lambda = 1.0_8 + this%particlesFraction = 1.0_8 else @@ -151,6 +151,13 @@ subroutine Species_setSpecie( this, speciesID) !! Adjust multiplicity this%multiplicity = this%multiplicity + 1 + !! Adjust eta + if(this%particles(1)%eta .ne. 0) then + this%eta = this%particles(1)%eta + this%lambda = this%particles(1)%eta + this%particlesFraction = 1.0_8/this%particles(1)%eta + end if + !! Adjust Occupation number this%ocupationNumber = this%ocupationNumber * this%particlesFraction diff --git a/src/core/Units.f90 b/src/core/Units.f90 index 3a51a81d..7f111dde 100644 --- a/src/core/Units.f90 +++ b/src/core/Units.f90 @@ -38,6 +38,7 @@ module Units_ real(8) , parameter :: DEBYES = 2.541764_8 real(8) , parameter :: ELECTRON_REST = 1.0_8 real(8) , parameter :: AMU = 1.0_8/1822.88850065855_8 + real(8) , parameter :: DALTON = 1822.88850065855_8 real(8) , parameter :: kg = 9.109382616D-31 real(8) , parameter :: DEGREES = 57.29577951_8 real(8) , parameter :: CM_NEG1 = 219476.0_8 diff --git a/src/ints/DirectIntegralManager.f90 b/src/ints/DirectIntegralManager.f90 index 4733500f..6d98d24c 100644 --- a/src/ints/DirectIntegralManager.f90 +++ b/src/ints/DirectIntegralManager.f90 @@ -32,7 +32,7 @@ module DirectIntegralManager_ use RysQuadrature_ use Matrix_ use Stopwatch_ - use ExternalPotential_ + use GTFPotential_ use String_ !# use RysQInts_ !! Please do not remove this line @@ -866,8 +866,8 @@ subroutine DirectIntegralManager_getExternalPotentialIntegrals(molSystem,species do i= 1, ExternalPotential_instance%ssize !if( trim(potential(i)%specie)==trim(interactNameSelected) ) then ! This does not work for UHF ! if ( String_findSubstring(trim( molSystem%species(speciesID)%name ), & - ! trim(String_getUpperCase(trim(ExternalPotential_instance%potentials(i)%specie)))) == 1 ) then - if ( trim( molSystem%species(speciesID)%symbol) == trim(String_getUpperCase(trim(ExternalPotential_instance%potentials(i)%specie))) ) then + ! trim(String_getUpperCase(trim(ExternalPotential_instance%potentials(i)%species)))) == 1 ) then + if ( trim( molSystem%species(speciesID)%symbol) == trim(String_getUpperCase(trim(ExternalPotential_instance%potentials(i)%species))) ) then potID=i exit end if diff --git a/src/ints/G12Integrals.f90 b/src/ints/G12Integrals.f90 index 19446252..1bfc7d4b 100644 --- a/src/ints/G12Integrals.f90 +++ b/src/ints/G12Integrals.f90 @@ -29,7 +29,7 @@ module G12Integrals_ use CONTROL_ use MolecularSystem_ use ContractedGaussian_ - use InterPotential_ + use GTFPotential_ implicit none #define contr(n,m) contractions(n)%contractions(m) @@ -253,7 +253,7 @@ subroutine G12Integrals_diskIntraSpecie(specieID) !Get potential ID do i=1, InterPotential_instance%ssize - if ( trim( MolecularSystem_instance%species(specieID)%symbol) == trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%specie))) .and. trim( MolecularSystem_instance%species(specieID)%symbol) == trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecie))) ) then + if ( trim( MolecularSystem_instance%species(specieID)%symbol) == trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%species))) .and. trim( MolecularSystem_instance%species(specieID)%symbol) == trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecies))) ) then potID=i exit end if @@ -918,14 +918,14 @@ subroutine G12Integrals_G12diskInterSpecie(nameOfSpecie, otherNameOfSpecie, spe !Get potential ID do i=1, InterPotential_instance%ssize if ( (trim(MolecularSystem_instance%species(specieID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%specie))) .and. & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%species))) .and. & trim(MolecularSystem_instance%species(otherSpecieID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecie)) ) & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecies)) ) & ) .or. & (trim( MolecularSystem_instance%species(otherSpecieID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%specie))) .and. & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%species))) .and. & trim( MolecularSystem_instance%species(specieID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecie)) ) & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecies)) ) & ) & ) then potID=i diff --git a/src/ints/IntegralManager.f90 b/src/ints/IntegralManager.f90 index 4c669083..ebba8eed 100644 --- a/src/ints/IntegralManager.f90 +++ b/src/ints/IntegralManager.f90 @@ -37,7 +37,7 @@ module IntegralManager_ use Matrix_ use CosmoCore_ use Stopwatch_ - use ExternalPotential_ + use GTFPotential_ use HarmonicIntegrals_ use FirstDerivativeIntegrals_ @@ -929,8 +929,8 @@ subroutine IntegralManager_writeThreeCenterIntegralsByProduct() do i= 1, ExternalPotential_instance%ssize !if( trim(potential(i)%specie)==trim(interactNameSelected) ) then ! This does not work for UHF ! if ( String_findSubstring(trim( MolecularSystem_instance%species(f)%name ), & - ! trim(String_getUpperCase(trim(ExternalPotential_instance%potentials(i)%specie)))) == 1 ) then - if ( trim( MolecularSystem_instance%species(f)%symbol) == trim(String_getUpperCase(trim(ExternalPotential_instance%potentials(i)%specie))) ) then + ! trim(String_getUpperCase(trim(ExternalPotential_instance%potentials(i)%species)))) == 1 ) then + if ( trim( MolecularSystem_instance%species(f)%symbol) == trim(String_getUpperCase(trim(ExternalPotential_instance%potentials(i)%species))) ) then potID=i exit end if diff --git a/src/ints/Libint2Interface.f90 b/src/ints/Libint2Interface.f90 index 23f5ffa1..78cc9d13 100644 --- a/src/ints/Libint2Interface.f90 +++ b/src/ints/Libint2Interface.f90 @@ -26,7 +26,7 @@ module Libint2Interface_ use, intrinsic :: iso_c_binding use MolecularSystem_ - use InterPotential_ + use GTFPotential_ use ContractedGaussian_ ! use Matrix_ @@ -824,9 +824,9 @@ subroutine Libint2Interface_computeG12Intraspecies_disk(speciesID) !Get potential ID do i=1, InterPotential_instance%ssize if ( trim( MolecularSystem_instance%species(speciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%specie))) .and. & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%species))) .and. & trim( MolecularSystem_instance%species(speciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecie))) ) then + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecies))) ) then potID=i exit end if @@ -883,14 +883,14 @@ subroutine Libint2Interface_computeG12Interspecies_disk(speciesID,otherSpeciesID !Get potential ID do i=1, InterPotential_instance%ssize if ( (trim(MolecularSystem_instance%species(speciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%specie))) .and. & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%species))) .and. & trim(MolecularSystem_instance%species(otherSpeciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecie)) ) & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecies)) ) & ) .or. & (trim( MolecularSystem_instance%species(otherSpeciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%specie))) .and. & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%species))) .and. & trim( MolecularSystem_instance%species(speciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecie)) ) & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecies)) ) & ) & ) then potID=i @@ -965,9 +965,9 @@ subroutine Libint2Interface_computeG12Intraspecies_direct(speciesID, density, tw !Get potential ID do i=1, InterPotential_instance%ssize if ( trim( molSys%species(speciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%specie))) .and. & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%species))) .and. & trim( molSys%species(speciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecie))) ) then + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecies))) ) then potID=i exit end if @@ -1027,14 +1027,14 @@ subroutine Libint2Interface_computeG12Interspecies_direct(speciesID,otherSpecies !Get potential ID do i=1, InterPotential_instance%ssize if ( (trim(molSys%species(speciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%specie))) .and. & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%species))) .and. & trim(molSys%species(otherSpeciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecie)) ) & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecies)) ) & ) .or. & (trim(molSys%species(otherSpeciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%specie))) .and. & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%species))) .and. & trim(molSys%species(speciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecie)) ) & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecies)) ) & ) & ) then potID=i @@ -1088,9 +1088,9 @@ subroutine Libint2Interface_compute2BodyIntraspecies_direct_all(speciesID, densi !Get potential ID do i=1, InterPotential_instance%ssize if ( trim( molSys%species(speciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%specie))) .and. & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%species))) .and. & trim( molSys%species(speciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecie))) ) then + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecies))) ) then potID=i exit end if @@ -1142,14 +1142,14 @@ subroutine Libint2Interface_compute2BodyInterspecies_direct_all(speciesID, other !Get potential ID do i=1, InterPotential_instance%ssize if ( (trim(molSys%species(speciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%specie))) .and. & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%species))) .and. & trim(molSys%species(otherSpeciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecie)) ) & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecies)) ) & ) .or. & (trim( molSys%species(otherSpeciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%specie))) .and. & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%species))) .and. & trim( molSys%species(speciesID)%symbol) == & - trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecie)) ) & + trim(String_getUpperCase(trim(InterPotential_instance%potentials(i)%otherSpecies)) ) & ) & ) then potID=i diff --git a/src/scf/DensityMatrixSCFGuess.f90 b/src/scf/DensityMatrixSCFGuess.f90 index a64a9a3a..f3a279a3 100644 --- a/src/scf/DensityMatrixSCFGuess.f90 +++ b/src/scf/DensityMatrixSCFGuess.f90 @@ -54,7 +54,7 @@ subroutine DensityMatrixSCFGuess_getGuess( speciesID, hcoreMatrix, transformatio type(MolecularSystem), pointer :: molSys type(Matrix) :: auxMatrix - character(30) :: nameOfSpecies + character(30) :: nameOfSpecies, symbolOfSpecies integer(8) :: orderOfMatrix, occupationNumber logical :: existPlain, existBinnary, readSuccess character(50) :: guessType @@ -71,6 +71,7 @@ subroutine DensityMatrixSCFGuess_getGuess( speciesID, hcoreMatrix, transformatio orderOfMatrix = MolecularSystem_getTotalnumberOfContractions(speciesID,molSys) nameOfSpecies = molSys%species(speciesID)%name + symbolOfSpecies = molSys%species(speciesID)%symbol occupationNumber = MolecularSystem_getOcupationNumber(speciesID,molSys) readSuccess=.false. @@ -79,11 +80,12 @@ subroutine DensityMatrixSCFGuess_getGuess( speciesID, hcoreMatrix, transformatio call Matrix_constructor(densityMatrix, int(orderOfMatrix,8), int(orderOfMatrix,8), 0.0_8 ) call Matrix_constructor(orbitals, int(orderOfMatrix,8), int(orderOfMatrix,8), 0.0_8 ) - + + readSuccess=.false. !!Verifica el archivo que contiene los coeficientes para una especie dada if ( CONTROL_instance%READ_FCHK ) then - call MolecularSystem_readFchk(trim(CONTROL_instance%INPUT_FILE)//trim(nameOfSpecies)//".fchk", orbitals, densityMatrix, nameOfSpecies ) - return + wfnFile=trim(CONTROL_instance%INPUT_FILE)//trim(symbolOfSpecies)//".fchk" + call MolecularSystem_readFchk(wfnFile, orbitals, densityMatrix, nameOfSpecies, readSuccess) else if ( CONTROL_instance%READ_COEFFICIENTS ) then wfnUnit = 30 @@ -107,10 +109,7 @@ subroutine DensityMatrixSCFGuess_getGuess( speciesID, hcoreMatrix, transformatio end if end if end if - - !check if the orbitals were read correctly - if(.not. allocated(orbitals%values) ) readSuccess=.false. - + if(readSuccess .and. printInfo ) print *, "Combination coefficients for ", trim(nameOfSpecies), " were read from ", trim(wfnFile) if(.not. readSuccess) then diff --git a/src/scf/OrbitalLocalizer.f90 b/src/scf/OrbitalLocalizer.f90 index 28e4ca88..006ab893 100644 --- a/src/scf/OrbitalLocalizer.f90 +++ b/src/scf/OrbitalLocalizer.f90 @@ -109,7 +109,7 @@ subroutine OrbitalLocalizer_erkaleLocal(speciesID,densityMatrix,fockMatrix,orbit type(Matrix) :: orbitalCoefficients type(Vector) :: orbitalEnergies - character(30) :: nameOfSpecies + character(30) :: nameOfSpecies, symbolOfSpecies integer :: statusSystem integer :: numberOfContractions @@ -118,11 +118,12 @@ subroutine OrbitalLocalizer_erkaleLocal(speciesID,densityMatrix,fockMatrix,orbit !! Convert lowdin fchk files to erkale chk files nameOfSpecies=MolecularSystem_getNameOfSpecies(speciesID) + symbolOfSpecies=MolecularSystem_getSymbolOfSpecies(speciesID) numberOfContractions = MolecularSystem_getTotalNumberOfContractions(speciesID) open(unit=30, file="erkale.read", status="replace", form="formatted") - write(30,*) "LoadFChk ", trim(CONTROL_instance%INPUT_FILE)//trim(nameOfSpecies)//".fchk" - write(30,*) "SaveChk ", trim(CONTROL_instance%INPUT_FILE)//trim(nameOfSpecies)//".chk" + write(30,*) "LoadFChk ", trim(CONTROL_instance%INPUT_FILE)//trim(symbolOfSpecies)//".fchk" + write(30,*) "SaveChk ", trim(CONTROL_instance%INPUT_FILE)//trim(symbolOfSpecies)//".chk" write(30,*) "Reorthonormalize true" close(30) @@ -132,8 +133,8 @@ subroutine OrbitalLocalizer_erkaleLocal(speciesID,densityMatrix,fockMatrix,orbit !! Localize orbitals open(unit=30, file="erkale.local", status="replace", form="formatted") - write(30,*) "LoadChk ", trim(CONTROL_instance%INPUT_FILE)//trim(nameOfSpecies)//".chk" - write(30,*) "SaveChk ", trim(CONTROL_instance%INPUT_FILE)//trim(nameOfSpecies)//".local.chk" + write(30,*) "LoadChk ", trim(CONTROL_instance%INPUT_FILE)//trim(symbolOfSpecies)//".chk" + write(30,*) "SaveChk ", trim(CONTROL_instance%INPUT_FILE)//trim(symbolOfSpecies)//".local.chk" write(30,*) "Method ", trim(CONTROL_instance%ERKALE_LOCALIZATION_METHOD) write(30,*) "Virtual false" write(30,*) "Maxiter 5000" @@ -149,15 +150,15 @@ subroutine OrbitalLocalizer_erkaleLocal(speciesID,densityMatrix,fockMatrix,orbit !!Convert erkale chk files to lowdin fchk files open(unit=30, file="erkale.write", status="replace", form="formatted") - write(30,*) "LoadChk ", trim(CONTROL_instance%INPUT_FILE)//trim(nameOfSpecies)//".local.chk" - write(30,*) "SaveFChk ", trim(CONTROL_instance%INPUT_FILE)//trim(nameOfSpecies)//".local.fchk" + write(30,*) "LoadChk ", trim(CONTROL_instance%INPUT_FILE)//trim(symbolOfSpecies)//".local.chk" + write(30,*) "SaveFChk ", trim(CONTROL_instance%INPUT_FILE)//trim(symbolOfSpecies)//".local.fchk" close(30) call system("erkale_fchkpt erkale.write") !! Read orbital coefficients from fchk files - call MolecularSystem_readFchk(trim(CONTROL_instance%INPUT_FILE)//trim(nameOfSpecies)//".local.fchk", orbitalCoefficients, densityMatrix, nameOfSpecies ) + call MolecularSystem_readFchk(trim(CONTROL_instance%INPUT_FILE)//trim(symbolOfSpecies)//".local.fchk", orbitalCoefficients, densityMatrix, nameOfSpecies ) orbitalEnergies%values=0.0 !! Molecular orbital fock operator expected value diff --git a/src/scf/SingleSCF.f90 b/src/scf/SingleSCF.f90 index 6cb85fd1..16c219fb 100644 --- a/src/scf/SingleSCF.f90 +++ b/src/scf/SingleSCF.f90 @@ -508,7 +508,7 @@ subroutine SingleSCF_readFrozenOrbitals(wfObject) if (CONTROL_instance%READ_FCHK) then call Matrix_constructor (auxiliaryMatrix, numberOfContractions, numberOfContractions) - call MolecularSystem_readFchk( trim(CONTROL_instance%INPUT_FILE)//trim(wfObject%name)//".fchk", & + call MolecularSystem_readFchk( trim(CONTROL_instance%INPUT_FILE)//trim(wfObject%molSys%species(wfObject%species)%symbol)//".fchk", & wfObject%waveFunctionCoefficients, auxiliaryMatrix, wfObject%name ) call Matrix_destructor(auxiliaryMatrix) diff --git a/test/CHDTHemu.massTest.lowdin b/test/CHDTHemu.massTest.lowdin new file mode 100644 index 00000000..85916399 --- /dev/null +++ b/test/CHDTHemu.massTest.lowdin @@ -0,0 +1,23 @@ + +GEOMETRY +e-(C) NAKAI-CC-PVTZ 0.0000000 0.0000000 0.0000000 +e-(H) NAKAI-CC-PVTZ 0.6283310 0.6283310 0.6283310 +e-(H) NAKAI-CC-PVTZ -0.6283310 -0.6283310 0.6283310 +e-(H) NAKAI-CC-PVTZ -0.6283310 0.6283310 -0.6283310 +e-(H) NAKAI-CC-PVTZ 0.6283310 -0.6283310 -0.6283310 +U- HEMU 0.6283310 0.6283310 0.6283310 +C_12 NAKAI-3-SP 0.0000000 0.0000000 0.0000000 +He_4 NAKAI-3-SP 0.6283310 0.6283310 0.6283310 +H_1 DZSPNB -0.6283310 -0.6283310 0.6283310 +H_2 DZSPNB -0.6283310 0.6283310 -0.6283310 +H_3 DZSPNB 0.6283310 -0.6283310 -0.6283310 +END GEOMETRY + +TASKS + method = "RHF" +END TASKS + +CONTROL + readCoefficients=F + removeTranslationalContamination=T +END CONTROL diff --git a/test/CHDTHemu.massTest.py b/test/CHDTHemu.massTest.py new file mode 100644 index 00000000..e1d5f5a4 --- /dev/null +++ b/test/CHDTHemu.massTest.py @@ -0,0 +1,84 @@ +#!/usr/bin/env python +from __future__ import print_function +import os +import sys +from colorstring import * + +if len(sys.argv)==2: + lowdinbin = sys.argv[1] +else: + lowdinbin = "lowdin2" + +testName = sys.argv[0][:-3] +inputName = testName + ".lowdin" +outputName = testName + ".out" +# Reference values and tolerance + +refValues = { + "Total mass" : [40383.2869, 5E-4], + "E- Kinetic energy" : [36.558549208752,1E-6], + "MUON Kinetic energy" : [4.823223401186,1E-6], + "C_12 Kinetic energy" : [0.020362665915,1E-6], + "HE_4 Kinetic energy" : [0.042798597334,1E-6], + "H_1 Kinetic energy" : [0.018810346458,1E-6], + "H_2 Kinetic energy" : [0.013654378944,1E-6], + "H_3 Kinetic energy" : [0.011128664979,1E-6], + "HF energy" : [-74.168958869716,1E-8] +} + +testValues = dict(refValues) #copy +for value in testValues: #reset + testValues[value] = 0 #reset + +# Run calculation + +status = os.system(lowdinbin + " -i " + inputName) + +if status: + print(testName + str_red(" ... NOT OK")) + sys.exit(1) + +output = open(outputName, "r") +outputRead = output.readlines() + +# Values +checkArray=[0,0,0,0] +for i in range(0,len(outputRead)): + line = outputRead[i] + if "TOTAL ENERGY =" in line: + testValues["HF energy"] = float(line.split()[3]) + if "MASS (m_e)" in line: + testValues["Total mass"] = float(line.split()[3]) + if "E- Kinetic energy" in line: + testValues["E- Kinetic energy"] = float(line.split()[4]) + if "MUON Kinetic energy" in line: + testValues["MUON Kinetic energy"] = float(line.split()[4]) + if "C_12 Kinetic energy" in line: + testValues["C_12 Kinetic energy"] = float(line.split()[4]) + if "HE_4 Kinetic energy" in line: + testValues["HE_4 Kinetic energy"] = float(line.split()[4]) + if "H_1 Kinetic energy" in line: + testValues["H_1 Kinetic energy"] = float(line.split()[4]) + if "H_2 Kinetic energy" in line: + testValues["H_2 Kinetic energy"] = float(line.split()[4]) + if "H_3 Kinetic energy" in line: + testValues["H_3 Kinetic energy"] = float(line.split()[4]) + +output.close() + +passTest = True + +for value in refValues: + diffValue = abs(refValues[value][0] - testValues[value]) + if ( diffValue <= refValues[value][1] ): + passTest = passTest * True + else : + passTest = passTest * False + print("%s %.8f %.8f %.2e" % ( value, refValues[value][0], testValues[value], diffValue)) + +if passTest : + print(testName + str_green(" ... OK")) +else: + print(testName + str_red(" ... NOT OK")) + sys.exit(1) + diff --git a/test/He2-C60potential-NOCI.lowdin b/test/He2-C60potential-NOCI.lowdin index c06d7379..4883911f 100644 --- a/test/He2-C60potential-NOCI.lowdin +++ b/test/He2-C60potential-NOCI.lowdin @@ -2,8 +2,8 @@ SYSTEM_DESCRIPTION='H' GEOMETRY N0 dirac 0.0 0.0 0.0 rotationPoint=1 -HEA3 HE2-1S 1.818207 0.0 -0.347193 rotateAround=1 -HEB3 HE2-1S -1.818207 0.0 0.347193 rotateAround=1 +HEA3 HE2-1S 1.818207 0.0 -0.347193 rotateAround=1 m=5494.8926 +HEB3 HE2-1S -1.818207 0.0 0.347193 rotateAround=1 m=5494.8926 END GEOMETRY TASKS diff --git a/test/He3-C60potential-NOCI.lowdin b/test/He3-C60potential-NOCI.lowdin index 16870c5d..a5efe877 100644 --- a/test/He3-C60potential-NOCI.lowdin +++ b/test/He3-C60potential-NOCI.lowdin @@ -2,9 +2,9 @@ SYSTEM_DESCRIPTION='H' GEOMETRY N0 dirac 0.0 0.0 0.0 rotationPoint=1 -HEA3 HE2-1S 2.159 0.0 0.0 rotateAround=1 -HEA3 HE2-1S 0.0 2.159 0.0 rotateAround=1 -HEB3 HE2-1S 0.0 0.0 2.159 rotateAround=1 +HEA3 HE2-1S 2.159 0.0 0.0 rotateAround=1 m=5494.8926 +HEA3 HE2-1S 0.0 2.159 0.0 rotateAround=1 m=5494.8926 +HEB3 HE2-1S 0.0 0.0 2.159 rotateAround=1 m=5494.8926 END GEOMETRY TASKS diff --git a/test/HemuH-CUSTOM_BASIS.lowdin b/test/HemuH-CUSTOM_BASIS.lowdin index 19b21087..27c3692a 100644 --- a/test/HemuH-CUSTOM_BASIS.lowdin +++ b/test/HemuH-CUSTOM_BASIS.lowdin @@ -2,9 +2,9 @@ GEOMETRY e-[H] cc-pvtz 0.0000 0.0000 0.00000 e-[H] CUSTOM_1 0.0000 0.0000 0.74144 - H_1 CUSTOM_2 0.0000 0.0000 0.00000 - U- CUSTOM_3 0.0000 0.0000 0.74144 - He_4 CUSTOM_3 0.0000 0.0000 0.74144 + H_1 CUSTOM_2 0.0000 0.0000 0.00000 m=1836.1527 + U- CUSTOM_3 0.0000 0.0000 0.74144 m=206.7683 + He_4 CUSTOM_3 0.0000 0.0000 0.74144 m=7349.6727 END GEOMETRY TASKS diff --git a/test/HemuH-CUSTOM_BASIS.py b/test/HemuH-CUSTOM_BASIS.py index 7d461f17..44c69e98 100644 --- a/test/HemuH-CUSTOM_BASIS.py +++ b/test/HemuH-CUSTOM_BASIS.py @@ -15,11 +15,11 @@ # Reference values and tolerance refValues = { -"HF energy" : [-343.383191892820,1E-8], -"U-HOMO" : [-371.890049816287,1E-1], -"H_1-HOMO" : [-1.019360160964,1E-4], -"He_4-HOMO" : [-652.366876763392,1E-1], -"e-HOMO" : [-0.585414450602,1E-4], +"HF energy" : [-343.383218426575,1E-8], +"U-HOMO" : [-371.889981903188,1E-1], +"H_1-HOMO" : [-1.019346186407,1E-4], +"He_4-HOMO" : [-652.365841581870,1E-1], +"e-HOMO" : [-0.585408097570,1E-4], } testValues = dict(refValues) #copy diff --git a/test/Ps2F2.RHF-customEta.lowdin b/test/Ps2F2.RHF-customEta.lowdin new file mode 100644 index 00000000..2484865d --- /dev/null +++ b/test/Ps2F2.RHF-customEta.lowdin @@ -0,0 +1,17 @@ +GEOMETRY + e-(F) aug-cc-pvdz 0.00 0.00 -1.33 addParticles=1 + e-(F) aug-cc-pvdz 0.00 0.00 1.33 addParticles=1 + e+ PSX-DZ 0.00 0.00 -1.33 eta=2 + e+ PSX-DZ 0.00 0.00 1.33 + F dirac 0.00 0.00 -1.33 + F dirac 0.00 0.00 1.33 +END GEOMETRY + +TASKS + method = "RHF" +END TASKS + +CONTROL + readCoefficients=F +END CONTROL + diff --git a/test/Ps2F2.RHF-customEta.py b/test/Ps2F2.RHF-customEta.py new file mode 100644 index 00000000..a812ee62 --- /dev/null +++ b/test/Ps2F2.RHF-customEta.py @@ -0,0 +1,68 @@ +#!/usr/bin/env python +from __future__ import print_function +import os +import sys +from colorstring import * + +if len(sys.argv)==2: + lowdinbin = sys.argv[1] +else: + lowdinbin = "lowdin2" + +testName = sys.argv[0][:-3] +inputName = testName + ".lowdin" +outputName = testName + ".out" + +# Reference values and tolerance +refValues = { + "HF energy" : [-199.213740198536,1E-8], + "eta e+" : [2.0,1E-8], + "occupation e+" : [1.0,1E-8] +} + +testValues = dict(refValues) #copy +for value in testValues: #reset + testValues[value] = 0 #reset + +# Run calculation + +status = os.system(lowdinbin + " -i " + inputName) + +if status: + print(testName + str_red(" ... NOT OK")) + sys.exit(1) + +output = open(outputName, "r") +outputRead = output.readlines() + +# Values +flagC=0 +for i in range(0,len(outputRead)): + line = outputRead[i] + if "TOTAL ENERGY =" in line: + testValues["HF energy"] = float(line.split()[3]) + if "CONSTANTS OF COUPLING" in line: + flagC=1 + if "E+" in line and flagC==1: + testValues["eta e+"] = float(line.split()[2]) + testValues["occupation e+"] = float(line.split()[4]) + flagC=0 +output.close() + +passTest = True + +for value in refValues: + diffValue = abs(refValues[value][0] - testValues[value]) + if ( diffValue <= refValues[value][1] ): + passTest = passTest * True + else : + passTest = passTest * False + print("%s %.8f %.8f %.2e" % ( value, refValues[value][0], testValues[value], diffValue)) + +if passTest : + print(testName + str_green(" ... OK")) +else: + print(testName + str_red(" ... NOT OK")) + sys.exit(1) + +output.close() diff --git a/test/TIP4P-dimer-singlet-NOCI.lowdin b/test/TIP4P-dimer-singlet-NOCI.lowdin index d886158d..213f2093 100644 --- a/test/TIP4P-dimer-singlet-NOCI.lowdin +++ b/test/TIP4P-dimer-singlet-NOCI.lowdin @@ -2,12 +2,12 @@ SYSTEM_DESCRIPTION='H' GEOMETRY N0 dirac 0 0 0 -X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 translationCenter=1 -Y0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 translationCenter=2 +X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 translationCenter=1 q=0.5564 m=1836.15267247 +Y0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 translationCenter=2 q=0.5564 m=1836.15267247 N0 dirac 0.0 0.0 6.0 -X1.1- dirac 0.0 0.0 6.292152 -X0.5+ dirac 0.0 1.430429 7.107157 -X0.5+ dirac 0.0 -1.430429 7.107157 +X1.1- dirac 0.0 0.0 6.292152 q=-1.1128 +X0.5+ dirac 0.0 1.430429 7.107157 q=0.5564 +X0.5+ dirac 0.0 -1.430429 7.107157 q=0.5564 END GEOMETRY TASKS @@ -36,3 +36,452 @@ INTERPOTENTIAL X0.5+ Y0.5+ VHH-CCSDT Y0.5+ Y0.5+ VHH-CCSDT END INTERPOTENTIAL + +EXTERPOTENTIAL + X0.5+ VOH-CCSDT + Y0.5+ VOH-CCSDT +END EXTERPOTENTIAL + +INTERPOTENTIAL + X0.5+ X0.5+ VHH-CCSDT + X0.5+ Y0.5+ VHH-CCSDT + Y0.5+ Y0.5+ VHH-CCSDT +END INTERPOTENTIAL + +BASIS H2O-1S1P1D +#optimized exponents +O-H-TIP X0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 + +O-H-TIP Y0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 +END BASIS + +POTENTIAL VOH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant RHH +O-X0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 + +O-Y0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 +END POTENTIAL + +POTENTIAL VHH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant ROH +O-X0.5+X0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 + +O-X0.5+Y0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 + +O-Y0.5+Y0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 +END POTENTIAL diff --git a/test/TIP4P-dimer-singlet-direct.lowdin b/test/TIP4P-dimer-singlet-direct.lowdin index 9e5e2eaa..3788cc16 100644 --- a/test/TIP4P-dimer-singlet-direct.lowdin +++ b/test/TIP4P-dimer-singlet-direct.lowdin @@ -1,11 +1,11 @@ GEOMETRY N0 dirac 0 0 0 -X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 -Y0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 +X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 q=0.5564 m=1836.15267247 +Y0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 q=0.5564 m=1836.15267247 N0 dirac 0.0 0.0 6.0 -X1.1- dirac 0.0 0.0 6.292152 -X0.5+ dirac 0.0 1.430429 7.107157 -X0.5+ dirac 0.0 -1.430429 7.107157 +X1.1- dirac 0.0 0.0 6.292152 q=-1.1128 +X0.5+ dirac 0.0 1.430429 7.107157 q=0.5564 +X0.5+ dirac 0.0 -1.430429 7.107157 q=0.5564 END GEOMETRY TASKS @@ -31,3 +31,440 @@ INTERPOTENTIAL Y0.5+ Y0.5+ VHH-CCSDT END INTERPOTENTIAL +BASIS H2O-1S1P1D +#optimized exponents +O-H-TIP X0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 + +O-H-TIP Y0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 +END BASIS + +POTENTIAL VOH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant RHH +O-X0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 + +O-Y0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 +END POTENTIAL + +POTENTIAL VHH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant ROH +O-X0.5+X0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 + +O-X0.5+Y0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 + +O-Y0.5+Y0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 +END POTENTIAL diff --git a/test/TIP4P-dimer-singlet-direct.py b/test/TIP4P-dimer-singlet-direct.py index f3b406ea..bc41eaa3 100644 --- a/test/TIP4P-dimer-singlet-direct.py +++ b/test/TIP4P-dimer-singlet-direct.py @@ -17,11 +17,11 @@ refValues = { "HF energy" : [-0.651377267238,1E-8], -"HA-TIP Ext Pot" : [0.005339817047,1E-4], -"HB-TIP Ext Pot" : [0.005265789260,1E-4], -"HA-TIP/HB-TIP Hartree" : [0.003393498016,1E-4], -"HA-TIP/Fixed interact." : [-0.025239485153,1E-4], -"HB-TIP/Fixed interact." : [-0.010156190638,1E-4] +"X0.5+ Ext Pot" : [0.005339817047,1E-4], +"Y0.5+ Ext Pot" : [0.005265789260,1E-4], +"X0.5+/Y0.5+ Hartree" : [0.003393498016,1E-4], +"X0.5+/Fixed interact." : [-0.025239485153,1E-4], +"Y0.5+/Fixed interact." : [-0.010156190638,1E-4] } testValues = dict(refValues) #copy @@ -44,16 +44,16 @@ line = outputRead[i] if "TOTAL ENERGY =" in line: testValues["HF energy"] = float(line.split()[3]) - if "HA-TIP Ext Pot" in line: - testValues["HA-TIP Ext Pot"] = float(line.split()[5]) - if "HB-TIP Ext Pot" in line: - testValues["HB-TIP Ext Pot"] = float(line.split()[5]) - if "HA-TIP/HB-TIP Hartree" in line: - testValues["HA-TIP/HB-TIP Hartree"] = float(line.split()[4]) - if "HA-TIP/Fixed interact." in line: - testValues["HA-TIP/Fixed interact."] = float(line.split()[4]) - if "HB-TIP/Fixed interact." in line: - testValues["HB-TIP/Fixed interact."] = float(line.split()[4]) + if "X0.5+ Ext Pot" in line: + testValues["X0.5+ Ext Pot"] = float(line.split()[5]) + if "Y0.5+ Ext Pot" in line: + testValues["Y0.5+ Ext Pot"] = float(line.split()[5]) + if "X0.5+/Y0.5+ Hartree" in line: + testValues["X0.5+/Y0.5+ Hartree"] = float(line.split()[4]) + if "X0.5+/Fixed interact." in line: + testValues["X0.5+/Fixed interact."] = float(line.split()[4]) + if "Y0.5+/Fixed interact." in line: + testValues["Y0.5+/Fixed interact."] = float(line.split()[4]) passTest = True diff --git a/test/TIP4P-dimer-singlet-memory.lowdin b/test/TIP4P-dimer-singlet-memory.lowdin index bc275bf3..1dc863f5 100644 --- a/test/TIP4P-dimer-singlet-memory.lowdin +++ b/test/TIP4P-dimer-singlet-memory.lowdin @@ -1,11 +1,11 @@ GEOMETRY N0 dirac 0 0 0 -X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 -Y0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 +X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 q=0.5564 m=1836.15267247 +Y0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 q=0.5564 m=1836.15267247 N0 dirac 0.0 0.0 6.0 -X1.1- dirac 0.0 0.0 6.292152 -X0.5+ dirac 0.0 1.430429 7.107157 -X0.5+ dirac 0.0 -1.430429 7.107157 +X1.1- dirac 0.0 0.0 6.292152 q=-1.1128 +X0.5+ dirac 0.0 1.430429 7.107157 q=0.5564 +X0.5+ dirac 0.0 -1.430429 7.107157 q=0.5564 END GEOMETRY TASKS @@ -31,3 +31,440 @@ INTERPOTENTIAL Y0.5+ Y0.5+ VHH-CCSDT END INTERPOTENTIAL +BASIS H2O-1S1P1D +#optimized exponents +O-H-TIP X0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 + +O-H-TIP Y0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 +END BASIS + +POTENTIAL VOH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant RHH +O-X0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 + +O-Y0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 +END POTENTIAL + +POTENTIAL VHH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant ROH +O-X0.5+X0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 + +O-X0.5+Y0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 + +O-Y0.5+Y0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 +END POTENTIAL diff --git a/test/TIP4P-dimer-singlet-memory.py b/test/TIP4P-dimer-singlet-memory.py index f3b406ea..bc41eaa3 100644 --- a/test/TIP4P-dimer-singlet-memory.py +++ b/test/TIP4P-dimer-singlet-memory.py @@ -17,11 +17,11 @@ refValues = { "HF energy" : [-0.651377267238,1E-8], -"HA-TIP Ext Pot" : [0.005339817047,1E-4], -"HB-TIP Ext Pot" : [0.005265789260,1E-4], -"HA-TIP/HB-TIP Hartree" : [0.003393498016,1E-4], -"HA-TIP/Fixed interact." : [-0.025239485153,1E-4], -"HB-TIP/Fixed interact." : [-0.010156190638,1E-4] +"X0.5+ Ext Pot" : [0.005339817047,1E-4], +"Y0.5+ Ext Pot" : [0.005265789260,1E-4], +"X0.5+/Y0.5+ Hartree" : [0.003393498016,1E-4], +"X0.5+/Fixed interact." : [-0.025239485153,1E-4], +"Y0.5+/Fixed interact." : [-0.010156190638,1E-4] } testValues = dict(refValues) #copy @@ -44,16 +44,16 @@ line = outputRead[i] if "TOTAL ENERGY =" in line: testValues["HF energy"] = float(line.split()[3]) - if "HA-TIP Ext Pot" in line: - testValues["HA-TIP Ext Pot"] = float(line.split()[5]) - if "HB-TIP Ext Pot" in line: - testValues["HB-TIP Ext Pot"] = float(line.split()[5]) - if "HA-TIP/HB-TIP Hartree" in line: - testValues["HA-TIP/HB-TIP Hartree"] = float(line.split()[4]) - if "HA-TIP/Fixed interact." in line: - testValues["HA-TIP/Fixed interact."] = float(line.split()[4]) - if "HB-TIP/Fixed interact." in line: - testValues["HB-TIP/Fixed interact."] = float(line.split()[4]) + if "X0.5+ Ext Pot" in line: + testValues["X0.5+ Ext Pot"] = float(line.split()[5]) + if "Y0.5+ Ext Pot" in line: + testValues["Y0.5+ Ext Pot"] = float(line.split()[5]) + if "X0.5+/Y0.5+ Hartree" in line: + testValues["X0.5+/Y0.5+ Hartree"] = float(line.split()[4]) + if "X0.5+/Fixed interact." in line: + testValues["X0.5+/Fixed interact."] = float(line.split()[4]) + if "Y0.5+/Fixed interact." in line: + testValues["Y0.5+/Fixed interact."] = float(line.split()[4]) passTest = True diff --git a/test/TIP4P-dimer-singlet.lowdin b/test/TIP4P-dimer-singlet.lowdin index 9d540025..7ffc1cfe 100644 --- a/test/TIP4P-dimer-singlet.lowdin +++ b/test/TIP4P-dimer-singlet.lowdin @@ -1,11 +1,11 @@ GEOMETRY N0 dirac 0 0 0 -X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 -Y0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 +X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 q=0.5564 m=1836.15267247 +Y0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 q=0.5564 m=1836.15267247 N0 dirac 0.0 0.0 6.0 -X1.1- dirac 0.0 0.0 6.292152 -X0.5+ dirac 0.0 1.430429 7.107157 -X0.5+ dirac 0.0 -1.430429 7.107157 +X1.1- dirac 0.0 0.0 6.292152 q=-1.1128 +X0.5+ dirac 0.0 1.430429 7.107157 q=0.5564 +X0.5+ dirac 0.0 -1.430429 7.107157 q=0.5564 END GEOMETRY TASKS @@ -31,3 +31,440 @@ INTERPOTENTIAL Y0.5+ Y0.5+ VHH-CCSDT END INTERPOTENTIAL +BASIS H2O-1S1P1D +#optimized exponents +O-H-TIP X0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 + +O-H-TIP Y0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 +END BASIS + +POTENTIAL VOH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant RHH +O-X0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 + +O-Y0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 +END POTENTIAL + +POTENTIAL VHH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant ROH +O-X0.5+X0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 + +O-X0.5+Y0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 + +O-Y0.5+Y0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 +END POTENTIAL diff --git a/test/TIP4P-dimer-singlet.py b/test/TIP4P-dimer-singlet.py index f3b406ea..bc41eaa3 100644 --- a/test/TIP4P-dimer-singlet.py +++ b/test/TIP4P-dimer-singlet.py @@ -17,11 +17,11 @@ refValues = { "HF energy" : [-0.651377267238,1E-8], -"HA-TIP Ext Pot" : [0.005339817047,1E-4], -"HB-TIP Ext Pot" : [0.005265789260,1E-4], -"HA-TIP/HB-TIP Hartree" : [0.003393498016,1E-4], -"HA-TIP/Fixed interact." : [-0.025239485153,1E-4], -"HB-TIP/Fixed interact." : [-0.010156190638,1E-4] +"X0.5+ Ext Pot" : [0.005339817047,1E-4], +"Y0.5+ Ext Pot" : [0.005265789260,1E-4], +"X0.5+/Y0.5+ Hartree" : [0.003393498016,1E-4], +"X0.5+/Fixed interact." : [-0.025239485153,1E-4], +"Y0.5+/Fixed interact." : [-0.010156190638,1E-4] } testValues = dict(refValues) #copy @@ -44,16 +44,16 @@ line = outputRead[i] if "TOTAL ENERGY =" in line: testValues["HF energy"] = float(line.split()[3]) - if "HA-TIP Ext Pot" in line: - testValues["HA-TIP Ext Pot"] = float(line.split()[5]) - if "HB-TIP Ext Pot" in line: - testValues["HB-TIP Ext Pot"] = float(line.split()[5]) - if "HA-TIP/HB-TIP Hartree" in line: - testValues["HA-TIP/HB-TIP Hartree"] = float(line.split()[4]) - if "HA-TIP/Fixed interact." in line: - testValues["HA-TIP/Fixed interact."] = float(line.split()[4]) - if "HB-TIP/Fixed interact." in line: - testValues["HB-TIP/Fixed interact."] = float(line.split()[4]) + if "X0.5+ Ext Pot" in line: + testValues["X0.5+ Ext Pot"] = float(line.split()[5]) + if "Y0.5+ Ext Pot" in line: + testValues["Y0.5+ Ext Pot"] = float(line.split()[5]) + if "X0.5+/Y0.5+ Hartree" in line: + testValues["X0.5+/Y0.5+ Hartree"] = float(line.split()[4]) + if "X0.5+/Fixed interact." in line: + testValues["X0.5+/Fixed interact."] = float(line.split()[4]) + if "Y0.5+/Fixed interact." in line: + testValues["Y0.5+/Fixed interact."] = float(line.split()[4]) passTest = True diff --git a/test/TIP4P-dimer-triplet-NOCI.lowdin b/test/TIP4P-dimer-triplet-NOCI.lowdin index b8411434..8a88242a 100644 --- a/test/TIP4P-dimer-triplet-NOCI.lowdin +++ b/test/TIP4P-dimer-triplet-NOCI.lowdin @@ -2,12 +2,12 @@ SYSTEM_DESCRIPTION='H' GEOMETRY N0 dirac 0 0 0 -X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 translationCenter=1 -X0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 translationCenter=2 +X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 translationCenter=1 q=0.5564 m=1836.15267247 +X0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 translationCenter=2 q=0.5564 m=1836.15267247 N0 dirac 0.0 0.0 6.0 -X1.1- dirac 0.0 0.0 6.292152 -X0.5+ dirac 0.0 1.430429 7.107157 -X0.5+ dirac 0.0 -1.430429 7.107157 +X1.1- dirac 0.0 0.0 6.292152 q=-1.1128 +X0.5+ dirac 0.0 1.430429 7.107157 q=0.5564 +X0.5+ dirac 0.0 -1.430429 7.107157 q=0.5564 END GEOMETRY TASKS @@ -33,3 +33,441 @@ END EXTERPOTENTIAL INTERPOTENTIAL X0.5+ X0.5+ VHH-CCSDT END INTERPOTENTIAL + +BASIS H2O-1S1P1D +#optimized exponents +O-H-TIP X0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 + +O-H-TIP Y0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 +END BASIS + +POTENTIAL VOH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant RHH +O-X0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 + +O-Y0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 +END POTENTIAL + +POTENTIAL VHH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant ROH +O-X0.5+X0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 + +O-X0.5+Y0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 + +O-Y0.5+Y0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 +END POTENTIAL diff --git a/test/TIP4P-dimer-triplet-direct.lowdin b/test/TIP4P-dimer-triplet-direct.lowdin index 15ab0abb..2db120e6 100644 --- a/test/TIP4P-dimer-triplet-direct.lowdin +++ b/test/TIP4P-dimer-triplet-direct.lowdin @@ -1,11 +1,11 @@ GEOMETRY N0 dirac 0 0 0 -X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 -X0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 +X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 q=0.5564 m=1836.15267247 +X0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 q=0.5564 m=1836.15267247 N0 dirac 0.0 0.0 6.0 -X1.1- dirac 0.0 0.0 6.292152 -X0.5+ dirac 0.0 1.430429 7.107157 -X0.5+ dirac 0.0 -1.430429 7.107157 +X1.1- dirac 0.0 0.0 6.292152 q=-1.1128 +X0.5+ dirac 0.0 1.430429 7.107157 q=0.5564 +X0.5+ dirac 0.0 -1.430429 7.107157 q=0.5564 END GEOMETRY TASKS @@ -28,3 +28,186 @@ INTERPOTENTIAL X0.5+ X0.5+ VHH-CCSDT END INTERPOTENTIAL + +BASIS H2O-1S1P1D +#optimized exponents +O-H-TIP X0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 +END BASIS + +POTENTIAL VOH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant RHH +O-X0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 +END POTENTIAL + +POTENTIAL VHH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant ROH +O-X0.5+X0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 +END POTENTIAL diff --git a/test/TIP4P-dimer-triplet-direct.py b/test/TIP4P-dimer-triplet-direct.py index 41f2ebc4..d13b500b 100644 --- a/test/TIP4P-dimer-triplet-direct.py +++ b/test/TIP4P-dimer-triplet-direct.py @@ -17,10 +17,10 @@ refValues = { "HF energy" : [-0.651377261423,1E-8], - "HA-TIP Ext Pot" : [0.010606635479,1E-4], - "HA-TIP/HA-TIP Hartree" : [1.217813835967,1E-4], - "HA-TIP Exchange" : [-1.214419735313,1E-4], - "HA-TIP/Fixed interact." : [-0.035396418562,1E-4] + "X0.5+ Ext Pot" : [0.010606635479,1E-4], + "X0.5+/X0.5+ Hartree" : [1.217813835967,1E-4], + "X0.5+ Exchange" : [-1.214419735313,1E-4], + "X0.5+/Fixed interact." : [-0.035396418562,1E-4] } testValues = dict(refValues) #copy @@ -43,14 +43,14 @@ line = outputRead[i] if "TOTAL ENERGY =" in line: testValues["HF energy"] = float(line.split()[3]) - if "HA-TIP Ext Pot" in line: - testValues["HA-TIP Ext Pot"] = float(line.split()[5]) - if "HA-TIP/HA-TIP Hartree" in line: - testValues["HA-TIP/HA-TIP Hartree"] = float(line.split()[4]) - if "HA-TIP Exchange" in line: - testValues["HA-TIP Exchange"] = float(line.split()[4]) - if "HA-TIP/Fixed interact." in line: - testValues["HA-TIP/Fixed interact."] = float(line.split()[4]) + if "X0.5+ Ext Pot" in line: + testValues["X0.5+ Ext Pot"] = float(line.split()[5]) + if "X0.5+/X0.5+ Hartree" in line: + testValues["X0.5+/X0.5+ Hartree"] = float(line.split()[4]) + if "X0.5+ Exchange" in line: + testValues["X0.5+ Exchange"] = float(line.split()[4]) + if "X0.5+/Fixed interact." in line: + testValues["X0.5+/Fixed interact."] = float(line.split()[4]) passTest = True diff --git a/test/TIP4P-dimer-triplet-memory.lowdin b/test/TIP4P-dimer-triplet-memory.lowdin index c37e1e48..5d36490e 100644 --- a/test/TIP4P-dimer-triplet-memory.lowdin +++ b/test/TIP4P-dimer-triplet-memory.lowdin @@ -1,11 +1,11 @@ GEOMETRY N0 dirac 0 0 0 -X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 -X0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 +X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 q=0.5564 m=1836.15267247 +X0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 q=0.5564 m=1836.15267247 N0 dirac 0.0 0.0 6.0 -X1.1- dirac 0.0 0.0 6.292152 -X0.5+ dirac 0.0 1.430429 7.107157 -X0.5+ dirac 0.0 -1.430429 7.107157 +X1.1- dirac 0.0 0.0 6.292152 q=-1.1128 +X0.5+ dirac 0.0 1.430429 7.107157 q=0.5564 +X0.5+ dirac 0.0 -1.430429 7.107157 q=0.5564 END GEOMETRY TASKS @@ -28,3 +28,194 @@ INTERPOTENTIAL X0.5+ X0.5+ VHH-CCSDT END INTERPOTENTIAL +EXTERPOTENTIAL + X0.5+ VOH-CCSDT +END EXTERPOTENTIAL + +INTERPOTENTIAL + X0.5+ X0.5+ VHH-CCSDT +END INTERPOTENTIAL + +BASIS H2O-1S1P1D +#optimized exponents +O-H-TIP X0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 +END BASIS + +POTENTIAL VOH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant RHH +O-X0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 + +END POTENTIAL + +POTENTIAL VHH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant ROH +O-X0.5+X0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 +END POTENTIAL diff --git a/test/TIP4P-dimer-triplet-memory.py b/test/TIP4P-dimer-triplet-memory.py index 41f2ebc4..d13b500b 100644 --- a/test/TIP4P-dimer-triplet-memory.py +++ b/test/TIP4P-dimer-triplet-memory.py @@ -17,10 +17,10 @@ refValues = { "HF energy" : [-0.651377261423,1E-8], - "HA-TIP Ext Pot" : [0.010606635479,1E-4], - "HA-TIP/HA-TIP Hartree" : [1.217813835967,1E-4], - "HA-TIP Exchange" : [-1.214419735313,1E-4], - "HA-TIP/Fixed interact." : [-0.035396418562,1E-4] + "X0.5+ Ext Pot" : [0.010606635479,1E-4], + "X0.5+/X0.5+ Hartree" : [1.217813835967,1E-4], + "X0.5+ Exchange" : [-1.214419735313,1E-4], + "X0.5+/Fixed interact." : [-0.035396418562,1E-4] } testValues = dict(refValues) #copy @@ -43,14 +43,14 @@ line = outputRead[i] if "TOTAL ENERGY =" in line: testValues["HF energy"] = float(line.split()[3]) - if "HA-TIP Ext Pot" in line: - testValues["HA-TIP Ext Pot"] = float(line.split()[5]) - if "HA-TIP/HA-TIP Hartree" in line: - testValues["HA-TIP/HA-TIP Hartree"] = float(line.split()[4]) - if "HA-TIP Exchange" in line: - testValues["HA-TIP Exchange"] = float(line.split()[4]) - if "HA-TIP/Fixed interact." in line: - testValues["HA-TIP/Fixed interact."] = float(line.split()[4]) + if "X0.5+ Ext Pot" in line: + testValues["X0.5+ Ext Pot"] = float(line.split()[5]) + if "X0.5+/X0.5+ Hartree" in line: + testValues["X0.5+/X0.5+ Hartree"] = float(line.split()[4]) + if "X0.5+ Exchange" in line: + testValues["X0.5+ Exchange"] = float(line.split()[4]) + if "X0.5+/Fixed interact." in line: + testValues["X0.5+/Fixed interact."] = float(line.split()[4]) passTest = True diff --git a/test/TIP4P-dimer-triplet.lowdin b/test/TIP4P-dimer-triplet.lowdin index 75996651..18fd85c8 100644 --- a/test/TIP4P-dimer-triplet.lowdin +++ b/test/TIP4P-dimer-triplet.lowdin @@ -1,11 +1,11 @@ GEOMETRY N0 dirac 0 0 0 -X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 -X0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 +X0.5+ H2O-1S1P1D -0.151399 0.000276 1.807200 q=0.5564 m=1836.15267247 +X0.5+ H2O-1S1P1D -1.715032 0.000843 -0.589570 q=0.5564 m=1836.15267247 N0 dirac 0.0 0.0 6.0 -X1.1- dirac 0.0 0.0 6.292152 -X0.5+ dirac 0.0 1.430429 7.107157 -X0.5+ dirac 0.0 -1.430429 7.107157 +X1.1- dirac 0.0 0.0 6.292152 q=-1.1128 +X0.5+ dirac 0.0 1.430429 7.107157 q=0.5564 +X0.5+ dirac 0.0 -1.430429 7.107157 q=0.5564 END GEOMETRY TASKS @@ -28,3 +28,186 @@ INTERPOTENTIAL X0.5+ X0.5+ VHH-CCSDT END INTERPOTENTIAL + +BASIS H2O-1S1P1D +#optimized exponents +O-H-TIP X0.5+ (1S) BASIS TYPE: 2 +3 +1 0 1 +14.509888498676842 1.0 +2 1 1 +6.885507269761004 1.0 +3 2 1 +9.023681376783887 1.0 +END BASIS + +POTENTIAL VOH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant RHH +O-X0.5+ +27 +1 0 +13.892336977 9.839539381 +0.0 0.0 0.0 +2 0 +10.000000000 15.916259785 +0.0 0.0 0.0 +3 0 +7.498942093 -18.199008237 +0.0 0.0 0.0 +4 0 +5.623413252 6.723316202 +0.0 0.0 0.0 +5 0 +4.216965034 4.671467422 +0.0 0.0 0.0 +6 0 +3.162277660 1.703771884 +0.0 0.0 0.0 +7 0 +2.371373706 0.297676262 +0.0 0.0 0.0 +8 0 +1.778279410 0.967600283 +0.0 0.0 0.0 +9 0 +1.333521432 1.034932150 +0.0 0.0 0.0 +10 0 +1.000000000 0.657621305 +0.0 0.0 0.0 +11 0 +0.749894209 -0.327381267 +0.0 0.0 0.0 +12 0 +0.562341325 0.221528399 +0.0 0.0 0.0 +13 0 +0.421696503 -0.255694044 +0.0 0.0 0.0 +14 0 +0.316227766 0.097853405 +0.0 0.0 0.0 +15 0 +0.237137371 0.744489008 +0.0 0.0 0.0 +16 0 +0.177827941 -0.750090202 +0.0 0.0 0.0 +17 0 +0.133352143 -0.309090719 +0.0 0.0 0.0 +18 0 +0.100000000 -0.998800172 +0.0 0.0 0.0 +19 0 +0.074989421 0.718795407 +0.0 0.0 0.0 +20 0 +0.056234133 0.228002288 +0.0 0.0 0.0 +21 0 +0.042169650 0.573744431 +0.0 0.0 0.0 +22 0 +0.031622777 -0.056969092 +0.0 0.0 0.0 +23 0 +0.023713737 -0.392241803 +0.0 0.0 0.0 +24 0 +0.017782794 -0.569609681 +0.0 0.0 0.0 +25 0 +0.013335214 0.759323484 +0.0 0.0 0.0 +26 0 +0.010000000 -0.194525832 +0.0 0.0 0.0 +27 0 +0.000000000 0.138912412 +0.0 0.0 0.0 +END POTENTIAL + +POTENTIAL VHH-CCSDT +#Fitted from CCSD(T)/def2-TZVPPD with constant ROH +O-X0.5+X0.5+ +26 +1 0 +10.000000000 1.313644267 +0.0 0.0 0.0 +2 0 +7.498942093 3.312208762 +0.0 0.0 0.0 +3 0 +5.623413252 -4.290599715 +0.0 0.0 0.0 +4 0 +4.216965034 -2.400490962 +0.0 0.0 0.0 +5 0 +3.162277660 6.160643793 +0.0 0.0 0.0 +6 0 +2.371373706 2.406172066 +0.0 0.0 0.0 +7 0 +1.778279410 -6.852617537 +0.0 0.0 0.0 +8 0 +1.333521432 -1.043562213 +0.0 0.0 0.0 +9 0 +1.000000000 7.556790649 +0.0 0.0 0.0 +10 0 +0.749894209 0.777037995 +0.0 0.0 0.0 +11 0 +0.562341325 -8.976316536 +0.0 0.0 0.0 +12 0 +0.421696503 0.575296040 +0.0 0.0 0.0 +13 0 +0.316227766 8.672468936 +0.0 0.0 0.0 +14 0 +0.237137371 -0.855388714 +0.0 0.0 0.0 +15 0 +0.177827941 -6.271012179 +0.0 0.0 0.0 +16 0 +0.133352143 1.152328020 +0.0 0.0 0.0 +17 0 +0.100000000 3.403807143 +0.0 0.0 0.0 +18 0 +0.074989421 -1.745300340 +0.0 0.0 0.0 +19 0 +0.056234133 -1.349326021 +0.0 0.0 0.0 +20 0 +0.042169650 1.605561812 +0.0 0.0 0.0 +21 0 +0.031622777 -0.267553445 +0.0 0.0 0.0 +22 0 +0.023713737 0.034566428 +0.0 0.0 0.0 +23 0 +0.017782794 -0.185356609 +0.0 0.0 0.0 +24 0 +0.013335214 -0.103128899 +0.0 0.0 0.0 +25 0 +0.010000000 0.128087663 +0.0 0.0 0.0 +26 0 +0.000000000 0.085213897 +0.0 0.0 0.0 +END POTENTIAL diff --git a/test/TIP4P-dimer-triplet.py b/test/TIP4P-dimer-triplet.py index 41f2ebc4..d13b500b 100644 --- a/test/TIP4P-dimer-triplet.py +++ b/test/TIP4P-dimer-triplet.py @@ -17,10 +17,10 @@ refValues = { "HF energy" : [-0.651377261423,1E-8], - "HA-TIP Ext Pot" : [0.010606635479,1E-4], - "HA-TIP/HA-TIP Hartree" : [1.217813835967,1E-4], - "HA-TIP Exchange" : [-1.214419735313,1E-4], - "HA-TIP/Fixed interact." : [-0.035396418562,1E-4] + "X0.5+ Ext Pot" : [0.010606635479,1E-4], + "X0.5+/X0.5+ Hartree" : [1.217813835967,1E-4], + "X0.5+ Exchange" : [-1.214419735313,1E-4], + "X0.5+/Fixed interact." : [-0.035396418562,1E-4] } testValues = dict(refValues) #copy @@ -43,14 +43,14 @@ line = outputRead[i] if "TOTAL ENERGY =" in line: testValues["HF energy"] = float(line.split()[3]) - if "HA-TIP Ext Pot" in line: - testValues["HA-TIP Ext Pot"] = float(line.split()[5]) - if "HA-TIP/HA-TIP Hartree" in line: - testValues["HA-TIP/HA-TIP Hartree"] = float(line.split()[4]) - if "HA-TIP Exchange" in line: - testValues["HA-TIP Exchange"] = float(line.split()[4]) - if "HA-TIP/Fixed interact." in line: - testValues["HA-TIP/Fixed interact."] = float(line.split()[4]) + if "X0.5+ Ext Pot" in line: + testValues["X0.5+ Ext Pot"] = float(line.split()[5]) + if "X0.5+/X0.5+ Hartree" in line: + testValues["X0.5+/X0.5+ Hartree"] = float(line.split()[4]) + if "X0.5+ Exchange" in line: + testValues["X0.5+ Exchange"] = float(line.split()[4]) + if "X0.5+/Fixed interact." in line: + testValues["X0.5+/Fixed interact."] = float(line.split()[4]) passTest = True