Commit b86654a9 authored by Adam Hospital's avatar Adam Hospital
Browse files

First commit

parents
Pipeline #5790 skipped
V = 0.2.6
EXEDIR =.
#
.PHONY: lib setup discrete
#
all: setup discrete
lib:
cd lib; make
discrete: lib
cd discrete; make
install: all
cp discrete/discrete $(EXEDIR)/discrete.exe
clean:
cd lib;make clean
cd discrete; make clean
rm $(EXEDIR)/discrete.exe
# GOdMD_different_natoms
cd dat/2lao_A-1lst_A-1/
../../discrete/discrete -i dmdtest.in -pdbin 2lao_A.REF.pdb -pdbtarg 1lst_A.pdb -p1 2lao_A.aln -p2 1lst_A.aln -trj 2lao_A-1lst_A_trajectory.crd -ener energy.dat -o 2lao_A-1lst_A_1.log
### Parameters description from ```discrete/paramSet.f```
- TSNAP: Frequency to output structures in the simulation. Low number gives more structures.
- TRECT: Frequency to apply the Maxwell-Demon algorithm
- TACT: Frequency to update the bonded distances.
- ACCEPTANCE: Fraction of sampled structures that would be accepted in Monte Carlo test. Low values give biased results. High values may lead to convergence problems. I recommend values between 0.5 and 0.8
- GOENER: Energy depth of the well kcal/mol
- SIGMAGO: Width of bonded bonds, in Armstrongs.
- WELLWIDTH: Width of non-bonded bonds, in Armstrongs.
- RGYRPERCENT: Maximum distance to consider non-bonded and changing interactions (double go wells) as a fraction of the average radius of gyration
- RCUT1: Maximum distance to consider non-bonded and non-changing interactions (single go wells)
- MINWENERWELL: Minimum energy interaction per well. Below this number the pseudo-Metadynamics algorithm will ignore this interaction.
- SIGMAFAKE: Bonded distance for non-consecutive particles (a gap in the structure)
\ No newline at end of file
DEBUG =
#DEBUG = -g -fbacktrace -fbounds-check
#DEBUG = -g -fbacktrace -fbounds-check -DDEBUG
DEBUGNOCHECK =
#DEBUGNOCHECK = -g -fbacktrace -DDEBUG
#PROFILE =
PROFILE = -pg
#
OPTIONS = -O3 -fimplicit-none -funroll-loops -mtune=generic -Wall -ftree-vectorizer-verbose=1 -ffast-math
LOPTIONS= -static-libgfortran
F77= /usr/local/bin/gfortran
CPP= cpp
LIBDIR = ../lib
CFLAGS = -c $(DEBUG) $(PROFILE) -I$(LIBDIR)
CFLAGSNODEBUG = -c $(PROFILE) -I$(LIBDIR)
CFLAGSNOCHECK = -c $(DEBUGNOCHECK) $(PROFILE) -I$(LIBDIR)
LFLAGS = $(PROFILE) -L$(LIBDIR)
ENDFLAG = -fconvert=little-endian
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
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
This diff is collapsed.
Molecules: 1 Residues: 238 Atoms: 238
Molecules: 1 Residues: 238 Atoms: 238
NATOM COMMON 238
=================================================
= =
= GOdMD =
= =
= P Sfriso, A. Emperador, JL. Gelpi, M.Orozco =
= =
= (c) 2013 =
=================================================
I N P U T F I L E S
====================
Settings : dmdtest.in
Initial PDB : 2lao_A.REF.pdb
Energies : energy.dat
Trajectory (PDB) : 2lao_A-1lst_A_trajectory.crd
Target PDB : 1lst_A.pdb
Calculation Log : 2lao_A-1lst_A_1.log
Table of same residu: 2lao_A.aln
Table of same residu: 1lst_A.aln
------------------------------------------------------------
| CALCULATION PARAMETERS |
------------------------------------------------------------
| Simulation settings |
------------------------------------------------------------
| Simulation Time (ps) (Nbloc x TSnap) | 5000.000 |
| Output structure (fs) | TSnap | 1000.000 |
| Re-scoring target (fs) | Trect | 100.000 |
| Update velocities (fs) | Tcorr | 100.000 |
| Update Lists, collision times (fs)| Tact | 5.000 |
| Min. accepted colision time (fs) | TMin | 0.00000010 |
| Temperature (K) | Temp | 36.08 |
------------------------------------------------------------
------------------------------------------------------------
| DIMS |
------------------------------------------------------------
|Acceptance |ACCEPTANCE| 0.76999998 |
|GO Energy (kcal) |GOENER | 0.25000000 |
|Changing distances well width(A) |sigmago | 0.10000000 |
|Tolerance between changes(A) |dtol | 0.30000001 |
|Not-Changing distances well width(A)|wellwidth | 0.13000000 |
|Dummy tolerance (A) |sigmaFake | 0.02000000 |
|Factor scaling NMA coupling |scalFactor| 0.30000001 |
|Cut-off in Radius Gyration | rgyrPercent| 0.50000000 |
|Cut-off not changing GO interact.(A)|rcut1 |100.00000000 |
|Minimum Energy of any well (kcal) | minEnerWell| 0.10000000 |
|Min time to eq. energy(ps) |Ener_evo_size | 120 |
|Number of ENM Eigenvectors |NEVECS | 5 |
|Fraction of energy to add at MD |DEGO_BASE | 0.03000000 |
|Minimum acceptable RMSD (A) |errorAcceptable| 1.25000000 |
|Number of frames to skip |NSKIP | 5 |
------------------------------------------------------------
| Other |
------------------------------------------------------------
| Random generator seed | 2835 |
| IDAB, IGOAB, IWR, ISOLV | 1 1 1 1 |
------------------------------------------------------------
PNS rgyr 18.161676487001049 CutOff 9.0808382435005246
UNO QUE SE MUEVE 6 5 3.8475082736698729 3.7694414300035253
UNO QUE SE MUEVE 13 12 3.7633583505571577 3.8434071459613461
UNO QUE SE MUEVE 73 72 3.7368136431476673 3.8148065644627778
UNO QUE SE MUEVE 112 111 3.7613276375573781 3.8663698270110842
UNO QUE SE MUEVE 116 115 3.7139506202096930 3.8323602804927526
UNO QUE SE MUEVE 122 121 3.7913098200750279 3.8832039056801846
UNO QUE SE MUEVE 142 141 3.7365821634840390 3.8529786020418011
UNO QUE SE MUEVE 158 157 3.7505540120924170 3.8595092803377535
UNO QUE SE MUEVE 161 160 3.7657241445026295 3.8456408888719804
This diff is collapsed.
MODEL 0
ATOM 1 CA ALA A1 16.538 89.143 48.164
ATOM 2 CA LEU A2 19.128 86.720 46.781
ATOM 3 CA PRO A3 22.331 86.237 48.900
ATOM 4 CA GLN A4 25.240 88.283 47.531
ATOM 5 CA THR A5 27.329 85.159 47.769
ATOM 6 CA VAL A6 25.631 81.739 48.242
ATOM 7 CA ARG A7 27.694 79.407 50.552
ATOM 8 CA ILE A8 26.992 75.834 49.546
ATOM 9 CA GLY A9 27.986 73.033 51.908
ATOM 10 CA THR A10 29.007 69.677 50.539
ATOM 11 CA ASP A11 30.636 66.544 51.948
ATOM 12 CA THR A12 33.646 65.863 49.689
ATOM 13 CA THR A13 34.021 62.181 50.371
ATOM 14 CA TYR A14 31.561 60.974 47.667
ATOM 15 CA ALA A15 32.984 60.558 44.199
ATOM 16 CA PRO A16 31.664 60.882 41.546
ATOM 17 CA PHE A17 29.373 63.524 43.090
ATOM 18 CA SER A18 31.658 65.498 45.437
ATOM 19 CA SER A19 35.381 64.983 46.000
ATOM 20 CA LYS A20 38.653 66.913 46.050
ATOM 21 CA ASP A21 41.296 66.550 43.452
ATOM 22 CA ALA A22 45.100 66.692 43.829
ATOM 23 CA LYS A23 44.976 70.518 43.839
ATOM 24 CA GLY A24 42.303 70.654 46.636
ATOM 25 CA GLU A25 39.609 71.837 44.090
ATOM 26 CA PHE A26 36.088 70.426 44.604
ATOM 27 CA ILE A27 35.214 68.189 41.669
ATOM 28 CA GLY A28 32.252 65.998 40.757
ATOM 29 CA PHE A29 28.781 65.976 39.305
ATOM 30 CA ASP A30 27.150 67.897 42.218
ATOM 31 CA ILE A31 29.899 70.506 41.870
ATOM 32 CA ASP A32 29.319 70.809 38.115
ATOM 33 CA LEU A33 25.558 71.101 38.551
ATOM 34 CA GLY A34 26.109 73.610 41.407
ATOM 35 CA ASN A35 28.528 75.827 39.521
ATOM 36 CA GLU A36 26.442 75.868 36.374
ATOM 37 CA MET A37 23.272 76.837 38.313
ATOM 38 CA CYS A38 25.206 79.706 40.069
ATOM 39 CA LYS A39 26.331 80.857 36.636
ATOM 40 CA ARG A40 22.763 81.018 35.294
ATOM 41 CA MET A41 21.514 82.568 38.581
ATOM 42 CA GLN A42 24.266 85.147 38.243
ATOM 43 CA VAL A43 25.443 84.941 41.881
ATOM 44 CA LYS A44 28.802 84.119 43.466
ATOM 45 CA CYS A45 28.752 80.693 45.173
ATOM 46 CA THR A46 31.372 79.414 47.530
ATOM 47 CA TRP A47 31.830 75.738 48.390
ATOM 48 CA VAL A48 32.252 74.839 52.024
ATOM 49 CA ALA A49 33.462 71.313 52.821
ATOM 50 CA SER A50 31.607 69.772 55.670
ATOM 51 CA ASP A51 30.715 66.478 57.187
ATOM 52 CA PHE A52 27.340 65.267 55.994
ATOM 53 CA ASP A53 25.769 65.153 59.431
ATOM 54 CA ALA A54 26.617 68.833 60.049
ATOM 55 CA LEU A55 25.061 70.051 56.828
ATOM 56 CA ILE A56 21.377 70.501 57.811
ATOM 57 CA PRO A57 22.522 71.978 61.256
ATOM 58 CA SER A58 24.753 74.406 59.407
ATOM 59 CA LEU A59 21.913 75.322 57.027
ATOM 60 CA LYS A60 19.447 75.962 59.815
ATOM 61 CA ALA A61 22.051 78.021 61.713
ATOM 62 CA LYS A62 22.742 80.061 58.519
ATOM 63 CA LYS A63 26.368 79.051 58.390
ATOM 64 CA ILE A64 25.676 77.996 54.763
ATOM 65 CA ASP A65 22.825 78.917 52.341
ATOM 66 CA ALA A66 22.389 75.657 50.400
ATOM 67 CA ILE A67 23.368 72.025 50.594
ATOM 68 CA ILE A68 24.363 70.166 47.397
ATOM 69 CA SER A69 25.654 66.856 48.595
ATOM 70 CA SER A70 23.502 64.044 47.042
ATOM 71 CA LEU A71 20.885 64.897 49.641
ATOM 72 CA SER A72 17.759 62.866 48.880
ATOM 73 CA ILE A73 14.379 64.443 49.109
ATOM 74 CA THR A74 12.768 62.420 51.844
ATOM 75 CA ASP A75 9.543 63.218 53.650
ATOM 76 CA LYS A76 11.371 63.507 56.979
ATOM 77 CA ARG A77 13.774 66.173 55.610
ATOM 78 CA GLN A78 10.871 68.034 53.935
ATOM 79 CA GLN A79 9.513 68.433 57.466
ATOM 80 CA GLU A80 12.676 70.477 58.206
ATOM 81 CA ILE A81 14.023 72.128 55.022
ATOM 82 CA ALA A82 13.005 72.937 51.490
ATOM 83 CA PHE A83 14.447 71.349 48.289
ATOM 84 CA SER A84 14.947 72.497 44.724
ATOM 85 CA ASP A85 13.542 70.363 41.969
ATOM 86 CA LYS A86 15.339 66.998 41.605
ATOM 87 CA LEU A87 18.941 67.083 40.400
ATOM 88 CA TYR A88 19.147 63.299 39.676
ATOM 89 CA ALA A89 17.698 60.071 41.006
ATOM 90 CA ALA A90 19.018 58.039 43.963
ATOM 91 CA ASP A91 17.301 54.647 43.457
ATOM 92 CA SER A92 18.347 51.100 44.224
ATOM 93 CA ARG A93 18.621 48.166 41.806
ATOM 94 CA LEU A94 19.280 44.424 41.926
ATOM 95 CA ILE A95 22.687 43.087 40.857
ATOM 96 CA ALA A96 22.938 39.461 39.715
ATOM 97 CA ALA A97 25.051 37.244 37.423
ATOM 98 CA LYS A98 25.157 38.364 33.737
ATOM 99 CA GLY A99 22.230 36.837 31.950
ATOM 100 CA SER A100 20.107 36.219 35.090
ATOM 101 CA PRO A101 16.330 36.280 34.409
ATOM 102 CA ILE A102 15.806 37.797 37.951
ATOM 103 CA GLN A103 13.484 40.786 38.191
ATOM 104 CA PRO A104 12.573 43.031 41.215
ATOM 105 CA THR A105 9.068 41.473 41.496
ATOM 106 CA LEU A106 7.745 39.266 44.188
ATOM 107 CA GLU A107 6.962 36.698 41.488
ATOM 108 CA SER A108 10.529 36.574 40.341
ATOM 109 CA LEU A109 12.191 36.704 43.780
ATOM 110 CA LYS A110 10.291 34.008 45.618
ATOM 111 CA GLY A111 12.801 31.304 46.171
ATOM 112 CA LYS A112 15.817 33.510 45.741
ATOM 113 CA HIS A113 18.274 34.731 48.404
ATOM 114 CA VAL A114 18.743 38.462 48.211
ATOM 115 CA GLY A115 21.486 40.124 50.210
ATOM 116 CA VAL A116 21.063 43.500 51.699
ATOM 117 CA LEU A117 23.367 45.675 53.860
ATOM 118 CA GLN A118 22.212 45.691 57.469
ATOM 119 CA GLY A119 20.632 48.932 58.614
ATOM 120 CA SER A120 20.367 50.553 55.181
ATOM 121 CA THR A121 17.311 51.918 53.314
ATOM 122 CA GLN A122 17.626 48.897 51.045
ATOM 123 CA GLU A123 17.001 46.452 53.999
ATOM 124 CA ALA A124 14.016 48.541 55.040
ATOM 125 CA TYR A125 12.534 48.320 51.537
ATOM 126 CA ALA A126 13.439 44.646 51.148
ATOM 127 CA ASN A127 12.103 43.739 54.600
ATOM 128 CA ASP A128 8.659 45.404 54.140
ATOM 129 CA ASN A129 8.132 44.605 50.482
ATOM 130 CA TRP A130 9.834 41.230 50.122
ATOM 131 CA ARG A131 10.906 39.435 53.301
ATOM 132 CA THR A132 7.389 39.706 54.868
ATOM 133 CA LYS A133 6.112 38.044 51.689
ATOM 134 CA GLY A134 8.519 35.088 51.885
ATOM 135 CA VAL A135 11.552 36.213 49.824
ATOM 136 CA ASP A136 14.678 35.117 51.797
ATOM 137 CA VAL A137 16.341 38.423 52.460
ATOM 138 CA VAL A138 19.736 38.176 54.222
ATOM 139 CA ALA A139 21.274 41.203 55.913
ATOM 140 CA TYR A140 25.043 41.639 55.942
ATOM 141 CA ALA A141 27.256 43.584 58.304
ATOM 142 CA ASN A142 29.547 44.440 55.479
ATOM 143 CA GLN A143 29.413 45.180 51.773
ATOM 144 CA ASP A144 32.507 43.036 51.038
ATOM 145 CA LEU A145 30.795 39.915 52.248
ATOM 146 CA ILE A146 27.763 40.653 50.050
ATOM 147 CA TYR A147 29.939 40.831 46.974
ATOM 148 CA SER A148 31.880 37.680 47.875
ATOM 149 CA ASP A149 28.693 35.731 48.558
ATOM 150 CA LEU A150 27.194 36.831 45.221
ATOM 151 CA THR A151 30.488 35.902 43.550
ATOM 152 CA ALA A152 30.472 32.571 45.394
ATOM 153 CA GLY A153 26.903 31.801 44.317
ATOM 154 CA ARG A154 25.625 32.016 47.897
ATOM 155 CA LEU A155 23.248 34.830 46.981
ATOM 156 CA ASP A156 21.046 35.119 43.932
ATOM 157 CA ALA A157 20.999 38.947 43.895
ATOM 158 CA ALA A158 21.914 41.928 45.979
ATOM 159 CA LEU A159 19.779 45.076 46.324
ATOM 160 CA GLN A160 22.109 48.053 46.292
ATOM 161 CA ASP A161 22.168 51.748 45.568
ATOM 162 CA GLU A162 22.406 51.949 41.727
ATOM 163 CA VAL A163 25.190 54.481 41.336
ATOM 164 CA ALA A164 27.158 52.891 44.167
ATOM 165 CA ALA A165 26.887 49.441 42.590
ATOM 166 CA SER A166 27.911 50.926 39.213
ATOM 167 CA GLU A167 30.942 52.792 40.349
ATOM 168 CA GLY A168 32.007 50.433 43.031
ATOM 169 CA PHE A169 31.551 46.875 41.777
ATOM 170 CA LEU A170 30.218 46.504 38.271
ATOM 171 CA LYS A171 33.123 48.702 37.148
ATOM 172 CA GLN A172 35.636 46.772 39.286
ATOM 173 CA PRO A173 37.371 43.573 37.981
ATOM 174 CA ALA A174 35.246 41.580 40.368
ATOM 175 CA GLY A175 31.984 42.764 38.774
ATOM 176 CA LYS A176 32.771 41.973 35.110
ATOM 177 CA GLU A177 30.509 38.909 35.071
ATOM 178 CA TYR A178 27.638 40.707 36.847
ATOM 179 CA ALA A179 24.985 43.174 35.835
ATOM 180 CA PHE A180 21.906 45.034 36.926
ATOM 181 CA ALA A181 18.974 42.680 36.903
ATOM 182 CA GLY A182 15.855 44.520 35.740
ATOM 183 CA PRO A183 15.058 48.245 36.131
ATOM 184 CA SER A184 15.781 50.331 39.242
ATOM 185 CA VAL A185 13.414 50.252 42.172
ATOM 186 CA LYS A186 11.921 53.602 42.899
ATOM 187 CA ASP A187 10.146 54.169 46.189
ATOM 188 CA LYS A188 9.931 57.714 47.714
ATOM 189 CA LYS A189 9.257 56.261 51.127
ATOM 190 CA TYR A190 12.621 54.478 51.381
ATOM 191 CA PHE A191 14.800 56.113 48.774
ATOM 192 CA GLY A192 13.306 59.610 48.612
ATOM 193 CA ASP A193 12.458 61.450 45.423
ATOM 194 CA GLY A194 15.857 61.989 43.872
ATOM 195 CA THR A195 18.331 64.402 45.413
ATOM 196 CA GLY A196 17.933 68.196 45.463
ATOM 197 CA VAL A197 19.600 71.411 46.622
ATOM 198 CA GLY A198 18.700 71.654 50.392
ATOM 199 CA LEU A199 17.538 75.181 51.277
ATOM 200 CA ARG A200 15.928 76.995 54.306
ATOM 201 CA LYS A201 12.128 76.703 54.125
ATOM 202 CA ASP A202 11.920 80.498 54.216
ATOM 203 CA ASP A203 14.223 81.168 51.320
ATOM 204 CA THR A204 11.747 81.053 48.544
ATOM 205 CA GLU A205 13.509 83.492 46.105
ATOM 206 CA LEU A206 16.673 81.363 46.383
ATOM 207 CA LYS A207 14.628 78.188 45.637
ATOM 208 CA ALA A208 13.025 79.933 42.648
ATOM 209 CA ALA A209 16.425 81.003 41.424
ATOM 210 CA PHE A210 17.846 77.401 41.689
ATOM 211 CA ASP A211 14.703 75.914 40.077
ATOM 212 CA LYS A212 14.789 78.361 37.145
ATOM 213 CA ALA A213 18.474 77.719 36.517
ATOM 214 CA LEU A214 18.084 73.938 36.703
ATOM 215 CA THR A215 15.205 73.928 34.082
ATOM 216 CA GLU A216 17.267 76.114 31.781
ATOM 217 CA LEU A217 20.387 73.986 31.838
ATOM 218 CA ARG A218 18.265 70.885 31.374
ATOM 219 CA GLN A219 16.422 72.395 28.373
ATOM 220 CA ASP A220 19.499 73.615 26.467
ATOM 221 CA GLY A221 21.263 70.224 26.900
ATOM 222 CA THR A222 24.010 71.299 29.355
ATOM 223 CA TYR A223 22.692 68.710 31.867
ATOM 224 CA ASP A224 23.357 65.858 29.427
ATOM 225 CA LYS A225 26.810 67.079 28.409
ATOM 226 CA MET A226 27.688 67.343 32.090
ATOM 227 CA ALA A227 26.275 63.917 33.012
ATOM 228 CA LYS A228 28.172 62.334 30.083
ATOM 229 CA LYS A229 31.498 62.822 31.938
ATOM 230 CA TYR A230 30.516 60.457 34.706
ATOM 231 CA PHE A231 27.569 58.362 33.769
ATOM 232 CA ASP A232 26.887 56.257 30.764
ATOM 233 CA PHE A233 23.249 55.923 31.878
ATOM 234 CA ASN A234 20.353 58.332 32.507
ATOM 235 CA VAL A235 21.192 59.487 36.033
ATOM 236 CA TYR A236 18.214 61.802 35.973
ATOM 237 CA GLY A237 15.674 58.929 35.996
ATOM 238 CA ASP A238 12.095 59.914 35.138
ENDMDL
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
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
&INPUT
TSNAP=1000
TEMP=300
seed=2835
&END
!
! File: ANM.f
! Author: psfriso
!
! Created on October 30, 2012, 12:45 PM
!
MODULE ANM
! Data dictionary and variable declaration
INTEGER, PARAMETER :: DBL=SELECTED_REAL_KIND(p=13)! Double precision real number definition processor independent
INTEGER, PARAMETER :: SGL=SELECTED_REAL_KIND(p=6) ! Single Precision real number definition processor independent
PRIVATE
PUBLIC ANM1
CONTAINS
!
SUBROUTINE ANM1(nres,coord,nevecs, evec,evals)
USE geometryDP
USE geometry
USE ls_rmsd
IMPLICIT NONE
INTEGER ::RES3,i,k,J,ioerr
INTEGER :: CUTOFF
INTEGER, ALLOCATABLE :: INDX(:)
REAL(SGL), DIMENSION(nevecs), INTENT(INOUT) :: evals
REAL(DBL), ALLOCATABLE :: DIS(:,:),G(:,:),HESS(:,:)
INTEGER ,intent(IN):: NRES ,nevecs ! Residue number
REAL(DBL), ALLOCATABLE :: W(:),V(:,:),WD(:),INVCONT(:,:)
REAL(SGL), PARAMETER :: KBT=0.597,EIGENCUT=1E-5
REAL(DBL), ALLOCATABLE :: rv1(:)
REAL, PARAMETER :: SL=3, IN=2.8
TYPE(pointDP), DIMENSION(nres), INTENT(IN) :: coord
TYPE(pointDP), INTENT(INOUT), DIMENSION(nevecs,nres) :: evec
REAL(DBL), dimension(nres, 3) :: r2
r2(:,1) =coord(:)%x
r2(:,2) =coord(:)%y
r2(:,3) =coord(:)%z
RES3=NRES*3
! WRITE (6,*)'Allocating arrays..'
ALLOCATE (V(RES3,RES3) ,W(RES3), WD(RES3), INDX(RES3), INVCONT(NRES,NRES)&
,DIS(NRES,NRES),G(NRES,NRES),HESS(RES3,RES3),rv1(RES3),STAT=ioerr)
call errorAllocmem (ioerr, "ANM ")
V=0.0_DBL
W=0.0_DBL
WD=0.0_DBL
INDX=0
INVCONT=0.0_DBL
! Setting cutoff
IF (NRES.LT.50) THEN
CUTOFF=8
ELSE
CUTOFF=INT(SL*LOG(REAL(NRES))-IN)
ENDIF
!
! WRITE (*,*) 'Computing C-alpha distances...'
CALL DISTANCES (NRES,R2,DIS)
! WRITE (*,*) 'Setting force constants...'
CALL SET_GAMMA(CUTOFF,NRES,DIS,G)
!
! WRITE (*,*) 'Now we calculate hessian...'
CALL HESSIAN(NRES,R2,DIS,G,HESS)
!
! write(*,*) 'Diagonalizing....'
CALL SVDCMP(RES3, HESS, RES3, RES3, W, V, rv1)
! write(*,*) 'Putting the evals in order....'
CALL INDEXX(RES3,W,INDX)
DO I=1,nevecs+6
k=1
IF(I.GT.6) THEN
DO j=1,3*nres,3
! write(*,*) "desde aqui" ,V(J, INDX(I))
evec(I-6,k)%x=V(j,INDX(I))
evec(I-6,k)%y=V(j+1,INDX(I))
evec(I-6,k)%z=V(j+2,INDX(I))
k=k+1
!
evals(I-6)=W(INDX(i))
!
END DO
END IF
ENDDO
! DO i=6,11
! WRITE(*,*)" INDX ", i, W(INDX(i))
! END DO
! WRITE(*,*) 'ANM Subroutine finished !!'
DEALLOCATE(V ,W, WD, INDX, INVCONT,DIS,G,HESS,rv1)
END SUBROUTINE ANM1
!********************************************************************************************
!* SUBROUTINES TO SET GAMMA NRES x NRES ARRAY ACCORDING TO MODE...
!********************************************************************************************
SUBROUTINE SET_GAMMA(CUTOFF,NRES,DIS,G)
INTEGER NRES,S,CUTOFF
REAL*4, PARAMETER :: Clin=10_SGL,Ckov=40_SGL,Ccart=6_SGL,Cseq=60_SGL
INTEGER, PARAMETER :: Slim=3,Ex=6
INTEGER:: K,J
REAL(DBL), INTENT(INOUT) :: G(NRES,NRES)
REAL(DBL),INTENT(IN) :: DIS(NRES,NRES)
! WRITE (*,*)'Allocating gamma array..'
G=0.0_DBL
! HESS=0.0
! WRITE (*,*) 'Computing with default values...'
!Here we build g(j,k) terms that are a combination between gamma parameter and kirchhoff elements(-1/0)
DO J=1,nres
DO K=1,nres
IF(J.NE.K) THEN
S=ABS(J-K)
IF(S.LE.Slim) THEN
G(J,K)=-(Cseq/(S**2))
ELSE
IF(DIS(J,K).LE.CUTOFF) THEN
G(J,K)=DBLE(-((Ccart/DIS(J,K))**Ex))
ELSE
G(J,K)=0_DBL
ENDIF
ENDIF
ENDIF
ENDDO
ENDDO
END SUBROUTINE SET_GAMMA
!********************************************************************************************
!* SUBROUTINES TO SET HESSIAN 3*NRES x 3*NRES ARRAY OF 2ND DERIVATIVES OF POTENTIAL
!********************************************************************************************
SUBROUTINE HESSIAN (NRES,R,DIS,G,HESS)
IMPLICIT NONE
INTEGER, INTENT(IN) :: NRES
INTEGER RES3
REAL(DBL) :: BX,BY,BZ,DIS2,GAMMA
REAL(DBL), INTENT(IN) :: R(NRES,NRES)
REAL(DBL), INTENT(IN) :: DIS(NRES,NRES),G(NRES,NRES)
REAL(DBL), INTENT(INOUT) :: HESS(3*NRES,3*NRES)
INTEGER :: J,K
RES3=3*NRES
! WRITE (*,*)'Allocating hessian array, dimension', RES3
HESS=0.0
! write(*,*) 'Calculating the hessian..'
DO J=1,NRES
DO K=1,NRES
bx=0.0
by=0.0
bz=0.0
dis2=0.0
gamma=0.0
IF(J.NE.K)then
BX=(r(J,1)-r(K,1))
BY=(r(J,2)-r(K,2))
BZ=(r(J,3)-r(K,3))
DIS2=DIS(J,K)**2
Gamma=G(J,K)
HESS(3*J-2,3*K-2)=+GAMMA*BX*BX/DIS2
HESS(3*J-2,3*K-1)=+GAMMA*BX*BY/DIS2
HESS(3*J-2,3*K)= +GAMMA*BX*BZ/DIS2
HESS(3*J-1,3*K-2)=+GAMMA*BY*BX/DIS2
HESS(3*J-1,3*K-1)=+GAMMA*BY*BY/DIS2
HESS(3*J-1,3*K)= +GAMMA*BY*BZ/DIS2
HESS(3*J,3*K-2)=+GAMMA*BZ*BX/DIS2
HESS(3*J,3*K-1)=+GAMMA*BZ*BY/DIS2
HESS(3*J,3*K)= +GAMMA*BZ*BZ/DIS2
!C SECOND: CREATION OF DIAGONAL Hii(Hii=-SUM(Hij/j#i) TERMS
HESS(3*J-2,3*J-2)=HESS(3*J-2,3*J-2)-GAMMA*BX*BX/DIS2
HESS(3*J-2,3*J-1)=HESS(3*J-2,3*J-1)-GAMMA*BX*BY/DIS2
HESS(3*J-2,3*J)=HESS(3*J-2,3*J) -GAMMA*BX*BZ/DIS2
HESS(3*J-1,3*J-2)=HESS(3*J-1,3*J-2)-GAMMA*BY*BX/DIS2
HESS(3*J-1,3*J-1)=HESS(3*J-1,3*J-1)-GAMMA*BY*BY/DIS2
HESS(3*J-1,3*J)=HESS(3*J-1,3*J) -GAMMA*BY*BZ/DIS2
HESS(3*J,3*J-2)=HESS(3*J,3*J-2) -GAMMA*BZ*BX/DIS2
HESS(3*J,3*J-1)=HESS(3*J,3*J-1) -GAMMA*BZ*BY/DIS2
HESS(3*J,3*J)=HESS(3*J,3*J) -GAMMA*BZ*BZ/DIS2
ENDIF
ENDDO
ENDDO
! WRITE (*,*) 'Hessian done!!'
! RETURN
END SUBROUTINE HESSIAN
!*******************************************************************************
!* SUBROUTINE DISTANCES - COMPUTES NRES x NRES ARRAY OF CALPHA-CALPHA DISTANCES
!*******************************************************************************
SUBROUTINE DISTANCES (NRES,R,DIS)
IMPLICIT NONE
INTEGER :: NRES,I,J
REAL(DBL) , INTENT(IN) :: R(NRES,3)
REAL(DBL), INTENT(INOUT), DIMENSION(NRES,NRES) :: Dis
REAL(DBL) :: x,y,z
! WRITE(*,*) 'Now we compute c-ALPHA distances...'
DO I=1,NRES
DO J=1,NRES
if(j.ne.i) then
X=(R(I,1)-R(J,1))**2
Y=(R(I,2)-R(J,2))**2
Z=(R(I,3)-R(J,3))**2
DIS(I,J)=DSQRT(X+Y+Z)
endif
ENDDO
ENDDO
END SUBROUTINE DISTANCES
!*******************************************************************************
!* SUBROUTINE INDEXX - ORDERS by eigenvalues
!*******************************************************************************
SUBROUTINE INDEXX(N,ARRIN,INDX)
INTEGER, INTENT(IN) :: N
INTEGER :: J,L,IR,INDXT,I
REAL(DBL), INTENT(IN), DIMENSION(N) :: ARRIN
INTEGER, INTENT(INOUT) :: INDX(N)
REAL(DBL) :: q
DO 11 J=1,N
INDX(J)=J
11 CONTINUE
IF(N.EQ.1)RETURN
L=N/2+1
IR=N
10 CONTINUE
IF(L.GT.1)THEN
L=L-1
INDXT=INDX(L)
Q=ARRIN(INDXT)
ELSE
INDXT=INDX(IR)
Q=ARRIN(INDXT)
INDX(IR)=INDX(1)
IR=IR-1
IF(IR.EQ.1)THEN
INDX(1)=INDXT
RETURN
ENDIF
ENDIF
I=L
J=L+L
20 IF(J.LE.IR)THEN
IF(J.LT.IR)THEN
IF(ARRIN(INDX(J)).LT.ARRIN(INDX(J+1)))J=J+1
ENDIF
IF(Q.LT.ARRIN(INDX(J)))THEN
INDX(I)=INDX(J)
I=J
J=J+J
ELSE
J=IR+1
ENDIF
GO TO 20
ENDIF
INDX(I)=INDXT
GO TO 10
END SUBROUTINE INDEXX
END MODULE ANM
File added
!
! File: DMD.f
! Author: gelpi
!
! Created on 14 de marzo de 2012, 16:54
!
subroutine dmdIntloop(nblist, r, v, nxm, tacact, tpart, ipart, TACT, TMIN, temps, iev, natom, toUpdate)
use intList
use geometryDP
integer, intent(IN) :: natom
real, intent(IN) :: nxm
real*8, intent(IN) :: TACT, TMIN
type(pointDP), intent(INOUT) :: r(natom), v(natom)
type(intpList), intent(INOUT) :: nblist(natom)
integer, intent(INOUT) :: iev
real*8, intent(INOUT) :: tacact, temps
logical toUpdate(natom)
!
real*8 tpart(natom)
real*8 tevent, tevent1
integer ipart(natom)
integer mem1, mem2, npair1, npair2, i,j
!
tacact = 0.
toUpdate = .true.
! evolucio temporal
!------------------------------------------------------------------------------
do while (tacact .lt. TACT)
! busca quina es la propera colisio
do i = 1, natom
if (toUpdate(i)) then
tpart(i) = 1.d15
do j = 1, nblist(i)%nats
if (nblist(i)%iData(j)%timp.lt.tpart(i)) then
tpart(i) = nblist(i)%iData(j)%timp
ipart(i) = j ! ipart recull index NO num d'atom
endif
toUpdate(i) = .false.
enddo
endif
enddo
tevent = 1.e15
call nextCol(mem1, mem2, npair1, npair2, tevent, ipart, tpart, natom, nblist)
!
tevent1 = tevent - temps
temps = tevent
tacact = tacact + tevent1
iev = iev + 1
! translacio i variacio dels temps
do i = 1, natom ! mantenim la versio inline degut a la barreja de tipus de real !!
r(i)%x = r(i)%x + v(i)%x * tevent1
r(i)%y = r(i)%y + v(i)%y * tevent1
r(i)%z = r(i)%z + v(i)%z * tevent1
enddo
call updateV(r, v, nblist(mem1)%iData(npair1)%deltak, nxm, nblist(mem1)%iData(npair1)%xsum, &
mem1, mem2, natom)
! Assegurem que trespassem la barrera en la direcció correcta
r(mem1)%x = r(mem1)%x + v(mem1)%x * tevent1/DBLE(100000.)
r(mem1)%y = r(mem1)%y + v(mem1)%y * tevent1/DBLE(100000.)
r(mem1)%z = r(mem1)%z + v(mem1)%z * tevent1/DBLE(100000.)
r(mem2)%x = r(mem2)%x + v(mem2)%x * tevent1/DBLE(100000.)
r(mem2)%y = r(mem2)%y + v(mem2)%y * tevent1/DBLE(100000.)
r(mem2)%z = r(mem2)%z + v(mem2)%z * tevent1/DBLE(100000.)
!
! ara actualitza els temps de colisio per a les dues particules que han xocat
! anulem la colisio per a que no es repeteixi
nblist(mem1)%iData(npair1)%timp = 1.d15
nblist(mem2)%iData(npair2)%timp = 1.d15
call updateXocPart(mem1, nblist, temps, r, v, TMIN, natom, toUpdate)
call updateXocPart(mem2, nblist, temps, r, v, TMIN, natom, toUpdate)
enddo
! deallocate (tpart, ipart)
! end do while (tacact.lt.tact)------------------------------------------------
end subroutine dmdIntloop
File added
!===============================================================================
subroutine readCommLine (files, unit_i, unit_o, NFILES)
USE commLine
integer, intent(IN) :: NFILES
integer, intent(OUT) :: unit_i, unit_o
type(commLineOption), intent(INOUT) :: files(NFILES)
!
call inputArgs(files)
unit_i = openFn(files, '-i')
if (fileName(files,'-o').ne.'log') then
unit_o = openFn(files, '-o')
else
unit_o = 6
endif
end subroutine readCommLine
!===============================================================================
subroutine programHeader(unit_o)
integer unit_o
#ifdef OMP
integer omp_get_max_threads
#endif
write (unit_o, *) "================================================="
write (unit_o, *) "= ="
write (unit_o, *) "= GOdMD ="
write (unit_o, *) "= ="
write (unit_o, *) "= P Sfriso, A. Emperador, JL. Gelpi, M.Orozco ="
write (unit_o, *) "= ="
write (unit_o, *) "= (c) 2013 ="
write (unit_o, *) "================================================="
write (unit_o, *)
#ifdef OMP
write (unit_o, '(" Running OpenMP version on ", i3, " processors")') omp_get_max_threads()
write (unit_o, *)
#endif
end subroutine programHeader
!===============================================================================
subroutine errorAllocmem (ioerr, text)
integer ioerr
character(len=*) text
if (ioerr.ne.0) then
write (0, '("Error allocating memory",a30)') text
stop 1
endif
end subroutine errorAllocmem
!===============================================================================
subroutine writeSnapshot (unit_traj, ibloc, r, atom, res, chain, rnum, natom)
use geometryDP
integer, intent(IN) :: natom, unit_traj, ibloc
type(pointDP), intent(IN) :: r(natom)
character(len=4), intent(IN) :: atom(natom), res(natom)
character(len=1), intent(IN) :: chain(natom)
character(len=5), intent(IN) :: rnum(natom)
integer i, j
write (unit_traj, '("MODEL",8X,I4)') ibloc
do i = 1,natom
j=i
write (unit_traj, '("ATOM",2X,I5,1X,A4,1X,A3,1X,A1,A5,3X,3F8.3)') i, atom(j), res(j), chain(j), rnum(j), r(i)
enddo
write (unit_traj, '("ENDMDL")')
end subroutine writeSnapshot
!===============================================================================
! subroutine writeRestartFiles(files, r, v, RSTV, natom, NFILES)
! use geometry
! use geometryDP
! use commline
! integer, intent(IN) :: natom, RSTV, NFILES
! type(pointDP), intent (IN) :: r(natom), v(natom)
! type(commLineOption), intent(IN) :: files(NFILES)
!
! type(point), allocatable :: rsp(:) ! Single precision version of r for binary I/O
! integer unit_rst, unit_rstv, i
!
! allocate (rsp(natom))
! do i = 1, natom
! rsp(i) = DPtoSP(r(i))
! enddo
! unit_rst = openFn(files, '-rst')
! write (unit_rst) natom
! write (unit_rst) (rsp(i), i = 1, natom)
! close(unit_rst)
!!
! if (RSTV .eq. 1) then
! unit_rstv = openFn(files, '-rstv')
! write (unit_rstv) natom
! write (unit_rstv) (v(i), i = 1, natom)
! close(unit_rstv)
! endif
! deallocate(rsp)
! !
!end subroutine writeRestartFiles
!===============================================================================
subroutine writeEnergies(unit_ener, unit_o, temps, epotgo, ekin0, natom,fractionDone,error,&
aux_acc,xbeta,saco,discardedEner,wrong,encero,WENER, WLOG,ierr)
! EN MAIN:
! call writeEnergies(unit_ener, unit_o, iev, ierr, temps,&
! epotgo, epot, ekin0, natom,fractionDone,error,aux_acc,xbeta ,.true.,.true.)
!
use constants
integer, intent(IN) :: unit_ener, unit_o, natom, saco,encero
INTEGER, INTENT(IN) :: wrong,ierr
real, intent(IN) :: epotgo, ekin0, xbeta,aux_acc,fractionDone
real*8,Intent(IN):: error
real*8, intent(IN) :: temps
real, INTENT(IN) :: discardedEner
logical, intent(IN) :: WENER, WLOG
real :: Kb= 8.314 ! J/(mol K)
character(len = 92), parameter :: &
fmtEner = '("Time ",f8.2, " Epot",f10.3, " Temp(K)",f8.2 ," Dener ",f12.2," Salto ",I4," Error ",I3)'
character(len = 105), parameter :: fmtLog = &
'(" Time(ps): ",f8.2, " Epot(kcal/mol): ",f10.3, " Done(%): ",f4.2," RMSd(A):",f8.3,f8.2,f5.2,X,I3,X,2I4)'
!
! FALTA Kb
if (WENER) &
write (unit_ener, fmtEner)temps * 1.e-3, epotgo, ekin0/1.5/natom/Kb ,discardedEner,saco,wrong
if (WLOG) &
write (unit_o, fmtLog) temps * 1.e-3, epotgo,fractionDone,error,xbeta,aux_acc,saco,encero,ierr
end subroutine writeEnergies
!===============================================================================
SUBROUTINE writePairList(unit_o,natom,contarBondedInt,stepPts)
use stepPotentials
!
IMPLICIT NONE
!
INTEGER, INTENT(IN) :: unit_o
INTEGER, INTENT(IN) :: natom
INTEGER, INTENT(IN) :: contarBondedInt
TYPE( stepPotInt), INTENT (IN) :: stepPts(natom,natom)
!
write (unit_o, '(" Initial Pairs list")')
write (unit_o, '(" ==================")')
write (unit_o, '(" Total: ",i9)') natom * (natom - 1) / 2
write (unit_o, '(" Bonded: ",i9)') contarBondedInt
write (unit_o, '(" Non Bonded: ",i9)') count(stepPts % active)
write (unit_o, '(" GO 1 SQ WELL ",i9)') count(stepPts % nstep .eq. 2 .and.stepPts%active.and.stepPts%tipInt==SSEC)
write (unit_o, '(" GO 2 SQ WELL ",i9)') count(stepPts % nstep .eq. 4 .and.stepPts%active)
write (unit_o, '(" Fake: ",i9)') count(stepPts % tipInt .eq. 3 .and.stepPts%active)
!
write (unit_o, *)
write (unit_o, '(" System setup completed")')
write (unit_o, *)
END SUBROUTINE writePairList
!=================================================================================
\ No newline at end of file
File added
include ../config.mk
#
OBJS = IO.o random.o energy.o colisio.o DMD.o dims_utils.o least_sq_fit.o ANM.o
MODS = paramSet.o Structure.o doubleDMD.o ls_rmsd.o
EXEFILE = discrete
discrete: ${MODS} $(OBJS) $(LIBDIR) main.o
$(F77) $(LFLAGS) $(ENDFLAG) -o $(EXEFILE) $(LOPTIONS) $(OBJS) $(MODS) main.o -l Discrete
main.o: vars.h
Structure.o: Structure.f
$(F77) -ffree-form $(OPTIONS) $(CFLAGSNOCHECK) Structure.f
.f.o: $<
$(F77) -ffree-form $(OPTIONS) $(CFLAGS) $<
.F.o: $<
$(F77) -ffree-form $(OPTIONS) $(CFLAGS) $<
clean:
rm *.o *.mod $(EXEFILE)
include ../config.mk
#
OBJS = IO.o random.o energy.o colisio.o DMD.o
MODS = paramSet.o
EXEFILE = discrete
discrete: ${MODS} $(OBJS) $(LIBDIR) main.o
$(F77) $(LFLAGS) $(ENDFLAG) -o $(EXEFILE) $(LOPTIONS) $(OBJS) $(MODS) main.o -lDiscrete
main.o: vars.h
.f.o: $<
$(F77) -ffree-form $(OPTIONS) $(CFLAGS) $<
.F.o: $<
$(F77) -ffree-form $(OPTIONS) $(CFLAGS) $<
clean:
rm *.o *.mod $(EXEFILE)
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment