anchor.diff:                      MPI_COMM_WORLD, stat, err)
   
                  if (rank == i) then
  -                  r_other = r
  +                  r_i = r
  +               else if (rank == nprocs - i)
  +                  r_np_i = r
                  end if
   
                  ! Broadcast i-th data
  -               call mpi_bcast(r_other, 3*N, MPI_REAL_XP, i, MPI_COMM_WORLD, err)
  +               call mpi_bcast(r_i, 3*N, MPI_REAL_XP, i, MPI_COMM_WORLD, err)
  +               ! Broadcast n-i-th data
  +               call mpi_bcast(r_np_i, 3*N, MPI_REAL_XP, nprocs - i, MPI_COMM_WORLD, err)
   
                  ! If i == rank, compute both diagonals
                  if (rank == i) then
  -                  r_right = r_other
  +                  r_right = r_i
                     call compute_force_diag(m, r, r, a)
                     call compute_force_diag(m, r_right, r_right, a_right)
                     a_comm = a
   
                     ! compute on right side but communicate nothing
                  else if (rank < i) then
  -                  call compute_force(m, r_other, r_right, a_right)
  +                  call compute_force(m, r_np_i, r_right, a_right)
                     a_comm = 0._xp
   
                     ! compute on left side and communicate interaction
                  else
  -                  call compute_force(m, r, r_other, a_comm)
  +                  call compute_force(m, r, r_i, a_comm)
                     a = a - a_comm
                  end if
   
  
positive.code: 
        call integrate(v, a, dt/2)  ! Compute v(t+dt/2)
        call integrate(r, v, dt)    ! Compute r(t+dt)
  
  
        select case(flag_compute_mpi)
           case(0)
  
negative.code:             stop "Unknown value of flag_compute_mpi"
        end select
        !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ FIX BELOW
        !--------------------------------
        ! reduce accelerations
        !--------------------------------
  