anchor.diff:  
   
   subroutine integrate(f, eu, ek, part, vel, npart, a, r_cut2, e_cut, dt, e_tot, temp_k)
  -integer                         :: npart
  -real(pr)                        :: r_cut2, a
  +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
  
positive.code:                 f(3*j - 1) = f(3*j - 1) - ff*dy
                  f(3*j)     = f(3*j)     - ff*dz
  
                  eu = eu + (4._pr*r6i*(r6i-1._pr) - e_cut)
              end if
          end do
  
negative.code: end subroutine vel_init
  
  
  subroutine force(part, npart, a, f, r_cut2, e_cut, eu, p)
  integer                         :: i, j
  integer                         :: npart
  real(pr)                        :: r2, r_cut2, a, ff, r2i, r6i
  real(pr)                        :: eu, e_cut, l, ln
  real(pr)                        :: dx, dy, dz. p
  real(pr), dimension(1:3*npart)  :: part, f
  
      f = 0
      p = 0
      eu = 0._pr
      ln = real(npart/4, pr)**(1./3.)
      l = ln * a
  