fortran 2003 and pointers

Question

I am writing a Fortran program. The program implements some numerical methods. Program speed is very important. I decided to get rid of dynamic arrays (whether it speeds up the program?), and faced with the following problem. I have 3d arrays (NXxNYxNZ = MAX elements) where I know MAX, but I dont know NX/NY/NZ ratio. It can be like this : 1x1xNZ or like this 2xNYx1 , etc. The solution I see - use pointers. Simplified 2D case:

program ptrtest
   parameter ( MAX = 50 )    ! I know this number on compile time.
   integer :: NX,NY          ! I will find this numbers on run time.
   real, target :: a(MAX)    ! Static Array
   real, pointer :: b(:,:)
   a = 5
   read(*,*) NX, NY          ! I get NX, NY (NX*NY == MAX)
   b (1:NX, 1:NY) => a       ! I can use b(:,:) <- this is my goal.
end program ptrtest

This example works, but i am afraid that such update slow down the real program where I am using even 5d arrays. Is it possible?

You don't show enough program to really comment (Are you working with a and b at the same time? Are you immediately passing b to a subroutine?), but in general performance will be worse when working with variables that have the TARGET or POINTER attribute because certain optimisations may be prevented. It will depend on specifics, but allocating an allocatable array of size (nx,ny) once is not likely to be that significant a performance hit relative to to something like the READ statement that populates nx and ny. — IanH, Dec 17 '14 at 23:58
If you need an array with run-time determined dimensions, I recommend using an allocatable array rather than a pointer array, unless you also need an additional feature of pointer. allocatables are safer. Also, its best to solve your problem, then identify slow portions, rather than guess at the start. Micro-optimizations rarely make much difference. Selecting an efficient algorithm is important. — M. S. B., Dec 18 '14 at 04:45
FWIW I concur with the other commentators. But if you know `MAX` at compile-time but not, until run-time, `nx`,`ny`, and `nz` you can always use the tried-and-tested (for which understand 'old-fashioned') approach of using a 1D array (of `MAX` elements) and roll-your-own index arithmetic. If considering this approach give careful thought to @M.S.B.'s observations on optimisation. — High Performance Mark, Dec 18 '14 at 08:25

score 3 · Accepted Answer · edited May 23 '17 at 10:26

You say that program speed is important, so old-style Fortran will give you the best:

allocate a 1D array of size MAX
pass is to a subroutine as well as NX,NY,and NZ. This will avoid the pointer indirection that will give you a loss of performance (unlike C, the Fortran standard assumes that subroutine array arguments don't overlap, see this post).

For example:

program noptrtest
   implicit none
   integer, parameter :: MAX = 50     ! I know this number on compile time.
   integer :: NX,NY          ! I will find this numbers on run time.
   real    :: a(MAX)         ! Static Array

   read(*,*) NX, NY          ! I get NX, NY (NX*NY == MAX)

   if (NX*NY /= MAX) then
     stop 'NX*NY /= MAX'
   endif
   call run(a,NX,NY)
end program noptrtest


subroutine run(b,NX,NY)
  integer :: NX,NY
  real    :: b(NX,NY)
  integer :: i,j

  do j=1,NY
    do i=1,NX
      ! use b(i,j) here
    enddo
  enddo
end

If performance really matters here are some other useful tricks:

Tell your compiler to align the arrays on a 32-byte boundary (with ifort, use the !DIR$ ATTRIBUTES ALIGN : 32 directive
Dimension physically your 2D array such that the size of leading dimension is a multiple of 32-bytes. For real*4, you need a multiple of 8 elements.
Tell the compiler that each column is properly aligned using the !DIR$ VECTOR ALIGNED directive

For example:

program noptrtest
   implicit none
   integer, parameter :: real_kind = 4 
   integer, parameter :: MAX = 50               ! I know this number on compile time.
   integer, parameter :: MAXDIM = MAX*(32/real_kind)    ! Max possible dimension required
   integer            :: NX,NY, NX_MOD         ! I will find this numbers on run time.
   real(real_kind)    :: a(MAXDIM)             ! Static Array

   !DIR$ ATTRIBUTES ALIGN : 32 :: a

   read(*,*) NX, NY          ! I get NX, NY (NX*NY == MAX)

   if (NX*NY /= MAX) then
     stop 'NX*NY /= MAX'
   endif

   if (mod(NX,real_kind) == 0) then
     NX_MOD = NX
   else
     NX_MOD = ((NX/real_kind)+1)*real_kind
   endif

   call run(a,NX_MOD,NX,NY)
end program noptrtest


subroutine run(b,NX_MOD,NX,NY)
  integer :: NX_MOD,NX,NY
  real    :: b(NX_MOD,NY)
  integer :: i,j

  do j=1,NY
    !DIR$ VECTOR ALIGNED
    do i=1,NX
    ! use b(i,j) here
    enddo
  enddo
end

Edit

This old-style Fortran trick avoids pointer aliasing.

References on pointer aliasing:

http://en.wikipedia.org/wiki/Aliasing_(computing) , section "Conflicts with optimization"
gcc C/C++ assume no pointer aliasing
Performance impact of -fno-strict-aliasing
http://scc.qibebt.cas.cn/docs/compiler/intel/11.1/Intel%20C++%20Compiler%20for%20Linux/main_cls/optaps/common/optaps_perf_run.htm

*so old-style Fortran will give you the best* Show us the numbers; without data this 'answer' is an opinion. — High Performance Mark, Dec 18 '14 at 13:02
It's not an opinion.Look at the -fno-alias option of C compilers for more explanations. As I said : "(unlike C, the Fortran standard assumes that subroutine array arguments don't overlap)". Pointer aliasing is a well known problem for compilers, and it was not present in Fortran before they introduced pointers. See — Anthony Scemama, Dec 18 '14 at 16:35

fortran 2003 and pointers

1 Answers1

Edit