Coarray nested loop

I am currently learning Coarray. Could anyone please show me how to translate the codes below to using Coarray?

  integer(8) :: n, m, i, j
  real(8) :: x(n, m), res
  
  res = 0d0

  !$omp parallel do collapse(2) reduction(+:res)
  do j = 1, m
    do i = 1, n
      res = res + x(i, j) 
    end do
  end do
1 Like

What have you tried so far?

For this specific problem (parallel sum), you don’t need coarrays. You can work with regular arrays and use Fortran 2018 collective subroutines. Here are a few hints:

  1. Declare array x such that each image has its portion of the domain. To do this you’ll need to make the array allocatable, and then allocate it on each image with start and end indices that depend on the values of this_image() and num_images(). Edit: start and end indices on each image don’t need to be different across images, the local arrays only need to have the same shape. And if you plan to work with coarrays, they will need to have the same extent (start and end indices) across images.
  2. Assign meaningful values to x on each image.
  3. Use the co_sum intrinsic subroutine to do the parallel reduction (sum). Careful with co_sum and other collective subroutines as they update the first argument in place on all images by default, or on a specified image if you use the result_image optional argument. See co_sum docs by GFortran.
  4. Now you have a section of your array summed across images, and you can use the local sum over that section to get the ultimate sum of the distributed array.
2 Likes

As stated by Milan, you could replace OpenMP directives by collective subroutines. I had started a project to explore coarrays and collective subroutines:

You can compare those two files:
https://github.com/vmagnin/exploring_coarrays/blob/main/pi_monte_carlo_openmp.f90
https://github.com/vmagnin/exploring_coarrays/blob/main/pi_monte_carlo_co_sum.f90

I have no time to progress further for the moment, but it can help you.

2 Likes

Like @milancurcic said you probably don’t need coarray here (or at least you don’t need an array explicitly declared with coindices). Here is an simple example following milan’s suggestions, but it’s for one dimensional array. For a 2D array it’s a bit more complicated as you have to judge if size(arr, dim=1) or size(arr, dim=2) is divisible by num_images() but the idea is the same.

module test_coarray
    use iso_fortran_env, only: i32 => int32, sp => real32
    implicit none

contains

    !> Add every element of an array to a value
    function co_add(val, arr) result(ret)
        real(sp), intent(in) :: val, arr(:)
        real(sp) :: ret
        integer(i32) :: njob, nrem, i

        ! Initialization
        ret = val 

        ! Number of jobs per image.
        njob = size(arr)/num_images() 

        ! For each image, sum part of the array.
        associate (x => this_image()*njob)
            ret = ret + sum([(arr(i), i=x - njob + 1, x)])
        end associate

        ! If size of the array is not divisible by num_images().
        nrem = mod(size(arr), num_images())
        if (nrem /= 0 .and. this_image() <= nrem) then
            ret = ret + arr(njob*num_images() + this_image())
        end if

        ! Wait for other images
        sync all 

        ! Result is stored on image 1.
        call co_sum(ret, result_image=1)
    end function co_add

end module test_coarray

program main

    use test_coarray
    implicit none

    real :: val

    val = co_add(0., [1., 2., 3., 4., 5., 6., 7., 8., 9., 10.])
    if (this_image() == 1) print *, val

end program main
1 Like