Prev: Integer return values in interfaces
Next: Cannot mex with gfortran on Snow Leopard (64bit)
From: James Van Buskirk on 16 Oct 2009 21:19
"Tobias Burnus" <burnus(a)net-b.de> wrote in message
news:4AD8DCF3.90905(a)net-b.de...

> On 16 October 2009 James Van Buskirk wrote:

>> > I read this as some sort of problem with the standard. At least it
>> > says that it lets you pass DCOS as an actual argument without having
>> > to declare the dummy argument as elemental, which would break codes
>> > dating back to at least f77.

> It is not completely clear to me why there is/would be a problem.

Thanks to Richard for decoding the above.

>> > It could be that whoever wrote the quoted paragraph assumed that for
>> > a dummy argument to be elemental was generally prohibited without
>> > actually checking, or that the prohibition existed in a previous
>> > draft but was dropped but not fixed in the quoted paragraph.

> Maybe - or I have simply missed to find the right spot.

>> > I can't see why the prohibition against passing a nonintrinsic
>> > elemental procedure as an actual argument exists anyhow.

> Me neither; I think the way elemental procedures are commonly
> implemented, it should just work - as ifort and NAG f95 show. (And
> gfortran would likely show if the error check were removed.) However,
> someone might have had a different implementation in mind
> (hypothetically or really existing) were allowing would have lead to
> problems.

>> > You know, there is a similar restriction regarding procedure
>> > pointers:

>> > abstract interface
>> > elemental function fun(x)
>> > [...]
>> > elemental function my_dcos(x)
>> > [...]
>> > procedure(fun), pointer :: f
>> > f => my_dcos
>> > write(*,*) f(x)

>> > Whoops, gfortran let this one slide!

> Actually, this program looks OK to me. Do you see anything in the
> standard which does not allow this one? (Note, the proc-pointer is not
> passed as actual argument but the procedure to which it is associated is
> just called.)

I was thinking in terms of N1601.pdf, C728:

"(R742) The proc-target shall not be a nonintrinsic elemental procedure."

Similar in wording to the actual argument case.

>> > OK, how about the other versions:

>> > abstract interface
>> > elemental function fun(x)
>> > [...]
>> > function my_dcos(x)
>> > [...]
>> > procedure(fun), pointer :: f
>> > f => my_dcos

>> > Now, we agree I think that gfortran should have caught that one as
>> > we shouldn't be able to change non-elemental procedures to
>> > syntactically elemental ones merely by pointing at them.

> I agree that it is invalid and that gfortran should have caught it. (As
> it is not a constraint, the compiler is not required to check this,
> however.)

You can always create situations where the compiler can't compare
interfaces, even if all procedures are module procedures, even in
the same module. That said, if the function invoked had one of
the things that don't make sense for an elemental function, the
results in all likelyhood wouldn't be pretty.

> From the F2003 standard (7.4.2.2 Procedure pointer assignment)
> "If proc-pointer-object has an explicit interface, its characteristics
> shall be the same as proc-target except that proc-target may be pure
> even if proc-pointer-object is not pure and proc-target may be an
> elemental intrinsic procedure even if proc-pointer-object is not
> elemental."

>> > The last version:

>> > abstract interface
>> > elemental function fun(x)

>> > procedure(fun), pointer :: f
>> > intrinsic dcos
>> > f => dcos

> Looks OK to me - assuming that "fun" has the same interface as "dcos"
> (which seems to be the case). Again, you do not use the proc-pointer as
> actual argument thus there should be no problem.

Yes, this I claim to be conforming. Note that in the passage from
7.4.2.2 that you quoted above, there is no question of any hint of
prohibiting this syntax. As a consequence we can construct an example
where a non-intrinsic elemental procedure can be passed to another
procedure. Sorry, no time for examples just now, so you will just
have to use your imagination, as also for the Cray pointer case.

--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end


From: Les Neilson on 19 Oct 2009 05:10

<analyst41(a)hotmail.com> wrote in message
news:4e146215-d595-42f7-af9e-faa2fd1fe4e5(a)a7g2000yqo.googlegroups.com...
> (1) Using Double Precision for latitude and longitude and using double
> precision for degrees to radians and radius of the earth (by the way,
> I have seen values that differ by about 4 miles)
>

Ignoring for the moment the ellipsoid problems.

Radius from center to sea level = R0
Radius from center to summit of Mt Everest = R0 + 5miles approximately

:-)
Les

From: jfh on 19 Oct 2009 17:29
On Oct 19, 10:10 pm, "Les Neilson" <l.neil...(a)nospam.co.uk> wrote:
> <analys...(a)hotmail.com> wrote in message
>
> news:4e146215-d595-42f7-af9e-faa2fd1fe4e5(a)a7g2000yqo.googlegroups.com...
>
> > (1) Using Double Precision for latitude and longitude and using double
> > precision for degrees to radians and radius of the earth (by the way,
> > I have seen values that differ by about 4 miles)
>
> Ignoring for the moment the ellipsoid problems.
>
> Radius from center to sea level = R0
> Radius from center to summit of Mt Everest = R0 + 5miles approximately
>
> :-)
> Les

But why ignore the ellipsoid problem? The OP seemed to want better
than 6-figure accuracy.
Assuming a spherical earth gives you 3-figure accuracy but not 4. (Try
calculating the
distance between 2 points on the equator that are 180 deg apart. The
shortest route goes
over a pole: the equatorial route is longer.)

The equatorial radius is 6378140 or 6378137 or 6378136 m according to
3 sources quoted in "Global Earth Physics: A Handbook of Physical
Constants" (T.J.Ahrens, ed., American Geophysical Union 1995), p.8,
and the flattening is 1/298.257 or 1/298.257222. A 4-metre
discrepancy in equatorial radius, not 4 miles. And anyway, shouldn't
high-accuracy geodesics be calculated on the geoid instead of any
ellipsoid?

BTW similar considerations apply to finding the direction from any
point P to Mecca. I do not know whether Muslims use the projection
onto the horizontal plane at P of the straight line through the Earth
to Mecca (ignoring general relativistic amendments to Euclidean
geometry?), or the direction in which the true geodesic on the surface
to Mecca leaves P, or whether they still assume that that geodesic is
the great circle calculated for a spherical earth. The differences are
large if P is close to the antipodes of Mecca, in French Polynesia.

-- John Harper
From: James Van Buskirk on 20 Oct 2009 01:36
"James Van Buskirk" <not_valid(a)comcast.net> wrote in message
news:hbb630$882$1(a)news.eternal-september.org...

> Yes, this I claim to be conforming. Note that in the passage from
> 7.4.2.2 that you quoted above, there is no question of any hint of
> prohibiting this syntax. As a consequence we can construct an example
> where a non-intrinsic elemental procedure can be passed to another
> procedure. Sorry, no time for examples just now, so you will just
> have to use your imagination, as also for the Cray pointer case.

Testing revealed that I was way overoptimistic about the possibility
of passing a non-intrinsic elemental procedure via a C_FUNPTR --
there are even more barriers here than in the procedure pointer case.

Well, that leaves the Cray pointer test:

C:\gfortran\clf\dcos>type cray.f90
module funcs
implicit none
integer, parameter :: dp = kind(1.0d0)
contains
! function my_dcos(x)
elemental function my_dcos(x)
integer, parameter :: dp = kind(1.0d0) !-- Or USE an appropriate module
real(dp), intent(in) :: x
real(dp) :: my_dcos
my_dcos = cos(x)
return
end function
end module funcs

program start
use funcs
implicit none
interface
elemental function fun(x)
import
implicit none
real(dp), intent(in) :: x
real(dp) fun
end function fun
end interface
pointer(f,fun)
real(dp) x(3)
intrinsic dcos
x = [1,2,3]
f = loc(my_dcos)
! f = loc(dcos)
write(*,*) fun(x)
end program start

C:\gfortran\clf\dcos>gfortran -fcray-pointer cray.f90 -ocray
cray.f90:30.11:

f = loc(my_dcos)
1
Error: ELEMENTAL non-INTRINSIC procedure 'my_dcos' is not allowed as an
actual a
rgument at (1)

I don't see why that should be a problem. After all, LOC is an
extension. Maybe that's just a blind spot I have in that I can see
the conflict if LOC was actually going to invoke MY_DCOS, but not if
it's just going to take its address. In that case LOC seems more like
an operator than a function, like operator(=>) (you can't overload
this one though, can you?) but gfortran is quite literal about LOC
being a function, not an operator... See what ifort does:

C:\gfortran\clf\dcos>ifort cray.f90
Intel(R) Fortran Compiler for Intel(R) EM64T-based applications, Version 9.1
Build 20061104
Copyright (C) 1985-2006 Intel Corporation. All rights reserved.

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

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

C:\gfortran\clf\dcos>cray
0.540302305868140 -0.416146836547142 -0.989992496600445

The 'B' version, with a non-elemental procedure:

C:\gfortran\clf\dcos>type cray.f90
module funcs
implicit none
integer, parameter :: dp = kind(1.0d0)
contains
function my_dcos(x)
! elemental function my_dcos(x)
integer, parameter :: dp = kind(1.0d0) !-- Or USE an appropriate module
real(dp), intent(in) :: x
real(dp) :: my_dcos
my_dcos = cos(x)
return
end function
end module funcs

program start
use funcs
implicit none
interface
elemental function fun(x)
import
implicit none
real(dp), intent(in) :: x
real(dp) fun
end function fun
end interface
pointer(f,fun)
real(dp) x(3)
intrinsic dcos
x = [1,2,3]
f = loc(my_dcos)
! f = loc(dcos)
write(*,*) fun(x)
end program start

C:\gfortran\clf\dcos>gfortran -fcray-pointer cray.f90 -ocray

C:\gfortran\clf\dcos>cray
0.54030230586813976 -0.41614683654714241 -0.98999249660044542

ifort likes this, too. In the general case I don't think the
compiler can tell that there is a mismatch, because we could have
computed the value of f arithmetically. In the intrinsic case:

C:\gfortran\clf\dcos>type cray.f90
module funcs
implicit none
integer, parameter :: dp = kind(1.0d0)
contains
function my_dcos(x)
! elemental function my_dcos(x)
integer, parameter :: dp = kind(1.0d0) !-- Or USE an appropriate module
real(dp), intent(in) :: x
real(dp) :: my_dcos
my_dcos = cos(x)
return
end function
end module funcs

program start
use funcs
implicit none
interface
elemental function fun(x)
import
implicit none
real(dp), intent(in) :: x
real(dp) fun
end function fun
end interface
pointer(f,fun)
real(dp) x(3)
intrinsic dcos
x = [1,2,3]
! f = loc(my_dcos)
f = loc(dcos)
write(*,*) fun(x)
end program start

C:\gfortran\clf\dcos>gfortran -fcray-pointer cray.f90 -ocray

C:\gfortran\clf\dcos>cray
0.54030230586813976 -0.41614683654714241 -0.98999249660044542

works as expected, as does ifort.

--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end


From: analyst41 on 20 Oct 2009 19:34
On Oct 19, 5:10 am, "Les Neilson" <l.neil...(a)nospam.co.uk> wrote:
> <analys...(a)hotmail.com> wrote in message
>
> news:4e146215-d595-42f7-af9e-faa2fd1fe4e5(a)a7g2000yqo.googlegroups.com...
>
> > (1) Using Double Precision for latitude and longitude and using double
> > precision for degrees to radians and radius of the earth (by the way,
> > I have seen values that differ by about 4 miles)
>
> Ignoring for the moment the ellipsoid problems.
>
> Radius from center to sea level = R0
> Radius from center to summit of Mt Everest = R0 + 5miles approximately
>
> :-)
> Les

Nothing so fancy. Every ZIP code has a latitude and longitude and the
idea is to calculate the distance from one to another.

I am using

FUNCTION ZIPDIST1(RLAT1,RLONG1,RLAT2,RLONG2)

DOUBLE PRECISION RLAT1,RLONG1,RLAT2,RLONG2,ACOSINPUT0,ACOSINPUT

ACOSINPUT0 = &
& SIN(RLAT1/57.2957795130823D0) * SIN(RLAT2/57.2957795130823D0) + COS
(RLAT1/57.2957795130823D0) &
& * COS(RLAT2/57.2957795130823D0) * COS(RLONG2/57.2957795130823D0 -
RLONG1/57.2957795130823D0)


ACOSINPUT = MIN(MAX( -1.0D0,ACOSINPUT0),1.0D0)

ZIPDIST1 = 3958.75586574D0 * acos(acosinput)

RETURN
END

All the Zipcodes will be in the US. I would say that one percent
accuracy would be more than sufficient - I just want to avoid any
known pitfalls with violent functions such as acos.
First  |  Prev  |  Next  |  Last
Pages: 1 2 3 4 5
Prev: Integer return values in interfaces
Next: Cannot mex with gfortran on Snow Leopard (64bit)