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