source/unres/src_MD-M/chainbuild.F

   1       subroutine chainbuild
   2 C
   3 C Build the virtual polypeptide chain. Side-chain centroids are moveable.
   4 C As of 2/17/95.
   5 C
   6       implicit real*8 (a-h,o-z)
   7       include 'DIMENSIONS'
   8       include 'COMMON.CHAIN'
   9       include 'COMMON.LOCAL'
  10       include 'COMMON.GEO'
  11       include 'COMMON.VAR'
  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
  18       lprn = .false.
  19       perbox=.false.
  20       fail=.false.
  21       call chainbuild_cart
  22       return
  23       end
  24 C#ifdef DEBUG
  25 C      if (perbox) then
  26 C first three CAs and SC(2).
  27 C
  28       subroutine chainbuild_extconf
  29 C
  30 C Build the virtual polypeptide chain. Side-chain centroids are moveable.
  31 C As of 2/17/95.
  32 C
  33       implicit real*8 (a-h,o-z)
  34       include 'DIMENSIONS'
  35       include 'COMMON.CHAIN'
  36       include 'COMMON.LOCAL'
  37       include 'COMMON.GEO'
  38       include 'COMMON.VAR'
  39       include 'COMMON.IOUNITS'
  40       include 'COMMON.NAMES'
  41       include 'COMMON.INTERACT'
  42       double precision e1(3),e2(3),e3(3)
  43       logical lprn /.false./,perbox,fail
  44
  45 c      write (iout,*) "Calling chainbuild_extconf"
  46       call orig_frame
  47 *
  48 * Build the alpha-carbon chain.
  49 *
  50       do i=4,nres
  51         call locate_next_res(i)
  52       enddo
  53 C
  54 C First and last SC must coincide with the corresponding CA.
  55 C
  56       do j=1,3
  57         dc(j,nres+1)=0.0D0
  58         dc_norm(j,nres+1)=0.0D0
  59         dc(j,nres+nres)=0.0D0
  60         dc_norm(j,nres+nres)=0.0D0
  61         c(j,nres+1)=c(j,1)
  62         c(j,nres+nres)=c(j,nres)
  63       enddo
  64 *
  65 * Temporary diagnosis
  66 *
  67       if (lprn) then
  68
  69       call cartprint
  70       write (iout,*) 'dc_norm'
  71       do i=1,nres
  72         write (iout,'(a3,1x,i3,3f10.5,5x,3f10.5)')
  73      &   restyp(itype(i)),i,(dc_norm(j,i),j=1,3),
  74      &   (dc_norm(j,i+nres),j=1,3)
  75       enddo
  76       write (iout,'(/a)') 'Recalculated internal coordinates'
  77       do i=2,nres-1
  78         do j=1,3
  79           c(j,maxres2)=0.5D0*(c(j,i-1)+c(j,i+1))
  80         enddo
  81         be=0.0D0
  82         if (i.gt.3) be=rad2deg*beta(i-3,i-2,i-1,i)
  83         be1=rad2deg*beta(nres+i,i,maxres2,i+1)
  84         alfai=0.0D0
  85         if (i.gt.2) alfai=rad2deg*alpha(i-2,i-1,i)
  86         write (iout,1212) restyp(itype(i)),i,dist(i-1,i),
  87      &  alfai,be,dist(nres+i,i),rad2deg*alpha(nres+i,i,maxres2),be1
  88       enddo
  89  1212 format (a3,'(',i3,')',2(f10.5,2f10.2))
  90
  91 C      endif
  92       endif
  93       return
  94       end
  95 C#endif
  96 c-------------------------------------------------------------------------
  97       subroutine orig_frame
  98 C
  99 C Define the origin and orientation of the coordinate system and locate
 100 C the first three atoms.
 101 C
 102       implicit real*8 (a-h,o-z)
 103       include 'DIMENSIONS'
 104       include 'COMMON.CHAIN'
 105       include 'COMMON.LOCAL'
 106       include 'COMMON.GEO'
 107       include 'COMMON.VAR'
 108       cost=dcos(theta(3))
 109       sint=dsin(theta(3))
 110       t(1,1,1)=-cost
 111       t(1,2,1)=-sint
 112       t(1,3,1)= 0.0D0
 113       t(2,1,1)=-sint
 114       t(2,2,1)= cost
 115       t(2,3,1)= 0.0D0
 116       t(3,1,1)= 0.0D0
 117       t(3,2,1)= 0.0D0
 118       t(3,3,1)= 1.0D0
 119       r(1,1,1)= 1.0D0
 120       r(1,2,1)= 0.0D0
 121       r(1,3,1)= 0.0D0
 122       r(2,1,1)= 0.0D0
 123       r(2,2,1)= 1.0D0
 124       r(2,3,1)= 0.0D0
 125       r(3,1,1)= 0.0D0
 126       r(3,2,1)= 0.0D0
 127       r(3,3,1)= 1.0D0
 128       do i=1,3
 129         do j=1,3
 130           rt(i,j,1)=t(i,j,1)
 131         enddo
 132       enddo
 133       do i=1,3
 134         do j=1,3
 135           prod(i,j,1)=0.0D0
 136           prod(i,j,2)=t(i,j,1)
 137         enddo
 138         prod(i,i,1)=1.0D0
 139       enddo
 140       c(1,1)=0.0D0
 141       c(2,1)=0.0D0
 142       c(3,1)=0.0D0
 143       c(1,2)=vbld(2)
 144       c(2,2)=0.0D0
 145       c(3,2)=0.0D0
 146       dc(1,0)=0.0d0
 147       dc(2,0)=0.0D0
 148       dc(3,0)=0.0D0
 149       dc(1,1)=vbld(2)
 150       dc(2,1)=0.0D0
 151       dc(3,1)=0.0D0
 152       dc_norm(1,0)=0.0D0
 153       dc_norm(2,0)=0.0D0
 154       dc_norm(3,0)=0.0D0
 155       dc_norm(1,1)=1.0D0
 156       dc_norm(2,1)=0.0D0
 157       dc_norm(3,1)=0.0D0
 158       do j=1,3
 159         dc_norm(j,2)=prod(j,1,2)
 160         dc(j,2)=vbld(3)*prod(j,1,2)
 161         c(j,3)=c(j,2)+dc(j,2)
 162       enddo
 163       call locate_side_chain(2)
 164       return
 165       end
 166 c-----------------------------------------------------------------------------
 167       subroutine locate_prev_res(i)
 168 C
 169 C Locate CA(i) and SC(i-1)
 170 C
 171       implicit real*8 (a-h,o-z)
 172       include 'DIMENSIONS'
 173       include 'COMMON.CHAIN'
 174       include 'COMMON.LOCAL'
 175       include 'COMMON.GEO'
 176       include 'COMMON.VAR'
 177       include 'COMMON.IOUNITS'
 178       include 'COMMON.NAMES'
 179       include 'COMMON.INTERACT'
 180 C
 181 C Define the rotation matrices corresponding to CA(i)
 182 C
 183 #ifdef OSF
 184       theti=theta(i+2)
 185       if (theti.ne.theti) theti=100.0
 186       phii=phi(i+3)
 187       if (phii.ne.phii) phii=180.0
 188 #else
 189       theti=theta(i+2)
 190       phii=phi(i+3)
 191 #endif
 192       cost=dcos(theti)
 193       sint=dsin(theti)
 194       cosphi=dcos(phii)
 195       sinphi=dsin(phii)
 196 * Define the matrices of the rotation about the virtual-bond valence angles
 197 * theta, T(i,j,k), virtual-bond dihedral angles gamma (miscalled PHI in this
 198 * program), R(i,j,k), and, the cumulative matrices of rotation RT
 199       t(1,1,i+2)=-cost
 200       t(1,2,i+2)=-sint
 201       t(1,3,i+2)= 0.0D0
 202       t(2,1,i+2)=-sint
 203       t(2,2,i+2)= cost
 204       t(2,3,i+2)= 0.0D0
 205       t(3,1,i+2)= 0.0D0
 206       t(3,2,i+2)= 0.0D0
 207       t(3,3,i+2)= 1.0D0
 208       r(1,1,i+2)= 1.0D0
 209       r(1,2,i+2)= 0.0D0
 210       r(1,3,i+2)= 0.0D0
 211       r(2,1,i+2)= 0.0D0
 212       r(2,2,i+2)=-cosphi
 213       r(2,3,i+2)= sinphi
 214       r(3,1,i+2)= 0.0D0
 215       r(3,2,i+2)= sinphi
 216       r(3,3,i+2)= cosphi
 217       rt(1,1,i+2)=-cost
 218       rt(1,2,i+2)=-sint
 219       rt(1,3,i+2)=0.0D0
 220       rt(2,1,i+2)=sint*cosphi
 221       rt(2,2,i+2)=-cost*cosphi
 222       rt(2,3,i+2)=sinphi
 223       rt(3,1,i+2)=-sint*sinphi
 224       rt(3,2,i+2)=cost*sinphi
 225       rt(3,3,i+2)=cosphi
 226       call matmult(prod(1,1,i+2),rt(1,1,i+2),prod(1,1,i+1))
 227       do j=1,3
 228         dc_norm(j,i)=-prod(j,1,i+1)
 229         dc(j,i)=-vbld(i)*prod(j,1,i+1)
 230         c(j,i)=c(j,i+1)+dc(j,i)
 231       enddo
 232 cd    print '(2i3,2(3f10.5,5x))', i-1,i,(dc(j,i+1),j=1,3),(c(j,i),j=1,3)
 233 C
 234 C Now calculate the coordinates of SC(i+1)
 235 C
 236       call locate_side_chain(i+1)
 237       return
 238       end
 239 c-----------------------------------------------------------------------------
 240       subroutine locate_next_res(i)
 241 C
 242 C Locate CA(i) and SC(i-1)
 243 C
 244       implicit real*8 (a-h,o-z)
 245       include 'DIMENSIONS'
 246       include 'COMMON.CHAIN'
 247       include 'COMMON.LOCAL'
 248       include 'COMMON.GEO'
 249       include 'COMMON.VAR'
 250       include 'COMMON.IOUNITS'
 251       include 'COMMON.NAMES'
 252       include 'COMMON.INTERACT'
 253 C
 254 C Define the rotation matrices corresponding to CA(i)
 255 C
 256 #ifdef OSF
 257       theti=theta(i)
 258       if (theti.ne.theti) theti=100.0
 259       phii=phi(i)
 260       if (phii.ne.phii) phii=180.0
 261 #else
 262       theti=theta(i)
 263       phii=phi(i)
 264 #endif
 265       cost=dcos(theti)
 266       sint=dsin(theti)
 267       cosphi=dcos(phii)
 268       sinphi=dsin(phii)
 269 * Define the matrices of the rotation about the virtual-bond valence angles
 270 * theta, T(i,j,k), virtual-bond dihedral angles gamma (miscalled PHI in this
 271 * program), R(i,j,k), and, the cumulative matrices of rotation RT
 272       t(1,1,i-2)=-cost
 273       t(1,2,i-2)=-sint
 274       t(1,3,i-2)= 0.0D0
 275       t(2,1,i-2)=-sint
 276       t(2,2,i-2)= cost
 277       t(2,3,i-2)= 0.0D0
 278       t(3,1,i-2)= 0.0D0
 279       t(3,2,i-2)= 0.0D0
 280       t(3,3,i-2)= 1.0D0
 281       r(1,1,i-2)= 1.0D0
 282       r(1,2,i-2)= 0.0D0
 283       r(1,3,i-2)= 0.0D0
 284       r(2,1,i-2)= 0.0D0
 285       r(2,2,i-2)=-cosphi
 286       r(2,3,i-2)= sinphi
 287       r(3,1,i-2)= 0.0D0
 288       r(3,2,i-2)= sinphi
 289       r(3,3,i-2)= cosphi
 290       rt(1,1,i-2)=-cost
 291       rt(1,2,i-2)=-sint
 292       rt(1,3,i-2)=0.0D0
 293       rt(2,1,i-2)=sint*cosphi
 294       rt(2,2,i-2)=-cost*cosphi
 295       rt(2,3,i-2)=sinphi
 296       rt(3,1,i-2)=-sint*sinphi
 297       rt(3,2,i-2)=cost*sinphi
 298       rt(3,3,i-2)=cosphi
 299       call matmult(prod(1,1,i-2),rt(1,1,i-2),prod(1,1,i-1))
 300       do j=1,3
 301         dc_norm(j,i-1)=prod(j,1,i-1)
 302         dc(j,i-1)=vbld(i)*prod(j,1,i-1)
 303         c(j,i)=c(j,i-1)+dc(j,i-1)
 304       enddo
 305 cd    print '(2i3,2(3f10.5,5x))', i-1,i,(dc(j,i-1),j=1,3),(c(j,i),j=1,3)
 306 C
 307 C Now calculate the coordinates of SC(i-1)
 308 C
 309       call locate_side_chain(i-1)
 310       return
 311       end
 312 c-----------------------------------------------------------------------------
 313       subroutine locate_side_chain(i)
 314 C
 315 C Locate the side-chain centroid i, 1 < i < NRES. Put in C(*,NRES+i).
 316 C
 317       implicit real*8 (a-h,o-z)
 318       include 'DIMENSIONS'
 319       include 'COMMON.CHAIN'
 320       include 'COMMON.LOCAL'
 321       include 'COMMON.GEO'
 322       include 'COMMON.VAR'
 323       include 'COMMON.IOUNITS'
 324       include 'COMMON.NAMES'
 325       include 'COMMON.INTERACT'
 326       dimension xx(3)
 327
 328 c      dsci=dsc(itype(i))
 329 c      dsci_inv=dsc_inv(itype(i))
 330       dsci=vbld(i+nres)
 331       dsci_inv=vbld_inv(i+nres)
 332 #ifdef OSF
 333       alphi=alph(i)
 334       omegi=omeg(i)
 335       if (alphi.ne.alphi) alphi=100.0
 336       if (omegi.ne.omegi) omegi=-100.0
 337 #else
 338       alphi=alph(i)
 339       omegi=omeg(i)
 340 #endif
 341       cosalphi=dcos(alphi)
 342       sinalphi=dsin(alphi)
 343       cosomegi=dcos(omegi)
 344       sinomegi=dsin(omegi)
 345       xp= dsci*cosalphi
 346       yp= dsci*sinalphi*cosomegi
 347       zp=-dsci*sinalphi*sinomegi
 348 * Now we have to rotate the coordinate system by 180-theta(i)/2 so as to get its
 349 * X-axis aligned with the vector DC(*,i)
 350       theta2=pi-0.5D0*theta(i+1)
 351       cost2=dcos(theta2)
 352       sint2=dsin(theta2)
 353       xx(1)= xp*cost2+yp*sint2
 354       xx(2)=-xp*sint2+yp*cost2
 355       xx(3)= zp
 356 cd    print '(a3,i3,3f10.5,5x,3f10.5)',restyp(itype(i)),i,
 357 cd   &   xp,yp,zp,(xx(k),k=1,3)
 358       do j=1,3
 359         xloc(j,i)=xx(j)
 360       enddo
 361 * Bring the SC vectors to the common coordinate system.
 362       xx(1)=xloc(1,i)
 363       xx(2)=xloc(2,i)*r(2,2,i-1)+xloc(3,i)*r(2,3,i-1)
 364       xx(3)=xloc(2,i)*r(3,2,i-1)+xloc(3,i)*r(3,3,i-1)
 365       do j=1,3
 366         xrot(j,i)=xx(j)
 367       enddo
 368       do j=1,3
 369         rj=0.0D0
 370         do k=1,3
 371           rj=rj+prod(j,k,i-1)*xx(k)
 372         enddo
 373         dc(j,nres+i)=rj
 374         dc_norm(j,nres+i)=rj*dsci_inv
 375         c(j,nres+i)=c(j,i)+rj
 376       enddo
 377       return
 378       end
 379 c------------------------------------------
 380       subroutine returnbox
 381       include 'DIMENSIONS'
 382 #ifdef MPI
 383       include 'mpif.h'
 384 #endif
 385       include 'COMMON.CONTROL'
 386       include 'COMMON.VAR'
 387       include 'COMMON.MD'
 388 #ifndef LANG0
 389       include 'COMMON.LANGEVIN'
 390 #else
 391       include 'COMMON.LANGEVIN.lang0'
 392 #endif
 393       include 'COMMON.CHAIN'
 394       include 'COMMON.DERIV'
 395       include 'COMMON.GEO'
 396       include 'COMMON.LOCAL'
 397       include 'COMMON.INTERACT'
 398       include 'COMMON.IOUNITS'
 399       include 'COMMON.NAMES'
 400       include 'COMMON.TIME1'
 401       include 'COMMON.REMD'
 402       include 'COMMON.SETUP'
 403       include 'COMMON.MUCA'
 404       include 'COMMON.HAIRPIN'
 405 C change suggested by Ana - begin
 406       integer allareout
 407 C change suggested by Ana - end
 408         j=1
 409         chain_beg=1
 410 C        do i=1,nres
 411 C       write(*,*) 'initial', i,j,c(j,i)
 412 C        enddo
 413 C change suggested by Ana - begin
 414         allareout=1
 415 C change suggested by Ana -end
 416         do i=1,nres-1
 417          if ((itype(i).eq.ntyp1).and.(itype(i+1).eq.ntyp1)) then
 418           chain_end=i
 419           if (allareout.eq.1) then
 420             ireturnval=int(c(j,i)/boxxsize)
 421             if (c(j,i).le.0) ireturnval=ireturnval-1
 422             do k=chain_beg,chain_end
 423               c(j,k)=c(j,k)-ireturnval*boxxsize
 424               c(j,k+nres)=c(j,k+nres)-ireturnval*boxxsize
 425             enddo
 426 C Suggested by Ana
 427             if (chain_beg.eq.1)
 428      &      dc_old(1,0)=dc_old(1,0)-ireturnval*boxxsize
 429 C Suggested by Ana -end
 430            endif
 431            chain_beg=i+1
 432            allareout=1
 433          else
 434           if (int(c(j,i)/boxxsize).eq.0) allareout=0
 435          endif
 436         enddo
 437          if (allareout.eq.1) then
 438             ireturnval=int(c(j,i)/boxxsize)
 439             if (c(j,i).le.0) ireturnval=ireturnval-1
 440             do k=chain_beg,nres
 441               c(j,k)=c(j,k)-ireturnval*boxxsize
 442               c(j,k+nres)=c(j,k+nres)-ireturnval*boxxsize
 443             enddo
 444           endif
 445 C NO JUMP
 446 C        do i=1,nres
 447 C        write(*,*) 'befor no jump', i,j,c(j,i)
 448 C        enddo
 449         nojumpval=0
 450         do i=2,nres
 451            if (itype(i).eq.ntyp1 .and. itype(i-1).eq.ntyp1) then
 452              difference=abs(c(j,i-1)-c(j,i))
 453 C             print *,'diff', difference
 454              if (difference.gt.boxxsize/2.0) then
 455                 if (c(j,i-1).gt.c(j,i)) then
 456                   nojumpval=1
 457                  else
 458                    nojumpval=-1
 459                  endif
 460               else
 461               nojumpval=0
 462               endif
 463               endif
 464               c(j,i)=c(j,i)+nojumpval*boxxsize
 465               c(j,i+nres)=c(j,i+nres)+nojumpval*boxxsize
 466          enddo
 467        nojumpval=0
 468         do i=2,nres
 469            if (itype(i).eq.ntyp1 .and. itype(i-1).eq.ntyp1) then
 470              difference=abs(c(j,i-1)-c(j,i))
 471              if (difference.gt.boxxsize/2.0) then
 472                 if (c(j,i-1).gt.c(j,i)) then
 473                   nojumpval=1
 474                  else
 475                    nojumpval=-1
 476                  endif
 477               else
 478               nojumpval=0
 479               endif
 480              endif
 481               c(j,i)=c(j,i)+nojumpval*boxxsize
 482               c(j,i+nres)=c(j,i+nres)+nojumpval*boxxsize
 483          enddo
 484
 485 C        do i=1,nres
 486 C        write(*,*) 'after no jump', i,j,c(j,i)
 487 C        enddo
 488
 489 C NOW Y dimension
 490 C suggesed by Ana begins
 491         allareout=1
 492 C suggested by Ana ends
 493         j=2
 494         chain_beg=1
 495         do i=1,nres-1
 496          if ((itype(i).eq.ntyp1).and.(itype(i+1).eq.ntyp1)) then
 497           chain_end=i
 498           if (allareout.eq.1) then
 499             ireturnval=int(c(j,i)/boxysize)
 500             if (c(j,i).le.0) ireturnval=ireturnval-1
 501             do k=chain_beg,chain_end
 502               c(j,k)=c(j,k)-ireturnval*boxysize
 503              c(j,k+nres)=c(j,k+nres)-ireturnval*boxysize
 504             enddo
 505 C Suggested by Ana
 506             if (chain_beg.eq.1)
 507      &      dc_old(1,0)=dc_old(1,0)-ireturnval*boxxsize
 508 C Suggested by Ana -end
 509            endif
 510            chain_beg=i+1
 511            allareout=1
 512          else
 513           if (int(c(j,i)/boxysize).eq.0) allareout=0
 514          endif
 515         enddo
 516          if (allareout.eq.1) then
 517             ireturnval=int(c(j,i)/boxysize)
 518             if (c(j,i).le.0) ireturnval=ireturnval-1
 519             do k=chain_beg,nres
 520               c(j,k)=c(j,k)-ireturnval*boxysize
 521               c(j,k+nres)=c(j,k+nres)-ireturnval*boxysize
 522             enddo
 523           endif
 524         nojumpval=0
 525         do i=2,nres
 526            if (itype(i).eq.ntyp1 .and. itype(i-1).eq.ntyp1) then
 527              difference=abs(c(j,i-1)-c(j,i))
 528              if (difference.gt.boxysize/2.0) then
 529                 if (c(j,i-1).gt.c(j,i)) then
 530                   nojumpval=1
 531                  else
 532                    nojumpval=-1
 533                  endif
 534              else
 535               nojumpval=0
 536               endif
 537            endif
 538               c(j,i)=c(j,i)+nojumpval*boxysize
 539               c(j,i+nres)=c(j,i+nres)+nojumpval*boxysize
 540          enddo
 541       nojumpval=0
 542         do i=2,nres
 543            if (itype(i).eq.ntyp1 .and. itype(i-1).eq.ntyp1) then
 544              difference=abs(c(j,i-1)-c(j,i))
 545              if (difference.gt.boxysize/2.0) then
 546                 if (c(j,i-1).gt.c(j,i)) then
 547                   nojumpval=1
 548                  else
 549                    nojumpval=-1
 550                  endif
 551               else
 552               nojumpval=0
 553               endif
 554             endif
 555               c(j,i)=c(j,i)+nojumpval*boxysize
 556               c(j,i+nres)=c(j,i+nres)+nojumpval*boxysize
 557          enddo
 558 C Now Z dimension
 559 C Suggested by Ana -begins
 560         allareout=1
 561 C Suggested by Ana -ends
 562        j=3
 563         chain_beg=1
 564         do i=1,nres-1
 565          if ((itype(i).eq.ntyp1).and.(itype(i+1).eq.ntyp1)) then
 566           chain_end=i
 567           if (allareout.eq.1) then
 568             ireturnval=int(c(j,i)/boxysize)
 569             if (c(j,i).le.0) ireturnval=ireturnval-1
 570             do k=chain_beg,chain_end
 571               c(j,k)=c(j,k)-ireturnval*boxzsize
 572               c(j,k+nres)=c(j,k+nres)-ireturnval*boxzsize
 573             enddo
 574 C Suggested by Ana
 575             if (chain_beg.eq.1)
 576      &      dc_old(1,0)=dc_old(1,0)-ireturnval*boxxsize
 577 C Suggested by Ana -end
 578            endif
 579            chain_beg=i+1
 580            allareout=1
 581          else
 582           if (int(c(j,i)/boxzsize).eq.0) allareout=0
 583          endif
 584         enddo
 585          if (allareout.eq.1) then
 586             ireturnval=int(c(j,i)/boxzsize)
 587             if (c(j,i).le.0) ireturnval=ireturnval-1
 588             do k=chain_beg,nres
 589               c(j,k)=c(j,k)-ireturnval*boxzsize
 590               c(j,k+nres)=c(j,k+nres)-ireturnval*boxzsize
 591             enddo
 592           endif
 593         nojumpval=0
 594         do i=2,nres
 595            if (itype(i).eq.ntyp1 .and. itype(i-1).eq.ntyp1) then
 596              difference=abs(c(j,i-1)-c(j,i))
 597              if (difference.gt.(boxzsize/2.0)) then
 598                 if (c(j,i-1).gt.c(j,i)) then
 599                   nojumpval=1
 600                  else
 601                    nojumpval=-1
 602                  endif
 603               else
 604               nojumpval=0
 605               endif
 606             endif
 607               c(j,i)=c(j,i)+nojumpval*boxzsize
 608               c(j,i+nres)=c(j,i+nres)+nojumpval*boxzsize
 609          enddo
 610        nojumpval=0
 611         do i=2,nres
 612            if (itype(i).eq.ntyp1 .and. itype(i-1).eq.ntyp1) then
 613              difference=abs(c(j,i-1)-c(j,i))
 614              if (difference.gt.boxzsize/2.0) then
 615                 if (c(j,i-1).gt.c(j,i)) then
 616                   nojumpval=1
 617                  else
 618                    nojumpval=-1
 619                  endif
 620               else
 621               nojumpval=0
 622               endif
 623             endif
 624              c(j,i)=c(j,i)+nojumpval*boxzsize
 625               c(j,i+nres)=c(j,i+nres)+nojumpval*boxzsize
 626          enddo
 627
 628         return
 629         end
 630