Home > Archive > Fortran > January 2006 > elemental assignment extension
You are viewing an archived Text-only version of the thread.
To view this thread in it's original format and/or if you want to reply to
this thread please [click here]
| Author |
elemental assignment extension
|
|
| inidos 2006-01-25, 7:57 am |
| Hi,
I'd like to ask the language experts about self defined
assignements. Executing the following code:
----------------------------
program test
use mod1
type(mytype) :: a(5), b(5)
do ii =3D 1, size(a)
call mytype_init(a(ii), (/ 1, 2, 3 /))
end do
b(:) =3D a(:)
print *, "ASSOCIATED a, b: ", associated(b(1)%array, a(1)%array)
end program test
----------------------------
module mod1
implicit none
type mytype
integer, pointer :: array(:)
end type mytype
interface assignment(=3D)
module procedure mytype_assign
end interface
contains
subroutine mytype_init(self, array)
type(mytype), intent(out) :: self
integer, intent(in) :: array(:)
allocate(self%array(size(array)))
self%array(:) =3D array(:)
end subroutine mytype_init
subroutine mytype_assign(left, right)
type(mytype), intent(in) :: right
type(mytype), intent(out) :: left
allocate(left%array(size(right%array)))
left%array(:) =3D right%array(:)
end subroutine mytype_assign
end module mod1
-----------------------------------------
Results in:
ASSOCIATED a, b: T
Specifying mytype_assign to be "elemental" results in
ASSOCIATED a, b: F
This means, in a case of an array assignement, my assignment will be
only called if its specified to be "elemental", otherwise the compiler
only copies the derived types by copying their fields, but without
creating new arrays. In the latter case, the pointers show then to the
same memory space as for the original instances of the derived types.
(This is bad, if I try later to destroy the elements of both array a
and b.)
This makes sense, but I did not find the appropriate specification in
the standard. Is this behaviour specified by the (F95) standard? (The
compilers I tried behaved similarly, so probably it is...) And what
would happen, if I would compile this with a Fortran 90 compiler? (In
Fortran 90 there are no elemental functions.) Would it use my
assignment, or the internal one?
Cheers,
B=E1lint
| |
| Michael Metcalf 2006-01-25, 7:10 pm |
|
"inidos" <aradi@phys.upb.de> wrote in message
news:1138193999.771111.293600@g49g2000cwa.googlegroups.com...
This means, in a case of an array assignement, my assignment will be
only called if its specified to be "elemental", otherwise the compiler
only copies the derived types by copying their fields, but without
creating new arrays. In the latter case, the pointers show then to the
same memory space as for the original instances of the derived types.
(This is bad, if I try later to destroy the elements of both array a
and b.)
You've correctly identified what's happening. If the subroutine is
elemental, it is used to define the assignment b(:) = a(:) on an
element-by-element basis. If it is not elemental, intrinsic assignment is
used and so no second allocation takes place. However, you can achieve the
same effect by replacing
b(:) = a(:) outside the loop by b(ii) = a(ii) inside it. This would work too
with an f90 compiler. You can't define an overload for an array assignment
in this case because it would imply a definition of an array of arrays,
which is not allowed in Fortran.
My question to you is why you are using pointers rather than allocatable
arrays? Many compilers support allocatable arrays as derived-type components
and that would seem to be what you're trying to do.
Regards,
Mike Metcalf
| |
| James Van Buskirk 2006-01-25, 7:10 pm |
| "Michael Metcalf" <michaelmetcalf@compuserve.com> wrote in message
news:6KLBf.4803$yE4.1137@news-wrt-01.rdc-nyc.rr.com...
> "inidos" <aradi@phys.upb.de> wrote in message
> news:1138193999.771111.293600@g49g2000cwa.googlegroups.com...
> This means, in a case of an array assignement, my assignment will be
> only called if its specified to be "elemental", otherwise the compiler
> only copies the derived types by copying their fields, but without
> creating new arrays. In the latter case, the pointers show then to the
> same memory space as for the original instances of the derived types.
> (This is bad, if I try later to destroy the elements of both array a
> and b.)
> You've correctly identified what's happening. If the subroutine is
> elemental, it is used to define the assignment b(:) = a(:) on an
> element-by-element basis.
I feel that this statement is a little imprecise. If it were taken
literally, the behavior of lf95, cvf, and ifort on
http://home.comcast.net/~kmbtib/For...elem_assign.f90
would be correct and that of g95 would be incorrect, whereas the
reverse is actually true. Assignment in Fortran of course happens
'as if' the right hand side were first copied to a temporary, then
then temporary is assigned, element by element if the assignment
is elemental, to the left hand side.
> If it is not elemental, intrinsic assignment is
> used and so no second allocation takes place. However, you can achieve the
> same effect by replacing
> b(:) = a(:) outside the loop by b(ii) = a(ii) inside it. This would work
too
> with an f90 compiler. You can't define an overload for an array assignment
> in this case because it would imply a definition of an array of arrays,
> which is not allowed in Fortran.
I don't see why you can't overload an array assignment. Isn't
the following code OK f90?
module mod1
implicit none
type mytype
integer, pointer :: array(:)
end type mytype
interface assignment(=)
module procedure mytype_assign_scalar, mytype_assign_vector
end interface
contains
subroutine mytype_init(self, array)
type(mytype), intent(out) :: self
integer, intent(in), optional :: array(:)
if(present(array)) then
allocate(self%array(size(array)))
self%array(:) = array(:)
else
nullify(self%array)
end if
end subroutine mytype_init
subroutine mytype_assign_scalar(left, right)
type(mytype), intent(in) :: right
type(mytype), intent(inout) :: left
if(associated(left%array)) deallocate(left%array)
allocate(left%array(size(right%array)))
left%array(:) = right%array(:)
end subroutine mytype_assign_scalar
subroutine mytype_assign_vector(left, right)
type(mytype), intent(in) :: right(:)
type(mytype), intent(inout) :: left(size(right))
integer ii
do ii = 1, size(right)
if(associated(left(ii)%array)) deallocate(left(ii)%array)
allocate(left(ii)%array(size(right(ii)%a
rray)))
left(ii)%array(:) = right(ii)%array(:)
end do
end subroutine mytype_assign_vector
end module mod1
!----------------------------
program test
use mod1
type(mytype) :: a(5), b(5)
do ii = 1, size(a)
call mytype_init(a(ii), (/ 1, 2, 3 /))
end do
do ii = 1, size(b)
call mytype_init(b(ii))
end do
b(:) = a(:)
print *, "ASSOCIATED a, b: ", associated(b(1)%array, a(1)%array)
end program test
--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end
| |
| Michael Metcalf 2006-01-25, 7:10 pm |
|
"James Van Buskirk" <not_valid@comcast.net> wrote in message
news:MI6dnesZO9OdBErenZ2dnUVZ_sydnZ2d@co
mcast.com...
> "Michael Metcalf" <michaelmetcalf@compuserve.com> wrote in message
> news:6KLBf.4803$yE4.1137@news-wrt-01.rdc-nyc.rr.com...
>
>
>
> I feel that this statement is a little imprecise. If it were taken
> literally, the behavior of lf95, cvf, and ifort on
>
> http://home.comcast.net/~kmbtib/For...elem_assign.f90
>
> would be correct and that of g95 would be incorrect, whereas the
> reverse is actually true. Assignment in Fortran of course happens
> 'as if' the right hand side were first copied to a temporary, then
> then temporary is assigned, element by element if the assignment
> is elemental, to the left hand side.
>
> too
>
> I don't see why you can't overload an array assignment. Isn't
> the following code OK f90?
>
> left(ii)%array(:) = right(ii)%array(:)
Yes, I agree. You avoid the array of array with this statement. I did it
otherwise.
Regards,
Mike Metcalf
| |
| inidos 2006-01-26, 7:58 am |
| > interface assignment(=3D)
> module procedure mytype_assign_scalar, mytype_assign_vector
> end interface
That's a nice solution, to provide an assignement extension for scalar
and
array, thank you. But, what happens, if you specify your scalar
assignement (mytype_assign_scalar) to be elemental while keeping also
the mytype_assign_vector routine? Which assignement would the compiler
(theoretically) have to take for the b(:) =3D a(:) assignement? The
elemental scalar or the non-elemantal vectorial one? Is there a
reference, where this "decision procedure" is properly explained?
Cheers,
B=E1lint
| |
| inidos 2006-01-26, 7:02 pm |
| Hi,
OK, the looped version should work also on a Fortan 90 compiler. But
what would an F90 compiler theoretically have to do with the original
example? Take the user defined assignment or take the intrinsic one?
The usage of allocatable arrays instead of pointers, as you suggested,
would be of course
very helpfull, since among others it would help me to avoid unintented
memory leaks (e.g. due to temporary copies, which can not be
explicitely destroyed). I wasn't aware of many memory leaks, until I
did not compile our code with g95, which let's you brutally face the
truth ;-).
However, I'm co-developing a simulation tool for atomistic modelling
and it's quite important to have a more or less portable source code.
Therefore, we restricted ourself to Fortran 95 constructs, without any
non-standardised extensions. Even so, we faced many compiler problems
(often compiler crashes!), when using sophisticated
(but standardised) Fortran 95 constructs, so I'm quite sceptical about
non-standard extensions (and about the future implementation of newer
Fortran standards by compiler vendors in general).
Cheers,
B=E1lint
| |
| Michael Metcalf 2006-01-26, 7:02 pm |
|
"inidos" <aradi@phys.upb.de> wrote in message
news:1138283981.935493.140740@g49g2000cwa.googlegroups.com...
Hi,
OK, the looped version should work also on a Fortan 90 compiler. But
what would an F90 compiler theoretically have to do with the original
example? Take the user defined assignment or take the intrinsic one?
The intrinsic one.
[snip] , so I'm quite sceptical about
non-standard extensions (and about the future implementation of newer
Fortran standards by compiler vendors in general).
I understand that.
Regards to Paderborn,
Mike Metcalf
| |
| Jan Vorbrüggen 2006-01-26, 7:02 pm |
| > However, I'm co-developing a simulation tool for atomistic modelling
> and it's quite important to have a more or less portable source code.
> Therefore, we restricted ourself to Fortran 95 constructs, without any
> non-standardised extensions.
The ALLOCATABLE TR, which Jugoslav mentioned, _is_ standardized, in fact
twice: once as a TR to F95 and now as a part of F03.
In all fairness, I have the impression that even relative newcomers such
as g95 now have a very good handle on the complexities of that TR, so you
should be as safe using its features as using other sophisticated, as you
say, F95 features. And the payoff of using it surely are substantial.
Jan
| |
|
| Well, the issue is that some compilers don't implement the TR at all.
While vendors of an f95 compiler might feel some obligation to fix bugs
wrt to the F95 standard, they might (for some time) happily ignore the
TR. On some machines you have no choice but work with the system
compiler. Therefore, at least in our project, we're also stuck using a
lot of pointers in our derived types.
Joost
| |
| Richard E Maine 2006-01-26, 7:03 pm |
| inidos <aradi@phys.upb.de> wrote:
> array, thank you. But, what happens, if you specify your scalar
> assignement (mytype_assign_scalar) to be elemental while keeping also
> the mytype_assign_vector routine?
You can't. That would be an illegal generic.
> Which assignement would the compiler
> (theoretically) have to take for the b(:) = a(:) assignement? The
> elemental scalar or the non-elemantal vectorial one? Is there a
> reference, where this "decision procedure" is properly explained?
For the most part, there is no decision procedure. The rules prohibit
you from making a generic where more than one specific could possibly
fit. However, we have a 2-level exception here.
Level 1. You can have a situation where both an elemental and a
nonelemental could fit. In that case, the non-elemental always "wins".
The idea is that the nonelemental version is more specific. Perhaps you
wanted to special-case one rank, but let the lemental handle other
ranks. There is an example in the standard on exactly that.
Level 2. I just noticed this one a few days ago and don't yet understand
the reason. I'm not sure it was intentional. I think it might be an
accidental holdover from before when users could write elemental
procedures, but... The level 1 exception doesn't "work" for the case of
having scalars be the special case. So you can't have an elemental and a
scalar version. You can have an elemental and any other rank, but not an
elemental and a scalar.
--
Richard Maine | Good judgment comes from experience;
email: my first.last at org.domain| experience comes from bad judgment.
org: nasa, domain: gov | -- Mark Twain
| |
| Keith Refson 2006-01-26, 7:03 pm |
| Joost wrote:
> Well, the issue is that some compilers don't implement the TR at all.
> While vendors of an f95 compiler might feel some obligation to fix bugs
> wrt to the F95 standard, they might (for some time) happily ignore the
> TR. On some machines you have no choice but work with the system
> compiler. Therefore, at least in our project, we're also stuck using a
> lot of pointers in our derived types.
We have been examining the same issue, for another project Joost may be familiar
with. In fact it turns out that the only compilers we found which don't implement
TR15581 at all are the SGI MIPSPro compilers, and Pathscale. (Portland seemed to
have a partial implementation which didn't include allocatable dummy arguments.)
This is a sufficiently small list, consisting of one legacy platform and another for
which a plethora of good alternatives exist that it's probably worth doing anyway.
Are there any others you know of?
Keith Refson
--
Dr Keith Refson,
Building R3
Rutherford Appleton Laboratory
Chilton
Didcot kr AT
Oxfordshire OX11 0QX isise D@T rl D.T ac D?T uk
| |
|
| > We have been examining the same issue, for another project Joost may be familiar
> with. In fact it turns out that the only compilers we found which don't implement
> TR15581 at all are the SGI MIPSPro compilers, and Pathscale. (Portland seemed to
> have a partial implementation which didn't include allocatable dummy arguments.)
>
> This is a sufficiently small list, consisting of one legacy platform and another for
> which a plethora of good alternatives exist that it's probably worth doing anyway.
These were the two I recalled as well, but e.g. did you test NEC SX5 or
cray XT3. In particular on the latter, pathscale might become available
and be a real alternative to the 'system compiler' (which ICEs out on
our current code after taking 1.5h to compile, grrr.). Provided one
avoids the TR.
Joost
| |
| James Van Buskirk 2006-01-26, 7:03 pm |
| "Richard E Maine" <nospam@see.signature> wrote in message
news:1h9rjup.4crym91nlmrtnN%nospam@see.signature...
> inidos <aradi@phys.upb.de> wrote:
[color=darkred]
> You can't. That would be an illegal generic.
> Level 1. You can have a situation where both an elemental and a
> nonelemental could fit. In that case, the non-elemental always "wins".
> The idea is that the nonelemental version is more specific. Perhaps you
> wanted to special-case one rank, but let the lemental handle other
> ranks. There is an example in the standard on exactly that.
I am surprised that you have the wrong answer followed by the
right answer here. The second paragraph directly contradicts
the first. Rare that you don't fix these things in your
ordinarily carefully worded and reviewed posts.
--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end
| |
| Richard E Maine 2006-01-26, 7:03 pm |
| James Van Buskirk <not_valid@comcast.net> wrote:
> "Richard E Maine" <nospam@see.signature> wrote in message
> news:1h9rjup.4crym91nlmrtnN%nospam@see.signature...
>
>
>
>
>
> I am surprised that you have the wrong answer followed by the
> right answer here. The second paragraph directly contradicts
> the first. Rare that you don't fix these things in your
> ordinarily carefully worded and reviewed posts.
Ah. I got the standard bits right (I think), but I too hastily read the
question. I was thinking that the question was about what I called the
level 2 exception, where there are specifics for elemental and scalar. I
now see that you are right - that he is talking about specifics for
elemental and rank 1.
--
Richard Maine | Good judgment comes from experience;
email: my first.last at org.domain| experience comes from bad judgment.
org: nasa, domain: gov | -- Mark Twain
|
|
|
|
|