Canonicalise WHERE with reduction intrinsics, which currently are placed into CodeBlocks

[sergisiso]

In ECMWF NEMO sbccpl.f90 this WHERE construct.

WHERE( picefr(:,:) > 1.e-10 ) ; zevap_ice(:,:,1) = frcv(jpr_ievp)%z3(:,:,1) * SUM( a_i_last_couple, dim=3 ) / picefr(:,:)
ELSEWHERE                     ; zevap_ice(:,:,1) = 0._wp
END WHERE

is converted to:

do widx2_2 = 1, SIZE(picefr, 2), 1
   do widx1_2 = 1, SIZE(picefr, 1), 1
      if (picefr(widx1_2,widx2_2) > 1.e-10) then                                                                                        1               
          zevap_ice(widx1_2,widx2_2,1) = frcv(jpr_ievp)%z3(widx1_2,widx2_2,1) * SUM(a_i_last_couple, dim=3) / picefr(widx1_2,widx2_2)
      else
          zevap_ice(widx1_2,widx2_2,1) = 0._wp
      end if
   enddo
enddo

Which produces the compiler error: The shapes of the array expressions do not conform. [ZEVAP_ICE]

I am not sure what the error is, but it is a nested WHERE (#1651 ?) with array elements ( #717 ?) and it uses 1, SIZE instead of L/UBound which could cause troubles with non-standard bounds.

[arporter]

I suspect it is the SUM(a_i_last_couple, dim=3) returning an array - essentially it reduces into a 2D (temporary) array which must then be indexed inot using widx1_2 and widx2_2. The only way for us to do this is to introduce the temporary array and replace the SUM(...) with it in the WHERE expression. Basically, at the moment if we see a SUM inside a WHERE we're going to have to use a CodeBlock (because to introduce the temporary will require type information - #1799).

[arporter]

Simon has also reported:

where(sum(v1(:), dim=2) > 0.0)   ! -> Issue 3: "do widx1 = 1, SIZE(v1, 1), 1" has incorrect bounds, it should be
                                 !             "do widx1 = lbound(v1, 1), rbound(v1, 1), 1"
   v1(:) = 1.0
end where

where(sum(v2(:,:), dim=2) > 0.0) ! -> Issue 4: the reduction "sum(v2(:,:), dim=2)" is ignored
  v3(:) = 1.0
end where

[arporter]

Looking at the intrinsics in fparser, I realise we are going to have similar problems with MAXVAL, MINVAL and PRODUCT.

[sergisiso]

As we talk about, these issues are also present in the arrayrange2loop transformation, although we correctly refuse them in the validate because they contain non-elemental intrinsics.

If we could convert

a(:) = b(:) * sum(c, dim=2 ) 

to

t = sum(c, dim=2)
do i = lbound(a,1), ubound(a,1)
    a(i) = b(i) * t(i)
end do

maybe we could reuse that transfomration to convert the wheres, but we need:

• Add support for querying the data type of a PSyIR expression. #1799
• Rationalise ArrayRange2Loop transformations. #2004

(trying to use the transformation for the where convertion can already be started by putting TransformationErrors to CodeBlocks - but it will have some limitations at least until #2004 )

[sergisiso]

@arporter Having a second thought on the previous comment, wouldn't the array to loop conversion from

a(:) = b(:) * sum(c, dim=2 ) 

be achieved without a temporary by doing

do i = lbound(a,1), ubound(a,1)
    a(i) = b(i) * sum(c(i,:))
end do

which it is also more memory efficient because skips the temporary store.

[arporter]

Ooh, that's clever.

[rupertford]

sum2code also does this I think but replaces the final sum as it is not sum2sum! So this functionality could be used here or extracted from here.

[sergisiso]

PR #2372 fixes some bugs on the WHERE parsing and adds the construct into a CodeBlock when there is a reduction that we can not safely canonicalise. This prevents the silent bug referenced at the top of this issue, but we leave the issue open until we can get rid of the CodeBlock.

[LonelyCat124]

PR #2883 handles the xfailing test for this issue, but I'm not sure it necessarily fixes the cases described here.