wham is prining energies; all "good" changes seems to be revoked improve in rmsd...
[unres.git] / source / unres / src_MD-M / energy_p_new_barrier.F
index 8a2d035..4a47935 100644 (file)
@@ -1107,7 +1107,7 @@ C------------------------------------------------------------------------
      &  escloc,wscloc,etors,wtor,etors_d,wtor_d,ehpb,wstrain,
      &  ecorr,wcorr,
      &  ecorr5,wcorr5,ecorr6,wcorr6,eel_loc,wel_loc,eello_turn3,wturn3,
-     &  eello_turn4,wturn4,eello_turn6,wturn6,esccor,wsccro,edihcnstr,
+     &  eello_turn4,wturn4,eello_turn6,wturn6,esccor,wsccor,edihcnstr,
      &  ethetacnstr,ebr*nss,Uconst,eliptran,wliptran,Eafmforc,
      &  etot
    10 format (/'Virtual-chain energies:'//
@@ -1780,6 +1780,7 @@ C lipbufthick is thickenes of lipid buffore
      &  +aa_aq(itypi,itypj)*(2.0d0-sslipi-sslipj)/2.0d0
       bb=bb_lip(itypi,itypj)*(sslipi+sslipj)/2.0d0
      &  +bb_aq(itypi,itypj)*(2.0d0-sslipi-sslipj)/2.0d0
+C      write(iout,*) "tu,", i,j,aa_lip(itypi,itypj),bb_lip(itypi,itypj)
 C      if (aa.ne.aa_aq(itypi,itypj)) write(63,'(2e10.5)')
 C     &(aa-aa_aq(itypi,itypj)),(bb-bb_aq(itypi,itypj))
 C      if (ssgradlipj.gt.0.0d0) print *,"??WTF??"
@@ -2056,6 +2057,7 @@ C lipbufthick is thickenes of lipid buffore
      &  +bb_aq(itypi,itypj)*(2.0d0-sslipi-sslipj)/2.0d0
 C      if (aa.ne.aa_aq(itypi,itypj)) write(63,'2e10.5') 
 C     &(aa-aa_aq(itypi,itypj)),(bb-bb_aq(itypi,itypj))
+C      write(iout,*) "tu,", i,j,aa,bb,aa_lip(itypi,itypj),sslipi,sslipj
       dist_init=(xj-xi)**2+(yj-yi)**2+(zj-zi)**2
       xj_safe=xj
       yj_safe=yj
@@ -3461,7 +3463,7 @@ c        end if
         num_cont_hb(i)=num_conti
       enddo
       do i=iturn4_start,iturn4_end
-        if (i.le.1) cycle
+        if (i.lt.1) cycle
         if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1
 C changes suggested by Ana to avoid out of bounds
 c     & .or.((i+5).gt.nres)
@@ -3591,7 +3593,7 @@ C I TU KURWA
         do j=ielstart(i),ielend(i)
 C          do j=16,17
 C          write (iout,*) i,j
-         if (j.le.1) cycle
+C         if (j.le.1) cycle
           if (itype(j).eq.ntyp1.or. itype(j+1).eq.ntyp1
 C changes suggested by Ana to avoid out of bounds
 c     & .or.((j+2).gt.nres)
@@ -3659,6 +3661,7 @@ C 13-go grudnia roku pamietnego...
       double precision unmat(3,3) /1.0d0,0.0d0,0.0d0,
      &                   0.0d0,1.0d0,0.0d0,
      &                   0.0d0,0.0d0,1.0d0/
+       integer xshift,yshift,zshift
 c          time00=MPI_Wtime()
 cd      write (iout,*) "eelecij",i,j
 c          ind=ind+1
@@ -3798,6 +3801,7 @@ cd     &      xmedi,ymedi,zmedi,xj,yj,zj
               write (iout,'(a6,2i5,0pf7.3,2i5,2e11.3)') 
      &'evdw1',i,j,evdwij
      &,iteli,itelj,aaa,evdw1
+              write (iout,*) sss
               write (iout,'(a6,2i5,0pf7.3,2f8.3)') 'ees',i,j,eesij,
      &fac_shield(i),fac_shield(j)
           endif
@@ -5848,6 +5852,7 @@ C       Checking if it involves dummy (NH3+ or COO-) group
          if (itype(i-1).eq.ntyp1 .or. itype(i).eq.ntyp1) then
 C YES   vbldpDUM is the equlibrium length of spring for Dummy atom
         diff = vbld(i)-vbldpDUM
+        if (energy_dec) write(iout,*) "dum_bond",i,diff 
          else
 C NO    vbldp0 is the equlibrium lenght of spring for peptide group
         diff = vbld(i)-vbldp0
@@ -5861,6 +5866,7 @@ C NO    vbldp0 is the equlibrium lenght of spring for peptide group
 c        write (iout,'(i5,3f10.5)') i,(gradb(j,i-1),j=1,3)
 c        endif
       enddo
+      
       estr=0.5d0*AKP*estr+estr1
 c
 c 09/18/07 AL: multimodal bond potential based on AM1 CA-SC PMF's included
@@ -8720,7 +8726,7 @@ c     & ' eij',eij,' eesij',ees0pij,ees0mij,' and ',k,l
 c     & ,' fcont ',ekl,' eeskl',ees0pkl,ees0mkl,' energy=',ekont*ees,
 c     & 'gradcorr_long'
 C Calculate the multi-body contribution to energy.
-c      ecorr=ecorr+ekont*ees
+C      ecorr=ecorr+ekont*ees
 C Calculate multi-body contributions to the gradient.
       coeffpees0pij=coeffp*ees0pij
       coeffmees0mij=coeffm*ees0mij
@@ -11120,7 +11126,7 @@ C       do i=1,1
         if (itype(i).eq.ntyp1) cycle
 
         positi=(mod(((c(3,i)+c(3,i+1))/2.0d0),boxzsize))
-        if (positi.le.0) positi=positi+boxzsize
+        if (positi.le.0.0) positi=positi+boxzsize
 C        print *,i
 C first for peptide groups
 c for each residue check if it is in lipid or lipid water border area