Re: Pointer-valued function to access inner components

On Feb 17, 5:56 am, Walter Spector <w6ws_xthiso...@xxxxxxxxxxxxx>
wrote:
Michael Metcalf wrote:

"Walter Spector" <w6ws_xthiso...@xxxxxxxxxxxxx> wrote in message
news:45D5C9C7.41820D64@xxxxxxxxxxxxxxxx
No, that points to a *copy* of the diagonal.

Agreed. Is this then what you meant (following Richard's hint)?

common/c/a...

Well, placing things in COMMON doesn't make them very dynamic.
But that is certainly one approach.

I was thinking more of the following:

program testdiag
implicit none

real, allocatable, target :: a(:,:)
real, pointer :: adiag(:)

integer :: i, j
integer :: n

write (*,'(a)', advance='no') 'size of matrix? '
read *, n

allocate (a(n,n))
a = 0

print *, 'A before:'
do, j=1, size (a,1)
write (*,'(100f5.2)') (a(j,i), i=1, size (a,2))
end do

adiag => getdiag (a, size (a,1))

adiag = 42.42

print *
print *, 'A after:'
do, j=1, size (a,1)
write (*,'(100f6.2)') (a(j,i), i=1, size (a,2))
end do

contains

function getdiag (t, tn)

! Return a F90 pointer which points to the diagonal of T.

integer, intent(in) :: tn
real, target :: t(tn+1,*)
real, pointer :: getdiag(:)

getdiag => t(1,:tn)

end function

end program

The above works fine, in practice, with a lot of different
compilers. I've used it with Intel, SGI IRIX, g95, gfortran,
and Salford compilers. However, it might be argued that it violates
the letter of the Standard... So one must be careful not to use it
in a situation where the compiler will potentially perform
a copyin/copyout on array A during the procedure call.

When copyin/copyout does get in the way, Cray pointers can come in
handy. Or, I presume, some of the new F2003 intrinsics. Time
will tell.

W.


using some more non-standard extensions and dirty tricks, one can even
implement a function
that can take the diagonal of an arbitrary array section (i.e., takes
assumed shape array as an argument). The following program works with
gfortran, g95 (with -Wno-globals), Intel and Pathscale:


module diagm
implicit none
contains
function diag(mat)
real,dimension(:,:),intent(in),target:: mat
real,dimension(:),pointer:: diag
integer:: n,inc
integer,parameter:: sizeof_real = 4
interface
subroutine dirty_pointer_1d(n,x,inc,ptr)
integer,intent(in):: n,inc
real:: x ! an intentional lie ...
real,dimension(:),pointer:: ptr
end subroutine
end interface
n = minval(shape(mat))
if (n == 0) then
nullify(diag)
else if (n == 1) then
diag => mat(1:1,1)
else
inc = (loc(mat(2,2)) - loc(mat(1,1))) / sizeof_real
call dirty_pointer_1d(n,mat(1,1),inc,diag)
end if
end function
end module
! this one has to be external
subroutine dirty_pointer_1d(n,x,inc,ptr)
integer,intent(in):: n,inc
real,target:: x(0:*)
real,dimension(:),pointer:: ptr
ptr => x(0:n*inc:inc)
end subroutine

program testp
use diagm
real,dimension(:,:),allocatable:: a
real,dimension(:),pointer:: dg
allocate(a(15,10))
a = 0
write(*,'(10F4.1)') transpose(a)
write(*,'(40("-"))')
read *
dg => diag(a)
dg = 1.1
write(*,'(10F4.1)') transpose(a)
write(*,'(40("-"))')
read *
dg => diag(a(4:12:2,1:6))
dg = 2.2
write(*,'(10F4.1)') transpose(a)
write(*,'(40("-"))')
read *
dg => diag(a(15:1:-1,:))
dg = 3.3
write(*,'(10F4.1)') transpose(a)
end program


.

Relevant Pages