Reputation: 379
Can I have a pointer array that points to an incontiguous section of a target array? In the following program:
program test
implicit none
integer :: i
integer, dimension(4), target :: orig
integer, dimension(:), pointer :: ref
orig = (/ (i , i = 1,4) /)
ref(1:2) => orig(1:2)
ref(3:3) => orig(4:4)
print *, 'orig = ', orig
print *, 'ref = ', ref
end program test
I want the pointer ref
to point to the incontiguous section (/1,2,4/)
of the target array orig
. The code compiles without a warning, however the output is not desired:
code output:
orig = 1 2 3 4
ref = 4
It appears as the compiler forces the pointer to point to a contiguous subsection even though I first point to indices 1:2
but compiler only forces pointing to 4:4
. As a result my pointer at the end only points to one entry.
I can pack the entries (/1,2,4/)
into another array and then use it contiguously but this would not be efficient Fortran coding for large data set. Do you have any idea how to make this working?
Upvotes: 3
Views: 183
Reputation: 379
Following "@IanH" suggestion, I implemented a wrapper which only has pointers to indices and data and everything is done IN PLACE so this should be expected to be more appropriate for large data. For small data set, I think "@Vladimir F" idea of vector indexing is better and more readable. The following code passed valgrind full memleak test and I tried for large data set (6GB double precision) and it doesn't copy data.
module any_stride
implicit none
private
! ------------------------------------------------------------
! data = pointer to the contiguous 1D double precision array
! indx = pointer to the non-uniform indices
! example: for a stride where:
! stride%data => (/ 1.0d0, 2.0d0, 3.0d0, 4.0d0 /)
! stride%indx => (/ 1, 4 /)
! is a wrapper around 1d data array with four elements
! that only elements 1 and 4 are visible to the outside world
!
type stride
private
real*8, dimension(:), pointer, public :: data => null()
integer, dimension(:), pointer, public :: indx => null()
contains
procedure :: assgn_from_scalar
procedure :: assgn_from_array
! ...
! ...
! ... add more overloaded procs if you want
generic :: assignment(=) => assgn_from_scalar &
, assgn_from_array
generic :: operator(+) => add_to_scalar &
, add_to_array
! ...
! ...
! ... add more overloaded procs if you want
end type stride
public :: stride
contains
subroutine assgn_from_scalar(this, rhs)
implicit none
class(stride), intent(inout) :: this
real*8, intent(in) :: rhs
! local vars
integer :: i
do i = 1, size(this%indx)
this%data(this%indx(i)) = rhs
end do
! done here
end subroutine assgn_from_scalar
subroutine assgn_from_array(this, rhs)
implicit none
class(stride), intent(inout) :: this
real*8, dimension(:), intent(in) :: rhs
! local vars
integer :: i
do i = 1, size(this%indx)
this%data(this%indx(i)) = rhs(i)
end do
! done here
end subroutine assgn_from_array
function add_to_scalar(lhs, rhs)
implicit none
class(stride), intent(in), target :: lhs
real*8, intent(in) :: rhs
class(stride), pointer :: add_to_scalar
! local vars
integer :: i
do i = 1, size(lhs%indx)
lhs%data(lhs%indx(i)) = lhs%data(lhs%indx(i)) + rhs
end do
add_to_scalar => lhs
! done here
end function add_to_scalar
function add_to_array(lhs, rhs)
implicit none
class(stride), intent(in), target :: lhs
real*8, dimension(:), intent(in) :: rhs
class(stride), pointer :: add_to_array
! local vars
integer :: i
do i = 1, size(lhs%indx)
lhs%data(lhs%indx(i)) = lhs%data(lhs%indx(i)) + rhs(i)
end do
add_to_array => lhs
! done here
end function add_to_array
end module any_stride
! a little tester program
!
! COMMENT and UNCOMMENT after this line
! when using the module in this file
!
program test
use any_stride
implicit none
! local vars
integer :: i
real*8, dimension(4), target :: data
real*8, dimension(2) :: rhs = (/ 3.0d0, 6.0d0 /)
integer, dimension(2), target :: indx
type(stride) :: astride
data = (/ 1.0d0, 2.0d0, 3.0d0, 4.0d0 /)
indx = (/ 1, 4 /)
astride = stride(data, indx)
astride = astride + rhs
print *, 'astride data= ', astride%data
! done here
end program test
the values 1 and 4 are only added to the first and the last entries of array data.
Upvotes: 0
Reputation: 21431
Your pointer assignments have what is known as bounds remapping. The "entire" pointer object is pointer assigned, but it is given the bounds indicated in the parentheses on the left hand side. This lets you have a lower bound other than one and have the rank of the pointer differ from the target. This feature is similar to the capabilities that you have with array argument association when the dummy array is explicit or assumed size.
You cannot have a pointer point to arbitrary elements of a target array - the elements in the target have to be "regular" - have a constant stride between them.
As an alternative to repacking you can obviously store those indices in an array and use the combination of the whole target and that array of indices to get the behaviour that you want, perhaps by wrapping a pointer to the whole target and the array of indices in a derived type with defined assignments.
Upvotes: 5