source/unres/src_MD-M-newcorr/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           write(iunit,'(a6,i4,1x,a3,i7,4x,a3,i7)')
  84      &         'SSBOND',i,'CYS',ihpb(i)-1-nres,
  85      &                    'CYS',jhpb(i)-1-nres
  86         enddo
  87       endif
  88
  89       iatom=0
  90       ichain=1
  91       ires=0
  92       do i=nnt,nct
  93         iti=itype(i)
  94         if (iti.eq.ntyp1) then
  95           ichain=ichain+1
  96           ires=0
  97           write (iunit,'(a)') 'TER'
  98         else
  99         ires=ires+1
 100         iatom=iatom+1
 101         ica(i)=iatom
 102         write (iunit,10) iatom,restyp(iti),chainid(ichain),
 103      &     ires,(c(j,i),j=1,3),vtot(i)
 104         if (iti.ne.10) then
 105           iatom=iatom+1
 106           write (iunit,20) iatom,restyp(iti),chainid(ichain),
 107      &      ires,(c(j,nres+i),j=1,3),
 108      &      vtot(i+nres)
 109         endif
 110         endif
 111       enddo
 112       write (iunit,'(a)') 'TER'
 113       do i=nnt,nct-1
 114         if (itype(i).eq.ntyp1) cycle
 115         if (itype(i).eq.10 .and. itype(i+1).ne.ntyp1) then
 116           write (iunit,30) ica(i),ica(i+1)
 117         else if (itype(i).ne.10 .and. itype(i+1).ne.ntyp1) then
 118           write (iunit,30) ica(i),ica(i+1),ica(i)+1
 119         else if (itype(i).ne.10 .and. itype(i+1).eq.ntyp1) then
 120           write (iunit,30) ica(i),ica(i)+1
 121         endif
 122       enddo
 123       if (itype(nct).ne.10) then
 124         write (iunit,30) ica(nct),ica(nct)+1
 125       endif
 126       do i=1,nss
 127         write (iunit,30) ica(ihpb(i)-nres)+1,ica(jhpb(i)-nres)+1
 128       enddo
 129       write (iunit,'(a6)') 'ENDMDL'
 130   10  FORMAT ('ATOM',I7,'  CA  ',A3,1X,A1,I4,4X,3F8.3,f15.3)
 131   20  FORMAT ('ATOM',I7,'  CB  ',A3,1X,A1,I4,4X,3F8.3,f15.3)
 132   30  FORMAT ('CONECT',8I5)
 133       return
 134       end
 135 c------------------------------------------------------------------------------
 136       subroutine MOL2out(etot,tytul)
 137 C Prints the Cartesian coordinates of the alpha-carbons in the Tripos mol2
 138 C format.
 139       implicit real*8 (a-h,o-z)
 140       include 'DIMENSIONS'
 141       include 'COMMON.CHAIN'
 142       include 'COMMON.INTERACT'
 143       include 'COMMON.NAMES'
 144       include 'COMMON.IOUNITS'
 145       include 'COMMON.HEADER'
 146       include 'COMMON.SBRIDGE'
 147       character*32 tytul,fd
 148       character*3 zahl
 149       character*6 res_num,pom,ucase
 150 #ifdef AIX
 151       call fdate_(fd)
 152 #elif (defined CRAY)
 153       call date(fd)
 154 #else
 155       call fdate(fd)
 156 #endif
 157       write (imol2,'(a)') '#'
 158       write (imol2,'(a)')
 159      & '#         Creating user name:           unres'
 160       write (imol2,'(2a)') '#         Creation time:                ',
 161      & fd
 162       write (imol2,'(/a)') '\@<TRIPOS>MOLECULE'
 163       write (imol2,'(a)') tytul
 164       write (imol2,'(5i5)') nct-nnt+1,nct-nnt+nss+1,nct-nnt+nss+1,0,0
 165       write (imol2,'(a)') 'SMALL'
 166       write (imol2,'(a)') 'USER_CHARGES'
 167       write (imol2,'(a)') '\@<TRIPOS>ATOM'
 168       do i=nnt,nct
 169         write (zahl,'(i3)') i
 170         pom=ucase(restyp(itype(i)))
 171         res_num = pom(:3)//zahl(2:)
 172         write (imol2,10) i,(c(j,i),j=1,3),i,res_num,0.0
 173       enddo
 174       write (imol2,'(a)') '\@<TRIPOS>BOND'
 175       do i=nnt,nct-1
 176         write (imol2,'(i5,2i6,i2)') i-nnt+1,i-nnt+1,i-nnt+2,1
 177       enddo
 178       do i=1,nss
 179         write (imol2,'(i5,2i6,i2)') nct-nnt+i,ihpb(i),jhpb(i),1
 180       enddo
 181       write (imol2,'(a)') '\@<TRIPOS>SUBSTRUCTURE'
 182       do i=nnt,nct
 183         write (zahl,'(i3)') i
 184         pom = ucase(restyp(itype(i)))
 185         res_num = pom(:3)//zahl(2:)
 186         write (imol2,30) i-nnt+1,res_num,i-nnt+1,0
 187       enddo
 188   10  FORMAT (I7,' CA      ',3F10.4,' C.3',I8,1X,A,F11.4,' ****')
 189   30  FORMAT (I7,1x,A,I14,' RESIDUE',I13,' ****  ****')
 190       return
 191       end
 192 c------------------------------------------------------------------------
 193       subroutine intout
 194       implicit real*8 (a-h,o-z)
 195       include 'DIMENSIONS'
 196       include 'COMMON.IOUNITS'
 197       include 'COMMON.CHAIN'
 198       include 'COMMON.VAR'
 199       include 'COMMON.LOCAL'
 200       include 'COMMON.INTERACT'
 201       include 'COMMON.NAMES'
 202       include 'COMMON.GEO'
 203       include 'COMMON.TORSION'
 204       write (iout,'(/a)') 'Geometry of the virtual chain.'
 205       write (iout,'(7a)') '  Res  ','         d','     Theta',
 206      & '       Phi','       Dsc','     Alpha','      Omega'
 207       do i=1,nres
 208         iti=itype(i)
 209         write (iout,'(a3,i4,6f10.3)') restyp(iti),i,vbld(i),
 210      &     rad2deg*theta(i),rad2deg*phi(i),vbld(nres+i),rad2deg*alph(i),
 211      &     rad2deg*omeg(i)
 212       enddo
 213       return
 214       end
 215 c---------------------------------------------------------------------------
 216       subroutine briefout(it,ener)
 217       implicit real*8 (a-h,o-z)
 218       include 'DIMENSIONS'
 219       include 'COMMON.IOUNITS'
 220       include 'COMMON.CHAIN'
 221       include 'COMMON.VAR'
 222       include 'COMMON.LOCAL'
 223       include 'COMMON.INTERACT'
 224       include 'COMMON.NAMES'
 225       include 'COMMON.GEO'
 226       include 'COMMON.SBRIDGE'
 227 c     print '(a,i5)',intname,igeom
 228 #if defined(AIX) || defined(PGI)
 229       open (igeom,file=intname,position='append')
 230 #else
 231       open (igeom,file=intname,access='append')
 232 #endif
 233       IF (NSS.LE.9) THEN
 234         WRITE (igeom,180) IT,ENER,NSS,(IHPB(I),JHPB(I),I=1,NSS)
 235       ELSE
 236         WRITE (igeom,180) IT,ENER,NSS,(IHPB(I),JHPB(I),I=1,9)
 237         WRITE (igeom,190) (IHPB(I),JHPB(I),I=10,NSS)
 238       ENDIF
 239 c     IF (nvar.gt.nphi) WRITE (igeom,200) (RAD2DEG*THETA(I),I=3,NRES)
 240       WRITE (igeom,200) (RAD2DEG*THETA(I),I=3,NRES)
 241       WRITE (igeom,200) (RAD2DEG*PHI(I),I=4,NRES)
 242 c     if (nvar.gt.nphi+ntheta) then
 243         write (igeom,200) (rad2deg*alph(i),i=2,nres-1)
 244         write (igeom,200) (rad2deg*omeg(i),i=2,nres-1)
 245 c     endif
 246       close(igeom)
 247   180 format (I5,F12.3,I2,9(1X,2I3))
 248   190 format (3X,11(1X,2I3))
 249   200 format (8F10.4)
 250       return
 251       end
 252 #ifdef WINIFL
 253       subroutine fdate(fd)
 254       character*32 fd
 255       write(fd,'(32x)')
 256       return
 257       end
 258 #endif
 259 c----------------------------------------------------------------
 260 #ifdef NOXDR
 261       subroutine cartout(time)
 262 #else
 263       subroutine cartoutx(time)
 264 #endif
 265       implicit real*8 (a-h,o-z)
 266       include 'DIMENSIONS'
 267       include 'COMMON.CHAIN'
 268       include 'COMMON.INTERACT'
 269       include 'COMMON.NAMES'
 270       include 'COMMON.IOUNITS'
 271       include 'COMMON.HEADER'
 272       include 'COMMON.SBRIDGE'
 273       include 'COMMON.DISTFIT'
 274       include 'COMMON.MD'
 275       double precision time
 276 #if defined(AIX) || defined(PGI)
 277       open(icart,file=cartname,position="append")
 278 #else
 279       open(icart,file=cartname,access="append")
 280 #endif
 281       write (icart,'(e15.8,2e15.5,f12.5,$)') time,potE,uconst,t_bath
 282       write (icart,'(i4,$)')
 283      &   nss,(ihpb(j),jhpb(j),j=1,nss)
 284        write (icart,'(i4,20f7.4)') nfrag+npair+3*nfrag_back,
 285      & (qfrag(i),i=1,nfrag),(qpair(i),i=1,npair),
 286      & (utheta(i),ugamma(i),uscdiff(i),i=1,nfrag_back)
 287       write (icart,'(8f10.5)')
 288      & ((c(k,j),k=1,3),j=1,nres),
 289      & ((c(k,j+nres),k=1,3),j=nnt,nct)
 290       close(icart)
 291       return
 292       end
 293 c-----------------------------------------------------------------
 294 #ifndef NOXDR
 295       subroutine cartout(time)
 296       implicit real*8 (a-h,o-z)
 297       include 'DIMENSIONS'
 298 #ifdef MPI
 299       include 'mpif.h'
 300       include 'COMMON.SETUP'
 301 #else
 302       parameter (me=0)
 303 #endif
 304       include 'COMMON.CHAIN'
 305       include 'COMMON.INTERACT'
 306       include 'COMMON.NAMES'
 307       include 'COMMON.IOUNITS'
 308       include 'COMMON.HEADER'
 309       include 'COMMON.SBRIDGE'
 310       include 'COMMON.DISTFIT'
 311       include 'COMMON.MD'
 312       double precision time
 313       integer iret,itmp
 314       real xcoord(3,maxres2+2),prec
 315
 316 #ifdef AIX
 317       call xdrfopen_(ixdrf,cartname, "a", iret)
 318       call xdrffloat_(ixdrf, real(time), iret)
 319       call xdrffloat_(ixdrf, real(potE), iret)
 320       call xdrffloat_(ixdrf, real(uconst), iret)
 321       call xdrffloat_(ixdrf, real(uconst_back), iret)
 322       call xdrffloat_(ixdrf, real(t_bath), iret)
 323       call xdrfint_(ixdrf, nss, iret)
 324       do j=1,nss
 325         call xdrfint_(ixdrf, ihpb(j), iret)
 326         call xdrfint_(ixdrf, jhpb(j), iret)
 327       enddo
 328       call xdrfint_(ixdrf, nfrag+npair+3*nfrag_back, iret)
 329       do i=1,nfrag
 330         call xdrffloat_(ixdrf, real(qfrag(i)), iret)
 331       enddo
 332       do i=1,npair
 333         call xdrffloat_(ixdrf, real(qpair(i)), iret)
 334       enddo
 335       do i=1,nfrag_back
 336         call xdrffloat_(ixdrf, real(utheta(i)), iret)
 337         call xdrffloat_(ixdrf, real(ugamma(i)), iret)
 338         call xdrffloat_(ixdrf, real(uscdiff(i)), iret)
 339       enddo
 340 #else
 341       call xdrfopen(ixdrf,cartname, "a", iret)
 342       call xdrffloat(ixdrf, real(time), iret)
 343       call xdrffloat(ixdrf, real(potE), iret)
 344       call xdrffloat(ixdrf, real(uconst), iret)
 345       call xdrffloat(ixdrf, real(uconst_back), iret)
 346       call xdrffloat(ixdrf, real(t_bath), iret)
 347       call xdrfint(ixdrf, nss, iret)
 348       do j=1,nss
 349         call xdrfint(ixdrf, ihpb(j), iret)
 350         call xdrfint(ixdrf, jhpb(j), iret)
 351       enddo
 352       call xdrfint(ixdrf, nfrag+npair+3*nfrag_back, iret)
 353       do i=1,nfrag
 354         call xdrffloat(ixdrf, real(qfrag(i)), iret)
 355       enddo
 356       do i=1,npair
 357         call xdrffloat(ixdrf, real(qpair(i)), iret)
 358       enddo
 359       do i=1,nfrag_back
 360         call xdrffloat(ixdrf, real(utheta(i)), iret)
 361         call xdrffloat(ixdrf, real(ugamma(i)), iret)
 362         call xdrffloat(ixdrf, real(uscdiff(i)), iret)
 363       enddo
 364 #endif
 365       prec=10000.0
 366       do i=1,nres
 367        do j=1,3
 368         xcoord(j,i)=c(j,i)
 369        enddo
 370       enddo
 371       do i=nnt,nct
 372        do j=1,3
 373         xcoord(j,nres+i-nnt+1)=c(j,i+nres)
 374        enddo
 375       enddo
 376
 377       itmp=nres+nct-nnt+1
 378 #ifdef AIX
 379       call xdrf3dfcoord_(ixdrf, xcoord, itmp, prec, iret)
 380       call xdrfclose_(ixdrf, iret)
 381 #else
 382       call xdrf3dfcoord(ixdrf, xcoord, itmp, prec, iret)
 383       call xdrfclose(ixdrf, iret)
 384 #endif
 385       return
 386       end
 387 #endif
 388 c-----------------------------------------------------------------
 389       subroutine statout(itime)
 390       implicit real*8 (a-h,o-z)
 391       include 'DIMENSIONS'
 392       include 'COMMON.CONTROL'
 393       include 'COMMON.CHAIN'
 394       include 'COMMON.INTERACT'
 395       include 'COMMON.NAMES'
 396       include 'COMMON.IOUNITS'
 397       include 'COMMON.HEADER'
 398       include 'COMMON.SBRIDGE'
 399       include 'COMMON.DISTFIT'
 400       include 'COMMON.MD'
 401       include 'COMMON.SETUP'
 402       integer itime
 403       double precision energia(0:n_ene)
 404       double precision gyrate
 405       external gyrate
 406       common /gucio/ cm
 407       character*256 line1,line2
 408       character*4 format1,format2
 409       character*30 format
 410 #ifdef AIX
 411       if(itime.eq.0) then
 412        open(istat,file=statname,position="append")
 413       endif
 414 #else
 415 #ifdef PGI
 416       open(istat,file=statname,position="append")
 417 #else
 418       open(istat,file=statname,access="append")
 419 #endif
 420 #endif
 421        if (refstr) then
 422          call rms_nac_nnc(rms,frac,frac_nn,co,.false.)
 423           write (line1,'(i10,f15.2,3f12.3,f7.2,4f6.3,3f12.3,i5,$)')
 424      &          itime,totT,EK,potE,totE,
 425      &          rms,frac,frac_nn,co,amax,kinetic_T,t_bath,gyrate(),me
 426           format1="a133"
 427         else
 428           write (line1,'(i10,f15.2,7f12.3,i5,$)')
 429      &           itime,totT,EK,potE,totE,
 430      &           amax,kinetic_T,t_bath,gyrate(),me
 431           format1="a114"
 432         endif
 433         if(usampl.and.totT.gt.eq_time) then
 434            write(line2,'(i5,2f9.4,300f7.4)') iset,uconst,uconst_back,
 435      &      (qfrag(ii1),ii1=1,nfrag),(qpair(ii2),ii2=1,npair),
 436      &      (utheta(i),ugamma(i),uscdiff(i),i=1,nfrag_back)
 437            write(format2,'(a1,i3.3)') "a",23+7*nfrag+7*npair
 438      &             +21*nfrag_back
 439         else
 440            format2="a001"
 441            line2=' '
 442         endif
 443         if (print_compon) then
 444           write(format,'(a1,a4,a1,a4,a10)') "(",format1,",",format2,
 445      &                                                     ",20f12.3)"
 446           write (istat,format) line1,line2,
 447      &      (potEcomp(print_order(i)),i=1,nprint_ene)
 448         else
 449           write(format,'(a1,a4,a1,a4,a1)') "(",format1,",",format2,")"
 450           write (istat,format) line1,line2
 451         endif
 452 #if defined(AIX)
 453         call flush(istat)
 454 #else
 455         close(istat)
 456 #endif
 457        return
 458       end
 459 c---------------------------------------------------------------
 460       double precision function gyrate()
 461       implicit real*8 (a-h,o-z)
 462       include 'DIMENSIONS'
 463       include 'COMMON.INTERACT'
 464       include 'COMMON.CHAIN'
 465       double precision cen(3),rg
 466
 467       do j=1,3
 468        cen(j)=0.0d0
 469       enddo
 470
 471       do i=nnt,nct
 472           do j=1,3
 473             cen(j)=cen(j)+c(j,i)
 474           enddo
 475       enddo
 476       do j=1,3
 477             cen(j)=cen(j)/dble(nct-nnt+1)
 478       enddo
 479       rg = 0.0d0
 480       do i = nnt, nct
 481         do j=1,3
 482          rg = rg + (c(j,i)-cen(j))**2
 483         enddo
 484       end do
 485       gyrate = sqrt(rg/dble(nct-nnt+1))
 486       return
 487       end
 488