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 #if defined(AIX) || defined(PGI)
 272       open(icart,file=cartname,position="append")
 273 #else
 274       open(icart,file=cartname,access="append")
 275 #endif
 276       write (icart,'(e15.8,2e15.5,f12.5,$)') time,potE,uconst,t_bath
 277       write (icart,'(i4,$)')
 278      &   nss,(ihpb(j),jhpb(j),j=1,nss)
 279        write (icart,'(i4,20f7.4)') nfrag+npair+3*nfrag_back,
 280      & (qfrag(i),i=1,nfrag),(qpair(i),i=1,npair),
 281      & (utheta(i),ugamma(i),uscdiff(i),i=1,nfrag_back)
 282       write (icart,'(8f10.5)')
 283      & ((c(k,j),k=1,3),j=1,nres),
 284      & ((c(k,j+nres),k=1,3),j=nnt,nct)
 285       close(icart)
 286       return
 287       end
 288 c-----------------------------------------------------------------
 289 #ifndef NOXDR
 290       subroutine cartout(time)
 291       implicit real*8 (a-h,o-z)
 292       include 'DIMENSIONS'
 293 #ifdef MPI
 294       include 'mpif.h'
 295       include 'COMMON.SETUP'
 296 #else
 297       parameter (me=0)
 298 #endif
 299       include 'COMMON.CHAIN'
 300       include 'COMMON.INTERACT'
 301       include 'COMMON.NAMES'
 302       include 'COMMON.IOUNITS'
 303       include 'COMMON.HEADER'
 304       include 'COMMON.SBRIDGE'
 305       include 'COMMON.DISTFIT'
 306       include 'COMMON.MD'
 307       double precision time
 308       integer iret,itmp
 309       real xcoord(3,maxres2+2),prec
 310
 311 #ifdef AIX
 312       call xdrfopen_(ixdrf,cartname, "a", iret)
 313       call xdrffloat_(ixdrf, real(time), iret)
 314       call xdrffloat_(ixdrf, real(potE), iret)
 315       call xdrffloat_(ixdrf, real(uconst), iret)
 316       call xdrffloat_(ixdrf, real(uconst_back), iret)
 317       call xdrffloat_(ixdrf, real(t_bath), iret)
 318       call xdrfint_(ixdrf, nss, iret)
 319       do j=1,nss
 320         call xdrfint_(ixdrf, ihpb(j), iret)
 321         call xdrfint_(ixdrf, jhpb(j), iret)
 322       enddo
 323       call xdrfint_(ixdrf, nfrag+npair+3*nfrag_back, iret)
 324       do i=1,nfrag
 325         call xdrffloat_(ixdrf, real(qfrag(i)), iret)
 326       enddo
 327       do i=1,npair
 328         call xdrffloat_(ixdrf, real(qpair(i)), iret)
 329       enddo
 330       do i=1,nfrag_back
 331         call xdrffloat_(ixdrf, real(utheta(i)), iret)
 332         call xdrffloat_(ixdrf, real(ugamma(i)), iret)
 333         call xdrffloat_(ixdrf, real(uscdiff(i)), iret)
 334       enddo
 335 #else
 336       call xdrfopen(ixdrf,cartname, "a", iret)
 337       call xdrffloat(ixdrf, real(time), iret)
 338       call xdrffloat(ixdrf, real(potE), iret)
 339       call xdrffloat(ixdrf, real(uconst), iret)
 340       call xdrffloat(ixdrf, real(uconst_back), iret)
 341       call xdrffloat(ixdrf, real(t_bath), iret)
 342       call xdrfint(ixdrf, nss, iret)
 343       do j=1,nss
 344         call xdrfint(ixdrf, ihpb(j), iret)
 345         call xdrfint(ixdrf, jhpb(j), iret)
 346       enddo
 347       call xdrfint(ixdrf, nfrag+npair+3*nfrag_back, iret)
 348       do i=1,nfrag
 349         call xdrffloat(ixdrf, real(qfrag(i)), iret)
 350       enddo
 351       do i=1,npair
 352         call xdrffloat(ixdrf, real(qpair(i)), iret)
 353       enddo
 354       do i=1,nfrag_back
 355         call xdrffloat(ixdrf, real(utheta(i)), iret)
 356         call xdrffloat(ixdrf, real(ugamma(i)), iret)
 357         call xdrffloat(ixdrf, real(uscdiff(i)), iret)
 358       enddo
 359 #endif
 360       prec=10000.0
 361       do i=1,nres
 362        do j=1,3
 363         xcoord(j,i)=c(j,i)
 364        enddo
 365       enddo
 366       do i=nnt,nct
 367        do j=1,3
 368         xcoord(j,nres+i-nnt+1)=c(j,i+nres)
 369        enddo
 370       enddo
 371
 372       itmp=nres+nct-nnt+1
 373 #ifdef AIX
 374       call xdrf3dfcoord_(ixdrf, xcoord, itmp, prec, iret)
 375       call xdrfclose_(ixdrf, iret)
 376 #else
 377       call xdrf3dfcoord(ixdrf, xcoord, itmp, prec, iret)
 378       call xdrfclose(ixdrf, iret)
 379 #endif
 380       return
 381       end
 382 #endif
 383 c-----------------------------------------------------------------
 384       subroutine statout(itime)
 385       implicit real*8 (a-h,o-z)
 386       include 'DIMENSIONS'
 387       include 'COMMON.CONTROL'
 388       include 'COMMON.CHAIN'
 389       include 'COMMON.INTERACT'
 390       include 'COMMON.NAMES'
 391       include 'COMMON.IOUNITS'
 392       include 'COMMON.HEADER'
 393       include 'COMMON.SBRIDGE'
 394       include 'COMMON.DISTFIT'
 395       include 'COMMON.MD'
 396       include 'COMMON.REMD'
 397       include 'COMMON.SETUP'
 398       integer itime
 399       double precision energia(0:n_ene)
 400       double precision gyrate
 401       external gyrate
 402       common /gucio/ cm
 403       character*256 line1,line2
 404       character*4 format1,format2
 405       character*30 format
 406 #ifdef AIX
 407       if(itime.eq.0) then
 408        open(istat,file=statname,position="append")
 409       endif
 410 #else
 411 #ifdef PGI
 412       open(istat,file=statname,position="append")
 413 #else
 414       open(istat,file=statname,access="append")
 415 #endif
 416 #endif
 417        if (refstr) then
 418          call rms_nac_nnc(rms,frac,frac_nn,co,.false.)
 419         if(tnp .or. tnp1 .or. tnh) then
 420         write (line1,'(i10,f15.2,3f12.3,f12.6,f7.2,4f6.3,3f12.3,i5,$)')
 421      &          itime,totT,EK,potE,totE,hhh,
 422      &          rms,frac,frac_nn,co,amax,kinetic_T,t_bath,gyrate(),me
 423           format1="a145"
 424         else
 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         endif
 430        else
 431         if(tnp .or. tnp1 .or. tnh) then
 432           write (line1,'(i10,f15.2,7f12.3,f12.6,i5,$)')
 433      &           itime,totT,EK,potE,totE,hhh,
 434      &           amax,kinetic_T,t_bath,gyrate(),me
 435           format1="a126"
 436         else
 437           write (line1,'(i10,f15.2,7f12.3,i5,$)')
 438      &           itime,totT,EK,potE,totE,
 439      &           amax,kinetic_T,t_bath,gyrate(),me
 440           format1="a114"
 441         endif
 442        endif
 443         if(usampl.and.totT.gt.eq_time) then
 444            write(line2,'(i5,2f9.4,300f7.4)') iset,uconst,uconst_back,
 445      &      (qfrag(ii1),ii1=1,nfrag),(qpair(ii2),ii2=1,npair),
 446      &      (utheta(i),ugamma(i),uscdiff(i),i=1,nfrag_back)
 447            write(format2,'(a1,i3.3)') "a",23+7*nfrag+7*npair
 448      &             +21*nfrag_back
 449         elseif(hremd.gt.0) then
 450            write(line2,'(i5)') iset
 451            format2="a005"
 452         else
 453            format2="a001"
 454            line2=' '
 455         endif
 456         if (print_compon) then
 457           if(itime.eq.0) then
 458            write(format,'(a1,a4,a1,a4,a10)') "(",format1,",",format2,
 459      &                                                     ",20a12)"
 460            write (istat,format) "#","",
 461      &      (ename(print_order(i)),i=1,nprint_ene)
 462           endif
 463           write(format,'(a1,a4,a1,a4,a10)') "(",format1,",",format2,
 464      &                                                     ",20f12.3)"
 465           write (istat,format) line1,line2,
 466      &      (potEcomp(print_order(i)),i=1,nprint_ene)
 467         else
 468           write(format,'(a1,a4,a1,a4,a1)') "(",format1,",",format2,")"
 469           write (istat,format) line1,line2
 470         endif
 471 #if defined(AIX)
 472         call flush(istat)
 473 #else
 474         close(istat)
 475 #endif
 476        return
 477       end
 478 c---------------------------------------------------------------
 479       double precision function gyrate()
 480       implicit real*8 (a-h,o-z)
 481       include 'DIMENSIONS'
 482       include 'COMMON.INTERACT'
 483       include 'COMMON.CHAIN'
 484       double precision cen(3),rg
 485
 486       do j=1,3
 487        cen(j)=0.0d0
 488       enddo
 489
 490       do i=nnt,nct
 491           do j=1,3
 492             cen(j)=cen(j)+c(j,i)
 493           enddo
 494       enddo
 495       do j=1,3
 496             cen(j)=cen(j)/dble(nct-nnt+1)
 497       enddo
 498       rg = 0.0d0
 499       do i = nnt, nct
 500         do j=1,3
 501          rg = rg + (c(j,i)-cen(j))**2
 502         enddo
 503       end do
 504       gyrate = sqrt(rg/dble(nct-nnt+1))
 505       return
 506       end
 507