homology energy correction and MREMD LANG restart
[unres.git] / source / unres / src-HCD-5D / MREMD.F
index 087b9be..38db8a8 100644 (file)
@@ -1,15 +1,25 @@
       subroutine MREMD
-      implicit real*8 (a-h,o-z)
+      implicit none
       include 'DIMENSIONS'
       include 'mpif.h'
       include 'COMMON.CONTROL'
       include 'COMMON.VAR'
       include 'COMMON.MD'
+#ifdef FIVEDIAG
+       include 'COMMON.LAGRANGE.5diag'
+#else
+       include 'COMMON.LAGRANGE'
+#endif
+      include 'COMMON.QRESTR'
 #ifndef LANG0
       include 'COMMON.LANGEVIN'
 #else
+#ifdef FIVEDIAG
+      include 'COMMON.LANGEVIN.lang0.5diag'
+#else
       include 'COMMON.LANGEVIN.lang0'
 #endif
+#endif
       include 'COMMON.CHAIN'
       include 'COMMON.DERIV'
       include 'COMMON.GEO'
       include 'COMMON.SETUP'
       include 'COMMON.MUCA'
       include 'COMMON.HAIRPIN'
-      integer ERRCODE
+      double precision time00,time01,time02,time03,time04,time05,
+     & time06,time07,time08,time001,tt0
+      double precision scalfac
+      integer i,j,k,il,il1,ii,iex,itmp,i_temp,i_mult,i_iset,i_mset,
+     & i_dir,i_temp1,i_mult1,i_mset1
+      integer ERRCODE,ierr,ierror
       double precision cm(3),L(3),vcm(3)
       double precision energia(0:n_ene)
       double precision remd_t_bath(maxprocs)
       external ilen
       character*50 tytul
       common /gucio/ cm
-      integer itime
+      integer itime,i_set_temp,itt,itime_master,irr,i_iset1
+      integer nharp,iharp(4,maxres/3)
 cold      integer nup(0:maxprocs),ndown(0:maxprocs)
       integer rep2i(0:maxprocs),ireqi(maxprocs)
       integer icache_all(maxprocs)
       integer status(MPI_STATUS_SIZE),statusi(MPI_STATUS_SIZE,maxprocs)
       logical synflag,end_of_run,file_exist /.false./,ovrtim
+      double precision t_bath_temp,delta,ene_iex_iex,ene_i_i,ene_iex_i,
+     & ene_i_iex,xxx,tmp,econstr_temp_iex,econstr_temp_i
+      integer iran_num
+      double precision ran_number
 
 cdeb      imin_itime_old=0
       ntwx_cache=0
@@ -155,6 +175,19 @@ cd           write (*,*) me," After broadcast: file_exist",file_exist
                enddo
               enddo
            endif 
+           stdfp=dsqrt(2*Rb*t_bath/d_time)
+           do i=1,ntyp
+             stdfsc(i)=dsqrt(2*Rb*t_bath/d_time)
+           enddo
+           if (lang.gt.0 .and. .not.surfarea) then
+             do i=nnt,nct-1
+               stdforcp(i)=stdfp*dsqrt(gamp)
+             enddo
+             do i=nnt,nct
+              if (itype(i).ne.ntyp1) stdforcsc(i)=stdfsc(iabs(itype(i)))
+     &                   *dsqrt(gamsc(iabs(itype(i))))
+             enddo
+           endif
        endif
 
        if(me.eq.king) then
@@ -217,6 +250,20 @@ cd           write (*,*) me," After broadcast: file_exist",file_exist
               write (iout,'(a6,i4,a1,100i4)') "ndowna",il,":",
      &                    (ndowna(i,il),i=1,ndowna(0,il))
              enddo
+             stdfp=dsqrt(2*Rb*t_bath/d_time)
+             do i=1,ntyp
+               stdfsc(i)=dsqrt(2*Rb*t_bath/d_time)
+             enddo
+             if (lang.gt.0 .and. .not.surfarea) then
+               do i=nnt,nct-1
+                 stdforcp(i)=stdfp*dsqrt(gamp)
+               enddo
+             do i=nnt,nct
+              if (itype(i).ne.ntyp1) stdforcsc(i)=stdfsc(iabs(itype(i)))
+     &                   *dsqrt(gamsc(iabs(itype(i))))
+               enddo
+             endif
+
         ELSE IF (.not.(rest.and.file_exist)) THEN
          do il=1,remd_m(1)
           ifirst(il)=il
@@ -1259,6 +1306,17 @@ co     &    " rescaling weights with temperature",t_bath
            stdfsc(i)=dsqrt(2*Rb*t_bath/d_time)
          enddo
 
+c Compute the standard deviations of stochastic forces for Langevin dynamics
+c if the friction coefficients do not depend on surface area
+         if (lang.gt.0 .and. .not.surfarea) then
+           do i=nnt,nct-1
+             stdforcp(i)=stdfp*dsqrt(gamp)
+           enddo
+           do i=nnt,nct
+             if (itype(i).ne.ntyp1) stdforcsc(i)=stdfsc(iabs(itype(i)))
+     &                   *dsqrt(gamsc(iabs(itype(i))))
+           enddo
+         endif
 cde         write(iout,*) 'REMD after',me,t_bath
            time08=MPI_WTIME()
            if (me.eq.king .or. .not. out1file) then
@@ -1306,10 +1364,17 @@ cd end
 
 c-----------------------------------------------------------------------
       subroutine write1rst(i_index)
-      implicit real*8 (a-h,o-z)
+      implicit none
       include 'DIMENSIONS'
       include 'mpif.h'
+      include 'COMMON.CONTROL'
       include 'COMMON.MD'
+#ifdef FIVEDIAG
+       include 'COMMON.LAGRANGE.5diag'
+#else
+       include 'COMMON.LAGRANGE'
+#endif
+      include 'COMMON.QRESTR'
       include 'COMMON.IOUNITS'
       include 'COMMON.REMD'
       include 'COMMON.SETUP'
@@ -1324,6 +1389,8 @@ c-----------------------------------------------------------------------
       integer*2 i_index
      &            (maxprocs/4,maxprocs/20,maxprocs/200,maxprocs/200)
        common /przechowalnia/ d_restart1,d_restart2
+      integer i,j,il1,il,ixdrf
+      integer ierr
 
        t5_restart1(1)=totT
        t5_restart1(2)=EK
@@ -1483,10 +1550,11 @@ c-----------------------------------------------------------------------
 
 
       subroutine write1traj
-      implicit real*8 (a-h,o-z)
+      implicit none
       include 'DIMENSIONS'
       include 'mpif.h'
       include 'COMMON.MD'
+      include 'COMMON.QRESTR'
       include 'COMMON.IOUNITS'
       include 'COMMON.REMD'
       include 'COMMON.SETUP'
@@ -1504,6 +1572,8 @@ c-----------------------------------------------------------------------
      &     p_uscdiff(100*maxprocs)
       real p_c(3,(maxres2+2)*maxprocs),r_c(3,maxres2+2)
       common /przechowalnia/ p_c
+      integer ii,i,il,j,ixdrf
+      integer ierr
 
       call mpi_bcast(ii_write,1,mpi_integer,
      &           king,CG_COMM,ierr)
@@ -1706,10 +1776,17 @@ c end debugging
 
 
       subroutine read1restart(i_index)
-      implicit real*8 (a-h,o-z)
+      implicit none
       include 'DIMENSIONS'
       include 'mpif.h'
+      include 'COMMON.CONTROL'
       include 'COMMON.MD'
+#ifdef FIVEDIAG
+       include 'COMMON.LAGRANGE.5diag'
+#else
+       include 'COMMON.LAGRANGE'
+#endif
+      include 'COMMON.QRESTR'
       include 'COMMON.IOUNITS'
       include 'COMMON.REMD'
       include 'COMMON.SETUP'
@@ -1721,6 +1798,8 @@ c end debugging
       integer*2 i_index
      &            (maxprocs/4,maxprocs/20,maxprocs/200,maxprocs/200)
       common /przechowalnia/ d_restart1
+      integer i,j,il,il1,ixdrf,iret,itmp
+      integer ierr
       write (*,*) "Processor",me," called read1restart"
 
          if(me.eq.king)then
@@ -1896,10 +1975,15 @@ c     &           CG_COMM,ierr)
         end
 
       subroutine read1restart_old
-      implicit real*8 (a-h,o-z)
+      implicit none
       include 'DIMENSIONS'
       include 'mpif.h'
       include 'COMMON.MD'
+#ifdef FIVEDIAG
+       include 'COMMON.LAGRANGE.5diag'
+#else
+       include 'COMMON.LAGRANGE'
+#endif
       include 'COMMON.IOUNITS'
       include 'COMMON.REMD'
       include 'COMMON.SETUP'
@@ -1909,6 +1993,8 @@ c     &           CG_COMM,ierr)
       real d_restart1(3,2*maxres*maxprocs),r_d(3,2*maxres),
      &                 t5_restart1(5)
       common /przechowalnia/ d_restart1
+      integer i,j,il,itmp
+      integer ierr
          if(me.eq.king)then
              open(irest2,file=mremd_rst_name,status='unknown')
              read (irest2,*) (i2rep(i),i=0,nodes-1)