mremd with homology constrains can start each replica from random model
authorCezary Czaplewski <czarek@chem.univ.gda.pl>
Thu, 28 Apr 2016 21:59:25 +0000 (23:59 +0200)
committerCezary Czaplewski <czarek@chem.univ.gda.pl>
Thu, 28 Apr 2016 21:59:25 +0000 (23:59 +0200)
add keyword START_FROM_MODELS, works only with READ2SIGMA
and without rest on MD keyword card

source/unres/src_MD/COMMON.CHAIN
source/unres/src_MD/COMMON.CONTROL
source/unres/src_MD/MD_A-MTS.F
source/unres/src_MD/readpdb.F
source/unres/src_MD/readrtns.F

index 6e19f8d..d4b60e4 100644 (file)
@@ -1,7 +1,7 @@
       integer nres,nsup,nstart_sup,nz_start,nz_end,iz_sc,
      &  nres0,nstart_seq
       double precision c,dc,dc_old,d_c_work,xloc,xrot,dc_norm,t,r,
-     & prod,rt,dc_work,cref,crefjlee,dc_norm2
+     & prod,rt,dc_work,cref,crefjlee,dc_norm2,chomo
       common /chain/ c(3,maxres2+2),dc(3,0:maxres2),dc_old(3,0:maxres2),
      & xloc(3,maxres),xrot(3,maxres),dc_norm(3,0:maxres2),
      & dc_norm2(3,0:maxres2),
@@ -11,3 +11,4 @@
       common /refstruct/ cref(3,maxres2+2),crefjlee(3,maxres2+2),
      & nsup,nstart_sup,nstart_seq
       common /from_zscore/ nz_start,nz_end,iz_sc
+      common /chomo_models/ chomo(3,maxres2+2,max_template)
\ No newline at end of file
index e512e53..4895b45 100644 (file)
@@ -7,13 +7,14 @@
      &                 sideadd,lsecondary,read_cart,unres_pdb,
      &                 vdisulf,searchsc,lmuca,dccart,extconf,out1file,
      &                 gnorm_check,gradout,split_ene,read2sigma
+     &                ,start_from_model
       common /cntrl/ modecalc,iscode,indpdb,indback,indphi,iranconf,
      & 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,
-     & homol_nset,read2sigma
+     & homol_nset,read2sigma,start_from_model
       common /homol/ waga_homology(maxprocs/20),
      & waga_dist, waga_angle, waga_theta, waga_d, dist_cut
 C... minim = .true. means DO minimization.
index c289095..a785ef0 100644 (file)
@@ -1827,6 +1827,30 @@ c Removing the velocity of the center of mass
        call flush(iout)
       endif
       if (.not.rest) then              
+         if (start_from_model) then
+          i_model=iran_num(1,constr_homology)
+          write (iout,*) 'starting from model ',i_model
+          do i=1,2*nres
+           do j=1,3
+            c(j,i)=chomo(j,i,i_model)
+           enddo
+          enddo
+
+          call int_from_cart(.true.,.false.)
+          call sc_loc_geom(.false.)
+          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
+         endif
          call chainbuild
          write (iout,*) "PREMINIM ",preminim
          if(iranconf.ne.0 .or. preminim) then
index 042a20f..eed985e 100644 (file)
@@ -742,6 +742,7 @@ C Copy the coordinates to reference coordinates
       do i=1,2*nres
         do j=1,3
           cref(j,i)=c(j,i)
+          chomo(j,i,k)=c(j,i)
         enddo
       enddo
 
index 74ba3d8..915e5f2 100644 (file)
@@ -2761,6 +2761,16 @@ c Alternative: reading from input
  
       call readi(controlcard,"HOMOL_NSET",homol_nset,1)       
       read2sigma=(index(controlcard,'READ2SIGMA').gt.0)
+      start_from_model=(index(controlcard,'START_FROM_MODELS').gt.0)
+      if(.not.read2sigma.and.start_from_model) then
+          write(iout,*) 'START_FROM_MODELS works only with READ2SIGMA'
+          start_from_model=.false.
+      endif
+      if(start_from_model) write(iout,*) 'START_FROM_MODELS is ON'
+      if(start_from_model .and. rest) then 
+        write(iout,*) 'START_FROM_MODELS is OFF'
+        write(iout,*) 'remove restart keyword from input'
+      endif
       if (homol_nset.gt.1)then
          call card_concat(controlcard)
          read(controlcard,*) (waga_homology(i),i=1,homol_nset)