After final fix by AL

diff --git a/source/unres/src_MD/COMMON.CONTROL b/source/unres/src_MD/COMMON.CONTROL
index 438efaf..8a3fbf5 100644 (file)
--- a/source/unres/src_MD/COMMON.CONTROL
+++ b/source/unres/src_MD/COMMON.CONTROL
@@ -2,17 +2,19 @@
      & inprint,i2ndstr,mucadyn,constr_dist,constr_homology,homol_nset
       real*8 waga_dist1,waga_angle1
       real*8 waga_dist, waga_angle, waga_theta, waga_d, dist_cut
-      logical minim,refstr,pdbref,outpdb,outmol2,overlapsc,energy_dec,
+      logical minim,refstr,pdbref,outpdb,outmol2,outx,overlapsc,
+     &                 energy_dec,
      &                 sideadd,lsecondary,read_cart,unres_pdb,
      &                 vdisulf,searchsc,lmuca,dccart,extconf,out1file,
      &                 gnorm_check,gradout,split_ene
       common /cntrl/ modecalc,iscode,indpdb,indback,indphi,iranconf,
-     & icheckgrad,minim,i2ndstr,refstr,pdbref,outpdb,outmol2,iprint,
+     & icheckgrad,minim,i2ndstr,refstr,pdbref,outpdb,outmol2,outx,
+     & iprint,
      & overlapsc,energy_dec,sideadd,lsecondary,read_cart,unres_pdb
      & ,vdisulf,searchsc,lmuca,dccart,mucadyn,extconf,out1file,
      & constr_dist,gnorm_check,gradout,split_ene,constr_homology,
-     & waga_dist, waga_angle, waga_theta, waga_d, dist_cut,
      & homol_nset
-     common /homol/ waga_dist1(maxprocs/20), waga_angle1(maxprocs/20)
+      common /homol/ waga_dist1(maxprocs/20), waga_angle1(maxprocs/20),
+     & waga_dist, waga_angle, waga_theta, waga_d, dist_cut
 C... minim = .true. means DO minimization.
 C... energy_dec = .true. means print energy decomposition matrix

diff --git a/source/unres/src_MD/COMMON.MD b/source/unres/src_MD/COMMON.MD
index b66ea2c..24fc1b3 100644 (file)
--- a/source/unres/src_MD/COMMON.MD
+++ b/source/unres/src_MD/COMMON.MD
@@ -52,7 +52,7 @@ c
       integer nresn,nyosh,nnos
       double precision glogs,qmass,vlogs,xlogs
       logical large,print_compon,tbf,rest,reset_moment,reset_vel,
-     & surfarea,rattle,usampl,mdpdb,RESPA,tnp,tnp1,tnh,xiresp
+     & surfarea,rattle,usampl,mdpdb,RESPA,tnp,tnp1,tnh,xiresp,preminim
       integer igmult_start,igmult_end,my_ng_count,ng_start,ng_counts,
      & nginv_start,nginv_counts,myginv_ng_count
       common /back_constr/ uconst_back,utheta,ugamma,uscdiff,
@@ -72,7 +72,7 @@ c
       common /mdpar/ v_ini,d_time,d_time0,scal_fric,
      & t_bath,tau_bath,dvmax,damax,n_timestep,mdpdb,
      & ntime_split,ntime_split0,maxtime_split,
-     & ntwx,ntwe,large,print_compon,tbf,rest,tnp,tnp1,tnh
+     & ntwx,ntwe,large,print_compon,tbf,rest,tnp,tnp1,tnh,preminim
       common /MDcalc/ totT,totE,potE,potEcomp,EK,amax,edriftmax,
      & kinetic_T
       common /lagrange/ d_t,d_t_old,d_t_new,d_t_work,

diff --git a/source/unres/src_MD/MD_A-MTS.F b/source/unres/src_MD/MD_A-MTS.F
index 95f174d..995b558 100644 (file)
--- a/source/unres/src_MD/MD_A-MTS.F
+++ b/source/unres/src_MD/MD_A-MTS.F
@@ -1828,7 +1828,8 @@ c Removing the velocity of the center of mass
       endif
       if (.not.rest) then              
          call chainbuild
-         if(iranconf.ne.0) then
+         write (iout,*) "PREMINIM ",preminim
+         if(iranconf.ne.0 .or. preminim) then
           if (overlapsc) then 
            print *, 'Calling OVERLAP_SC'
            call overlap_sc(fail)
@@ -1850,8 +1851,12 @@ c Removing the velocity of the center of mass
            call minimize(etot,varia,iretcode,nfun)
            call var_to_geom(nvar,varia)
           endif
-          if(me.eq.king.or..not.out1file)
-     &       write(iout,*) 'SUMSL return code is',iretcode,' eval ',nfun
+          if(me.eq.king.or..not.out1file) then
+             write(iout,*) "Minimized energy is",etot
+             write(iout,*) 'SUMSL return code is',iretcode,' eval ',nfun
+             call etotal(potEcomp)
+             call enerprint(potEcomp)
+          endif
          endif
       endif      
       call chainbuild_cart

diff --git a/source/unres/src_MD/Makefile_MPICH_ifort b/source/unres/src_MD/Makefile_MPICH_ifort
index fb95aa7..8a07ff9 100644 (file)
--- a/source/unres/src_MD/Makefile_MPICH_ifort
+++ b/source/unres/src_MD/Makefile_MPICH_ifort
@@ -102,13 +102,13 @@ readpdb.o : readpdb.F
        ${FC} ${FFLAGS2} ${CPPFLAGS} readpdb.F
 
 sumsld.o : sumsld.f
-       ${FC} ${FFLAGS} ${CPPFLAGS} sumsld.f
+       ${FC} ${FFLAGS2} ${CPPFLAGS} sumsld.f
         
 cored.o : cored.f
-       ${FC} ${FFLAGS1} ${CPPFLAGS} cored.f
+       ${FC} ${FFLAGS2} ${CPPFLAGS} cored.f
  
 rmdd.o : rmdd.f
-       ${FC} ${FFLAGS} ${CPPFLAGS} rmdd.f
+       ${FC} ${FFLAGS2} ${CPPFLAGS} rmdd.f
 
 energy_p_new_barrier.o : energy_p_new_barrier.F
        ${FC} ${FFLAGSE} ${CPPFLAGS} energy_p_new_barrier.F

diff --git a/source/unres/src_MD/cinfo.f b/source/unres/src_MD/cinfo.f
index e00e4af..81455fd 100644 (file)
--- a/source/unres/src_MD/cinfo.f
+++ b/source/unres/src_MD/cinfo.f
@@ -1,10 +1,10 @@
 C DO NOT EDIT THIS FILE - IT HAS BEEN GENERATED BY COMPINFO.C
-C 0 40376 2
+C 0 40376 39
       subroutine cinfo
       include 'COMMON.IOUNITS'
       write(iout,*)'++++ Compile info ++++'
-      write(iout,*)'Version 0.40376 build 2'
-      write(iout,*)'compiled Wed Feb 25 14:06:14 2015'
+      write(iout,*)'Version 0.40376 build 39'
+      write(iout,*)'compiled Tue Mar  3 12:57:02 2015'
       write(iout,*)'compiled by felipe@piasek4'
       write(iout,*)'OS name:    Linux '
       write(iout,*)'OS release: 3.2.0-70-generic '

diff --git a/source/unres/src_MD/convert.f b/source/unres/src_MD/convert.f
index dc0cccd..9cb5741 100644 (file)
--- a/source/unres/src_MD/convert.f
+++ b/source/unres/src_MD/convert.f
@@ -15,7 +15,8 @@ C
       include 'COMMON.GEO'
       include 'COMMON.CHAIN'
       double precision x(n)
-cd    print *,'nres',nres,' nphi',nphi,' ntheta',ntheta,' nvar',nvar
+c      print *,'nres',nres,' nphi',nphi,' ntheta',ntheta,' nvar',nvar
+c      print *,"n",n
       do i=4,nres
         x(i-3)=phi(i)
 cd      print *,i,i-3,phi(i)
@@ -28,9 +29,9 @@ cd      print *,i,i-2+nphi,theta(i)
       if (n.eq.nphi+ntheta) return
       do i=2,nres-1
        if (ialph(i,1).gt.0) then
+c          print *,i,ialph(i,1),ialph(i,1)+nside,alph(i),omeg(i)
          x(ialph(i,1))=alph(i)
          x(ialph(i,1)+nside)=omeg(i)
-cd        print *,i,ialph(i,1),ialph(i,1)+nside,alph(i),omeg(i)
         endif
       enddo      
       return

diff --git a/source/unres/src_MD/geomout.F b/source/unres/src_MD/geomout.F
index ac2c759..77d342e 100644 (file)
--- a/source/unres/src_MD/geomout.F
+++ b/source/unres/src_MD/geomout.F
@@ -268,6 +268,7 @@ c----------------------------------------------------------------
       include 'COMMON.DISTFIT'
       include 'COMMON.MD'
       double precision time
+      write (iout,*) "cartout: cartname ",cartname
 #if defined(AIX) || defined(PGI)
       open(icart,file=cartname,position="append")
 #else

diff --git a/source/unres/src_MD/minimize_p.F b/source/unres/src_MD/minimize_p.F
index c7922c7..9192bc1 100644 (file)
--- a/source/unres/src_MD/minimize_p.F
+++ b/source/unres/src_MD/minimize_p.F
@@ -483,6 +483,7 @@ c     v(25)=4.0D0
         endif
       enddo
 
+      call flush(iout)
       call sumsl(k,d,x,func_dc,grad_dc,iv,liv,lv,v,idum,rdum,fdum)      
 
       k=0

diff --git a/source/unres/src_MD/readrtns.F b/source/unres/src_MD/readrtns.F
index 8c3938e..8a68708 100644 (file)
--- a/source/unres/src_MD/readrtns.F
+++ b/source/unres/src_MD/readrtns.F
@@ -8,6 +8,7 @@
       include 'COMMON.CONTROL'
       include 'COMMON.SBRIDGE'
       include 'COMMON.IOUNITS'
+      include 'COMMON.CHAIN'
       logical file_exist
 C Read force-field parameters except weights
       call parmread
@@ -130,6 +131,7 @@ C Set up the time limit (caution! The time must be input in minutes!)
       sideadd=(index(controlcard,'SIDEADD').gt.0)
       energy_dec=(index(controlcard,'ENERGY_DEC').gt.0)
       outpdb=(index(controlcard,'PDBOUT').gt.0)
+      outx=(index(controlcard,'XOUT').gt.0)
       outmol2=(index(controlcard,'MOL2OUT').gt.0)
       pdbref=(index(controlcard,'PDBREF').gt.0)
       refstr=pdbref .or. (index(controlcard,'REFSTR').gt.0)
@@ -386,6 +388,11 @@ C
       large = index(controlcard,"LARGE").gt.0
       print_compon = index(controlcard,"PRINT_COMPON").gt.0
       rattle = index(controlcard,"RATTLE").gt.0
+      preminim = index(controlcard,"PREMINIM").gt.0
+      if (preminim) then
+        dccart=(index(controlcard,'CART').gt.0)
+        call read_minim
+      endif
 c  if performing umbrella sampling, fragments constrained are read from the fragment file 
       nset=0
       if(usampl) then
@@ -426,6 +433,10 @@ c  if performing umbrella sampling, fragments constrained are read from the frag
        write (iout,'(a60,i10)') "Frequency of coordinate output:",ntwx
        if (rattle) write (iout,'(a60)') 
      &  "Rattle algorithm used to constrain the virtual bonds"
+       if (preminim .or. iranconf.gt.0) then
+         write (iout,'(a60)')
+     &      "Initial structure will be energy-minimized" 
+       endif
       endif
       reset_fricmat=1000
       if (lang.gt.0) then
@@ -928,6 +939,12 @@ c        print *,'Begin reading pdb data'
             crefjlee(j,i)=c(j,i)
           enddo
         enddo
+#ifdef DEBUG
+        do i=1,nres
+          write (iout,'(i5,3f8.3,5x,3f8.3)') i,(crefjlee(j,i),j=1,3),
+     &      (crefjlee(j,i+nres),j=1,3)
+        enddo
+#endif
 c        print *,'Finished reading pdb data'
         if(me.eq.king.or..not.out1file)
      &   write (iout,'(a,i3,a,i3)')'nsup=',nsup,
@@ -1179,6 +1196,36 @@ c        write (iout,*) "constr_dist",constr_dist,nstart_sup,nsup
             enddo
           enddo
         endif
+#ifdef DEBUG
+        write (iout,*) "Array C"
+        do i=1,nres
+          write (iout,'(i5,3f8.3,5x,3f8.3)') i,(c(j,i),j=1,3),
+     &      (c(j,i+nres),j=1,3)
+        enddo
+        write (iout,*) "Array Cref"
+        do i=1,nres
+          write (iout,'(i5,3f8.3,5x,3f8.3)') i,(cref(j,i),j=1,3),
+     &      (cref(j,i+nres),j=1,3)
+        enddo
+#endif
+       call int_from_cart1(.false.)
+       call sc_loc_geom(.false.)
+       do i=1,nres
+         thetaref(i)=theta(i)
+         phiref(i)=phi(i)
+       enddo
+       do i=1,nres-1
+         do j=1,3
+           dc(j,i)=c(j,i+1)-c(j,i)
+           dc_norm(j,i)=dc(j,i)*vbld_inv(i+1)
+         enddo
+       enddo
+       do i=2,nres-1
+         do j=1,3
+           dc(j,i+nres)=c(j,i+nres)-c(j,i)
+           dc_norm(j,i+nres)=dc(j,i+nres)*vbld_inv(i+nres)
+         enddo
+       enddo
       else
         homol_nset=0
       endif
@@ -1196,6 +1243,8 @@ C initial geometry.
           if(me.eq.king.or..not.out1file .and.fg_rank.eq.0)
      &     write (iout,'(a)') 'Initial geometry will be read in.'
           if (read_cart) then
+            read (inp,*) time,potE,uconst,t_bath,
+     &       nss,(ihpb(j),jhpb(j),j=1,nss), nn, (qfrag(i),i=1,nn)
             read(inp,'(8f10.5)',end=36,err=36)
      &       ((c(l,k),l=1,3),k=1,nres),
      &       ((c(l,k+nres),l=1,3),k=nnt,nct)
@@ -1214,7 +1263,7 @@ C initial geometry.
                 enddo
               endif
             enddo
-            return
+c            return
           else
             call read_angles(inp,*36)
           endif
@@ -1298,6 +1347,8 @@ C Generate distance constraints, if the PDB structure is to be regularized.
       if (nthread.gt.0) then
         call read_threadbase
       endif
+      write (iout,*) "READRTNS: Calling setup_var"
+      call flush(iout)
       call setup_var
       if (me.eq.king .or. .not. out1file)
      & call intout
@@ -2276,6 +2327,8 @@ c      print *,"Processor",myrank," fg_rank",fg_rank
      &  //'.pdb'
       mol2name=prefix(:lenpre)//'_'//pot(:lenpot)//
      &  liczba(:ilen(liczba))//'.mol2'
+      cartname=prefix(:lenpre)//'_'//pot(:lenpot)//
+     &  liczba(:ilen(liczba))//'.x'
       statname=prefix(:lenpre)//'_'//pot(:lenpot)//
      &  liczba(:ilen(liczba))//'.stat'
       if (lentmp.gt.0)
@@ -2294,6 +2347,7 @@ c      print *,"Processor",myrank," fg_rank",fg_rank
       intname=prefix(:lenpre)//'_'//pot(:lenpot)//'.int'
       pdbname=prefix(:lenpre)//'_'//pot(:lenpot)//'.pdb'
       mol2name=prefix(:lenpre)//'_'//pot(:lenpot)//'.mol2'
+      cartname=prefix(:lenpre)//'_'//pot(:lenpot)//'.x'
       statname=prefix(:lenpre)//'_'//pot(:lenpot)//'.stat'
       if (lentmp.gt.0)
      &  call copy_to_tmp(pref_orig(:ilen(pref_orig))//'_'//pot(:lenpot)

diff --git a/source/unres/src_MD/unres.F b/source/unres/src_MD/unres.F
index 632374b..6d9de55 100644 (file)
--- a/source/unres/src_MD/unres.F
+++ b/source/unres/src_MD/unres.F
@@ -219,7 +219,6 @@ c---------------------------------------------------------------------------
       time_ene=tcpu()-time00
 #endif
       write (iout,*) "Time for energy evaluation",time_ene
-      print *,"after etotal"
       etota = energy(0)
       etot =etota
       call enerprint(energy(0))
@@ -228,18 +227,14 @@ c---------------------------------------------------------------------------
       if (minim) then
 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 
 
         if (dccart) then
-          print *, 'Calling MINIM_DC'
 #ifdef MPI
           time1=MPI_WTIME()
 #else
@@ -252,7 +247,6 @@ crc overlap test
             call chainbuild
           endif
           call geom_to_var(nvar,varia)
-          print *,'Calling MINIMIZE.'
 #ifdef MPI
           time1=MPI_WTIME()
 #else
@@ -285,6 +279,12 @@ crc overlap test
         if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
         call briefout(0,etot)
       endif
+      potE=etot
+#ifdef NOXDR
+       if (outx) call cartout(0.0d0)
+#else
+       if (outx) call cartoutx(0.0d0)
+#endif
       if (outpdb) call pdbout(etot,titel(:32),ipdb)
       if (outmol2) call mol2out(etot,titel(:32))
       return
@@ -323,6 +323,12 @@ c---------------------------------------------------------------------------
       call enerprint(energy(0))
       call intout
       call briefout(0,etot)
+      potE=etot
+#ifdef NOXDR
+      if (outx) call cartout(0.0d0)
+#else
+      if (outx) call cartoutx(0.0d0)
+#endif
       if (outpdb) call pdbout(etot,titel(:32),ipdb)
       if (outmol2) call mol2out(etot,titel(:32))
       if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)