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