source: bldmol.f

Last change on this file was 38d77eb, checked in by baerbaer <baerbaer@…>, 14 years ago

Redirected standard out to logString.

SMMP produced a lot of log messages. This became an issue when run in massively
parallel environments. I replaced all writes to standard out to a write to logString.
The next step is to pass this string to a function that writes the messages to a log
file according to their importance and the chosen log level.

git-svn-id: svn+ssh://svn.berlios.de/svnroot/repos/smmp/trunk@34 26dc1dd8-5c4e-0410-9ffe-d298b4865968

  • Property mode set to 100644
File size: 9.7 KB
Line 
1! **************************************************************
2!
3! This file contains the subroutines: bldmol, fnd3ba,eyring,
4! setsys,setgbl
5!
6! Copyright 2003-2005 Frank Eisenmenger, U.H.E. Hansmann,
7! Shura Hayryan, Chin-Ku
8! Copyright 2007 Frank Eisenmenger, U.H.E. Hansmann,
9! Jan H. Meinke, Sandipan Mohanty
10!
11! **************************************************************
12
13 subroutine bldmol(nml)
14
15! .................................................
16! PURPOSE: calculate coordinates for molecule 'nml'
17!
18! OUTPUT: xat,yat,zat,xbaat,ybaat,zbaat,xtoat,ytoat,
19! ztoat (via 'eyring')
20!
21! 1st backbone atom of 1st residue of 'nml':
22!
23! - it's position: from 'gbpr(1-3,nml)'
24! - it's axes: from 'setgbl' according to
25! global angles 'gbpr(4-5,nml)'
26!
27! CALLS: eyring, fnd3ba,setgbl,setsys
28! .................................................
29
30 include 'INCL.H'
31 ! Subroutine arguments
32 integer nml
33
34 integer i, i1, i2, i3, ifirs, jow, j, jj
35 double precision x1, x2, x3, y1, y2, y3, z1, z2, z3
36
37 double precision xg, zg
38 dimension xg(3),zg(3)
39
40
41 call fnd3ba(nml,i1,i2,i3)
42! ------------------------------ first 3 bb atoms of 'nml'
43 ixrfpt(1,nml)=i1
44 ixrfpt(2,nml)=i2
45 ixrfpt(3,nml)=i3
46
47! ------------------------------ position of 1st bb atom
48 xat(i1) = gbpr(1,nml)
49 yat(i1) = gbpr(2,nml)
50 zat(i1) = gbpr(3,nml)
51
52 rfpt(1,nml)=xat(i1)
53 rfpt(2,nml)=yat(i1)
54 rfpt(3,nml)=zat(i1)
55
56 call setgbl(nml,i1,i2,i3, xg, zg)
57
58 xbaat(i1) = zg(1)
59 ybaat(i1) = zg(2)
60 zbaat(i1) = zg(3)
61
62 xtoat(i1) = xg(1)
63 ytoat(i1) = xg(2)
64 ztoat(i1) = xg(3)
65
66 ifirs=irsml1(nml)
67
68 do i=ifirs,irsml2(nml)
69 if (i.eq.ifirs) then ! not construct
70 jj=iatrs1(i)+1 ! 1st bb atom again
71 else
72 jj=iatrs1(i)
73 endif
74 do j=jj,iatrs2(i)
75 jow=iowat(j)
76 call eyring(j,jow)
77 enddo
78 enddo
79
80 call setsys(i1,i2,i3, x1,x2,x3,y1,y2,y3,z1,z2,z3)
81
82 xrfax(1,nml)=x1 ! J
83 xrfax(2,nml)=x2
84 xrfax(3,nml)=x3
85 yrfax(1,nml)=-y1 ! K
86 yrfax(2,nml)=-y2
87 yrfax(3,nml)=-y3
88 zrfax(1,nml)=-z1 ! L
89 zrfax(2,nml)=-z2
90 zrfax(3,nml)=-z3
91
92 return
93 end
94! ***********************************
95 subroutine fnd3ba(nml,i1,i2,i3)
96
97! .................................................
98! PURPOSE: return indices 'i1,i2,i3' of
99! first 3 backbone atoms in molecule 'nml'
100!
101! CALLS: fndbrn
102! .................................................
103
104 include 'INCL.H'
105! arguments
106 integer nml, i1, i2, i3
107
108 integer ibd, i, irs, ix
109 dimension ibd(4)
110 logical bb
111
112
113 irs=irsml1(nml)
114
115! --------------------------- 1st bb atom
116 i1=iatrs1(irs)
117
118 call fndbrn(nml,irs,i1,i,ix,i2,bb)
119
120! --------------------------- 2nd bb atom
121 i2=i+1
122
123! ------------------------ check for ring
124
125 ibd(1)=iowat(i1)
126 ibd(2)=ibdat(1,i1)
127 ibd(3)=ibdat(2,i1)
128 ibd(4)=ibdat(3,i1)
129
130 ix = 0
131 do i=1,nbdat(i1)+1
132
133 if (iowat(ibd(i)).ne.i1) then
134
135 if (ix.ne.0) then
136 write (logString, '(2a,i3)')
137 & ' fnd3ba> Can handle only simple ring at 1st',
138 & ' atom of molecule #',nml
139 stop
140 endif
141
142 ix=ibd(i)
143 endif
144 enddo
145
146! --------------------------- 3rd bb atom
147
148 ix=ixatrs(irs)
149
150 if (i2.eq.ix) then
151
152 if (irs.lt.irsml2(nml)) then
153 i3 = iatrs1(irs+1)
154 return
155 endif
156
157 else
158
159 do i = ix,i2+1,-1
160 if (iowat(i).eq.i2) then
161 i3=i
162 return
163 endif
164 enddo
165
166 endif
167
168 write (logString, '(4a,i4,a,i4)')
169 & ' fnd3ba> Cannot find backbone atom following ',nmat(i2),
170 & ' of residue ',seq(irs),irs,' in molecule #',nml
171 stop
172
173 end
174! ***************************
175 subroutine eyring(i,ia)
176
177! .........................................................
178! PURPOSE: calculate cartesian coordinates of atom 'i'
179! using EYRING's algorithm
180! INPUT: i - index of atom to be constructed
181! for 'i': blat,csbaat,snbaat,cstoat,sntoat
182! ia- index of atom from which 'i' is to be built
183! OUTPUT: for 'i': xat,yat,zat,xbaat,ybaat,zbaat,xtoat,ytoat,ztoat
184!
185! CALLS: none
186! .........................................................
187
188 include 'INCL.H'
189! arguments
190 integer i, ia
191 double precision ct, st, ca, sa, bl, x1, x2, x3, z1, z2, z3
192 double precision y1, y2, y3, h2, h3, dx, dz
193
194 ct=cstoat(i)
195 st=sntoat(i)
196 ca=csbaat(i)
197 sa=snbaat(i)
198 bl=blat(i)
199
200 x1=xtoat(ia)
201 x2=ytoat(ia)
202 x3=ztoat(ia)
203
204 z1=xbaat(ia)
205 z2=ybaat(ia)
206 z3=zbaat(ia)
207
208 y1=z2*x3-z3*x2
209 y2=z3*x1-z1*x3
210 y3=z1*x2-z2*x1
211
212 h2=-sa*ct
213 h3=-sa*st
214
215 x1=-ca*x1+h2*y1+h3*z1
216 x2=-ca*x2+h2*y2+h3*z2
217 x3=-ca*x3+h2*y3+h3*z3
218
219 dx=one/sqrt(x1*x1+x2*x2+x3*x3)
220 x1=x1*dx
221 x2=x2*dx
222 x3=x3*dx
223
224 xtoat(i)=x1
225 ytoat(i)=x2
226 ztoat(i)=x3
227
228 z1=-st*y1+ct*z1
229 z2=-st*y2+ct*z2
230 z3=-st*y3+ct*z3
231
232 dz=one/sqrt(z1*z1+z2*z2+z3*z3)
233 xbaat(i)=z1*dz
234 ybaat(i)=z2*dz
235 zbaat(i)=z3*dz
236
237 xat(i)=xat(ia)+x1*bl
238 yat(i)=yat(ia)+x2*bl
239 zat(i)=zat(ia)+x3*bl
240
241 return
242 end
243! ***********************************************************
244 subroutine setsys(i1,i2,i3, x1,x2,x3,y1,y2,y3,z1,z2,z3)
245
246! ..........................................................
247! PURPOSE: calculate axes X,Y,Z of right-handed orthogonal
248! system given by three atom positions R1, R2, R3
249!
250! X = (R2-R1)/ |( )|
251! Z = {X x (R2-R3)} / |{ }|
252! Y = Z x X
253!
254! INPUT: i1, i2, i3 - indices of three atoms
255! OUTPUT: x1,x2,x3 |
256! y1,y2,y3 | -direction cosines of X,Y,Z
257! z1,z2,z3 |
258!
259! CALLS: none
260! ...................................................
261
262
263 include 'INCL.H'
264! arguments
265 integer i1, i2, i3
266 double precision x1, x2, x3, y1, y2, y3, z1, z2, z3, h1, h2, h3
267 double precision dz, dx
268
269 h1=xat(i2)
270 h2=yat(i2)
271 h3=zat(i2)
272
273 x1=h1-xat(i1)
274 x2=h2-yat(i1)
275 x3=h3-zat(i1)
276
277 y1=h1-xat(i3)
278 y2=h2-yat(i3)
279 y3=h3-zat(i3)
280
281 z1=x2*y3-x3*y2
282 z2=x3*y1-x1*y3
283 z3=x1*y2-x2*y1
284
285 dx=one/sqrt(x1*x1+x2*x2+x3*x3)
286 x1=x1*dx
287 x2=x2*dx
288 x3=x3*dx
289
290 dz=one/sqrt(z1*z1+z2*z2+z3*z3)
291 z1=z1*dz
292 z2=z2*dz
293 z3=z3*dz
294
295 y1=z2*x3-z3*x2
296 y2=z3*x1-z1*x3
297 y3=z1*x2-z2*x1
298
299 return
300 end
301
302! *****************************************
303 subroutine setgbl(nml,i1,i2,i3,xg,zg)
304
305! ...................................................
306!
307! PURPOSE: 1. Obtain global axes (J,K,L)
308! related to x(1 0 0),y(0 1 0),z(0 0 1)
309! by 3 rotations (gbl. parameters #4-#6):
310!
311! - round z by angle alpha
312! - round x' by a. beta
313! - round y" by a. gamma
314!
315! 2. Return x-axis (xg) & z-axis (zg)
316! for atom #1 in order to orientate J
317! along the bond from backbone atom #1
318! to bb.a. #2 and L according to the
319! cross product [ bond(#1->#2) x
320! bond(#2->#3) ] when using Eyring's
321! algorithm to get the coordinates
322!
323! CALLS: none
324! ..............................................
325
326 include 'INCL.H'
327! arguments
328 integer nml, i1, i2, i3
329 double precision xg, zg, ag
330 dimension xg(3),zg(3),ag(3,3)
331
332 integer i
333 double precision ca, sa, cb,sb, cg, sg, d, ct2, ct3, dx, dz
334 double precision x1, x2, x3, y1, y2, y3, z1, z2, z3, st2, sa2, st3
335
336 ca = cos(gbpr(4,nml)) ! alpha
337 sa = sin(gbpr(4,nml))
338 cb = cos(gbpr(5,nml)) ! beta
339 sb = sin(gbpr(5,nml))
340 cg = cos(gbpr(6,nml)) ! gamma
341 sg = sin(gbpr(6,nml))
342
343! ----------------------------- J
344 x1 = ca*cg - sa*sb*sg
345 x2 = sa*cg + ca*sb*sg
346 x3 = -cb*sg
347
348 d = sqrt(x1**2+x2**2+x3**2)
349 ag(1,1) = x1/d
350 ag(2,1) = x2/d
351 ag(3,1) = x3/d
352! ----------------------------- K
353 y1 = -sa*cb
354 y2 = ca*cb
355 y3 = sb
356
357 d = sqrt(y1**2+y2**2+y3**2)
358 ag(1,2) = y1/d
359 ag(2,2) = y2/d
360 ag(3,2) = y3/d
361! ----------------------------- L
362 z1 = ca*sg + sa*sb*cg
363 z2 = sa*sg - ca*sb*cg
364 z3 = cb*cg
365
366 d = sqrt(z1**2+z2**2+z3**2)
367 ag(1,3) = z1/d
368 ag(2,3) = z2/d
369 ag(3,3) = z3/d
370
371! ------------------------------------ X1
372 ct2 = cstoat(i2)
373 st2 = sntoat(i2)
374 sa2 = snbaat(i2)
375
376 x1 = -csbaat(i2)
377 x2 = -sa2*ct2
378 x3 = -sa2*st2
379
380 dx = sqrt(x1**2+x2**2+x3**2)
381 x1 = x1/dx
382 x2 = x2/dx
383 x3 = x3/dx
384! ------------------------------------- Z1
385 st3 = sntoat(i3)
386 ct3 = cstoat(i3)
387
388 z1 = -st3*(st2*x3+ct2*x2)
389 z2 = st3*ct2*x1 + ct3*st2
390 z3 = st3*st2*x1 - ct3*ct2
391
392 dz = sqrt(z1**2+z2**2+z3**2)
393 z1 = z1/dz
394 z2 = z2/dz
395 z3 = z3/dz
396! ------------------------------------- Y1
397 y1 = z2 * x3 - z3 * x2
398 y3 = z1 * x2 - z2 * x1 ! do not need y2
399
400! ----------------------------- into global system
401
402 xg(1) = ag(1,1)*x1 + ag(1,2)*y1 + ag(1,3)*z1
403 xg(2) = ag(2,1)*x1 + ag(2,2)*y1 + ag(2,3)*z1
404 xg(3) = ag(3,1)*x1 + ag(3,2)*y1 + ag(3,3)*z1
405
406 dx = sqrt(xg(1)**2+xg(2)**2+xg(3)**2)
407
408 zg(1) = ag(1,1)*x3 + ag(1,2)*y3 + ag(1,3)*z3
409 zg(2) = ag(2,1)*x3 + ag(2,2)*y3 + ag(2,3)*z3
410 zg(3) = ag(3,1)*x3 + ag(3,2)*y3 + ag(3,3)*z3
411
412 dz = sqrt(zg(1)**2+zg(2)**2+zg(3)**2)
413
414 do i=1,3
415 xg(i) = xg(i)/dx
416 zg(i) = zg(i)/dz
417 enddo
418
419 return
420 end
421
Note: See TracBrowser for help on using the repository browser.