source/unres/src-HCD-5D/geomout.F

   1       subroutine pdbout(etot,tytul,iunit)
   2       implicit real*8 (a-h,o-z)
   3       include 'DIMENSIONS'
   4       include 'COMMON.CHAIN'
   5       include 'COMMON.INTERACT'
   6       include 'COMMON.NAMES'
   7       include 'COMMON.IOUNITS'
   8       include 'COMMON.HEADER'
   9       include 'COMMON.SBRIDGE'
  10       include 'COMMON.DISTFIT'
  11       include 'COMMON.MD'
  12       character*50 tytul
  13       character*1 chainid(10) /'A','B','C','D','E','F','G','H','I','J'/
  14       dimension ica(maxres)
  15       write (iunit,'(3a,1pe15.5)') 'REMARK ',tytul,' ENERGY ',etot
  16 cmodel      write (iunit,'(a5,i6)') 'MODEL',1
  17       if (nhfrag.gt.0) then
  18        do j=1,nhfrag
  19         iti=itype(hfrag(1,j))
  20         itj=itype(hfrag(2,j))
  21         if (j.lt.10) then
  22            write (iunit,'(a5,i5,1x,a1,i1,2x,a3,i7,2x,a3,i7,i3,t76,i5)')
  23      &           'HELIX',j,'H',j,
  24      &           restyp(iti),hfrag(1,j)-1,
  25      &           restyp(itj),hfrag(2,j)-1,1,hfrag(2,j)-hfrag(1,j)
  26         else
  27              write (iunit,'(a5,i5,1x,a1,i2,1x,a3,i7,2x,a3,i7,i3)')
  28      &           'HELIX',j,'H',j,
  29      &           restyp(iti),hfrag(1,j)-1,
  30      &           restyp(itj),hfrag(2,j)-1,1,hfrag(2,j)-hfrag(1,j)
  31         endif
  32        enddo
  33       endif
  34
  35       if (nbfrag.gt.0) then
  36
  37        do j=1,nbfrag
  38
  39         iti=itype(bfrag(1,j))
  40         itj=itype(bfrag(2,j)-1)
  41
  42         write (iunit,'(a5,i5,1x,a1,i1,i3,1x,a3,i6,2x,a3,i6,i3)')
  43      &           'SHEET',1,'B',j,2,
  44      &           restyp(iti),bfrag(1,j)-1,
  45      &           restyp(itj),bfrag(2,j)-2,0
  46
  47         if (bfrag(3,j).gt.bfrag(4,j)) then
  48
  49          itk=itype(bfrag(3,j))
  50          itl=itype(bfrag(4,j)+1)
  51
  52          write (iunit,'(a5,i5,1x,a1,i1,i3,1x,a3,i6,2x,a3,i6,i3,
  53      &              2x,a1,2x,a3,i6,3x,a1,2x,a3,i6)')
  54      &           'SHEET',2,'B',j,2,
  55      &           restyp(itl),bfrag(4,j),
  56      &           restyp(itk),bfrag(3,j)-1,-1,
  57      &           "N",restyp(itk),bfrag(3,j)-1,
  58      &           "O",restyp(iti),bfrag(1,j)-1
  59
  60         else
  61
  62          itk=itype(bfrag(3,j))
  63          itl=itype(bfrag(4,j)-1)
  64
  65
  66         write (iunit,'(a5,i5,1x,a1,i1,i3,1x,a3,i6,2x,a3,i6,i3,
  67      &              2x,a1,2x,a3,i6,3x,a1,2x,a3,i6)')
  68      &           'SHEET',2,'B',j,2,
  69      &           restyp(itk),bfrag(3,j)-1,
  70      &           restyp(itl),bfrag(4,j)-2,1,
  71      &           "N",restyp(itk),bfrag(3,j)-1,
  72      &           "O",restyp(iti),bfrag(1,j)-1
  73
  74
  75
  76         endif
  77
  78        enddo
  79       endif
  80
  81       if (nss.gt.0) then
  82         do i=1,nss
  83          if (dyn_ss) then
  84           write(iunit,'(a6,i4,1x,a3,i7,4x,a3,i7)')
  85      &         'SSBOND',i,'CYS',idssb(i)-nnt+1,
  86      &                    'CYS',jdssb(i)-nnt+1
  87          else
  88           write(iunit,'(a6,i4,1x,a3,i7,4x,a3,i7)')
  89      &         'SSBOND',i,'CYS',ihpb(i)-nnt+1-nres,
  90      &                    'CYS',jhpb(i)-nnt+1-nres
  91          endif
  92         enddo
  93       endif
  94
  95       iatom=0
  96       ichain=1
  97       ires=0
  98       do i=nnt,nct
  99         iti=itype(i)
 100         if ((iti.eq.ntyp1).and.((itype(i+1)).eq.ntyp1)) then
 101           ichain=ichain+1
 102           ires=0
 103           write (iunit,'(a)') 'TER'
 104         else
 105         ires=ires+1
 106         iatom=iatom+1
 107         ica(i)=iatom
 108         if (iti.ne.ntyp1) then
 109         write (iunit,10) iatom,restyp(iti),chainid(ichain),
 110      &     ires,(c(j,i),j=1,3),vtot(i)
 111         if (iti.ne.10) then
 112           iatom=iatom+1
 113           write (iunit,20) iatom,restyp(iti),chainid(ichain),
 114      &      ires,(c(j,nres+i),j=1,3),
 115      &      vtot(i+nres)
 116          endif
 117         endif
 118         endif
 119       enddo
 120       write (iunit,'(a)') 'TER'
 121       do i=nnt,nct-1
 122         if (itype(i).eq.ntyp1) cycle
 123         if (itype(i).eq.10 .and. itype(i+1).ne.ntyp1) then
 124           write (iunit,30) ica(i),ica(i+1)
 125         else if (itype(i).ne.10 .and. itype(i+1).ne.ntyp1) then
 126           write (iunit,30) ica(i),ica(i+1),ica(i)+1
 127         else if (itype(i).ne.10 .and. itype(i+1).eq.ntyp1) then
 128           write (iunit,30) ica(i),ica(i)+1
 129         endif
 130       enddo
 131       if (itype(nct).ne.10) then
 132         write (iunit,30) ica(nct),ica(nct)+1
 133       endif
 134       do i=1,nss
 135        if (dyn_ss) then
 136         write (iunit,30) ica(idssb(i))+1,ica(jdssb(i))+1
 137        else
 138         write (iunit,30) ica(ihpb(i)-nres)+1,ica(jhpb(i)-nres)+1
 139        endif
 140       enddo
 141       write (iunit,'(a6)') 'ENDMDL'
 142   10  FORMAT ('ATOM',I7,'  CA  ',A3,1X,A1,I4,4X,3F8.3,f15.3)
 143   20  FORMAT ('ATOM',I7,'  CB  ',A3,1X,A1,I4,4X,3F8.3,f15.3)
 144   30  FORMAT ('CONECT',8I5)
 145       return
 146       end
 147 c------------------------------------------------------------------------------
 148       subroutine MOL2out(etot,tytul)
 149 C Prints the Cartesian coordinates of the alpha-carbons in the Tripos mol2
 150 C format.
 151       implicit real*8 (a-h,o-z)
 152       include 'DIMENSIONS'
 153       include 'COMMON.CHAIN'
 154       include 'COMMON.INTERACT'
 155       include 'COMMON.NAMES'
 156       include 'COMMON.IOUNITS'
 157       include 'COMMON.HEADER'
 158       include 'COMMON.SBRIDGE'
 159       character*32 tytul,fd
 160       character*3 zahl
 161       character*6 res_num,pom,ucase
 162 #ifdef AIX
 163       call fdate_(fd)
 164 #elif (defined CRAY)
 165       call date(fd)
 166 #else
 167       call fdate(fd)
 168 #endif
 169       write (imol2,'(a)') '#'
 170       write (imol2,'(a)')
 171      & '#         Creating user name:           unres'
 172       write (imol2,'(2a)') '#         Creation time:                ',
 173      & fd
 174       write (imol2,'(/a)') '\@<TRIPOS>MOLECULE'
 175       write (imol2,'(a)') tytul
 176       write (imol2,'(5i5)') nct-nnt+1,nct-nnt+nss+1,nct-nnt+nss+1,0,0
 177       write (imol2,'(a)') 'SMALL'
 178       write (imol2,'(a)') 'USER_CHARGES'
 179       write (imol2,'(a)') '\@<TRIPOS>ATOM'
 180       do i=nnt,nct
 181         write (zahl,'(i3)') i
 182         pom=ucase(restyp(itype(i)))
 183         res_num = pom(:3)//zahl(2:)
 184         write (imol2,10) i,(c(j,i),j=1,3),i,res_num,0.0
 185       enddo
 186       write (imol2,'(a)') '\@<TRIPOS>BOND'
 187       do i=nnt,nct-1
 188         write (imol2,'(i5,2i6,i2)') i-nnt+1,i-nnt+1,i-nnt+2,1
 189       enddo
 190       do i=1,nss
 191         write (imol2,'(i5,2i6,i2)') nct-nnt+i,ihpb(i),jhpb(i),1
 192       enddo
 193       write (imol2,'(a)') '\@<TRIPOS>SUBSTRUCTURE'
 194       do i=nnt,nct
 195         write (zahl,'(i3)') i
 196         pom = ucase(restyp(itype(i)))
 197         res_num = pom(:3)//zahl(2:)
 198         write (imol2,30) i-nnt+1,res_num,i-nnt+1,0
 199       enddo
 200   10  FORMAT (I7,' CA      ',3F10.4,' C.3',I8,1X,A,F11.4,' ****')
 201   30  FORMAT (I7,1x,A,I14,' RESIDUE',I13,' ****  ****')
 202       return
 203       end
 204 c------------------------------------------------------------------------
 205       subroutine intout
 206       implicit real*8 (a-h,o-z)
 207       include 'DIMENSIONS'
 208       include 'COMMON.IOUNITS'
 209       include 'COMMON.CHAIN'
 210       include 'COMMON.VAR'
 211       include 'COMMON.LOCAL'
 212       include 'COMMON.INTERACT'
 213       include 'COMMON.NAMES'
 214       include 'COMMON.GEO'
 215       include 'COMMON.TORSION'
 216       write (iout,'(/a)') 'Geometry of the virtual chain.'
 217       write (iout,'(7a)') '  Res  ','         d','     Theta',
 218      & '       Phi','       Dsc','     Alpha','      Omega'
 219       do i=1,nres
 220         iti=itype(i)
 221         write (iout,'(a3,i4,6f10.3)') restyp(iti),i,vbld(i),
 222      &     rad2deg*theta(i),rad2deg*phi(i),vbld(nres+i),rad2deg*alph(i),
 223      &     rad2deg*omeg(i)
 224       enddo
 225       return
 226       end
 227 c---------------------------------------------------------------------------
 228       subroutine briefout(it,ener)
 229       implicit real*8 (a-h,o-z)
 230       include 'DIMENSIONS'
 231       include 'COMMON.IOUNITS'
 232       include 'COMMON.CHAIN'
 233       include 'COMMON.VAR'
 234       include 'COMMON.LOCAL'
 235       include 'COMMON.INTERACT'
 236       include 'COMMON.NAMES'
 237       include 'COMMON.GEO'
 238       include 'COMMON.SBRIDGE'
 239 c     print '(a,i5)',intname,igeom
 240 #if defined(AIX) || defined(PGI) || defined(CRAY)
 241       open (igeom,file=intname,position='append')
 242 #else
 243       open (igeom,file=intname,access='append')
 244 #endif
 245       IF (NSS.LE.9) THEN
 246         WRITE (igeom,180) IT,ENER,NSS,(IHPB(I),JHPB(I),I=1,NSS)
 247       ELSE
 248         WRITE (igeom,180) IT,ENER,NSS,(IHPB(I),JHPB(I),I=1,9)
 249         WRITE (igeom,190) (IHPB(I),JHPB(I),I=10,NSS)
 250       ENDIF
 251 c     IF (nvar.gt.nphi) WRITE (igeom,200) (RAD2DEG*THETA(I),I=3,NRES)
 252       WRITE (igeom,200) (RAD2DEG*THETA(I),I=3,NRES)
 253       WRITE (igeom,200) (RAD2DEG*PHI(I),I=4,NRES)
 254 c     if (nvar.gt.nphi+ntheta) then
 255         write (igeom,200) (rad2deg*alph(i),i=2,nres-1)
 256         write (igeom,200) (rad2deg*omeg(i),i=2,nres-1)
 257 c     endif
 258       close(igeom)
 259   180 format (I5,F12.3,I2,9(1X,2I3))
 260   190 format (3X,11(1X,2I3))
 261   200 format (8F10.4)
 262       return
 263       end
 264 #ifdef WINIFL
 265       subroutine fdate(fd)
 266       character*32 fd
 267       write(fd,'(32x)')
 268       return
 269       end
 270 #endif
 271 c----------------------------------------------------------------
 272 #ifdef NOXDR
 273       subroutine cartout(time)
 274 #else
 275       subroutine cartoutx(time)
 276 #endif
 277       implicit real*8 (a-h,o-z)
 278       include 'DIMENSIONS'
 279       include 'COMMON.CHAIN'
 280       include 'COMMON.INTERACT'
 281       include 'COMMON.NAMES'
 282       include 'COMMON.IOUNITS'
 283       include 'COMMON.HEADER'
 284       include 'COMMON.SBRIDGE'
 285       include 'COMMON.DISTFIT'
 286       include 'COMMON.MD'
 287       double precision time
 288 #if defined(AIX) || defined(PGI) || defined(CRAY)
 289       open(icart,file=cartname,position="append")
 290 #else
 291       open(icart,file=cartname,access="append")
 292 #endif
 293       write (icart,'(e15.8,2e15.5,f12.5,$)') time,potE,uconst,t_bath
 294       if (dyn_ss) then
 295        write (icart,'(i4,$)')
 296      &   nss,(idssb(j)+nres,jdssb(j)+nres,j=1,nss)
 297       else
 298        write (icart,'(i4,$)')
 299      &   nss,(ihpb(j),jhpb(j),j=1,nss)
 300        endif
 301        write (icart,'(i4,20f7.4)') nfrag+npair+3*nfrag_back,
 302      & (qfrag(i),i=1,nfrag),(qpair(i),i=1,npair),
 303      & (utheta(i),ugamma(i),uscdiff(i),i=1,nfrag_back)
 304       write (icart,'(8f10.5)')
 305      & ((c(k,j),k=1,3),j=1,nres),
 306      & ((c(k,j+nres),k=1,3),j=nnt,nct)
 307       close(icart)
 308       return
 309       end
 310 c-----------------------------------------------------------------
 311 #ifndef NOXDR
 312       subroutine cartout(time)
 313       implicit real*8 (a-h,o-z)
 314       include 'DIMENSIONS'
 315 #ifdef MPI
 316       include 'mpif.h'
 317       include 'COMMON.SETUP'
 318 #else
 319       parameter (me=0)
 320 #endif
 321       include 'COMMON.CHAIN'
 322       include 'COMMON.INTERACT'
 323       include 'COMMON.NAMES'
 324       include 'COMMON.IOUNITS'
 325       include 'COMMON.HEADER'
 326       include 'COMMON.SBRIDGE'
 327       include 'COMMON.DISTFIT'
 328       include 'COMMON.MD'
 329       double precision time
 330       integer iret,itmp
 331       real xcoord(3,maxres2+2),prec
 332
 333 #ifdef AIX
 334       call xdrfopen_(ixdrf,cartname, "a", iret)
 335       call xdrffloat_(ixdrf, real(time), iret)
 336       call xdrffloat_(ixdrf, real(potE), iret)
 337       call xdrffloat_(ixdrf, real(uconst), iret)
 338       call xdrffloat_(ixdrf, real(uconst_back), iret)
 339       call xdrffloat_(ixdrf, real(t_bath), iret)
 340       call xdrfint_(ixdrf, nss, iret)
 341       do j=1,nss
 342        if (dyn_ss) then
 343         call xdrfint_(ixdrf, idssb(j)+nres, iret)
 344         call xdrfint_(ixdrf, jdssb(j)+nres, iret)
 345        else
 346         call xdrfint_(ixdrf, ihpb(j), iret)
 347         call xdrfint_(ixdrf, jhpb(j), iret)
 348        endif
 349       enddo
 350       call xdrfint_(ixdrf, nfrag+npair+3*nfrag_back, iret)
 351       do i=1,nfrag
 352         call xdrffloat_(ixdrf, real(qfrag(i)), iret)
 353       enddo
 354       do i=1,npair
 355         call xdrffloat_(ixdrf, real(qpair(i)), iret)
 356       enddo
 357       do i=1,nfrag_back
 358         call xdrffloat_(ixdrf, real(utheta(i)), iret)
 359         call xdrffloat_(ixdrf, real(ugamma(i)), iret)
 360         call xdrffloat_(ixdrf, real(uscdiff(i)), iret)
 361       enddo
 362 #else
 363       call xdrfopen(ixdrf,cartname, "a", iret)
 364 c      write (iout,*) "Writing conformation: time",time," potE",potE,
 365 c     & " uconst",uconst," uconst_back",uconst_back," t_bath",t_bath,
 366 c     & " nss",nss
 367       call xdrffloat(ixdrf, real(time), iret)
 368       call xdrffloat(ixdrf, real(potE), iret)
 369       call xdrffloat(ixdrf, real(uconst), iret)
 370       call xdrffloat(ixdrf, real(uconst_back), iret)
 371       call xdrffloat(ixdrf, real(t_bath), iret)
 372       call xdrfint(ixdrf, nss, iret)
 373       do j=1,nss
 374        if (dyn_ss) then
 375         call xdrfint(ixdrf, idssb(j)+nres, iret)
 376         call xdrfint(ixdrf, jdssb(j)+nres, iret)
 377        else
 378         call xdrfint(ixdrf, ihpb(j), iret)
 379         call xdrfint(ixdrf, jhpb(j), iret)
 380        endif
 381       enddo
 382       call xdrfint(ixdrf, nfrag+npair+3*nfrag_back, iret)
 383       do i=1,nfrag
 384         call xdrffloat(ixdrf, real(qfrag(i)), iret)
 385       enddo
 386       do i=1,npair
 387         call xdrffloat(ixdrf, real(qpair(i)), iret)
 388       enddo
 389       do i=1,nfrag_back
 390         call xdrffloat(ixdrf, real(utheta(i)), iret)
 391         call xdrffloat(ixdrf, real(ugamma(i)), iret)
 392         call xdrffloat(ixdrf, real(uscdiff(i)), iret)
 393       enddo
 394 #endif
 395       prec=10000.0
 396       do i=1,nres
 397        do j=1,3
 398         xcoord(j,i)=c(j,i)
 399        enddo
 400       enddo
 401       do i=nnt,nct
 402        do j=1,3
 403         xcoord(j,nres+i-nnt+1)=c(j,i+nres)
 404        enddo
 405       enddo
 406
 407       itmp=nres+nct-nnt+1
 408 #ifdef AIX
 409       call xdrf3dfcoord_(ixdrf, xcoord, itmp, prec, iret)
 410       call xdrfclose_(ixdrf, iret)
 411 #else
 412       call xdrf3dfcoord(ixdrf, xcoord, itmp, prec, iret)
 413       call xdrfclose(ixdrf, iret)
 414 #endif
 415       return
 416       end
 417 #endif
 418 c-----------------------------------------------------------------
 419       subroutine statout(itime)
 420       implicit real*8 (a-h,o-z)
 421       include 'DIMENSIONS'
 422       include 'COMMON.CONTROL'
 423       include 'COMMON.CHAIN'
 424       include 'COMMON.INTERACT'
 425       include 'COMMON.NAMES'
 426       include 'COMMON.IOUNITS'
 427       include 'COMMON.HEADER'
 428       include 'COMMON.SBRIDGE'
 429       include 'COMMON.DISTFIT'
 430       include 'COMMON.MD'
 431       include 'COMMON.REMD'
 432       include 'COMMON.SETUP'
 433       integer itime
 434       double precision energia(0:n_ene)
 435       double precision gyrate
 436       external gyrate
 437       common /gucio/ cm
 438       character*256 line1,line2
 439       character*4 format1,format2
 440       character*30 format
 441 #ifdef AIX
 442       if(itime.eq.0) then
 443        open(istat,file=statname,position="append")
 444       endif
 445 #else
 446 #if defined(PGI) || defined(CRAY)
 447       open(istat,file=statname,position="append")
 448 #else
 449       open(istat,file=statname,access="append")
 450 #endif
 451 #endif
 452        if (AFMlog.gt.0) then
 453          if (refstr) then
 454          call rms_nac_nnc(rms,frac,frac_nn,co,.false.)
 455           write (line1,'(i10,f15.2,3f12.3,f7.2,2f6.3,4f12.3,i5,$)')
 456      &          itime,totT,EK,potE+potEcomp(27),totE+potEcomp(27),
 457      &          rms,frac,frac_nn,kinetic_T,t_bath,gyrate(),
 458      &          potEcomp(23),me
 459           format1="a133"
 460          else
 461 C          print *,'A CHUJ',potEcomp(23)
 462           write (line1,'(i10,f15.2,7f12.3,i5,$)')
 463      &           itime,totT,EK,potE+potEcomp(27),totE+potEcomp(27),
 464      &           kinetic_T,t_bath,gyrate(),
 465      &           potEcomp(23),me
 466           format1="a114"
 467          endif
 468        else if (selfguide.gt.0) then
 469        distance=0.0
 470        do j=1,3
 471        distance=distance+(c(j,afmend)-c(j,afmbeg))**2
 472        enddo
 473        distance=dsqrt(distance)
 474         if (refstr) then
 475          call rms_nac_nnc(rms,frac,frac_nn,co,.false.)
 476           write (line1,'(i10,f15.2,3f12.3,f7.2,2f6.3,f12.3,f10.1,2f8.2,
 477      &    f9.3,i5,$)')
 478      &          itime,totT,EK,potE+potEcomp(27),totE+potEcomp(27),
 479      &          rms,frac,frac_nn,kinetic_T,t_bath,gyrate(),
 480      &          distance,potEcomp(23),me
 481           format1="a133"
 482 C          print *,"CHUJOWO"
 483         else
 484 C          print *,'A CHUJ',potEcomp(23)
 485           write (line1,'(i10,f15.2,8f12.3,i5,$)')
 486      &           itime,totT,EK,potE+potEcomp(27),totE+potEcomp(27),
 487      &           kinetic_T,t_bath,gyrate(),
 488      &           distance,potEcomp(23),me
 489           format1="a114"
 490         endif
 491        else
 492         if (refstr) then
 493          call rms_nac_nnc(rms,frac,frac_nn,co,.false.)
 494           write (line1,'(i10,f15.2,3f12.3,f7.2,4f6.3,3f12.3,i5,$)')
 495      &          itime,totT,EK,potE+potEcomp(27),totE+potEcomp(27),
 496      &          rms,frac,frac_nn,co,amax,kinetic_T,t_bath,gyrate(),me
 497           format1="a133"
 498         else
 499           write (line1,'(i10,f15.2,7f12.3,i5,$)')
 500      &           itime,totT,EK,potE+potEcomp(27),totE+potEcomp(27),
 501      &           amax,kinetic_T,t_bath,gyrate(),me
 502           format1="a114"
 503         endif
 504        endif
 505        if(usampl.and.totT.gt.eq_time) then
 506            if (loc_qlike) then
 507            write(line2,'(i5,2f9.4,300f7.4)') iset,uconst,uconst_back,
 508      &      (qfrag(ii1),ii1=1,nfrag),(qpair(ii2),ii2=1,npair),
 509      &      (utheta(i),ugamma(i),uscdiff(i),i=1,nfrag_back),
 510      &      ((qloc(j,i),j=1,3),i=1,nfrag_back)
 511            write(format2,'(a1,i3.3)') "a",23+7*nfrag+7*npair
 512      &             +42*nfrag_back
 513            else
 514            write(line2,'(i5,2f9.4,300f7.4)') iset,uconst,uconst_back,
 515      &      (qfrag(ii1),ii1=1,nfrag),(qpair(ii2),ii2=1,npair),
 516      &      (utheta(i),ugamma(i),uscdiff(i),i=1,nfrag_back)
 517            write(format2,'(a1,i3.3)') "a",23+7*nfrag+7*npair
 518      &             +21*nfrag_back
 519            endif
 520        else
 521            format2="a001"
 522            line2=' '
 523        endif
 524        if (print_compon) then
 525           if(itime.eq.0) then
 526            write(format,'(a1,a4,a1,a4,a10)') "(",format1,",",format2,
 527      &                                                     ",100a12)"
 528            write (istat,format) "#"," ",
 529      &      (ename(print_order(i)),i=1,nprint_ene)
 530           endif
 531           write(format,'(a1,a4,a1,a4,a10)') "(",format1,",",format2,
 532      &                                                     ",100f12.3)"
 533           write (istat,format) line1,line2,
 534      &      (potEcomp(print_order(i)),i=1,nprint_ene)
 535        else
 536           write(format,'(a1,a4,a1,a4,a1)') "(",format1,",",format2,")"
 537           write (istat,format) line1,line2
 538        endif
 539 #if defined(AIX)
 540         call flush(istat)
 541 #else
 542         close(istat)
 543 #endif
 544        return
 545       end
 546 c---------------------------------------------------------------
 547       double precision function gyrate()
 548       implicit real*8 (a-h,o-z)
 549       include 'DIMENSIONS'
 550       include 'COMMON.INTERACT'
 551       include 'COMMON.CHAIN'
 552       double precision cen(3),rg
 553
 554       do j=1,3
 555        cen(j)=0.0d0
 556       enddo
 557
 558       ii=0
 559       do i=nnt,nct
 560         if (itype(i).eq.ntyp1) cycle
 561         ii=ii+1
 562         do j=1,3
 563           cen(j)=cen(j)+c(j,i)
 564         enddo
 565       enddo
 566       do j=1,3
 567         cen(j)=cen(j)/dble(ii)
 568       enddo
 569       rg = 0.0d0
 570       do i = nnt, nct
 571         if (itype(i).eq.ntyp1) cycle
 572         do j=1,3
 573          rg = rg + (c(j,i)-cen(j))**2
 574         enddo
 575       end do
 576       gyrate = dsqrt(rg/dble(ii))
 577       return
 578       end