source/wham/src-M/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 C      write (iout,*) "tu2", nres,nsup
 159       noverlap=nres
 160       if (nres.gt.nsup+nnt-1) noverlap=nsup+nnt-1
 161 c      write (iout,*) "tu3,",noverlap
 162       do i=1,symetr
 163        nperm=nperm*i
 164       enddo
 165       do kkk=1,nperm
 166        write (iout,*) "kkk",kkk
 167        nnsup=0
 168        do i=1,noverlap
 169         if (itype(i).ne.ntyp1) then
 170           nnsup=nnsup+1
 171           do j=1,3
 172             cc(j,nnsup)=c(j,i)
 173             ccref(j,nnsup)=cref(j,i,kkk)
 174           enddo
 175               write (iout,'(2i5,3f10.5,5x,3f10.5)')
 176      &        i,nnsup,(cc(j,nnsup),j=1,3),(ccref(j,nnsup),j=1,3)
 177         endif
 178        enddo
 179
 180        call fitsq(rms,cc(1,1),ccref(1,1),nnsup,przes,obrot,non_conv)
 181        if (non_conv) then
 182         print *,'Error: FITSQ non-convergent, jcon',jcon,i
 183         rms=1.0d10
 184        else if (rms.lt.-1.0d-6) then
 185         print *,'Error: rms^2 = ',rms,jcon,i
 186         rms = 1.0d10
 187        else if (rms.ge.1.0d-6 .and. rms.lt.0) then
 188         rms=0.0d0
 189        endif
 190        if (rms.le.rminrms) rminrms=rms
 191 c       write (iout,*) "kkk",kkk," rmsnat",rms , rminrms
 192       enddo
 193       rmsnat = dsqrt(rminrms)
 194 C      write (iout,*)  "analysys",rmsnat, anatemp
 195 C      liczenie rmsdla pojedynczego lancucha
 196       return
 197       end
 198 c-----------------------------------------------------------------------------
 199       double precision function gyrate(jcon)
 200       implicit real*8 (a-h,o-z)
 201       include 'DIMENSIONS'
 202       include 'COMMON.INTERACT'
 203       include 'COMMON.CHAIN'
 204       double precision cen(3),rg
 205
 206       do j=1,3
 207        cen(j)=0.0d0
 208       enddo
 209
 210       do i=nnt,nct
 211           do j=1,3
 212             cen(j)=cen(j)+c(j,i)
 213           enddo
 214       enddo
 215       do j=1,3
 216             cen(j)=cen(j)/dble(nct-nnt+1)
 217       enddo
 218       rg = 0.0d0
 219       do i = nnt, nct
 220         do j=1,3
 221          rg = rg + (c(j,i)-cen(j))**2
 222         enddo
 223       end do
 224       gyrate = dsqrt(rg/dble(nct-nnt+1))
 225       return
 226       end