*******> update.14

 Author: Tyler Luchko

 Date: March 21, 2019

 Programs: rism1d

 Description: fixes rism1d test cases for changes in MDL files
    * added ions94 directory for old ion MDL files (some of which are used for tests)
    * updated test files to use the new MDL files for cSPC/E and cTIP3P
    * updated test files to find the old ions94 MDLs in their new location
 --------------------------------------------------------------------------------
 AmberTools/test/rism1d/spc-hnc/Run.spc-lj          |   2 +-
 AmberTools/test/rism1d/spc-hnc/spc-lj.therm.save   |   2 +-
 AmberTools/test/rism1d/spc-kh/Run.spc              |   2 +-
 AmberTools/test/rism1d/spc-kh/spc.therm.save       |   2 +-
 AmberTools/test/rism1d/spc-mv0/Run.spc             |   2 +-
 AmberTools/test/rism1d/spc-mv0/spc.therm.save      |   2 +-
 AmberTools/test/rism1d/spc-polyt/Run.spc           |   2 +-
 AmberTools/test/rism1d/spc-polyt/spc.therm.save    |   2 +-
 AmberTools/test/rism1d/spc-psen/Run.spc-na-kh      |   4 +-
 AmberTools/test/rism1d/spc-psen/Run.spc-nacl-3     |   6 +-
 .../test/rism1d/spc-psen/spc-nacl-3.therm.save     |   2 +-
 AmberTools/test/rism1d/tip3p-kh/Run.tip3p          |   2 +-
 AmberTools/test/rism1d/tip3p-kh/tip3p.therm.save   |   2 +-
 dat/rism1d/mdl/ions94/CIO.mdl                      |  31 ++
 dat/rism1d/mdl/ions94/Cl-.mdl                      |  31 ++
 dat/rism1d/mdl/ions94/Cs+.mdl                      |  31 ++
 dat/rism1d/mdl/ions94/IB.mdl                       |  31 ++
 dat/rism1d/mdl/ions94/K+.mdl                       |  31 ++
 dat/rism1d/mdl/ions94/Li+.mdl                      |  31 ++
 dat/rism1d/mdl/ions94/MG2.mdl                      |  31 ++
 dat/rism1d/mdl/ions94/Na+.mdl                      |  31 ++
 dat/rism1d/mdl/ions94/Rb+.mdl                      |  31 ++
 dat/rism1d/model/prepdat.pl                        | 472 ---------------------
 23 files changed, 295 insertions(+), 488 deletions(-)

diff --git AmberTools/test/rism1d/spc-hnc/Run.spc-lj AmberTools/test/rism1d/spc-hnc/Run.spc-lj
index ca81fd559d..7752c1a0fa 100755
--- AmberTools/test/rism1d/spc-hnc/Run.spc-lj
+++ AmberTools/test/rism1d/spc-hnc/Run.spc-lj
@@ -21,7 +21,7 @@ NSP=2
 !DENSITY=55.5,
 !corresponds very closely to 0.0333 1/A3
 DENSITY=55.296d0,
-MODEL="../../../../dat/rism1d/model/SPC.mdl"
+MODEL="../../../../dat/rism1d/mdl/cSPCE.mdl"
 /
 &SPECIES
 DENSITY=0
diff --git AmberTools/test/rism1d/spc-hnc/spc-lj.therm.save AmberTools/test/rism1d/spc-hnc/spc-lj.therm.save
index 856228b1af..c464f82e66 100644
--- AmberTools/test/rism1d/spc-hnc/spc-lj.therm.save
+++ AmberTools/test/rism1d/spc-hnc/spc-lj.therm.save
@@ -5,7 +5,7 @@ Pressure_(Free_Energy)            Pfe         [MPa]           4.94946637E+002
 Excess_free_energy                FE          [kcal/mol/A^3] -2.29996704E-001
 
 #Species properties               
-#Description                      Variable    Units           SPC              LJ              
+#Description                      Variable    Units           cSPC              LJ              
 Excess_chemical_potential_PR      EXCHEMsp_PR [kcal/mol]     -3.37896237E+000  5.93901592E+000 
 Partial_molar_volume              PMV         [A^-3]          3.00300236E+001  4.13870618E+001 
 
diff --git AmberTools/test/rism1d/spc-kh/Run.spc AmberTools/test/rism1d/spc-kh/Run.spc
index 9acaeae4dd..27e43a54c6 100755
--- AmberTools/test/rism1d/spc-kh/Run.spc
+++ AmberTools/test/rism1d/spc-kh/Run.spc
@@ -21,7 +21,7 @@ NSP=1
 &SPECIES
 !corresponds very closely to 0.0333 1/A3
 DENSITY=55.296d0,
-MODEL="../../../../dat/rism1d/model/SPC.mdl"
+MODEL="../../../../dat/rism1d/mdl/cSPCE.mdl"
 /
 EOF
 
diff --git AmberTools/test/rism1d/spc-kh/spc.therm.save AmberTools/test/rism1d/spc-kh/spc.therm.save
index a0e566b3cd..91bbbe2a63 100644
--- AmberTools/test/rism1d/spc-kh/spc.therm.save
+++ AmberTools/test/rism1d/spc-kh/spc.therm.save
@@ -6,7 +6,7 @@ Pressure_(Free_Energy)            Pfe         [MPa]           5.85603370E+002
 Excess_free_energy                FE          [kcal/mol/A^3] -2.24797798E-001
 
 #Species properties               
-#Description                      Variable    Units           SPC             
+#Description                      Variable    Units           cSPC             
 Excess_chemical_potential_PR      EXCHEMsp_PR [kcal/mol]     -2.79190338E+000 
 Excess_chemical_potential_SC      EXCHEMsp_SC [kcal/mol]     -6.30982097E+000 
 Excess_chemical_potential_SM      EXCHEMsp_SM [kcal/mol]     -4.81172182E+000 
diff --git AmberTools/test/rism1d/spc-mv0/Run.spc AmberTools/test/rism1d/spc-mv0/Run.spc
index 8cac2e96d9..c062f6eb4d 100755
--- AmberTools/test/rism1d/spc-mv0/Run.spc
+++ AmberTools/test/rism1d/spc-mv0/Run.spc
@@ -20,7 +20,7 @@ NSP=1
 !DENSITY=55.5,
 !corresponds very closely to 0.0333 1/A3
 DENSITY=55.296d0,
-MODEL="../../../../dat/rism1d/model/SPC.mdl"
+MODEL="../../../../dat/rism1d/mdl/cSPCE.mdl"
 /
 EOF
 
diff --git AmberTools/test/rism1d/spc-mv0/spc.therm.save AmberTools/test/rism1d/spc-mv0/spc.therm.save
index 35d628d0fc..1266c40526 100644
--- AmberTools/test/rism1d/spc-mv0/spc.therm.save
+++ AmberTools/test/rism1d/spc-mv0/spc.therm.save
@@ -3,7 +3,7 @@
 Compressibility                   xi          [1/MPa]         3.13478444E-004
 
 #Species properties               
-#Description                      Variable    Units           SPC             
+#Description                      Variable    Units           cSPC             
 Partial_molar_volume              PMV         [A^-3]          3.00300236E+001 
 
 #Site properties                  
diff --git AmberTools/test/rism1d/spc-polyt/Run.spc AmberTools/test/rism1d/spc-polyt/Run.spc
index 3372b7f226..e9adc022d0 100755
--- AmberTools/test/rism1d/spc-polyt/Run.spc
+++ AmberTools/test/rism1d/spc-polyt/Run.spc
@@ -20,7 +20,7 @@ NSP=1
 &SPECIES
 !corresponds very closely to 0.0333 1/A3
 DENSITY=55.296d0,
-MODEL="../../../../dat/rism1d/model/SPC.mdl"
+MODEL="../../../../dat/rism1d/mdl/cSPCE.mdl"
 /
 EOF
 
diff --git AmberTools/test/rism1d/spc-polyt/spc.therm.save AmberTools/test/rism1d/spc-polyt/spc.therm.save
index c4a3425497..04ea2cc247 100644
--- AmberTools/test/rism1d/spc-polyt/spc.therm.save
+++ AmberTools/test/rism1d/spc-polyt/spc.therm.save
@@ -6,7 +6,7 @@ Pressure_(Virial)               Pvir      [kPa]        2.31065822E+003
 Excess_free_energy              FE        [kcal/mol]  -1.03963153E+003
 
 #Species properties             
-#Description                    Variable  Units        SPC             
+#Description                    Variable  Units        cSPC             
 Excess_chemical_potential       EXCHEMsp  [kcal/mol]   2.74658196E-310 
 Partial_molar_volume            PMV       [A^-3]       3.00300236E+001 
 
diff --git AmberTools/test/rism1d/spc-psen/Run.spc-na-kh AmberTools/test/rism1d/spc-psen/Run.spc-na-kh
index 99a99dd28a..f92b676257 100755
--- AmberTools/test/rism1d/spc-psen/Run.spc-na-kh
+++ AmberTools/test/rism1d/spc-psen/Run.spc-na-kh
@@ -20,11 +20,11 @@ NSP=2
 !DENSITY=55.5,
 !corresponds very closely to 0.0333 1/A3
 DENSITY=55.296d0,
-MODEL="../../../../dat/rism1d/model/SPC.mdl"
+MODEL="../../../../dat/rism1d/mdl/cSPCE.mdl"
 /
 &SPECIES
 DENSITY=2.35d0,
-MODEL="../../../../dat/rism1d/model/Na+.mdl"
+MODEL="../../../../dat/rism1d/mdl/ions94/Na+.mdl"
 /
 EOF
 
diff --git AmberTools/test/rism1d/spc-psen/Run.spc-nacl-3 AmberTools/test/rism1d/spc-psen/Run.spc-nacl-3
index e0ea0e42b4..6c85d13186 100755
--- AmberTools/test/rism1d/spc-psen/Run.spc-nacl-3
+++ AmberTools/test/rism1d/spc-psen/Run.spc-nacl-3
@@ -20,15 +20,15 @@ NSP=3
 !DENSITY=55.5,
 !corresponds very closely to 0.0333 1/A3
 DENSITY=55.296d0,
-MODEL="../../../../dat/rism1d/model/SPC.mdl"
+MODEL="../../../../dat/rism1d/mdl/cSPCE.mdl"
 /
 &SPECIES
 DENSITY=.1d0,
-MODEL="../../../../dat/rism1d/model/Na+.mdl"
+MODEL="../../../../dat/rism1d/mdl/ions94/Na+.mdl"
 /
 &SPECIES
 DENSITY=.1d0,
-MODEL="../../../../dat/rism1d/model/Cl-.mdl"
+MODEL="../../../../dat/rism1d/mdl/ions94/Cl-.mdl"
 /
 EOF
 
diff --git AmberTools/test/rism1d/spc-psen/spc-nacl-3.therm.save AmberTools/test/rism1d/spc-psen/spc-nacl-3.therm.save
index 8b690f36b6..dc29b5d8ab 100644
--- AmberTools/test/rism1d/spc-psen/spc-nacl-3.therm.save
+++ AmberTools/test/rism1d/spc-psen/spc-nacl-3.therm.save
@@ -6,7 +6,7 @@ Pressure_(Free_Energy)            Pfe         [MPa]           5.29103932E+002
 Excess_free_energy                FE          [kcal/mol/A^3] -2.37705468E-001
 
 #Species properties               
-#Description                      Variable    Units           SPC              Na+              Cl-             
+#Description                      Variable    Units           cSPC              Na+              Cl-             
 Excess_chemical_potential_PR      EXCHEMsp_PR [kcal/mol]     -2.94840611E+000 -8.02432203E+001 -8.06716307E+001 
 Excess_chemical_potential_SC      EXCHEMsp_SC [kcal/mol]     -6.84667208E+000 -7.89150847E+001 -7.93435001E+001 
 Excess_chemical_potential_SM      EXCHEMsp_SM [kcal/mol]     -5.19507070E+000 -7.55619280E+001 -7.62537016E+001 
diff --git AmberTools/test/rism1d/tip3p-kh/Run.tip3p AmberTools/test/rism1d/tip3p-kh/Run.tip3p
index 3b59136383..206ada50fe 100755
--- AmberTools/test/rism1d/tip3p-kh/Run.tip3p
+++ AmberTools/test/rism1d/tip3p-kh/Run.tip3p
@@ -18,7 +18,7 @@ NSP=1
 /
 &SPECIES
 DENSITY=55.296,
-MODEL="../../../../dat/rism1d/model/TP3.mdl"
+MODEL="../../../../dat/rism1d/mdl/cTIP3P.mdl"
 /
 EOF
 
diff --git AmberTools/test/rism1d/tip3p-kh/tip3p.therm.save AmberTools/test/rism1d/tip3p-kh/tip3p.therm.save
index b8f9ef8ff5..2968863954 100644
--- AmberTools/test/rism1d/tip3p-kh/tip3p.therm.save
+++ AmberTools/test/rism1d/tip3p-kh/tip3p.therm.save
@@ -6,7 +6,7 @@ Pressure_(Free_Energy)            Pfe         [MPa]           5.75349142E+002
 Excess_free_energy                FE          [kcal/mol/A^3] -2.06489638E-001
 
 #Species properties               
-#Description                      Variable    Units           TP3             
+#Description                      Variable    Units           cTIP             
 Excess_chemical_potential_PR      EXCHEMsp_PR [kcal/mol]     -2.14794896E+000 
 Solvation_energy                  ESOLVsp     [kcal/mol]     -1.06935504E+001 
 -Temperature*solvation_entropy_PR -TSsp_PR    [kcal/mol]      8.54560148E+000 
diff --git dat/rism1d/mdl/ions94/CIO.mdl dat/rism1d/mdl/ions94/CIO.mdl
new file mode 100644
index 0000000000..59f2091b8c
--- /dev/null
+++ dat/rism1d/mdl/ions94/CIO.mdl
@@ -0,0 +1,31 @@
+%VERSION  VERSION_STAMP = V0001.000  DATE = 05/12/09  17:44:15
+%FLAG TITLE
+%FORMAT(20a4)
+ CIO
+%FLAG POINTERS
+%FORMAT(10I8)
+       1       1
+%FLAG ATMTYP
+%FORMAT(10I8)
+       1
+%FLAG ATMNAME
+%FORMAT(20a4)
+Na+ 
+%FLAG MASS
+%FORMAT(5e16.8)
+  2.29900000e+01
+%FLAG CHG
+%FORMAT(5e16.8)
+  1.82223000e+01
+%FLAG LJEPSILON
+%FORMAT(5e16.8)
+  2.77000001e-03
+%FLAG LJSIGMA
+%FORMAT(5e16.8)
+  1.86800000e+00
+%FLAG MULTI
+%FORMAT(10I8)
+       1
+%FLAG COORD
+%FORMAT(5e16.8)
+  0.00000000e+00  0.00000000e+00  0.00000000e+00
diff --git dat/rism1d/mdl/ions94/Cl-.mdl dat/rism1d/mdl/ions94/Cl-.mdl
new file mode 100644
index 0000000000..c60530a3b8
--- /dev/null
+++ dat/rism1d/mdl/ions94/Cl-.mdl
@@ -0,0 +1,31 @@
+%VERSION  VERSION_STAMP = V0001.000  DATE = 05/12/09  17:44:16
+%FLAG TITLE
+%FORMAT(20a4)
+ Cl-
+%FLAG POINTERS
+%FORMAT(10I8)
+       1       1
+%FLAG ATMTYP
+%FORMAT(10I8)
+       1
+%FLAG ATMNAME
+%FORMAT(20a4)
+Cl- 
+%FLAG MASS
+%FORMAT(5e16.8)
+  3.54500000e+01
+%FLAG CHG
+%FORMAT(5e16.8)
+ -1.82223000e+01
+%FLAG LJEPSILON
+%FORMAT(5e16.8)
+  1.00000000e-01
+%FLAG LJSIGMA
+%FORMAT(5e16.8)
+  2.47000000e+00
+%FLAG MULTI
+%FORMAT(10I8)
+       1
+%FLAG COORD
+%FORMAT(5e16.8)
+  0.00000000e+00  0.00000000e+00  0.00000000e+00
diff --git dat/rism1d/mdl/ions94/Cs+.mdl dat/rism1d/mdl/ions94/Cs+.mdl
new file mode 100644
index 0000000000..745ec9555b
--- /dev/null
+++ dat/rism1d/mdl/ions94/Cs+.mdl
@@ -0,0 +1,31 @@
+%VERSION  VERSION_STAMP = V0001.000  DATE = 05/12/09  17:44:17
+%FLAG TITLE
+%FORMAT(20a4)
+ Cs+
+%FLAG POINTERS
+%FORMAT(10I8)
+       1       1
+%FLAG ATMTYP
+%FORMAT(10I8)
+       1
+%FLAG ATMNAME
+%FORMAT(20a4)
+Cs+ 
+%FLAG MASS
+%FORMAT(5e16.8)
+  1.32910000e+02
+%FLAG CHG
+%FORMAT(5e16.8)
+  1.82223000e+01
+%FLAG LJEPSILON
+%FORMAT(5e16.8)
+  8.05999997e-05
+%FLAG LJSIGMA
+%FORMAT(5e16.8)
+  3.39500000e+00
+%FLAG MULTI
+%FORMAT(10I8)
+       1
+%FLAG COORD
+%FORMAT(5e16.8)
+  0.00000000e+00  0.00000000e+00  0.00000000e+00
diff --git dat/rism1d/mdl/ions94/IB.mdl dat/rism1d/mdl/ions94/IB.mdl
new file mode 100644
index 0000000000..398d657480
--- /dev/null
+++ dat/rism1d/mdl/ions94/IB.mdl
@@ -0,0 +1,31 @@
+%VERSION  VERSION_STAMP = V0001.000  DATE = 05/12/09  17:44:17
+%FLAG TITLE
+%FORMAT(20a4)
+  IB
+%FLAG POINTERS
+%FORMAT(10I8)
+       1       1
+%FLAG ATMTYP
+%FORMAT(10I8)
+       1
+%FLAG ATMNAME
+%FORMAT(20a4)
+IB  
+%FLAG MASS
+%FORMAT(5e16.8)
+  1.31000000e+02
+%FLAG CHG
+%FORMAT(5e16.8)
+  1.82223000e+01
+%FLAG LJEPSILON
+%FORMAT(5e16.8)
+  1.00000000e-01
+%FLAG LJSIGMA
+%FORMAT(5e16.8)
+  5.00000000e+00
+%FLAG MULTI
+%FORMAT(10I8)
+       1
+%FLAG COORD
+%FORMAT(5e16.8)
+  0.00000000e+00  0.00000000e+00  0.00000000e+00
diff --git dat/rism1d/mdl/ions94/K+.mdl dat/rism1d/mdl/ions94/K+.mdl
new file mode 100644
index 0000000000..75f0aec0e1
--- /dev/null
+++ dat/rism1d/mdl/ions94/K+.mdl
@@ -0,0 +1,31 @@
+%VERSION  VERSION_STAMP = V0001.000  DATE = 06/12/09  12:22:01
+%FLAG TITLE
+%FORMAT(20a4)
+  K+
+%FLAG POINTERS
+%FORMAT(10I8)
+       1       1
+%FLAG ATMTYP
+%FORMAT(10I8)
+       1
+%FLAG ATMNAME
+%FORMAT(20a4)
+K+  
+%FLAG MASS
+%FORMAT(5e16.8)
+  3.91000000e+01
+%FLAG CHG
+%FORMAT(5e16.8)
+  1.82223000e+01
+%FLAG LJEPSILON
+%FORMAT(5e16.8)
+  3.28000003e-04
+%FLAG LJSIGMA
+%FORMAT(5e16.8)
+  2.65800000e+00
+%FLAG MULTI
+%FORMAT(10I8)
+       1
+%FLAG COORD
+%FORMAT(5e16.8)
+  0.00000000e+00  0.00000000e+00  0.00000000e+00
diff --git dat/rism1d/mdl/ions94/Li+.mdl dat/rism1d/mdl/ions94/Li+.mdl
new file mode 100644
index 0000000000..7672c6bde5
--- /dev/null
+++ dat/rism1d/mdl/ions94/Li+.mdl
@@ -0,0 +1,31 @@
+%VERSION  VERSION_STAMP = V0001.000  DATE = 05/12/09  17:44:16
+%FLAG TITLE
+%FORMAT(20a4)
+ Li+
+%FLAG POINTERS
+%FORMAT(10I8)
+       1       1
+%FLAG ATMTYP
+%FORMAT(10I8)
+       1
+%FLAG ATMNAME
+%FORMAT(20a4)
+Li+ 
+%FLAG MASS
+%FORMAT(5e16.8)
+  6.94000000e+00
+%FLAG CHG
+%FORMAT(5e16.8)
+  1.82223000e+01
+%FLAG LJEPSILON
+%FORMAT(5e16.8)
+  1.83000000e-02
+%FLAG LJSIGMA
+%FORMAT(5e16.8)
+  1.13700000e+00
+%FLAG MULTI
+%FORMAT(10I8)
+       1
+%FLAG COORD
+%FORMAT(5e16.8)
+  0.00000000e+00  0.00000000e+00  0.00000000e+00
diff --git dat/rism1d/mdl/ions94/MG2.mdl dat/rism1d/mdl/ions94/MG2.mdl
new file mode 100644
index 0000000000..a0d18a712a
--- /dev/null
+++ dat/rism1d/mdl/ions94/MG2.mdl
@@ -0,0 +1,31 @@
+%VERSION  VERSION_STAMP = V0001.000  DATE = 06/12/09  12:22:00
+%FLAG TITLE
+%FORMAT(20a4)
+ MG2
+%FLAG POINTERS
+%FORMAT(10I8)
+       1       1
+%FLAG ATMTYP
+%FORMAT(10I8)
+       1
+%FLAG ATMNAME
+%FORMAT(20a4)
+MG  
+%FLAG MASS
+%FORMAT(5e16.8)
+  2.43050000e+01
+%FLAG CHG
+%FORMAT(5e16.8)
+  3.64446000e+01
+%FLAG LJEPSILON
+%FORMAT(5e16.8)
+  8.94700001e-01
+%FLAG LJSIGMA
+%FORMAT(5e16.8)
+  7.92600000e-01
+%FLAG MULTI
+%FORMAT(10I8)
+       1
+%FLAG COORD
+%FORMAT(5e16.8)
+  0.00000000e+00  0.00000000e+00  0.00000000e+00
diff --git dat/rism1d/mdl/ions94/Na+.mdl dat/rism1d/mdl/ions94/Na+.mdl
new file mode 100644
index 0000000000..36066df670
--- /dev/null
+++ dat/rism1d/mdl/ions94/Na+.mdl
@@ -0,0 +1,31 @@
+%VERSION  VERSION_STAMP = V0001.000  DATE = 05/12/09  17:44:16
+%FLAG TITLE
+%FORMAT(20a4)
+ Na+
+%FLAG POINTERS
+%FORMAT(10I8)
+       1       1
+%FLAG ATMTYP
+%FORMAT(10I8)
+       1
+%FLAG ATMNAME
+%FORMAT(20a4)
+Na+ 
+%FLAG MASS
+%FORMAT(5e16.8)
+  2.29900000e+01
+%FLAG CHG
+%FORMAT(5e16.8)
+  1.82223000e+01
+%FLAG LJEPSILON
+%FORMAT(5e16.8)
+  2.77000001e-03
+%FLAG LJSIGMA
+%FORMAT(5e16.8)
+  1.86800000e+00
+%FLAG MULTI
+%FORMAT(10I8)
+       1
+%FLAG COORD
+%FORMAT(5e16.8)
+  0.00000000e+00  0.00000000e+00  0.00000000e+00
diff --git dat/rism1d/mdl/ions94/Rb+.mdl dat/rism1d/mdl/ions94/Rb+.mdl
new file mode 100644
index 0000000000..f83bf2fe72
--- /dev/null
+++ dat/rism1d/mdl/ions94/Rb+.mdl
@@ -0,0 +1,31 @@
+%VERSION  VERSION_STAMP = V0001.000  DATE = 05/12/09  17:44:17
+%FLAG TITLE
+%FORMAT(20a4)
+ Rb+
+%FLAG POINTERS
+%FORMAT(10I8)
+       1       1
+%FLAG ATMTYP
+%FORMAT(10I8)
+       1
+%FLAG ATMNAME
+%FORMAT(20a4)
+Rb+ 
+%FLAG MASS
+%FORMAT(5e16.8)
+  8.54700000e+01
+%FLAG CHG
+%FORMAT(5e16.8)
+  1.82223000e+01
+%FLAG LJEPSILON
+%FORMAT(5e16.8)
+  1.70000000e-04
+%FLAG LJSIGMA
+%FORMAT(5e16.8)
+  2.95600000e+00
+%FLAG MULTI
+%FORMAT(10I8)
+       1
+%FLAG COORD
+%FORMAT(5e16.8)
+  0.00000000e+00  0.00000000e+00  0.00000000e+00
diff --git dat/rism1d/model/prepdat.pl dat/rism1d/model/prepdat.pl
deleted file mode 100644
index 0b18f9750c..0000000000
--- dat/rism1d/model/prepdat.pl
+++ /dev/null
@@ -1,472 +0,0 @@
-#!/usr/bin/perl -w
-use POSIX;
-use strict;
-
-#  Note: this script illustrates one way to create mdl files for RISM1D.
-#  It was written by Tyler Luchko back in 2009, before the Amber
-#  community updated the way ions are handled (with parameters specific
-#  to each water model.)  Hence, the actual files created here were
-#  deleted in 2018, since they do not contain any currently recommended
-#  parameters.  More up to date versions are in the ../mdl directory.
-#  We are keeping this script here to help users to wish to create their
-#  own mdl files.
-
-# This script creates parameter files to be read by RISM1D for all
-# predefined solvent and ion parameters in the Amber force field (FF99SB is
-# used).  This is done in two parts for each molecule.  First, a PRMTOP and
-# CRD file are created in TLEAP. These two resulting files are parsed and
-# the relavent data placed into a single file following the PRMTOP format
-# and can be read by NXTSEC.
-
-# detect atoms with the same sigma, epsilon and charge and assign the
-# correct multiplicity value.  This makes for a small number of solvent
-# atoms and makes 3D-RISM more efficient
-my $multiplicity=1;
-
-#use the Amber PARM7 format. This is really the only supported format in RISM1D.
-#The other format is namelist but this is difficult to read for variable 
-#   sized arrays
-my $ambParm=1;
-
-#solvent models
-my %solvent=(
-	  "TP3"=>["WAT","TIP3P"],
-	  "SPC"=>["WAT","SPCE"],
-	  "SPF"=>["WAT","SPCFw"],
-	  "PL3"=>["WAT","POL3"],
-	  "TP4"=>["WAT","TIP4P"],
-	  "T4E"=>["WAT","TIP4PEw"],
-	  "DC4"=>["WAT","DC4"],
-	  "TP5"=>["WAT","TIP5P"],
-	  "Li+"=>["ION","Li+"],
-	  "Na+"=>["ION","Na+"],
-	  "K+"=>["ION","K+"],
-	  "Rb+"=>["ION","Rb+"],
-	  "Cs+"=>["ION","Cs+"],
-	  "Cl-"=>["ION","Cl-"],
-	  "MG2"=>["ION","MG"],
-	  "IB"=>["ION","IB"],
-	  "CIO"=>["ION","CIO"],
-	  "CL3"=>["OTH","CHCL3"],
-	  "MOH"=>["OTH","MEOH"],
-	  "NMA"=>["OTH","NMA"]);
-
-#output directory
-my $datdir=".";
-
-#proceed through each solvent
-foreach my $solvent (keys (%solvent)) {
-    print "solvent $solvent\n";
-    tleap($solvent,@{$solvent{$solvent}},$datdir);
-    buildParm($solvent,$datdir);
-}
-
-exit;
-
-#creates a single residue system in TLEAP for the given molecule.
-#$solvent - residue name
-#$type    - i.e. WAT, ION, OTH...
-#$alt     - frcmod suffix (typically not needed for IONs)
-#$datdir  - destination for resulting PARM7 and RST7 files
-sub tleap {
-    my ($solvent,$type,$alt,$datdir) = @_;
-
-    my $tleap = "tleap";
-    my $inp = "$datdir/$solvent.inp";
-    my $parm = "$datdir/$solvent.parm7";
-    my $rst = "$datdir/$solvent.rst7";
-    my $log = "$datdir/$solvent.log";
-    open(INP,">$inp") || die "ERROR: could not open >$inp:$!\n";
-    print INP "source leaprc.ff99SB\n";
-    if($type eq "WAT"){
-	print INP "WAT = $solvent\n";
-	if($solvent ne "TP3"){
-	    print INP "loadAmberParams frcmod.".(lc($alt))."\n";
-	}
-    }elsif($type eq "OTH"){
-	print INP "loadAmberParams frcmod.".(lc($alt))."\n"
-	    ."loadAmberPrep ".lc($alt).".in\n";
-    }elsif($type eq "ION"){
-	
-    }else{
-	die "ERROR: unknown solvent type: $solvent $type\n";
-    }
-    print INP "sol = sequence { $solvent }\n";
-    print INP "saveAmberParm sol $parm $rst\n"
-	."quit\n";
-    close INP;
-
-    my $cmd="$tleap -f $inp >& $log";
-    !system($cmd) || die "ERROR: could not run $cmd:$!\n";
-
-    unlink $inp;
-    unlink $log;
-    unlink "leap.log";
-}
-
-#create the MDL file for the molecule.  A PARM7 and RST7 file should exist in 
-#$datdir.  These are read in and processed for output in the MDL file and then 
-#deleted.
-#$solvent - solvent molecule name.  Used as the root name for all files
-#$datdir  - output location for all work
-sub buildParm {
-     my ($solvent,$datdir) = @_;
-     local $_;
-     my $rst="$solvent.rst7";
-     my $parm="$solvent.parm7";
-     my $mdl="$solvent.mdl";
-     #all the parameters we will write out
-     my (@coord, @name, @charge, @mass, @typeNum, @typeName, @lja, @ljb, 
-	 @radius, @ljIndex, @ljsigma, @ljepsilon,@multi);
-     my $i=0;
-     #get coordinates
-     open(RST,"<$rst") || die "ERROR: could not open <$rst:$!\n";
-     <RST>;
-     <RST>;
-     while(<RST>){
-	 
-	 my @temp=split();
-	 push(@coord, [@temp[0,1,2]]);	 
-	 $i++;
-	 if($#temp==5){
-	     push(@coord, [@temp[3,4,5]]);
-	     $i++;	     
-	 }
-     }
-     close RST;
-
-     #read all parameters into a hash
-     my %top = readPRMTOP($parm);
-
-     #compute the well depth and radius for each atom.
-     my $nType = sqrt($#{$top{NONBONDED_PARM_INDEX}}+1);
-     for(my $i=0;$i<=$#{$top{ATOM_TYPE_INDEX}};$i++){
-	 my $index = ${$top{NONBONDED_PARM_INDEX}}[$nType*(${$top{ATOM_TYPE_INDEX}}[$i]-1) + ${$top{ATOM_TYPE_INDEX}}[$i]-1]-1;
-	 if($index < 0 || ${$top{LENNARD_JONES_BCOEF}}[$index] == 0 || ${$top{LENNARD_JONES_ACOEF}}[$index] == 0){
-#	     push(@ljepsilon, 0.046);
-#	     push(@ljsigma, .2*pow(2,1./6.));
-	     #if the radius is zero, assign a sero for further processing
-	     push(@ljepsilon, 0);
-	     push(@ljsigma, 0);
-	     next;
-	 }
-	 push(@ljepsilon, ${$top{LENNARD_JONES_BCOEF}}[$index]*${$top{LENNARD_JONES_BCOEF}}[$index]/(4*${$top{LENNARD_JONES_ACOEF}}[$index]));
-	 push(@ljsigma, pow(2*${$top{LENNARD_JONES_ACOEF}}[$index]/(${$top{LENNARD_JONES_BCOEF}}[$index]),1./6.)/2.);
-     }
-
-     #use the coincident radius for atoms that have none defined
-     for(my $index = 0; $index <= $#ljsigma; $index++){
-	 !$ljsigma[$index] || next;
-	 my @parent = getParents(\%top,$index);
-	 my @bondlength;
-	 for my $parent (@parent){
-	     push(@bondlength,getBondLength(\%top,$index,$parent));
-	 }
-	 ($ljsigma[$index],$ljepsilon[$index]) = calcCoincidentRadius(\@parent,\@bondlength,\@ljsigma,\@ljepsilon);
-	 #Radius must be greater than 0
-	 $ljsigma[$index] >0 || warn "Negative coincident radius for atom ".($index+1).
-	     "${$top{ATOM_NAME}}[$index]\n";
-     }
-     
-     if($multiplicity){
-	 removeMultiples(\@coord, \@{$top{ATOM_NAME}}, \@{$top{CHARGE}}, \@{$top{MASS}}, \@{$top{ATOM_TYPE_INDEX}}, \@{$top{AMBER_ATOM_TYPE}}, \@{$top{RADII}}, \@{$top{NONBONDED_PARM_INDEX}}, \@ljsigma, \@ljepsilon, \@multi);
-     }
-
-     #write the output in PARM7 or namelist format
-     open(MDL,">$solvent.mdl") || die "ERROR: could not open >$solvent.mdl:$!\n";
-     if($ambParm){
-	 writeAmbHeader(*MDL);
-	 writeAmbParm(*MDL,"20a4","TITLE",[$solvent]);
-     }else{
-	 print MDL "&SOLVENT\n";
-     }
-
-     #write array sizes
-     if($ambParm){
-	 my @sizes = ($#coord+1,$#{$top{ATOM_TYPE_INDEX}}+1);
-	 writeAmbParm(*MDL,"10I8","POINTERS",\@sizes);
-     }
-
-     #write atom types
-     if($ambParm){
-	 writeAmbParm(*MDL,"10I8","ATMTYP",\@{$top{ATOM_TYPE_INDEX}});
-     }else{
-	 writeParam(*MDL,"!ATOM TYPES","ATMTYP",\@{$top{ATOM_TYPE_INDEX}});
-     }
-
-     #write atom names
-     if($ambParm){
-	 writeAmbParm(*MDL,"20a4","ATMNAME",\@{$top{ATOM_NAME}});
-     }else{
-	 writeParam(*MDL,"!ATOM NAMES","ATMNAME",\@{$top{ATOM_NAME}});
-     }
-
-     #write masses
-     if($ambParm){
-	 writeAmbParm(*MDL,"5e16.8","MASS",\@{$top{MASS}});
-     }else{
-	 writeParam(*MDL,"!MASS","MASS",\@{$top{MASS}});
-     }
-
-     #write charges
-     if($ambParm){
-	 writeAmbParm(*MDL,"5e16.8","CHG",\@{$top{CHARGE}});
-     }else{
-	 writeParam(*MDL,"!CHARGES","CHG",\@{$top{CHARGE}});
-     }
-
-     #write LJ epsilon parameters
-     if($ambParm){
-	 writeAmbParm(*MDL,"5e16.8","LJEPSILON",\@ljepsilon);
-     }else{
-	 writeParam(*MDL,"!LENNARD-JONES EPSILON PARAMETER","LJEPSILON",\@ljepsilon);
-     }
-
-     #write LJ sigma parameters
-     if($ambParm){
-	 writeAmbParm(*MDL,"5e16.8","LJSIGMA",\@ljsigma);
-     }else{
-	 writeParam(*MDL,"!LENNARD-JONES SIGMA PARAMETER","LJSIGMA",\@ljsigma);
-     }
-     #write multiplicity
-     if($multiplicity){
-	 if($ambParm){
-	     writeAmbParm(*MDL,"10I8","MULTI",\@multi);
-	 }else{
-	     writeParam(*MDL,"!ATOM TYPE MULTIPLICITY","MULTI",\@multi);
-	 }
-     }
-     
-     
-     #write coordinates
-     if($ambParm){
-	 my @flatcoord;
-	 for(my $i=0; $i<=$#coord; $i++){
-	     for(my $j=0; $j<=$#{$coord[$i]}; $j++){
-		 push(@flatcoord,${$coord[$i]}[$j]);
-	     }
-	 }
-	 writeAmbParm(*MDL,"5e16.8","COORD",\@flatcoord);
-     }else{
-	 print MDL "/\n";
-	 for(my $i=0;$i<=$#coord;$i++){
-	     print MDL "@{$coord[$i]}\n";
-	 }
-     }
-     close MDL;
-     unlink "$solvent.parm7";
-     unlink "$solvent.rst7";
-}
-
-#looks for atoms with identical charge and Lennard-Jones parameters.  
-#If found, the later is deleted and the multiplicity of the remainder is incremented.
-sub removeMultiples {
-    my ($coord, $name, $charge, $mass, $typeNum, $typeName, $radius, $ljIndex, $ljsigma, $ljepsilon, $multi) = @_;
-    for(my $i=0; $i<=$#{$name}; $i++){
-	@$multi[$i]=1;
-    }
-    for(my $i=0; $i<=$#{$name}; $i++){
-	for(my $j=$i+1; $j<=$#{$name}; $j++){
-	    if(@$ljsigma[$i] == @$ljsigma[$j] && @$ljepsilon[$i] == @$ljepsilon[$j] &&
-	       @$charge[$i] == @$charge[$j]){
-		my $coordi=0;
-		for(my $k=0;$k<=$i;$k++){
-		    $coordi+=@$multi[$k];
-		}
-		my $coordj=0;
-		for(my $k=0;$k<=$j;$k++){
-		    $coordj+=@$multi[$k];
-		}
-		$coordj--;
-		(@$coord[$coordi],@$coord[$coordj])=(@$coord[$coordj],@$coord[$coordi]);
-		splice(@$name,$j,1);
-		splice(@$charge,$j,1);
-		splice(@$mass,$j,1);
-		splice(@$typeNum,$j,1);
-		splice(@$typeName,$j,1);
-		splice(@$radius,$j,1);
-		splice(@$ljIndex,$j,1);
-		splice(@$ljsigma,$j,1);
-		splice(@$ljepsilon,$j,1);
-		splice(@$multi,$j,1);
-		@$multi[$i]++;
-		$j--;
-	    }
-	}
-    }
-}
-
-#reads in the given PARM7 file into a hash using the %FLAG for each section as 
-#the key.  Each key then contains a 1D array of values.  The entire hash is returned.
-sub readPRMTOP {
-  my $prmtop = shift;
-  my %prmHash;
-  open(PRMTOP,"<$prmtop") || die "ERROR: could not open $prmtop:$!\n";
-  my ($flag,$format,$size);
-  while(<PRMTOP>){
-    if(/^\%VERSION/){
-      
-    }elsif(/^\%FLAG/){
-      chomp;
-      m/\%FLAG ([^\s]+)/;
-      $flag = $1;
-    }elsif(/^\%FORMAT/){
-      chomp;
-      m/\%FORMAT\(([^\s]+)\)/;
-      $format = $1;
-      $format =~ m/[A-Za-z](\d+)/;
-      $size = $1;
-    }else{
-      my @data = split(/(.{$size})/);
-      for(my $idata=1;$idata <=$#data;$idata+=2){
-	push(@{$prmHash{$flag}},$data[$idata]);
-      }
-    }
-  }
-  close PRMTOP;
-  return %prmHash;
-}
-
-#writes a parameter section using the namelist format
-#$mdl     - filehandle
-#$comment - comment to preceed
-#$varName - name of variable to be written
-#$data    - 1D array of data values
-sub writeParam {
-    my ($mdl,$comment,$varName,$data)=@_;
-    print $mdl "$comment\n";
-    print $mdl "$varName = ";
-    foreach my $data (@$data){
-	print $mdl "$data, ";
-    }
-    print $mdl "\n";
-}
-
-#writes the PARM7 format header
-#$mdl - filehandle
-sub writeAmbHeader {
-    my ($mdl) = @_;
-    my ($sec, $min, $hr, $day, $mon, $year) = localtime;
-    $mon++;
-    $year = sprintf("%02d", $year % 100);
-    print $mdl "%";
-    printf $mdl "VERSION  VERSION_STAMP = V0001.000  DATE = %02d/%02d/%02d  %02d:%02d:%02d\n"
-	,$mon,$day,$year,$hr,$min,$sec;
-}
-
-#writes a parameter section for the PARM7 format
-#$mdl     - filehandle
-#$format  - string conversion format for reading/writing
-#$varname - name of variable to be written
-#$data    - 1D array of data to be written
-sub writeAmbParm {
-    my ($mdl,$format,$varName,$data)=@_;
-    print $mdl "\%FLAG $varName\n";
-    print $mdl "\%FORMAT($format)\n";
-    $format =~ m/(\d+)([A-Za-z]+)([\d\.]+)/;
-    my ($col,$type,$precision)=($1,$2,$3);
-    $type = lc($type);
-    $type=~s/a/s/;
-    for(my $i=0; $i<=$#{$data}; $i+=$col){
-	for(my $j=0; $i+$j<=$#{$data} && $j < $col; $j++){
-	    printf $mdl "%$precision$type",@$data[$i+$j];
-	}
-	print $mdl "\n";
-    }
-
-}
-
-#gets 'parents' or other atoms bonded to $index.  return as an array of indices 
-#with 0 offset
-#$top   - PARM7 hash
-#$index - atom index to get the bonded neighbours of
-sub getParents {
-    my $top = shift;
-    my $index = shift;
-    my @parent;
-    #find all of the entrys in BONDS_WITHOUT_HYDROGEN and
-    #BONDS_INC_HYDROGEN that contain this index
-    for(my $ibond=0; $ibond <= $#{$top->{BONDS_INC_HYDROGEN}}; $ibond+=3){
-	if(${$top->{BONDS_INC_HYDROGEN}}[$ibond] == $index*3){
-	    push(@parent,${$top->{BONDS_INC_HYDROGEN}}[$ibond+1]/3);
-	}
-	if(${$top->{BONDS_INC_HYDROGEN}}[$ibond+1] == $index*3){
-	    push(@parent,${$top->{BONDS_INC_HYDROGEN}}[$ibond]/3);
-	}
-    }
-    for(my $ibond=0; $ibond <= $#{$top->{BONDS_WITHOUT_HYDROGEN}}; $ibond+=3){
-	if(${$top->{BONDS_WITHOUT_HYDROGEN}}[$ibond] == $index*3){
-	    push(@parent,${$top->{BONDS_WITHOUT_HYDROGEN}}[$ibond+1]/3);
-	}
-	if(${$top->{BONDS_WITHOUT_HYDROGEN}}[$ibond+1] == $index*3){
-	    push(@parent,${$top->{BONDS_WITHOUT_HYDROGEN}}[$ibond]/3);
-	}
-    }
-    return @parent;
-}
-
-#given two bonded atoms, return the equilibrium bond length
-#$top - PARM7 hash
-#$A   - atom 1 index
-#$B   - atom 2 index
-sub getBondLength {
-    my $top = shift;
-    my $A = shift;
-    my $B = shift;
-    for(my $ibond=0; $ibond <= $#{$top->{BONDS_INC_HYDROGEN}}; $ibond+=3){
-	if((${$top->{BONDS_INC_HYDROGEN}}[$ibond] == $A*3 && 
-	   ${$top->{BONDS_INC_HYDROGEN}}[$ibond+1] == $B*3) ||
-	   (${$top->{BONDS_INC_HYDROGEN}}[$ibond] == $B*3 && 
-	   ${$top->{BONDS_INC_HYDROGEN}}[$ibond+1] == $A*3)){
-	    return ${$top->{BOND_EQUIL_VALUE}}[${$top->{BONDS_INC_HYDROGEN}}[$ibond+2]-1];
-	}
-    }
-    for(my $ibond=0; $ibond <= $#{$top->{BONDS_WITHOUT_HYDROGEN}}; $ibond+=3){
-	if((${$top->{BONDS_WITHOUT_HYDROGEN}}[$ibond] == $A*3 && 
-	   ${$top->{BONDS_WITHOUT_HYDROGEN}}[$ibond+1] == $B*3) ||
-	   (${$top->{BONDS_WITHOUT_HYDROGEN}}[$ibond] == $B*3 && 
-	   ${$top->{BONDS_WITHOUT_HYDROGEN}}[$ibond+1] == $A*3)){
-	    return ${$top->{BOND_EQUIL_VALUE}}[${$top->{BONDS_WITHOUT_HYDROGEN}}[$ibond+2]-1];
-	}
-    }
-    #atoms apparently are not bonded
-    warn "bond not found $A $B\n";
-}
-
-#given a number of atoms and there bondlengths to the atom in questions, calculates
-#the coincident radius.  This is the radius such that LJ potential of our enclosed
-#atom is 0 at the same radial distance from the enclosing atom along the bond.  
-#I.e. the enclosing radius - the bond length.  This can be a negative number if the
-#target atom is not actually enclosed.  This sub routine picks the largest coincident
-#radius and returns it along with a epsilon value one tenth that of the enclosing
-#atoms.
-#$parents     - indices of atoms bonded to the target atom
-#$bondlengths - equilibrium bond lengths between parent atoms and the target
-#$ljsigma     - list of all sigma values so far determined
-#$ljepsilon   - list of all epsilon values so far determined
-sub calcCoincidentRadius {
-    my $parents = shift;
-    my $bondlengths = shift;
-    my $ljsigma = shift;
-    my $ljepsilon = shift;
-    
-    my @radius;
-    for(my $i = 0; $i<= $#{$parents}; $i++){
-	push(@radius, ((@$ljsigma[@$parents[$i]])/pow(2,1./6.) -@$bondlengths[$i])*pow(2,1./6.));
-    }
-    
-    my $index = maxIndex(@radius);
-    return ($radius[$index],0.1*@$ljepsilon[@$parents[$index]]);
-}
-
-#returns the index of the maximum value in an array
-sub maxIndex {
-    my @a = @_;
-    my $index = 0;
-    for(my $i = 1; $i<=$#a; $i++){
-	if($a[$i] > $a[$index]){
-	    $index=$i;
-	}
-    }
-    return $index;
-}