source/wham/src-M-Pawel/rmscalc.f

   1       double precision function rmscalc(ishif,i,j,jcon,lprn)
   2       implicit real*8 (a-h,o-z)
   3       include 'DIMENSIONS'
   4       include 'DIMENSIONS.ZSCOPT'
   5       include 'DIMENSIONS.COMPAR'
   6       include 'COMMON.IOUNITS'
   7       include 'COMMON.COMPAR'
   8       include 'COMMON.CHAIN'
   9       include 'COMMON.INTERACT'
  10       include 'COMMON.VAR'
  11       include 'COMMON.CONTROL'
  12       double precision przes(3),obrot(3,3)
  13       double precision creff(3,maxres2),cc(3,maxres2)
  14       logical iadded(maxres)
  15       integer inumber(2,maxres)
  16       common /ccc/ creff,cc,iadded,inumber
  17       logical lprn
  18       logical non_conv
  19       integer ishif,i,j
  20       if (lprn) then
  21         write (iout,*) "i",i," j",j," jcont",jcon," ishif",ishif
  22         write (iout,*) "npiece",npiece(j,i)
  23         call flush(iout)
  24       endif
  25 c      write (iout,*) "symetr",symetr
  26 c      call flush(iout)
  27       nperm=1
  28       do idup=1,symetr
  29       nperm=nperm*idup
  30       enddo
  31 c      write (iout,*) "nperm",nperm
  32 c      call flush(iout)
  33       do kkk=1,nperm
  34       idup=0
  35       do l=1,nres
  36         iadded(l)=.false.
  37       enddo
  38 c      write (iout,*) "kkk",kkk
  39 c      call flush(iout)
  40       do k=1,npiece(j,i)
  41         if (i.eq.1) then
  42           if (lprn) then
  43             write (iout,*) "Level 1: j=",j,"k=",k," adding fragment",
  44      &         ifrag(1,k,j),ifrag(2,k,j)
  45             call flush(iout)
  46           endif
  47           call cprep(ifrag(1,k,j),ifrag(2,k,j),ishif,idup,kkk)
  48 c          write (iout,*) "Exit cprep"
  49 c          call flush(iout)
  50 c          write (iout,*) "ii=",ii
  51         else
  52           kk = ipiece(k,j,i)
  53 c          write (iout,*) "kk",kk," npiece",npiece(kk,1)
  54           do l=1,npiece(kk,1)
  55             if (lprn) then
  56               write (iout,*) "Level",i,": j=",j,"k=",k," kk=",kk,
  57      &          " l=",l," adding fragment",
  58      &          ifrag(1,l,kk),ifrag(2,l,kk)
  59               call flush(iout)
  60             endif
  61             call cprep(ifrag(1,l,kk),ifrag(2,l,kk),ishif,idup,kkk)
  62 c            write (iout,*) "After cprep"
  63 c            call flush(iout)
  64           enddo
  65         endif
  66       enddo
  67       enddo
  68       if (lprn) then
  69         write (iout,*) "tuszukaj"
  70         do kkk=1,nperm
  71           do k=1,idup
  72             write(iout,'(5i4,2(3f10.5,5x))') i,j,k,inumber(1,k),
  73      &        inumber(2,k),(creff(l,k),l=1,3),(cc(l,k),l=1,3)
  74           enddo
  75         enddo
  76         call flush(iout)
  77       endif
  78       rminrms=1.0d10
  79       do kkk=1,nperm
  80       call fitsq(rms,cc(1,1),creff(1,1),idup,przes,obrot,non_conv)
  81       if (non_conv) then
  82         print *,'Error: FITSQ non-convergent, jcon',jcon,i
  83         rms = 1.0d10
  84       else if (rms.lt.-1.0d-6) then
  85         print *,'Error: rms^2 = ',rms,jcon,i
  86         rms = 1.0d10
  87       else if (rms.ge.1.0d-6 .and. rms.lt.0) then
  88         rms = 0.0d0
  89       endif
  90 c      write (iout,*) "rmsmin", rminrms, "rms", rms
  91       if (rms.le.rminrms) rminrms=rms
  92       enddo
  93       rmscalc = dsqrt(rminrms)
  94 c      write (iout, *) "analysys", rmscalc,anatemp
  95       return
  96       end
  97 c-------------------------------------------------------------------------
  98       subroutine cprep(if1,if2,ishif,idup,kwa)
  99       implicit real*8 (a-h,o-z)
 100       include 'DIMENSIONS'
 101       include 'DIMENSIONS.ZSCOPT'
 102       include 'DIMENSIONS.COMPAR'
 103       include 'COMMON.CONTROL'
 104       include 'COMMON.IOUNITS'
 105       include 'COMMON.COMPAR'
 106       include 'COMMON.CHAIN'
 107       include 'COMMON.INTERACT'
 108       include 'COMMON.VAR'
 109       double precision przes(3),obrot(3,3)
 110       double precision creff(3,maxres2),cc(3,maxres2)
 111       logical iadded(maxres)
 112       integer inumber(2,maxres),iistrart,kwa,blar
 113       common /ccc/ creff,cc,iadded,inumber
 114 c      write (iout,*) "Calling cprep symetr",symetr," kwa",kwa
 115       nperm=1
 116       do blar=1,symetr
 117       nperm=nperm*blar
 118       enddo
 119 c      write (iout,*) "nperm",nperm
 120       kkk=kwa
 121 c      ii=0
 122       do l=if1,if2
 123 c        write (iout,*) "l",l," iadded",iadded(l)
 124 c        call flush(iout)
 125         if (l+ishif.gt.1 .and. l+ishif.le.nres .and. .not.iadded(l))
 126      &  then
 127           idup=idup+1
 128           iadded(l)=.true.
 129           inumber(1,idup)=l
 130           inumber(2,idup)=l+ishif
 131           do m=1,3
 132             creff(m,idup)=cref(m,l,kkk)
 133             cc(m,idup)=c(m,l+ishif)
 134           enddo
 135         endif
 136       enddo
 137       return
 138       end
 139 c-------------------------------------------------------------------------
 140       double precision function rmsnat(jcon)
 141       implicit real*8 (a-h,o-z)
 142       include 'DIMENSIONS'
 143       include 'DIMENSIONS.ZSCOPT'
 144       include 'DIMENSIONS.COMPAR'
 145       include 'COMMON.IOUNITS'
 146       include 'COMMON.COMPAR'
 147       include 'COMMON.CHAIN'
 148       include 'COMMON.INTERACT'
 149       include 'COMMON.VAR'
 150       include 'COMMON.CONTROL'
 151       double precision przes(3),obrot(3,3),cc(3,2*maxres),
 152      &   ccref(3,2*maxres)
 153       logical non_conv
 154       integer ishif,i,j,resprzesun
 155       rminrms=10.0d10
 156       rmsminsing=10d10
 157       nperm=1
 158       do i=1,symetr
 159        nperm=nperm*i
 160       enddo
 161       do kkk=1,nperm
 162        nnsup=0
 163        do i=1,nres
 164         if (itype(i).ne.ntyp1) then
 165           nnsup=nnsup+1
 166           do j=1,3
 167             cc(j,nnsup)=c(j,i)
 168             ccref(j,nnsup)=cref(j,i,kkk)
 169           enddo
 170         endif
 171        enddo
 172        call fitsq(rms,cc(1,1),ccref(1,1),nnsup,przes,obrot,non_conv)
 173        if (non_conv) then
 174         print *,'Error: FITSQ non-convergent, jcon',jcon,i
 175         rms=1.0d10
 176        else if (rms.lt.-1.0d-6) then
 177         print *,'Error: rms^2 = ',rms,jcon,i
 178         rms = 1.0d10
 179        else if (rms.ge.1.0d-6 .and. rms.lt.0) then
 180         rms=0.0d0
 181        endif
 182        if (rms.le.rminrms) rminrms=rms
 183 c       write (iout,*) "kkk",kkk," rmsnat",rms , rminrms
 184       enddo
 185       rmsnat = dsqrt(rminrms)
 186 C      write (iout,*)  "analysys",rmsnat, anatemp
 187 C      liczenie rmsdla pojedynczego lancucha
 188       return
 189       end
 190 c-----------------------------------------------------------------------------
 191       double precision function gyrate(jcon)
 192       implicit real*8 (a-h,o-z)
 193       include 'DIMENSIONS'
 194       include 'COMMON.INTERACT'
 195       include 'COMMON.CHAIN'
 196       double precision cen(3),rg
 197
 198       do j=1,3
 199        cen(j)=0.0d0
 200       enddo
 201
 202       do i=nnt,nct
 203           do j=1,3
 204             cen(j)=cen(j)+c(j,i)
 205           enddo
 206       enddo
 207       do j=1,3
 208             cen(j)=cen(j)/dble(nct-nnt+1)
 209       enddo
 210       rg = 0.0d0
 211       do i = nnt, nct
 212         do j=1,3
 213          rg = rg + (c(j,i)-cen(j))**2
 214         enddo
 215       end do
 216       gyrate = dsqrt(rg/dble(nct-nnt+1))
 217       return
 218       end