anchor.diff:  
      implicit none
   
  -   real(kind = xp) :: m(npoints)                 ! Masses of the particles
  -   real(kind = xp) :: r(3,npoints)               ! Positions of the particles (3-dim vectors)
  -   real(kind = xp) :: v(3,npoints)               ! Speeds of the particles (3-dim vectors)
  -   real(kind = xp) :: a(3,npoints)               ! Acceleration of the particles (3-dim vectors)
  +   real(kind = xp), allocatable :: m(:)                 ! Masses of the particles
  +   real(kind = xp), allocatable :: r(:, :)               ! Positions of the particles (3-dim vectors)
  +   real(kind = xp), allocatable :: v(:, :)               ! Speeds of the particles (3-dim vectors)
  +   real(kind = xp), allocatable :: a(:, : )               ! Acceleration of the particles (3-dim vectors)
      real(kind = xp), allocatable :: a_comm(:, :), a_right(:, :)
      real(kind = xp), allocatable :: r_i(:, :), r_np_i(:, :), r_right(:, :)
  -   real(kind = xp) :: a_reduced(3,npoints)
  +   real(kind = xp), allocatable :: a_reduced(:, :)
      real(kind = xp) :: Ec                         ! Total kinetic energy
      real(kind = xp) :: Ep                         ! Total potential energy
      real(kind = xp) :: E                          ! Total energy
  
positive.code:    N = ceiling(1._xp * npoints / nprocs)
  
     !---------------------------------------------
     ! Read initial positions
     !---------------------------------------------
     open(newunit=un, file='initial_conditions.dat', status="old")
  
negative.code:    iend = npoints
     select case (flag_compute_force)
        case(0)
           call compute_force(m, r, istart, iend, r, istart, iend, a)
        case(1)
           call compute_force_omp(m, r, istart, iend, r, istart, iend, a)
        case(2)
           call compute_force_omp_nn_1(m, r, istart, iend, r, istart, iend, a)
        case default
           stop "Unknown value of flag_compute_force"
     end select
  