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