1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
|
#!/usr/bin/perl
# usage: make_gromos_nb.pl gromos_atoms > gromos_nb.itp
# this script generates a GROMOS96 nonbonded forcefield file for GROMACS,
# with as input the file gromos_atoms, which contains records of the form
#: 1 O 0.04756 0.8611E-3 1.125E-3 0.0
#: #CS6 CS12 parameters LJ14PAIR
#: 0.04756 0.8611E-3
#: 1 1 2 2 2 2 2 2 2 2 1 1 1 1 1 1 1 1 1 1
#: 2 2 2 2 2 2 2 1 1 1 1 1 2 2 1 1 1 1 1 1
#: 1 1 1
#: #---
# taken directly from ifp43a1.dat. For a description of what the numbers
# mean, see the GROMOS96 manual and the fie ifp43a1.dat
# set the input separator to #--- so we read one atom entry at a time
# make sure the records are actually separated by #---\n and not by say #--\n!!
# don't put a separator at the end of the file, only between records
$/= "#---\n";
# start arrays with 1 instead of 0
$[=1;
$i=1;
# read in all entries (43 in ifp43a1.dat). This will BREAK for other input!
while(<>) {
@lines = split('\n',$_);
($number[$i],$name[$i],$c6[$i],$c12_1[$i],$c12_2[$i],$c12_3[$i]) =
split(' ',$lines[1]);
($c6_14[$i],$c12_14[$i]) = split(' ',$lines[3]);
$combination[$i] = $lines[4] . $lines[5] . $lines[6];
# one type is called P,SI, the same LJ parameters for both P and SI
# treat P,SI different: create a 44th type for Si, just a copy of Si,P
# and rename P,SI to P
if ($name[$i] =~ /P,SI/) {
$number[44] = 44; $name[44] = "SI"; $c6[44] = $c6[$i];
$c12_1[44] = $c12_1[$i]; $c12_2[44] = $c12_2[$i];
$c12_3[44] = $c12_3[$i]; $combination[44] = $combination[$i];
$name[$i] = "P";
$P = $i;
}
$i++;
}
# now $number[$i] has the atom nr, $name[$i] the name, $c6 the C6 LJ parameter
# $c12_1[$i], $c12_2[$i], $c12_3[$i] the three C12 parameters and
# $combination[$i] the matrix with 43 elements that tells which C12 parameter
# you need in combination with each of the 44 atomtypes. $i runs from 1 to 44,
# one for each atom type. This goes nicely wrong because of the Stupid SI,P
# entry so we have to give an extra element 44 with the same value as that of
# P
# start printing to the output file: header for gromos-nb.itp
print "[ atomtypes ]\n";
print ";name mass charge ptype c6 c12\n";
# print out the atomtypes, plus the C6 and C12 for interactions with themselves
# the masses and charges are set to 0, the masses should come from gromos.atp
for ($j=1;$j<=44;$j++) {
# lookup the C12 with itself
@c12types = split(' ',$combination[$j]);
$c12types[44] = $c12types[$P]; # SI has the same as P
if ($c12types[$j] == 1)
{$c12 = $c12_1[$j];}
elsif ($c12types[$j] == 2)
{$c12 = $c12_2[$j];}
elsif ($c12types[$j] == 3)
{$c12 = $c12_3[$j];}
else {die "Error [atomtypes]: c12 type is not 1,2,3:j=$j,c12=$c12\n";}
print sprintf("%5s 0.000 0.000 A %10.8g %10.8g\n",
$name[$j],$c6[$j]*$c6[$j],$c12*$c12);
}
# Now make the LJ matrix. It's ok to do some double work in shell scripts,
# trust me.
print "\n";
print "[ nonbond_params ]\n";
print "; i j func c6 c12\n";
for ($j=1;$j<=44;$j++) {
for ($k=1;$k<$j;$k++) {
# lookup the C12 of j for k
@c12types_j = split(' ',$combination[$j]);
$c12types_j[44] = $c12types_j[$P]; # SI has the same as P
if ($c12types_j[$k] == 1)
{$c12_j = $c12_1[$j];}
elsif ($c12types_j[$k] == 2)
{$c12_j = $c12_2[$j];}
elsif ($c12types_j[$k] == 3)
{$c12_j = $c12_3[$j];}
else {
die "Error [nonbond-params] j=$j,k=$k: c12 type is not one of 1,2,3\n";
}
# lookup the C12 of k for j
@c12types_k = split(' ',$combination[$k]);
$c12types_k[44] = $c12types_k[$P]; # SI has the same as P
if ($c12types_k[$j] == 1)
{$c12_k = $c12_1[$k];}
elsif ($c12types_k[$j] == 2)
{$c12_k = $c12_2[$k];}
elsif ($c12types_k[$j] == 3)
{$c12_k = $c12_3[$k];}
else {
die "Error [nonbond-params] j=$j,k=$k: c12 type is not one of 1,2,3\n";
}
print sprintf("%8s %8s 1 %10.8g %10.8g\n", $name[$j], $name[$k],
$c6[$j]*$c6[$k], $c12_j*$c12_k);
}
}
# Now do the same for the 1-4 interactions
print "\n";
print "[ pairtypes ]\n";
print "; i j func c6 c12\n";
for ($j=1;$j<=44;$j++) {
for ($k=1;$k<=$j;$k++) {
print sprintf("%8s %8s 1 %10.8g %10.8g\n", $name[$j], $name[$k],
$c6_14[$j]*$c6_14[$k], $c12_14[$j]*$c12_14[$k]);
}
}
|
