nanophere working
[unres.git] / source / unres / src_MD-M / energy_p_new_barrier.F
index d44e78f..c037ff3 100644 (file)
@@ -12659,23 +12659,37 @@ C for vectorization reasons we will sumup at the end to avoid depenence of previ
        sc_aa_tube=sc_aa_tube_par(iti)
        sc_bb_tube=sc_bb_tube_par(iti)
        enetube(i+nres)=sc_aa_tube/rdiff6**2.0d0+sc_bb_tube/rdiff6
+C       enetube(i+nres)=0.0d0
 C pep_aa_tube and pep_bb_tube are precomputed values A=4eps*sigma^12 B=4eps*sigma^6
 C now we calculate gradient
        fac=-12.0d0*sc_aa_tube/rdiff6**2.0d0/rdiff-
      &       6.0d0*sc_bb_tube/rdiff6/rdiff
+C       fac=0.0
 C now direction of gg_tube vector
 C Now cavity term E=a(x+bsqrt(x)+c)/(1+dx^12)
-         if (acavtub(iti).eq.0.0d0) go to 667
+         if (acavtub(iti).eq.0.0d0) then
+C go to 667
+         enecavtube(i)=0.0
+         faccav=0.0
+         else
          denominator=(1.0+dcavtub(iti)*rdiff6*rdiff6)
          enecavtube(i)=
-     &   acavtub(iti)*(rdiff+bcavtub(iti)*sqrt(rdiff)+ccavtub(iti))
+     &   (bcavtub(iti)*rdiff+acavtub(iti)*sqrt(rdiff)+ccavtub(iti))
      &   /denominator
-         faccav=(acavtub(iti)*(1.0+bcavtub(iti)/2.0/sqrt(rdiff))
-     &   *denominator-acavtub(iti)*(rdiff+bcavtub(iti)*sqrt(rdiff)
-     &   +ccavtub(iti))*rdiff6**2.0d0/rdiff*dcavtub(iti))
+C         enecavtube(i)=0.0
+         faccav=((bcavtub(iti)*1.0d0+acavtub(iti)/2.0d0/sqrt(rdiff))
+     &   *denominator-(bcavtub(iti)*rdiff+acavtub(iti)*sqrt(rdiff)
+     &   +ccavtub(iti))*rdiff6**2.0d0/rdiff*dcavtub(iti)*12.0d0)
      &   /denominator**2.0d0
+C         faccav=0.0
          fac=fac+faccav
- 667     continue
+C 667     continue
+         endif
+         print *,"TUT",i,iti,rdiff,rdiff6,acavtub(iti),denominator,
+     &   enecavtube(i),faccav
+         print *,"licz=",
+     & (bcavtub(iti)*rdiff+acavtub(iti)*sqrt(rdiff)+ccavtub(iti))
+         print *,"finene=",enetube(i+nres)+enecavtube(i)
          do j=1,3
           gg_tube_SC(j,i)=gg_tube_SC(j,i)+vectube(j)*fac
           gg_tube(j,i-1)=gg_tube(j,i-1)+vectube(j)*fac
@@ -12690,7 +12704,7 @@ C        if (acavtub(iti).eq.0.0) cycle
 
 
         do i=itube_start,itube_end
-          Etube=Etube+enetube(i)+enetube(i+nres)
+          Etube=Etube+enetube(i)+enetube(i+nres)+enecavtube(i)
         enddo
 C        print *,"ETUBE", etube
         return