anchor.diff:                  f(3*j - 1) = f(3*j - 1) - ff*dy
                   f(3*j)     = f(3*j)     - ff*dz
   
  +                p = p + ff*(dx+dy+dz)
  +
                   eu = eu + (4._pr*r6i*(r6i-1._pr) - e_cut)
               end if
           end do
  
positive.code: 
      ! nuevas fuerzas
      f_old = f
      call force(x_new, npart, a, f, r_cut2, e_cut, eu)
  
      !nuevas velocidades
      v_new  = vel + dt * (f + f_old)/2._pr
      sumvv2 = sum(v_new**2)
  
      ! calculo temperatura y energia
      temp_k = sumvv2/real(3*npart, pr)
      ek = 0.5_pr * sumvv2
      e_tot  = eu + ek
  
  
negative.code: 
  
  subroutine integrate(f, eu, ek, part, vel, npart, a, r_cut2, e_cut, dt, e_tot, temp_k)
  integer                         :: npart, ii, jj
  real(pr)                        :: r_cut2, a, p
  real(pr)                        :: eu, ek, e_cut, temp_k, e_tot
  real(pr)                        :: sumvv, sumvv2, dt
  real(pr), dimension(1:3*npart)  :: f_old, x_new, v_new
  