Intel® Fortran Compiler
Build applications that can scale for the future with optimized code designed for Intel® Xeon® and compatible processors.

application of an elemental function to assumed-rank array

Aleksey_K_
Beginner
1,153 Views
Hi folks, I wrote some really beautiful function to compare arrays of arbitrary rank using Intel ifort 19.3. The code below compiles and runs perfectly: it gives right results for whatever integer arrays I throw into the function. However I could not find any reference that it is (current or any future) Fortran standard conforming. PURE LOGICAL FUNCTION array_is_equal(arr1,arr2) RESULT(is_equal) INTEGER, DIMENSION(..), INTENT(in) :: arr1 INTEGER, DIMENSION(..), INTENT(in) :: arr2 is_equal = .FALSE. IF (rank(arr1)==rank(arr2)) THEN IF (all(shape(arr1)==shape(arr2))) THEN is_equal = all(is_equal_val(arr1,arr2)) END IF END IF CONTAINS ELEMENTAL LOGICAL FUNCTION is_equal_val(v1,v2) INTEGER, INTENT(in) :: v1 !< First value. INTEGER, INTENT(in) :: v2 !< Second value. is_equal_val = v1==v2 END FUNCTION is_equal_val END FUNCTION array_is_equal The code uses so-called assumed-rank array so that the function can receive arrays of any rank for comparison. All I could find out is that an intrinsic inquiry function can be used with an assumed-rank array, nothing about a possibility to apply an elemental function element-wise to such an array. Therefore I wonder if it is a bug or a feature in Ifort 19.3. I really enjoy this code working correctly and hope a lot that it is a feature. It would be really sad if this function would stop working with some future release of Intel Fortran due to some fixed incompatibility with the standard. I would like to know if the function above conforms with some Fortran standard so that I can confidently use it in production code. --Aleksey
0 Kudos
21 Replies
Steve_Lionel
Honored Contributor III
1,029 Views

I agree with you - this is not allowed by the standard and the compiler should complain. I sort of understand what's happening here as the expression doesn't depend on the rank of the arrays, just the total number of elements, but that's a slippery slope. I'd encourage you to report this as a bug through the Intel Online Service Center, as I'm sure other compilers would not allow this.

0 Kudos
FortranFan
Honored Contributor II
1,029 Views

@Aleksey K.,

Line 8 of your code,

is_equal = all(is_equal_val(arr1,arr2))

appears non-conformant with the Fortran standard with Intel Fortran compiler in violation with regard to the constraint

15 C838 An assumed-rank variable name shall not appear in a designator or expression except as an actual
16      argument that corresponds to a dummy argument that is assumed-rank, the argument of the function
17      C_LOC or C_SIZEOF from the intrinsic module ISO_C_BINDING (18.2), the first dummy argument
18      of an intrinsic inquiry function, or the selector of a SELECT RANK statement

as described in Fortran 2018.

I personally think this reflects a terribly missed opportunity in Fortran 2018 to offer some valuable facilities with rank genericity in the language on its own merit rather than in situations corresponding to interoperability with a companion C processor.  It's quite a shame really since it will be year 2035 or later when the next set of compiler implementations might become viable with a friendlier future revision of the standard and when you can feel safe about pushing anything like your "beautiful" code into production i.e., if at all the language succeeds in offering better generics!  Until then, ELEMENTAL procedures and/or generic interfaces (i.e., overloaded methods) are your standard-conforming options...

As suggested, you may want to file a bug report with Intel.

 

0 Kudos
jimdempseyatthecove
Honored Contributor III
1,029 Views

I think one could inelegantly use C_F_POINTER on the CLOC of the array after obtaining the size as an arguments to the is_equal_val.

I think it a pity that the language has assumed rank, but that you cannot functionally use it in the language (you can only pass it on to something else).

Jim Dempsey

0 Kudos
Steve_Lionel
Honored Contributor III
1,029 Views

Fortran 2018 has SELECT RANK (not yet in ifort) that allows you full access, though it can still be cumbersome to deal with, especially if you have more than one assumed-rank argument. There are various proposals for the next Fortran under the general term of "generics/templates" that could make this easier. The fundamental problem is that the language is designed to always know what the rank of a data object is. As you say, there are ways around it using C interoperability.

0 Kudos
jimdempseyatthecove
Honored Contributor III
1,029 Views

I should caution on my post #4, one could have a Fortran pointer with strided ranks pointing into an array, and then this pointer is passed into a subroutine that uses assumed rank. IOW the input arg(..) does not represent a contiguous chunk of memory. Thus the C_F_POINTER on the CLOC of the first element of the input arg (strided pointer) would not represent the correct memory association of the entire argument.

IMHO Due to the descriptor being passed by reference (or a temporary descriptor being generated then passed by reference) then I see no reasonable reason why the compiler (and Language) should not support treating the array as a whole: e.g. A = B + C, or passed to an elemental routine. Additionally, one should be able to use ASSOCIATE (with limitations) to construct in effect a pointer/reference to the data blob with desired rank (some conformity test at runtime may be necessary).

Jim Dempsey

0 Kudos
Steve_Lionel
Honored Contributor III
1,029 Views

If one allowed an assumed-rank array as the argument to an elemental, the result of the call would necessarily be the same shape, but completely unknown at compile time, including rank. What could you do with it?

0 Kudos
FortranFan
Honored Contributor II
1,029 Views

Steve Lionel (Ret.) (Blackbelt) wrote:

Fortran 2018 has SELECT RANK (not yet in ifort) that allows you full access, though it can still be cumbersome to deal with, especially if you have more than one assumed-rank argument. ..

It's with the cumbersome aspect of SELECT RANK I didn't suggest the option; in many use cases, my thinking at present is generic interfaces (overloaded methods) with perhaps INCLUDE files will be the "lesser of the evil".

Anyways, the Fortran standard states:

rank.png

so as and when Intel Fortran implementation finds SELECT RANK facility, the following code (caveat: it's been tested using only my 'buggy' mental Fortran compiler!!) might prove standard-conforming and performant also - it will be up to OP then to determine whether it's "beautiful" to consider pushing to production!!!

   pure logical function array_is_equal(arr1,arr2) result(is_equal)

      integer, dimension(..), intent(in) :: arr1
      integer, dimension(..), intent(in) :: arr2

      ! Local variables
      integer :: sz_arr

      is_equal = .false.
      if ( rank(arr1)==rank(arr2) ) then
         if ( all(shape(arr1)==shape(arr2)) ) then
            sz_arr = size(arr1)
            srnk1: select rank ( a1 => arr1 )
               rank (*)
                  srnk2: select rank ( a2 => arr2 )
                     rank (*)
                        ! Note: instruction below employs 'associating entities'
                        is_equal = all(is_equal_val(a1(1:sz_arr),a2(1:sz_arr)))
                  end select srnk2
            end select srnk1
         end if
      end if

   contains

      elemental logical function is_equal_val(v1,v2)
         integer, intent(in) :: v1 !< First value.
         integer, intent(in) :: v2 !< Second value.
         is_equal_val = (v1 == v2)
      end function is_equal_val

   end function array_is_equal

 

0 Kudos
Steve_Lionel
Honored Contributor III
1,029 Views

Ah, I had forgotten about the RANK DEFAULT case. That will indeed be useful.

0 Kudos
jimdempseyatthecove
Honored Contributor III
1,029 Views

Steve Lionel (Ret.) (Blackbelt) wrote:

If one allowed an assumed-rank array as the argument to an elemental, the result of the call would necessarily be the same shape, but completely unknown at compile time, including rank. What could you do with it?

The compiler knows that the assumed rank argument is an array descriptor. With that available, the compiler can generate code that can use the descriptor information to direct the code flow through a general procedure to process the array. A although  not optimal as it could with knowing the rank (and extents of each, and stride), it still permits one to code with assumed rank (should they be inclined to do so).

Please note, that interpreters handle this. IOW completely unknown at compile time (as there is no compile time).

Please note that I am not requesting this be implemented. Once implemented, the noobs would tend to use it, then complain that the code runs slowly.

Jim Dempsey

0 Kudos
Steve_Lionel
Honored Contributor III
1,029 Views

Elemental calls aren't handled that way. The compiler generates one or more loops, calls the scalar routine many times, and builds the result one element at a time. My point above was that once you have this array of unknown rank, how do you get it into something else you can use?

0 Kudos
jimdempseyatthecove
Honored Contributor III
1,029 Views

The compiler knows the argument is an array descriptor. And in there it has everything it needs to know (# ranks, dim of each rank, and possibly stride of each rank). It is fully capable (I should say the compiler writers) could easily peel out elements of the arbitrary array descriptor. You would likely implement this as two sections of code: a) when you can determine all ranks have stride of 1 then treat it (blob of data) as a rank-1 array, or b) treat the list of ranks as a tumbler Where you iterate walking through rank like you do days on a calendar. I will see if I can write something up.

Jim Dempsey

 

0 Kudos
kondratyev__aleksey
1,029 Views

FortranFan wrote:

Quote:

Steve Lionel (Ret.) (Blackbelt) wrote:

 

Fortran 2018 has SELECT RANK (not yet in ifort) that allows you full access, though it can still be cumbersome to deal with, especially if you have more than one assumed-rank argument. ..

      is_equal = .false.
      if ( rank(arr1)==rank(arr2) ) then
         if ( all(shape(arr1)==shape(arr2)) ) then
            sz_arr = size(arr1)
            srnk1: select rank ( a1 => arr1 )
               rank (*)
                  srnk2: select rank ( a2 => arr2 )
                     rank (*)
                        ! Note: instruction below employs 'associating entities'
                        is_equal = all(is_equal_val(a1(1:sz_arr),a2(1:sz_arr)))
                  end select srnk2
            end select srnk1
         end if
      end if

Thank you for the inspiration! It looks like it would indeed be a portable way to achieve exactly the same effect as my original code intends to do (compare integer arrays of any rank without a separate case for each rank). The drawback would be the SELECT RANK and related stuff but it would be a small price for compliance with Fortran standard. We just have to wait a bit until most mainstream Fortran compilers would be able to handle SELECT RANK construct.

It is actually curious how RANK(*) would generally turn a multi-dimensional array with possibly different strides along dimensions into a single dimensional one. Maybe by some reshape-like functionality. Unfortunately it would likely require to create a copy of the original array in some cases. A fine implementation could use the same array but with reinterpreted shape as one dimensional, e.g. if the array is contiguous (an extreme case).

Anyway, the implementation of RANK(*) seems like a serious challenge to compiler developers.

--Aleksey

0 Kudos
kondratyev__aleksey
1,029 Views

jimdempseyatthecove wrote:

The compiler knows the argument is an array descriptor. And in there it has everything it needs to know (# ranks, dim of each rank, and possibly stride of each rank). It is fully capable (I should say the compiler writers) could easily peel out elements of the arbitrary array descriptor. You would likely implement this as two sections of code: a) when you can determine all ranks have stride of 1 then treat it (blob of data) as a rank-1 array, or b) treat the list of ranks as a tumbler Where you iterate walking through rank like you do days on a calendar. I will see if I can write something up.

Jim Dempsey

 

Can you write it in Fortran? AFAIK there is no way to access the descriptor of an array (including strides) from the Fortran side. I can imagine a C function implementing your proposal but not a Fortran one. It is already a very good thing that such C code can be portable since F2018 has finally provided a standard for data structures to represent array descriptors (section 18.4 in F2018 draft). However it would be truly awesome to avoid any C code and write this functionality in Fortran alone.

--Aleksey

0 Kudos
Steve_Lionel
Honored Contributor III
1,029 Views

[Deleted]

0 Kudos
FortranFan
Honored Contributor II
1,029 Views

FortranFan wrote:

.. It's with the cumbersome aspect of SELECT RANK I didn't suggest the option ..

.. as and when Intel Fortran implementation finds SELECT RANK facility, the following code (caveat: it's been tested using only my 'buggy' mental Fortran compiler!!) might prove standard-conforming ..

As I stated above, SELECT RANK is cumbersome and as to whether it's to be preferred over generic interfaces might be more of a matter of taste.

More importantly, I now notice the verbiage I missed out earlier in the Fortran standard with "A RANK ( * ) statement matches the selector if the selector is argument associated with an assumed-size array" and which indicates the code I showed in Quote #8 is faulty.  Something like the variant below might be the conforming version:

program p
   
   implicit none
   
   blk1: block
      integer, allocatable :: foo(:), bar(:)
      foo = [ 1, 2, 3 ]
      bar = [ 1, 1, 1 ]
      print *, "Block 1: Rank 1 with same shape, different values"
      print *, "foo = bar? ", array_is_equal(foo, bar)
      print *, "Expected is F.", new_line("")
   end block blk1
   
   blk2: block
      integer, allocatable :: foo(:,:), bar(:,:)
      integer :: i
      foo = reshape( [( i, i = 1, 2*3 )], shape=[ 2, 3 ] )
      bar = foo
      print *, "Block 2: Rank 2 with same shape and values"
      print *, "foo = bar? ", array_is_equal(foo, bar)
      print *, "Expected is T.", new_line("")
   end block blk2
   
   blk3: block
      integer, allocatable :: foo(:,:,:), bar(:,:,:)
      integer :: i
      foo = reshape( [( i, i= 1,1*2*3 )], shape=[ 1, 2, 3 ] )
      bar = reshape( [( i, i = 1, 2*3*4 )], shape=[ 2, 3, 4 ] )
      print *, "Block 3: Rank 3 with different shape"
      print *, "foo = bar? ", array_is_equal(foo, bar)
      print *, "Expected is F.", new_line("")
   end block blk3
   
   blk4: block
      integer :: foo, bar
      foo = 42 ; bar = foo
      print *, "Block 4: Rank 0 with same values"
      print *, "foo = bar? ", array_is_equal(foo, bar)
      print *, "Expected is T", new_line("")
   end block blk4
   
   stop

contains

   logical function array_is_equal(arr1,arr2) result(is_equal)

      integer, dimension(..), intent(in) :: arr1
      integer, dimension(..), intent(in) :: arr2

      is_equal = .false.
      if ( rank(arr1)==rank(arr2) ) then
         if ( all(shape(arr1)==shape(arr2)) ) then
            srnk1: select rank ( a1 => arr1 )
               rank ( 0 )
                  select rank ( a2 => arr2 )
                     rank ( 0 )
                        ! Note: instruction below employs 'associating entities'
                        is_equal = is_equal_val(a1,a2)
                  end select
               rank ( 1 )
                  select rank ( a2 => arr2 )
                     rank ( 1 )
                        ! Note: instruction below employs 'associating entities'
                        is_equal = all( is_equal_val(a1,a2) )
                  end select
               rank ( 2 )
                  select rank ( a2 => arr2 )
                     rank ( 2 )
                        ! Note: instruction below employs 'associating entities'
                        is_equal = all( is_equal_val(a1,a2) )
                  end select
               rank ( 3 )
                  select rank ( a2 => arr2 )
                     rank ( 3 )
                        ! Note: instruction below employs 'associating entities'
                        is_equal = all( is_equal_val(a1,a2) )
                  end select
            end select srnk1
         end if
      end if

   end function array_is_equal

   elemental logical function is_equal_val(v1,v2)
      integer, intent(in) :: v1 !< First value.
      integer, intent(in) :: v2 !< Second value.
      is_equal_val = (v1 == v2)
   end function is_equal_val
   
end

The new version of Intel Fortran Compiler 19.1 Pre-Release Beta works as I expect with the above code - see below:

C:\Temp>ifort /standard-semantics /warn:all /stand:f18 p.f90
Intel(R) Visual Fortran Intel(R) 64 Compiler for applications running on Intel(R) 64, Version 19.1.0.056 Pre-Release Beta Build 20190321
Copyright (C) 1985-2019 Intel Corporation.  All rights reserved.

ifort: NOTE: The Beta evaluation period for this product ends on 9-oct-2019 UTC.
Microsoft (R) Incremental Linker Version 14.16.27027.1
Copyright (C) Microsoft Corporation.  All rights reserved.

-out:p.exe
-subsystem:console
p.obj

C:\Temp>p.exe
 Block 1: Rank 1 with same shape, different values
 foo = bar?  F
 Expected is F.

 Block 2: Rank 2 with same shape and values
 foo = bar?  T
 Expected is T.

 Block 3: Rank 3 with different shape
 foo = bar?  F
 Expected is F.

 Block 4: Rank 0 with same values
 foo = bar?  T
 Expected is T


C:\Temp>

 

0 Kudos
Aleksey_K_
Beginner
1,029 Views

FortranFan wrote:

As I stated above, SELECT RANK is cumbersome and as to whether it's to be preferred over generic interfaces might be more of a matter of taste.

Well, generality of any form tends to indulge sloth: I would  prefer to write the compaison function once, not to add a new piece of code every time an array of a higher rank (more than three in your sample code) is to be handled.

FortranFan wrote:

More importantly, I now notice the verbiage I missed out earlier in the Fortran standard with "A RANK ( * ) statement matches the selector if the selector is argument associated with an assumed-size array" and which indicates the code I showed in Quote #8 is faulty.  Something like the variant below might be the conforming version:

Uh, this is really cumbersome. A more viable alternative is desperately needed. Possibly a small function in C that receives the array descriptors CFI_cdesc_t and works with them to run through arrays would do the job in a small piece of code to be written once and would never be extended no matter what arrays are thrown to the function. AFAIK, F2018 guarantees portability of CFI_cdesc_t and related declarations from ISO_Fortran_binding.h header. This way would probably have an advantage over SELECT RANK for maintenance of the source. If there was any chance to write the same in Fortran, I would gladly go for it but there does not seem to be any at the moment.

Surely it would still be important that the C code runs fast but it is an extensive and very different topic: How to introduce vectorization opportunities in the C code and how well would this opportunities be exploited by current C compilers.

At least the solution with a small C function should never become as bad as the code in some very old project that I have seen: The whole thing was written in F77, surely without any possibilities for dynamic memory. However the developers found out that they can use malloc() and free() functions from C if they link their code with C run-time. The whole mixture of F77 with calls to C heap management functions looked weird to say the least, mostly the tricks to turn the result of a call to malloc() into a F77 array ;)

--Aleksey

0 Kudos
FortranFan
Honored Contributor II
1,029 Views

Aleksey K. wrote:

 .. Possibly a small function in C ..

..

Well, if one is willing to bring in a function in C which conceivably will bring in the requirements of interoperable types then, one can short-circuit any such pursuit by employing interoperable types to being with and using intrinsic procedures from ISO_C_BINDING module to setup pointers without copying any data and achieve the desired result.  See below, I think it conforms to the standard:

program p

   use, intrinsic :: iso_c_binding, only : c_int
   
   implicit none
   
   blk1: block
      integer(c_int), allocatable :: foo(:), bar(:)
      foo = [ 1, 2, 3 ]
      bar = [ 1, 1, 1 ]
      print *, "Block 1: Rank 1 with same shape, different values"
      print *, "foo = bar? ", array_is_equal(foo, bar)
      print *, "Expected is F.", new_line("")
   end block blk1
   
   blk2: block
      integer(c_int), allocatable :: foo(:,:), bar(:,:)
      integer :: i
      foo = reshape( [( i, i = 1, 2*3 )], shape=[ 2, 3 ] )
      bar = foo
      print *, "Block 2: Rank 2 with same shape and values"
      print *, "foo = bar? ", array_is_equal(foo, bar)
      print *, "Expected is T.", new_line("")
   end block blk2
   
   blk3: block
      integer(c_int), allocatable :: foo(:,:,:), bar(:,:,:)
      integer :: i
      foo = reshape( [( i, i= 1,1*2*3 )], shape=[ 1, 2, 3 ] )
      bar = reshape( [( i, i = 1, 2*3*4 )], shape=[ 2, 3, 4 ] )
      print *, "Block 3: Rank 3 with different shape"
      print *, "foo = bar? ", array_is_equal(foo, bar)
      print *, "Expected is F.", new_line("")
   end block blk3
   
   blk4: block
      integer(c_int) :: foo, bar
      foo = 42 ; bar = foo
      print *, "Block 4: Rank 0 with same values"
      print *, "foo = bar? ", array_is_equal(foo, bar)
      print *, "Expected is T", new_line("")
   end block blk4
   
   stop

contains

   logical function array_is_equal(arr1,arr2) result(is_equal)

      integer(c_int), intent(in), target :: arr1(..)
      integer(c_int), intent(in), target :: arr2(..)

      ! Local variables
      integer :: sz_arr

      is_equal = .false.
      if ( rank(arr1)==rank(arr2) ) then
         if ( all(shape(arr1)==shape(arr2)) ) then
            sz_arr = size(arr1)
            blk: block
               use, intrinsic :: iso_c_binding, only : c_loc, c_f_pointer
               integer(c_int), pointer :: a1(:) => null()
               integer(c_int), pointer :: a2(:) => null()
               call c_f_pointer( cptr=c_loc(arr1), fptr=a1, shape=[ sz_arr ] )
               call c_f_pointer( cptr=c_loc(arr2), fptr=a2, shape=[ sz_arr ] )
               is_equal = all( is_equal_val(a1, a2) )
               a1 => null()
               a2 => null()
            end block blk
         end if
      end if

   end function array_is_equal

   elemental logical function is_equal_val(v1,v2)
      integer, intent(in) :: v1 !< First value.
      integer, intent(in) :: v2 !< Second value.
      is_equal_val = (v1 == v2)
   end function is_equal_val
   
end

C:\Temp>ifort /standard-semantics /warn:all x.f90
Intel(R) Visual Fortran Intel(R) 64 Compiler for applications running on Intel(R) 64, Version 19.0.3.203 Build 20190206
Copyright (C) 1985-2019 Intel Corporation.  All rights reserved.

Microsoft (R) Incremental Linker Version 14.16.27027.1
Copyright (C) Microsoft Corporation.  All rights reserved.

-out:x.exe
-subsystem:console
x.obj

C:\Temp>x.exe
 Block 1: Rank 1 with same shape, different values
 foo = bar?  F
 Expected is F.

 Block 2: Rank 2 with same shape and values
 foo = bar?  T
 Expected is T.

 Block 3: Rank 3 with different shape
 foo = bar?  F
 Expected is F.

 Block 4: Rank 0 with same values
 foo = bar?  T
 Expected is T


C:\Temp>

 

0 Kudos
kondratyev__aleksey
1,029 Views

FortranFan wrote:

Quote:

Well, if one is willing to bring in a function in C which conceivably will bring in the requirements of interoperable types then, one can short-circuit any such pursuit by employing interoperable types to being with and using intrinsic procedures from ISO_C_BINDING module to setup pointers without copying any data and achieve the desired result.  See below, I think it conforms to the standard:

Your code is fine. It looks compliant with the standard although I  am no expert in Fortran standard (just a software developer;). It will work perfectly with contiguous arrays. However the function is_equal_val will compare a lot of junk once the arrays arr1 or arr2 have "non-trivial" strides, i.e. there are memory gaps between neighboring elements in arrays. I do not see another way to handle such a situation but to read the strides from array descriptor CFI_cdesc_t, namely it "sm" field of "dim[]" field of CFI_cdesc_t, all from ISO_Fortran_binding.h

typedef struct CFI_dim_t_ {
    CFI_index_t     extent;
    CFI_index_t     sm;
    CFI_index_t     lower_bound;
} CFI_dim_t;

typedef struct CFI_cdesc_t_ {
    void * base_addr;
    size_t elem_len;
    int version;
    CFI_attribute_t attribute;
    CFI_rank_t rank;
    CFI_type_t type;
    // Begin Intel-specific fields - these are for
    // Intel compiler/library use only and should 
    // not be modified
    . . . . .
    /*  Intel-specific fields */
    CFI_dim_t dim[];
} CFI_cdesc_t;

Note the stride in sm is in bytes. I do not think there is any sane way to make use of such kind of data in Fortran. I can imagine a weird way to access sm from Fortran by wrapping the descriptor types (which is IMHO a terrible idea anyway)

TYPE, BIND(C) :: F_CFI_dim_t
 . . .
END TYPE F_CFI_CFI_dim_t

TYPE, BIND(C) :: F_CFI_cdesc_t
 . . .
END TYPE F_CFI_cdesc_t

but it would require a very ugly and presumably inefficient code in Fortran to use the values of "sm" (given in bytes!) for different dimensions to calculate the memory address of each element for comparison.

Interestingly, array stride is the only array characteristic in Fortran without an intrinsic function to access it. Probably it is just because there is nothing meaningful to do with a stride in bytes on the Fortran side. At least without terrible hacks that use transfer, c_f_pointer and things like that. Note that transfer, c_f_pointer, etc are perfectly fine functions, still the use of them with stride information to decode the layout of an array manually (i.e. not by compiler generated code) would probably be evil if not to say more ;)

Imho, C is a proper language for such things. Therefore, it still looks like a viable option to compare arrays by a small C function which

  1. Receives a pair of arrays by their descriptor pointers CFI_cdesc_t*.
  2. Gathers all needed info about the layout of array elements by traversing the descriptors, especially the dim[] fields.
  3. Accesses both arrays element-wise using the layout info collected from array descriptors.
  4. Does the comparison element-wise.

It would be challenging (imho still possible) to this equip such a function in C with an ample of vectorization opportunities.

Or, the ideal case: somebody with an influence on Fortran standard board would read this branch and propose an intrinsic that does the four steps above. This way Fortran developers would have to write the code to iterate over arrays of an arbitrary dimension and any stride along each dimension, not me.

--Aleksey

0 Kudos
andrew_4619
Honored Contributor II
1,029 Views

Would I be correct in assuming that   all( is_equal_val(a1, a2) ) would always cause the evaluation of the equality test for all elements of the array? If so then in the general case it is going to be very inefficient with large arrays as we lose the possibility of early exit.  

0 Kudos
Aleksey_K_
Beginner
937 Views

andrew_4619 wrote:

Would I be correct in assuming that   all( is_equal_val(a1, a2) ) would always cause the evaluation of the equality test for all elements of the array? If so then in the general case it is going to be very inefficient with large arrays as we lose the possibility of early exit.  

There can be no assumptions about the number or the order of calls of is_equal_val. That is why is_equal_val must be elemental. It is up to implementation of intrinsic functions like "all" in a specific Fortran compiler. The compiler may choose to call is_equal_val in any order and any number of times. AFAIK, we can only hope that compiler developers are very smart people and would employ vectorization along with some kind of short circuit evaluation of booleans inside an intrinsic like "all" or "any" and without the necessity to build the whole logical array as the argument for the intrinsic, at least with optimization enabled. The good news is that modern Fortran compilers are very capable in this respect.

Note that intrinsic functions are kind of built in the language. Therefore the compiler knows precisely what the intrinsic does and can exploit this knowledge for a truly impressive optimization.

I did not run performance tests. Therefore I can not say for sure but I an willing to bet that a variant with all() intrinsic would generally perform better than a DO loop equivalent with explicit short circuit evaluation of conjunction (logical .AND.).

It is also interesting that all() intrinsic is allowed inside OMP WORKSHARE construct. I wonder if such a use of openmp

!$OMP PARALLEL WORKSHARE
result = all(is_equal_val(a1,a2))
!$OMP END PARALLEL WORKSHARE

would ever give any sensible performance gain unless is_equal_val is computationally expensive or some very exotic platform is used.

--Aleksey

0 Kudos
Reply