Merge branch 'lipid' of mmka.chem.univ.gda.pl:unres into lipid
[unres.git] / source / unres / src_Eshel / unres.F
diff --git a/source/unres/src_Eshel/unres.F b/source/unres/src_Eshel/unres.F
new file mode 100644 (file)
index 0000000..b039983
--- /dev/null
@@ -0,0 +1,210 @@
+CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
+C                                                                              C
+C                                U N R E S                                     C
+C                                                                              C
+C Version for use in parameter optimization on fixed set of decoys             C
+C                                                                              C
+C A. Liwo 12/14/13                                                             C
+C                                                                              C
+CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
+      implicit real*8 (a-h,o-z)
+      include 'DIMENSIONS'
+      include 'COMMON.TIME1'
+      include 'COMMON.INTERACT'
+      include 'COMMON.NAMES'
+      include 'COMMON.GEO'
+      include 'COMMON.HEADER'
+      include 'COMMON.CONTROL'
+      include 'COMMON.CONTACTS'
+      include 'COMMON.CHAIN'
+      include 'COMMON.VAR'
+      include 'COMMON.IOUNITS'
+      include 'COMMON.FFIELD'
+      include 'COMMON.REMD'
+      include 'COMMON.MD'
+      include 'COMMON.SBRIDGE'
+      double precision hrtime,mintime,sectime
+      double precision energia(0:n_ene)
+
+C Initializing variables, setting constants, translation of residue types
+C to numeric types, etc. MUST be called.
+      call init_task
+C Read force field parameters and job setup data
+      call readrtns
+C
+C
+#ifdef MPI
+With fine-graining energy, slaves compute only their energy chunks.
+      if (fg_rank.gt.0) then
+C Fine-grain slaves just do energy and gradient components.
+        call ergastulum ! slave workhouse in Latin
+      else
+#endif
+        if (indpdb.gt.0) then
+          do ifile=1,npdbfile
+            call pdbinput(ifile) 
+            call exec_eeval_or_minim(energia(0))
+            call statout(ifile,energia(0))
+          enddo
+        else
+          call exec_eeval_or_minim(energia(0))
+          call statout(1,energia(0))
+        endif
+#ifdef MPI
+      endif
+C Finish task.
+      if (fg_rank.eq.0) call finish_task
+#ifdef TIMING
+       call print_detailed_timing
+#endif
+      call MPI_Finalize(ierr)
+      stop 'Bye Bye...'
+#else
+      call dajczas(tcpu(),hrtime,mintime,sectime)
+      stop '********** Program terminated normally.'
+#endif
+      end
+c--------------------------------------------------------------------------
+      subroutine exec_eeval_or_minim(energy)
+      implicit real*8 (a-h,o-z)
+      include 'DIMENSIONS'
+#ifdef MPI
+      include 'mpif.h'
+#endif
+      include 'COMMON.SETUP'
+      include 'COMMON.TIME1'
+      include 'COMMON.INTERACT'
+      include 'COMMON.NAMES'
+      include 'COMMON.GEO'
+      include 'COMMON.HEADER'
+      include 'COMMON.CONTROL'
+      include 'COMMON.CONTACTS'
+      include 'COMMON.CHAIN'
+      include 'COMMON.VAR'
+      include 'COMMON.IOUNITS'
+      include 'COMMON.FFIELD'
+      include 'COMMON.REMD'
+      include 'COMMON.MD'
+      include 'COMMON.SBRIDGE'
+      include 'COMMON.ECOMPON'
+      common /srutu/ icall
+      double precision energy(0:n_ene)
+      double precision energy_long(0:n_ene),energy_short(0:n_ene)
+      double precision varia(maxvar)
+      if (indpdb.eq.0) call chainbuild
+#ifdef MPI
+      time00=MPI_Wtime()
+#else
+      time00=tcpu()
+#endif
+      do i=nnt+1,nct-1
+        do j=1,3
+          dc(j,i)=c(j,i+1)-c(j,i)
+        enddo
+      enddo
+      if (nnt.gt.1) then
+        do j=1,3
+          dc(j,1)=dc(j,3)
+          c(j,1)=c(j,2)-dc(j,1)
+        enddo
+      endif
+      if (nct.lt.nres) then
+        do j=1,3
+          dc(j,nres-1)=dc(j,nres-3)
+          c(j,nres)=c(j,nres)+dc(j,nres-1)
+        enddo
+      endif 
+      call cartprint
+      call int_from_cart(.true.,.false.)
+      write (iout,*) "Time for energy evaluation",time_ene
+
+crc overlap test
+      if (overlapsc) then 
+        print *, 'Calling OVERLAP_SC'
+        call overlap_sc(fail)
+      endif 
+
+      if (searchsc) then 
+        call sc_move(2,nres-1,10,1d10,nft_sc,etot)
+        print *,'SC_move',nft_sc,etot
+        write(iout,*) 'SC_move',nft_sc,etot
+      endif 
+
+      call etotal(energy(0))
+#ifdef MPI
+      time_ene=MPI_Wtime()-time00
+#else
+      time_ene=tcpu()-time00
+#endif
+      if (regular) then
+        call gen_dist_constr
+        call regularize(nct-nnt+1,etot,rms,cref(1,nnt),iretcode)
+        call etotal(energy(0))
+        energy(0)=energy(0)-energy(14)
+        etot=energy(0)
+        call enerprint(energy(0))
+        call intout
+        if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
+        write (iout,'(a,i3)') 'SUMSL return code:',iretcode
+      else if (minim) then
+        print *, 'Calling MINIM_DC'
+#ifdef MPI
+        time1=MPI_WTIME()
+#else
+        time1=tcpu()
+#endif
+        call minim_dc(etot,iretcode,nfun)
+        print *,'SUMSL return code is',iretcode,' eval ',nfun
+#ifdef MPI
+        evals=nfun/(MPI_WTIME()-time1)
+#else
+        evals=nfun/(tcpu()-time1)
+#endif
+        print *,'# eval/s',evals
+        print *,'refstr=',refstr
+        call etotal(energy(0))
+        etot = energy(0)
+        call enerprint(energy(0))
+
+        call intout
+        call briefout(0,etot)
+        if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
+        write (iout,'(a,i3)') 'SUMSL return code:',iretcode
+        write (iout,'(a,i20)') '# of energy evaluations:',nfun+1
+        write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals
+      else
+        if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
+      endif
+      call enerprint(energy)
+#define DEBUG
+#ifdef DEBUG
+c Check the correctness of energy components
+      evdw=0.0d0
+      etor=0.0d0
+      esccor=0.0d0
+      do i=1,ntyp
+        do j=1,ntyp
+          evdw=evdw+vdwcompon(i,j)
+          etor=etor+torcompon(i,j)
+          esccor=esccor+sccorcompon(i,j)
+        enddo
+      enddo
+      ebe=0.0d0
+      esc=0.0d0
+      evdw2=0.0d0
+      etord=0.0d0
+      do i=1,ntyp
+        ebe=ebe+becompon(i)
+        esc=esc+sccompon(i)
+        evdw2=evdw2+vdw2compon(i)
+        etord=etord+tordcompon(i) 
+      enddo
+      write (iout,'(a8,1pe15.6)') "EVDW",evdw,"EVDW2",EVDW2,"EBEND",ebe,
+     &  "ESCLOC",esc,"ETOR",etor,"ETORD",etord,"ESCCOR",esccor
+#endif
+#undef DEBUG
+      if (outpdb) call pdbout(etot,titel(:32),ipdb)
+      if (outmol2) call mol2out(etot,titel(:32))
+      return
+      end
+c---------------------------------------------------------------------------