3 C Build the virtual polypeptide chain. Side-chain centroids are moveable.
6 implicit real*8 (a-h,o-z)
12 include 'COMMON.IOUNITS'
13 include 'COMMON.NAMES'
14 include 'COMMON.INTERACT'
15 double precision e1(3),e2(3),e3(3)
16 logical lprn,perbox,fail
17 C Set lprn=.true. for debugging
26 C first three CAs and SC(2).
28 subroutine chainbuild_extconf
30 C Build the virtual polypeptide chain. Side-chain centroids are moveable.
33 implicit real*8 (a-h,o-z)
35 include 'COMMON.CHAIN'
36 include 'COMMON.LOCAL'
39 include 'COMMON.IOUNITS'
40 include 'COMMON.NAMES'
41 include 'COMMON.INTERACT'
42 double precision e1(3),e2(3),e3(3)
43 logical lprn,perbox,fail
46 c write (iout,*) "Calling chainbuild_extconf"
49 * Build the alpha-carbon chain.
52 call locate_next_res(i)
55 C First and last SC must coincide with the corresponding CA.
59 dc_norm(j,nres+1)=0.0D0
61 dc_norm(j,nres+nres)=0.0D0
63 c(j,nres+nres)=c(j,nres)
71 write (iout,'(/a)') 'Recalculated internal coordinates'
74 c(j,maxres2)=0.5D0*(c(j,i-1)+c(j,i+1))
77 if (i.gt.3) be=rad2deg*beta(i-3,i-2,i-1,i)
78 be1=rad2deg*beta(nres+i,i,maxres2,i+1)
80 if (i.gt.2) alfai=rad2deg*alpha(i-2,i-1,i)
81 write (iout,1212) restyp(itype(i)),i,dist(i-1,i),
82 & alfai,be,dist(nres+i,i),rad2deg*alpha(nres+i,i,maxres2),be1
84 1212 format (a3,'(',i3,')',2(f10.5,2f10.2))
91 c-------------------------------------------------------------------------
94 C Define the origin and orientation of the coordinate system and locate
95 C the first three atoms.
97 implicit real*8 (a-h,o-z)
99 include 'COMMON.CHAIN'
100 include 'COMMON.LOCAL'
154 dc_norm(j,2)=prod(j,1,2)
155 dc(j,2)=vbld(3)*prod(j,1,2)
156 c(j,3)=c(j,2)+dc(j,2)
158 call locate_side_chain(2)
161 c-----------------------------------------------------------------------------
162 subroutine locate_next_res(i)
164 C Locate CA(i) and SC(i-1)
166 implicit real*8 (a-h,o-z)
168 include 'COMMON.CHAIN'
169 include 'COMMON.LOCAL'
172 include 'COMMON.IOUNITS'
173 include 'COMMON.NAMES'
174 include 'COMMON.INTERACT'
176 C Define the rotation matrices corresponding to CA(i)
180 if (theti.ne.theti) theti=100.0
182 if (phii.ne.phii) phii=180.0
191 * Define the matrices of the rotation about the virtual-bond valence angles
192 * theta, T(i,j,k), virtual-bond dihedral angles gamma (miscalled PHI in this
193 * program), R(i,j,k), and, the cumulative matrices of rotation RT
215 rt(2,1,i-2)=sint*cosphi
216 rt(2,2,i-2)=-cost*cosphi
218 rt(3,1,i-2)=-sint*sinphi
219 rt(3,2,i-2)=cost*sinphi
221 call matmult(prod(1,1,i-2),rt(1,1,i-2),prod(1,1,i-1))
223 dc_norm(j,i-1)=prod(j,1,i-1)
224 dc(j,i-1)=vbld(i)*prod(j,1,i-1)
225 c(j,i)=c(j,i-1)+dc(j,i-1)
227 cd print '(2i3,2(3f10.5,5x))', i-1,i,(dc(j,i-1),j=1,3),(c(j,i),j=1,3)
229 C Now calculate the coordinates of SC(i-1)
231 call locate_side_chain(i-1)
234 c-----------------------------------------------------------------------------
235 subroutine locate_side_chain(i)
237 C Locate the side-chain centroid i, 1 < i < NRES. Put in C(*,NRES+i).
239 implicit real*8 (a-h,o-z)
241 include 'COMMON.CHAIN'
242 include 'COMMON.LOCAL'
245 include 'COMMON.IOUNITS'
246 include 'COMMON.NAMES'
247 include 'COMMON.INTERACT'
251 c dsci_inv=dsc_inv(itype(i))
253 dsci_inv=vbld_inv(i+nres)
257 if (alphi.ne.alphi) alphi=100.0
258 if (omegi.ne.omegi) omegi=-100.0
268 yp= dsci*sinalphi*cosomegi
269 zp=-dsci*sinalphi*sinomegi
270 * Now we have to rotate the coordinate system by 180-theta(i)/2 so as to get its
271 * X-axis aligned with the vector DC(*,i)
272 theta2=pi-0.5D0*theta(i+1)
275 xx(1)= xp*cost2+yp*sint2
276 xx(2)=-xp*sint2+yp*cost2
278 cd print '(a3,i3,3f10.5,5x,3f10.5)',restyp(itype(i)),i,
279 cd & xp,yp,zp,(xx(k),k=1,3)
283 * Bring the SC vectors to the common coordinate system.
285 xx(2)=xloc(2,i)*r(2,2,i-1)+xloc(3,i)*r(2,3,i-1)
286 xx(3)=xloc(2,i)*r(3,2,i-1)+xloc(3,i)*r(3,3,i-1)
293 rj=rj+prod(j,k,i-1)*xx(k)
296 dc_norm(j,nres+i)=rj*dsci_inv
297 c(j,nres+i)=c(j,i)+rj
301 c------------------------------------------
307 include 'COMMON.CONTROL'
311 include 'COMMON.LANGEVIN'
313 include 'COMMON.LANGEVIN.lang0'
315 include 'COMMON.CHAIN'
316 include 'COMMON.DERIV'
318 include 'COMMON.LOCAL'
319 include 'COMMON.INTERACT'
320 include 'COMMON.IOUNITS'
321 include 'COMMON.NAMES'
322 include 'COMMON.TIME1'
323 include 'COMMON.REMD'
324 include 'COMMON.SETUP'
325 include 'COMMON.MUCA'
326 include 'COMMON.HAIRPIN'
327 C change suggested by Ana - begin
329 C change suggested by Ana - end
333 C write(*,*) 'initial', i,j,c(j,i)
335 C change suggested by Ana - begin
337 C change suggested by Ana -end
339 if ((itype(i).eq.ntyp1).and.(itype(i+1).eq.ntyp1)) then
341 if (allareout.eq.1) then
342 ireturnval=int(c(j,i)/boxxsize)
343 if (c(j,i).le.0) ireturnval=ireturnval-1
344 do k=chain_beg,chain_end
345 c(j,k)=c(j,k)-ireturnval*boxxsize
346 c(j,k+nres)=c(j,k+nres)-ireturnval*boxxsize
350 & dc_old(1,0)=dc_old(1,0)-ireturnval*boxxsize
351 C Suggested by Ana -end
356 if (int(c(j,i)/boxxsize).eq.0) allareout=0
359 if (allareout.eq.1) then
360 ireturnval=int(c(j,i)/boxxsize)
361 if (c(j,i).le.0) ireturnval=ireturnval-1
363 c(j,k)=c(j,k)-ireturnval*boxxsize
364 c(j,k+nres)=c(j,k+nres)-ireturnval*boxxsize
369 C write(*,*) 'befor no jump', i,j,c(j,i)
373 if (itype(i).eq.ntyp1 .and. itype(i-1).eq.ntyp1) then
374 difference=abs(c(j,i-1)-c(j,i))
375 C print *,'diff', difference
376 if (difference.gt.boxxsize/2.0) then
377 if (c(j,i-1).gt.c(j,i)) then
386 c(j,i)=c(j,i)+nojumpval*boxxsize
387 c(j,i+nres)=c(j,i+nres)+nojumpval*boxxsize
391 if (itype(i).eq.ntyp1 .and. itype(i-1).eq.ntyp1) then
392 difference=abs(c(j,i-1)-c(j,i))
393 if (difference.gt.boxxsize/2.0) then
394 if (c(j,i-1).gt.c(j,i)) then
403 c(j,i)=c(j,i)+nojumpval*boxxsize
404 c(j,i+nres)=c(j,i+nres)+nojumpval*boxxsize
408 C write(*,*) 'after no jump', i,j,c(j,i)
412 C suggesed by Ana begins
414 C suggested by Ana ends
418 if ((itype(i).eq.ntyp1).and.(itype(i+1).eq.ntyp1)) then
420 if (allareout.eq.1) then
421 ireturnval=int(c(j,i)/boxysize)
422 if (c(j,i).le.0) ireturnval=ireturnval-1
423 do k=chain_beg,chain_end
424 c(j,k)=c(j,k)-ireturnval*boxysize
425 c(j,k+nres)=c(j,k+nres)-ireturnval*boxysize
429 & dc_old(1,0)=dc_old(1,0)-ireturnval*boxxsize
430 C Suggested by Ana -end
435 if (int(c(j,i)/boxysize).eq.0) allareout=0
438 if (allareout.eq.1) then
439 ireturnval=int(c(j,i)/boxysize)
440 if (c(j,i).le.0) ireturnval=ireturnval-1
442 c(j,k)=c(j,k)-ireturnval*boxysize
443 c(j,k+nres)=c(j,k+nres)-ireturnval*boxysize
448 if (itype(i).eq.ntyp1 .and. itype(i-1).eq.ntyp1) then
449 difference=abs(c(j,i-1)-c(j,i))
450 if (difference.gt.boxysize/2.0) then
451 if (c(j,i-1).gt.c(j,i)) then
460 c(j,i)=c(j,i)+nojumpval*boxysize
461 c(j,i+nres)=c(j,i+nres)+nojumpval*boxysize
465 if (itype(i).eq.ntyp1 .and. itype(i-1).eq.ntyp1) then
466 difference=abs(c(j,i-1)-c(j,i))
467 if (difference.gt.boxysize/2.0) then
468 if (c(j,i-1).gt.c(j,i)) then
477 c(j,i)=c(j,i)+nojumpval*boxysize
478 c(j,i+nres)=c(j,i+nres)+nojumpval*boxysize
481 C Suggested by Ana -begins
483 C Suggested by Ana -ends
487 if ((itype(i).eq.ntyp1).and.(itype(i+1).eq.ntyp1)) then
489 if (allareout.eq.1) then
490 ireturnval=int(c(j,i)/boxysize)
491 if (c(j,i).le.0) ireturnval=ireturnval-1
492 do k=chain_beg,chain_end
493 c(j,k)=c(j,k)-ireturnval*boxzsize
494 c(j,k+nres)=c(j,k+nres)-ireturnval*boxzsize
498 & dc_old(1,0)=dc_old(1,0)-ireturnval*boxxsize
499 C Suggested by Ana -end
504 if (int(c(j,i)/boxzsize).eq.0) allareout=0
507 if (allareout.eq.1) then
508 ireturnval=int(c(j,i)/boxzsize)
509 if (c(j,i).le.0) ireturnval=ireturnval-1
511 c(j,k)=c(j,k)-ireturnval*boxzsize
512 c(j,k+nres)=c(j,k+nres)-ireturnval*boxzsize
517 if (itype(i).eq.ntyp1 .and. itype(i-1).eq.ntyp1) then
518 difference=abs(c(j,i-1)-c(j,i))
519 if (difference.gt.(boxzsize/2.0)) then
520 if (c(j,i-1).gt.c(j,i)) then
529 c(j,i)=c(j,i)+nojumpval*boxzsize
530 c(j,i+nres)=c(j,i+nres)+nojumpval*boxzsize
534 if (itype(i).eq.ntyp1 .and. itype(i-1).eq.ntyp1) then
535 difference=abs(c(j,i-1)-c(j,i))
536 if (difference.gt.boxzsize/2.0) then
537 if (c(j,i-1).gt.c(j,i)) then
546 c(j,i)=c(j,i)+nojumpval*boxzsize
547 c(j,i+nres)=c(j,i+nres)+nojumpval*boxzsize