SIZE with 64bit gFortran

[campbel...@gmail.com]

I thought that SIZE with 64 bit gFortran returned integer*8, but the following does not work when I tested. I am using gcc version 4.9.2 (x86_64-posix-seh-rev2, Built by MinGW-W64 project), which I thought would work.
Any ideas as to what is wrong ?
Makes me worry that gFortran does not fully support arrays with more than 2^31 elements.

module important_variables
integer*8 ipos
real*8, allocatable, dimension(:) :: Profile
end module important_variables

use important_variables

integer*4 status
integer*8 n1, n2
!
n1 = 2905
n2 = 1000000
ipos = n1*n2
write (*,*) 'ipos =',ipos
!
allocate ( profile(ipos), stat=status )
write (*,*) '[Profile] allocated; status =',status,': size =',size(profile)
!
n2 = size(profile)
write (*,*) n2
end

response is:
ipos = 2905000000
[Profile] allocated; status = 0 : size = -1389967296
-1389967296


(if your paging file does not support 22gb, try real*4 Profile)

[Arjen Markus]

Op maandag 28 maart 2016 07:30:53 UTC+2 schreef campbel...@gmail.com:

Try this instead:

write (*,*) '[Profile] allocated; status =',status,': size =',size(profile,kind=kind(n1))

Regards,

Arjen

[campbel...@gmail.com]

Arjen,

Thanks very much for the advice. It fixed the problem.
Other compilers I use default to integer(8) for 64-bit, but I suspect gFortran may be standard conforming (?)

John

[Clive Page]

You say that as though it is a disappointment to you - I hope not.

I think gfortran does conform to the Standard, but all Fortran Standards
have been very carefully drafted to not specify any particular number of
bits for any data type. In fact your code, using things like INTEGER*8
and REAL*8,is not standard-conforming, though it may work with many
compilers. It is a bit of a nuisance but there are really only two
ways of getting the behaviour you want in a completely portable and
standard-conforming way:

(a) invoking the intrinsic functions like SELECTED_REAL_KIND and
SELECTED_INT_KIND, or

(b) using the intrinsic module ISO_FORTRAN_ENV and then using e.g.
INTEGER(INT64) or REAL(REAL64).


--
Clive Page

[campbel...@gmail.com]

Clive,

I don't want to start a discussion on this, but for the mix of compilers I use, INTEGER*8 is the most portable and concise way I have to document the code I write. ( for portability Fortran 90 should have given guidance on the values for kind, rather than the needless verbosity that was preached. We should have seen that coming! )
Back in 80's and 70's especially I transferred many codes between hardware, where the precision of variables was vague. Having some verbose definition in an include file or module (that was often not included with the dump) was an annoyance. "Was the code developed on a Vax, Cyber or some English variation" was always the unknown, which meant the byte precision was uncertain. REAL*8 gave a lot of certainty, without requiring a search of files, hoping to find the answer. It has certainly become the go-to kind. Only once in 40 years have I written code with flexible precision, which did not give a useful solution. Typically the algorithm requires a certain precision, being 4, 8, 10 or 16 bytes and changing rarely achieves anything meaningful. They are typically the only alternatives we have.

Anyway, I presume SIZE should default to default integer kind, while with 64_bit, Intel and Silverfrost maintain integer*4 default for integers but provide integer*8 for size and memory address intrinsic, as that is more useful.

Portability is a tricky beast.

[Tim Prince]

I don't find documented any deviation from standard for SIZE in Intel
Fortran. In fact, my ifort has the same problem with your posted
example, regardless of whether I change your *8 notation to standard (8).
c_loc has a value return of type c_ptr, so it should always make the
adjustment you appear to want, under USE iso_c_binding. Intel Fortran
doesn't check standard usage of c_loc even under the option -stand, so
it's possible to write code which "works" in ifort but not gfortran.

[robin....@gmail.com]

On Monday, March 28, 2016 at 10:32:10 PM UTC+11, campbel...@gmail.com wrote:
> Clive,
>
> I don't want to start a discussion on this, but for the mix of compilers
> I use, INTEGER*8 is the most portable and concise way I have to document
> the code I write.

It's not standard, and not portable.

> ( for portability Fortran 90 should have given guidance on the values
> for kind, rather than the needless verbosity that was preached.

What's difficult about
dp = kind (1.0d0)
integer (kind = dp) :: x
?

[Paul van Delst]

On 03/28/16 07:32, campbel...@gmail.com wrote:
> Clive,
>

> I don't want to start a discussion on this,...

Nah, mate, you can't preface an only slightly camouflaged "the new
method is bad and the old non-standard method is better" and then state
you don't want to start (or presumably participate in) a discussion!

Hoo boy. :o)

Your comments immediately put you at a disadvantage because you are
advocating non-standard syntax before stating "Portability is a tricky
beast".

Darn tootin' - it's even trickier when you use code that is, by
definition, non-portable.

The fact that the *8/*4/*whatever syntax works for you and what you need
to do (or needed to do way back when) is great, but your statements
approach the pontification for which you exhibited disdain in your post:
"the needless verbosity that was preached". Really? (I tried to verbose
my sentence up, too, but I ended up deleting the pleonastic tautological
circularities. ha ha).

O.k., I definitely need more coffee. Or maybe less?

cheers,

paulv



p.s. For the record, I would *much* rather have a module like so:

module type_kinds
implicit none
private
public :: byte, short, long
public :: single, double

integer, parameter :: byte = selected_int_kind(1)
integer, parameter :: short = selected_int_kind(4)
integer, parameter :: long = selected_int_kind(8)

integer, parameter :: single = selected_real_kind(6)
integer, parameter :: double = selected_real_kind(15)
end module type_kinds

that I can then use anywhere in any of my code so I can do stuff like

use type_kinds
integer(byte) :: b
real(single) :: x
real(double) :: y

etc..

Is it more work up front? Yes.

Is it portable? Yes.

(Bonus question) Is it self documenting? Yes.

[FortranFan]

On Monday, March 28, 2016 at 8:37:55 AM UTC-4, Tim Prince wrote:

> ..

> >
> I don't find documented any deviation from standard for SIZE in Intel
> Fortran. In fact, my ifort has the same problem with your posted
> example, regardless of whether I change your *8 notation to standard (8).

Standard clearly states in Section 13.7.156 (document WD 1539-1, J3/10-007r1, 24th November 2010)

-- begin snip from standard --
..
18 13.7.156 SIZE (ARRAY [, DIM, KIND])
19 1 Description. Size of an array or one extent.
20 2 Class. Inquiry function.
21 3 Arguments.
22 ARRAY shall be an array of any type. It shall not be an unallocated
allocatable variable or a pointer that
23 is not associated. If ARRAY is an assumed-size array, DIM shall be
present with a value less than
24 the rank of ARRAY.
25 DIM (optional) shall be an integer scalar with a value in the range
1 <= DIM <= n, where n is the rank of ARRAY.
26 KIND (optional) shall be a scalar integer constant expression.
27 4 Result Characteristics. Integer scalar. If KIND is present, the kind
type parameter is that specified by the
28 value of KIND; otherwise the kind type parameter is that of default
integer type.
..

-- end snip --

From what I notice too with current Intel Fortran version (16, update 2) with standard-semantics compiler option, it conforms to the standard.

As to the following,
On Monday, March 28, 2016 at 6:01:39 AM UTC-4, campbel...@gmail.com wrote:
> .. Other compilers I use default to integer(8) for 64-bit ..

and

On Monday, March 28, 2016 at 7:32:10 AM UTC-4, campbel...@gmail.com wrote:
> .. while with 64_bit, Intel and Silverfrost maintain integer*4 default for integers but provide integer*8 for size and memory address intrinsic, as that is more useful. ..


Not exactly sure what OP is trying to convey here; the statements appear misleading and incorrect, at least with Intel Fortran.

[Richard Maine]

<campbel...@gmail.com> wrote:

> Anyway, I presume SIZE should default to default integer kind

Yes, the standard is quite explicit and clear on that.

It is a little unfortunate that 64-bit compilers often don't make
default integer (and real) 64-bits. Some compilers do; others don't. The
standard allows either choice.

Probably the main reason that man compilers keep default real and
integer at 32 bits even for 64-bit address spaces is that there is an
awful lot of code written with the non-portable assumption that the
default sizes are 32 bits. That assumption is not justified by the
standard, but lots of code makes it anyway. Many vendors have decided
that it is more important to keep that old code working than to
facilitate new code that uses large arrays. You can, of course, still
use large arrays, but it is a little more awkward because you have to
explicitly specify the kind in several places, as in how Arjen pointed
out with the SIZE intrinsic.

I find it minorly ironic (but far from unprecedented) to read about an
awkwardness that essentially arises from vendors prioritizing support of
code that uses nonstandard assumptions, while I also read the same
author advocating for a nonstandard coding syntax. But no, I don't
really want to discuss that again either; this one paragraph is (more
than) enough on that subject from me in this thread.

--
Richard Maine
email: last name at domain . net
dimnain: summer-triangle

[Ron Shepard]

On 3/28/16 12:30 AM, campbel...@gmail.com wrote:
> I thought that SIZE with 64 bit gFortran returned integer*8, but the following does not work when I tested.

I'm pretty sure that the SIZE() intrinsic must return a default integer
kind by default. SIZE() has the option to change the kind of its
returned value if desired, as do many of the intrinsics in fortran.

It would be possible for the compiler to treat 8-byte integers as the
default, but that would also require treating 8-byte logical and reals
as the default and 16-byte reals as DOUBLE PRECISION in order to be
standard conforming. But most byte-addressable machines don't do that,
except possibly with compiler options, choosing instead to keep 4-byte
integers as the default with either 32-bit or 64-bit addressing. This is
also consistent with the default behavior with other languages, such as
C, when switching from 32-bit to 64-bit addressing.

This switchover from 32-bit to 64-bit addressing occurred over a decade
ago for most of us, and a lot of fortran code has been written in the
meantime, so I doubt that you are going to get any of the major
compilers to switch their default behavior now.

Instead, I recommend that you place all of these KIND definitions in a
single module, and USE that module everywhere in your code. This
localizes the machine dependency to a single place in your code in case
you ever need to change things. If you use the SELECTED_REAL_KIND() (and
so on) intrinsics appropriately, then you might not ever even need to
make any changes to this module as you move to new hardware and compiler
combinations. Like you, I also wrote code in a wide variety of computers
in the 70s, 80s, and 90s (with default integers with 16, 24, 32, 36, and
64 bits, and with corresponding REAL and DOUBLE PRECSIION data sizes),
so this is an important feature of modern fortran to me.

$.02 -Ron Shepard

[Ron Shepard]

On 3/28/16 10:06 AM, Ron Shepard wrote:
> Like you, I also wrote code in a wide variety of computers in the 70s,
> 80s, and 90s (with default integers with 16, 24, 32, 36, and 64 bits,
> and with corresponding REAL and DOUBLE PRECSIION data sizes), so this is
> an important feature of modern fortran to me.

Oops, I should mention 60-bit (CDC) too.

$.02 -Ron Shepard

[Tim Prince]

You could include an integer type based on c_ptr if, like John, you want
one which switches automatically between 32- and 64-bit target modes,
for those compilers where that is relevant.
Up to now, there have been platforms which run a 64-bit OS but support
much more efficient 32-bit array indexing. I don't know whether next
year's Intel CPUs will correct that hardware level deficiency in 64-bit
support. Even if the hardware changes so that a change in default would
be satisfactory, inertia will prevent changes in compiler defaults.
For example, gfortran i386 32-bit still defaults to generation of code
for CPUs which went out of production 15 years ago. People still get
confused over the various levels of double precision evaluation of
single precision expressions.
5 years ago, Intel hardware fixed the deficiency in support of IEEE754
gradual underflow mode, but Intel compilers still set abrupt underflow
as a default, even when specifically targeted to those newer CPUs. That
hardware fix seemed motivated by the requirements of gcc and gfortran.

[Tim Prince]

Reading further about AVX512, there is support for 52-bit integer fma,
so there is a case where 64-bit integers will remain less efficient than
smaller ones (aside from inevitable memory bandwidth questions).

[Clive Page]

On 28/03/2016 12:32, campbel...@gmail.com wrote:
> Portability is a tricky beast.

I guess we can all agree on that, if on nothing else.

The discussion so far hasn't brought out that there are two different
circumstances in which you want to specify precision:

(a) To ensure that your calculations have enough bits that no
significant information is lost - for this the use of
SELECTED_REAL_KIND(N) is fine, but it only guarantees you that the data
type selected will have enough bits to store at least N decimal digits,
it doesn't, in theory at least, specify a particular number of bits. you
might get more than you need, but for computations that doesn't matter.

(b) To read or write a binary (unformatted) file where you know that the
values have a particular number of bits and you want variables with
exactly that number of bits, no more and no fewer. Here the
introduction in Fortran 2008 of the ISO_FORTRAN_ENV module allowing
INTEGER(int64) etc. is valuable and I've used it quite a bit.

In my experience, both of these are portable, but I guess there are
compilers around which so far support only (a) and not (b).

--
Clive Page

[kargl]

campbel...@gmail.com wrote:

> I thought that SIZE with 64 bit gFortran returned integer*8, but the following does not work
> when I tested. I am using gcc version 4.9.2 (x86_64-posix-seh-rev2,
Built by MinGW-W64 project), which I thought would work.
> Any ideas as to what is wrong ?

You did not read the documentation that came with your compiler.
You did not read the Fortran standard.

--
steve

[steve kargl]

Richard Maine wrote:

> <campbel...@gmail.com> wrote:
>
>> Anyway, I presume SIZE should default to default integer kind
>
> Yes, the standard is quite explicit and clear on that.
>
> It is a little unfortunate that 64-bit compilers often don't make
> default integer (and real) 64-bits. Some compilers do; others don't. The
> standard allows either choice.
>
> Probably the main reason that man compilers keep default real and
> integer at 32 bits even for 64-bit address spaces is that there is an
> awful lot of code written with the non-portable assumption that the
> default sizes are 32 bits. That assumption is not justified by the
> standard, but lots of code makes it anyway. Many vendors have decided
> that it is more important to keep that old code working than to
> facilitate new code that uses large arrays. You can, of course, still
> use large arrays, but it is a little more awkward because you have to
> explicitly specify the kind in several places, as in how Arjen pointed
> out with the SIZE intrinsic.
>

I suspect the main reason is the consequence of storage association.
A 64-bit REAL requires a 128-bit DOUBLE PRECSION. Common hardware
does not support a 128-bit hardware floating point; so, tall the DOUBLE
PRECISION with use a software implementation (which is slowwwwwww).
I suppose one could enforce the 128-bit storage association, and then
do computations in whatever extended precision is supported by the
hardware (e.g., Intel 80-bit extended double precision).

--
steve

[Richard Maine]

steve kargl <s...@REMOVEtroutmask.apl.washington.edu> wrote:

> Richard Maine wrote:

> > Probably the main reason that man compilers keep default real and
> > integer at 32 bits even for 64-bit address spaces is that there is an
> > awful lot of code written with the non-portable assumption that the

> > default sizes are 32 bits....

> I suspect the main reason is the consequence of storage association.
> A 64-bit REAL requires a 128-bit DOUBLE PRECSION. Common hardware
> does not support a 128-bit hardware floating point; so, tall the DOUBLE
> PRECISION with use a software implementation (which is slowwwwwww).
> I suppose one could enforce the 128-bit storage association, and then
> do computations in whatever extended precision is supported by the
> hardware (e.g., Intel 80-bit extended double precision).

Yes, it would be perfectly standard conforming for double to use 80-bit
reals stored in 128-bits for storage association purposes. I've used
compilers that did exactly that. I don't see that as probably having
been a big driver against making default real 64 bits. Granted I could
be wrong, but that would not be my guess as to a significant reason.

When the kind argument got added to all the intrinsics like SIZE, I
recall thinking that it was not worth the bother because surely
compilers that wanted to support large address spaces would just make 64
bits the default. But I was obviously wrong on that. :-)

[Ian Harvey]

On 2016-03-29 4:00 AM, Clive Page wrote:
...

> (b) To read or write a binary (unformatted) file where you know that the
> values have a particular number of bits and you want variables with
> exactly that number of bits, no more and no fewer. Here the
> introduction in Fortran 2008 of the ISO_FORTRAN_ENV module allowing
> INTEGER(int64) etc. is valuable and I've used it quite a bit.

A subtle point, which may not be relevant on commonly available Fortran
processors, is that INT64, REAL64, etc, relate to the storage size of an
array element in memory, which is not the same thing as the number of
bits required to store the object in an unformatted file (and it says
nothing about how those bits are to be interpreted to determine the
value of the object), and it is not the same thing as the "value" of the
thing having a certain number of bits.

[FX]

> (b) using the intrinsic module ISO_FORTRAN_ENV and then using e.g.
> INTEGER(INT64) or REAL(REAL64).

If the intent is to use a integer kind that directly maps the processor's
memory access (32-bit or 64-bit on respective hardware), then
ISO_C_BINDING's C_SIZE_T or C_INTPTR_T is the best choice.

$ cat b.f90
use, intrinsic :: iso_c_binding
print *, c_size_t
end
$ gfortran -m32 b.f90 && ./a.out
4
$ gfortran -m64 b.f90 && ./a.out
8

--
FX

[campbel...@gmail.com]

How did it become acceptable for gFortran to store real*10 in 16 bytes ?

[campbel...@gmail.com]

On Tuesday, March 29, 2016 at 12:06:43 AM UTC+11, robin....@gmail.com wrote:
>
> It's not standard, and not portable.
>
> > ( for portability Fortran 90 should have given guidance on the values
> > for kind, rather than the needless verbosity that was preached.
>
> What's difficult about
> dp = kind (1.0d0)
> integer (kind = dp) :: x
> ?
>

Robin,

We may have had similar experiences with different Fortran compilers in the past, but if real*8 is not portable, then that is a compiler which I would not use. How many compilers are in use in this millennium and fail this portability test ?

As for what's difficult about "integer (kind = dp) :: x" lots !!

A concise way of using 64 bit integer constants also has problems; dp = kind(??)

[herrman...@gmail.com]

On Monday, March 28, 2016 at 3:25:40 PM UTC-7, campbel...@gmail.com wrote:

(snip)

> How did it become acceptable for gFortran to store real*10 in 16 bytes ?

Alignment is always complicated, but I believe that on many processors fetch/store of temporary real is faster if appropriately aligned.

The 8087, with a 16 bit data bus, runs faster with data on even addresses.
The 80486, with a 32 bit data bus, prefers four byte alignment.

Much software written in the 80486 era only supplied four byte alignment, such as for the return value of C's malloc().

Starting with the Pentium, there is a significant speed penalty for data not eight byte aligned.
It took some time for all the 80486 era software to disappear.

Since one often wants speed with floating point, I suspect that compilers try to align REAL*10 for appropriate speed.

[herrman...@gmail.com]

On Monday, March 28, 2016 at 3:50:21 PM UTC-7, campbel...@gmail.com wrote:

(snip)

> We may have had similar experiences with different Fortran compilers in the past,
> but if real*8 is not portable, then that is a compiler which I would not use.
> How many compilers are in use in this millennium and fail this portability test ?

> As for what's difficult about "integer (kind = dp) :: x" lots !!

> A concise way of using 64 bit integer constants also has problems; dp = kind(??)

The REAL*8 form is not in any ANSI or ISO standard, but came through from IBM
and Fortran IV.

In its manuals, IBM carefully documents its extensions to Fortran 66, but to
stay competitive, others followed IBM on this one.

VAX/VMS even includes REAL*16 and COMPLEX*32, as IBM used in
the H extended compiler.

So, REAL*8 is not standard, but reasonably portable.

The other option, REAL(8) is standard, but not portable. The standard
gives no meaning to specific numerical values for KIND.

KIND(1.D0) is often 8, but even on byte addressable machines
with a 64 bit double precision data type, it isn't required to be 8.

I often use the REAL*8 form in quick programs meant to test out
features, and might even post them to the group.

INTEGER*8 is a little different. I don't know that IBM ever supported
a 64 bit integer type for its Fortran IV compilers. I don't remember
it being in VS Fortran. DEC might have it in VAX/VMS, though.

[campbel...@gmail.com]

This original post may have gone in another direction, but for the original example:

n2 = size(profile)

To get a value of -1389967296 is not acceptable performance, what ever the lame excuse.

Fortran must be useable and not just standard conforming. I hope there are some who understand this.

Ever used y = sin (x) and got an unacceptable answer

[steve kargl]

gfortran is apart of GCC. See the gcc.info documentation for the
-m96bit-long-double and -m128bit-long-double options.

--
steve

[steve kargl]

campbel...@gmail.com wrote:

> This original post may have gone in another direction, but for the original example:
>
> n2 = size(profile)
>
> To get a value of -1389967296 is not acceptable performance, what ever the lame excuse.
>

Your original program is not standard conforming (and I'm not talking about the
de facto standard REAL* form). The standard is very clear on this point. From
F2003, 13.7:

A program is prohibited from invoking an intrinsic procedure
under circumstances where a value to be returned in a subroutine
argument or function result is outside the range of values
representable by objects of the specified type and type parameters,
unless the intrinsic module IEEE ARITHMETIC (section 14) is
accessible and there is support for an infinite or a NaN result,
as appropriate.

The size of profile exceeds HUGE(1_4), so the processor is free to do
whatever it wants. In this case, you appear to be getting the well
known wrap around semantics of 2-complement integers. If you have
an array whose size exceeds HUGE(1_4), then use the optional
argument SIZE(profile, 8)

--
steve

[Richard Maine]

<campbel...@gmail.com> wrote:

> How did it become acceptable for gFortran to store real*10 in 16 bytes ?

To the extent that your "acceptable" means "acceptable to the standard",
i.e. standard conforming, the standard has allowed things like that from
day one. Storage alignment requirements are distinct from precision
requirements. Many examples of mismatches between the two have existed.

As noted elsewhere, I have used multiple compilers that stored 80-bit
real values with 128-bit alignment. No, gFortran is not the only one.

A simillar case is that old Microsoft compilers used to store 16-bit
integers with 32-bit alognment. The 32-bit alignment was so that they
matched the storage used for reals as required by the standard, and the
16 bits was because that was faster for some processors of the day (yes,
it was a while ago). If I recall correctly, details of this could be
controlled by compiler switches.

To the extent that your "acceptable" means something more along the
lines of useful to users, I'd suggest that the existance of hardware
that supports 8-bit reals is pretty solid evidence that there is use for
such things. There is a lot of such hardware from multiple vendors. Many
applications need something a bit more than 64 bits, but don't need a
full 128 bits.

As to why those 80-bit reals are often aligned on 128-bit boundaries,
Glenn mentioned that is is sometimes more efficient to use that
alignment, and there is even some hardware that literally cannot deal
with data not aligned on at least 4-byte boundaries.

Plus, the Fortran standard requires the 128-bit alignment for doubles if
single is 64 bits. Yes, there are vendors who do pay attention to the
standard. And in the other direction, the standard tends to be
influenced by needs of vendors and users. The standard requires double
to have exactly twice the storage of single, but only requires that its
precision be "greater"; how much greater is deliberately left
unspecified.

[campbel...@gmail.com]

The standard is an ass !
Thank goodness we didn't get:
y = sin (x, kind=selected_real_kind(15) )

[Richard Maine]

<campbel...@gmail.com> wrote:

> This original post may have gone in another direction, but for the
> original example:
>
> n2 = size(profile)
>
> To get a value of -1389967296 is not acceptable performance, what ever the
> lame excuse.
>
> Fortran must be useable and not just standard conforming. I hope there are
> some who understand this.

I seriously doubt that this subthread will end up going anywhere useful.
I probably shouldn't answer at all, but I'll go with this one post. I
won't argue farther about it. (Really.)

Many things get wrong answers if you don't follow the rules (aka the
standard). Heck, your sin(x) example can get wrong answers even if you
do follow the rules. (Try it with very large values of x; some compilers
will try to cope, others probably give garbage, and the core problem is
that the best answer just isn't very well defined for finite precision.)
For almost every thing that one person thinks is a wrong answer, the
proposed "fix" usually makes other things wrong.

If you think that default integer should be 64 bits on machines with a
64-bit address space, well, as noted before, I actually tend to agree
with that. And then your n2=size(profile) would get the expected answer.
But then there would be a lot of people whose code would get wrong
answers because they had assumed defaults were 32 bits. You might call
that a lame excuse, but I suspect they would have equally hostile
replies. I think I'll jump out of the middle of that one, as I'm not one
of those people. But I rather doubt that you are going to make much
progress with that argument.

I'm not sure what other fix you might have in mind, but all the ones I
can think of would break various code that was formerly perfectly
standard conforming. For example, if you want to do something like have
SIZE return something other than default integer, then that breaks
existing things in fairly "obvious" ways. (No, I don't feel it worth
bothering to spell out some for people who might not think it so
obvious, because I have less than zero interest in arguing about how
important those cases might or might not be - the cases exist.) I
suspect the people whose code that would break would consider your
reason for breaking it "lame".

I would tend to agree that it would be nice if integer overflow
generated an exception. The Fortran standard certainly allows for things
like that, but it doesn't demand it. Why that doesn't tend to happen on
most machines is a topic outside of just the Fortran language.

The standard is actually there in order to help make things useable. It
certanly has compromises (and just plain imperfections), but that is the
core of its purpose. But if you don't follow it, then it isn't going to
help much in any case. (Nor is kvetching here.)

[Ron Shepard]

On 3/28/16 7:30 PM, campbel...@gmail.com wrote:
> The standard is an ass !
> Thank goodness we didn't get:
> y = sin (x, kind=selected_real_kind(15) )

I think you are making much ado over nothing regarding the SIZE()
intrinsic, which is perfectly understandable and perhaps even intuitive.
But if you think otherwise, then I recommend that you avoid using the
CMPLX() intrinsic function.

$.02 -Ron Shepard :-)

[herrman...@gmail.com]

On Monday, March 28, 2016 at 5:10:57 PM UTC-7, Richard Maine wrote:

> <campbe...@gmail.com> wrote:

> > How did it become acceptable for gFortran to store real*10 in 16 bytes ?

> To the extent that your "acceptable" means "acceptable to the standard",
> i.e. standard conforming, the standard has allowed things like that from
> day one. Storage alignment requirements are distinct from precision
> requirements. Many examples of mismatches between the two have existed.

Day one gave 36 bit REAL and 16 bit INTEGER in 36 bit words on the 704,
so yes not all bits were used.

> As noted elsewhere, I have used multiple compilers that stored 80-bit
> real values with 128-bit alignment. No, gFortran is not the only one.

> A simillar case is that old Microsoft compilers used to store 16-bit
> integers with 32-bit alognment. The 32-bit alignment was so that they
> matched the storage used for reals as required by the standard, and the
> 16 bits was because that was faster for some processors of the day (yes,
> it was a while ago). If I recall correctly, details of this could be
> controlled by compiler switches.

Pretty much all compilers I knew for 16 bit minicomputers, such
as the PDP-11, used 16 bits for INTEGER, and 32 bits for REAL.
With 64K (or often only 16K or 32K) you don't waste bits.

The compilers likely have an option to generate 32 bit integers,
and use only 16 of the bits. I never knew anyone to use the option.

The DEC PDP-11 compilers have INTEGER*2 and INTEGER*4, both
holding, as well as I remember, 16 bit integers.

> To the extent that your "acceptable" means something more along the
> lines of useful to users, I'd suggest that the existance of hardware
> that supports 8-bit reals is pretty solid evidence that there is use for
> such things. There is a lot of such hardware from multiple vendors. Many
> applications need something a bit more than 64 bits, but don't need a
> full 128 bits.

> As to why those 80-bit reals are often aligned on 128-bit boundaries,
> Glenn mentioned that is is sometimes more efficient to use that
> alignment, and there is even some hardware that literally cannot deal
> with data not aligned on at least 4-byte boundaries.

Pretty common for RISC processors.

In the Fortran 66 days, it was believed that COMMON didn't
allow padding. Compilers generated unaligned variables,
and then fixed up the alignment exception at run-time.

Fortunately, that has changed.

> Plus, the Fortran standard requires the 128-bit alignment for doubles if
> single is 64 bits. Yes, there are vendors who do pay attention to the
> standard. And in the other direction, the standard tends to be
> influenced by needs of vendors and users. The standard requires double
> to have exactly twice the storage of single, but only requires that its
> precision be "greater"; how much greater is deliberately left
> unspecified.

DOUBLE PRECISION is required to be twice as large, though I don't
believe that the standard requires a specific alignment for it.
(They could be 128 bits big, but only need 64 bit alignment.)

I was always disappointed that 128 bit floating point hardware,
and software to support it, was so rare.

Early in my Fortran programming (maybe my second program),
I tried to compute square roots with QSQRT. It seems that only
the more expensive H extended compiler supported them.

I had the manual, but no access to the compiler.

IBM supported REAL*16 back to the 360/85, and on all
S/370 and successor models. (Except that divide, DXR,
was done in software. IBM did studies to show that it was
rare enough that software was fast enough. In ESA/390
years, hardware DXR was finally added.)

The VAX architecture includes H-float, but it is an optional extra
charge microcode option on most processors.
The 11/730, low end machine, has as standard.

Numerically, as calculations get larger, there is
more need for intermediates with extra precision.

The 8087 temporary real came from the IEEE standard attempt to allow
for extra precision for intermediate values.

128 bit floating point is an optional feature of RISC V.
It will be interesting to see how many implement it.

[herrman...@gmail.com]

On Monday, March 28, 2016 at 5:30:52 PM UTC-7, campbel...@gmail.com wrote:

(snip)

> The standard is an ass !
> Thank goodness we didn't get:
> y = sin (x, kind=selected_real_kind(15) )

Don't forget CMPLX, though.

[James Van Buskirk]

wrote in message
news:19ad0ba4-f7bc-4188...@googlegroups.com...

> Fortran must be useable and not just standard conforming.
> I hope there are some who understand this.

> Ever used y = sin (x) and got an unacceptable answer

Yeah, y = sin(x) is quite controversial.

CASE 1:

x = 6.022d23
y = sin(x)

Now, some folks would like y = NaN because it effectively has negative
significant figures. Others think that x should be reduced during
execution of the SIN function modulo the exact value of 2*pi, no
matter how long it takes. Others think that reduction via the machine
representation of 2*pi in double precision should take place, and still
others might count on reduction with a special 66-bit approximation to
2*pi like the x87 FPU does. In the last two instances probably the user
hopes that there is a damping factor that wipes out these results.

CASE 2:

x = atan(1.0d0)+(0,1)*nearest(1.5*log(2.0d0)+log(huge(1.0d0)),-1.0d0)
y = sin(x)

ifort gives y = (1.797693134862125E+308,1.797693134862125E+308),
gfortran says y = ( Infinity, Infinity)

Do you want to guard against such cases that are close to overflow or
not take the time for such checks on every invocation of the SIN
function?

[Ron Shepard]

On 3/28/16 6:42 PM, campbel...@gmail.com wrote:
> Ever used y = sin (x) and got an unacceptable answer

For an appropriate definition of unacceptable, yes. Consider for example
how to do argument reduction for large values of x -- no matter how it
is done, there are situations in which the choice is the wrong way for
the application.

I have also gotten what might be an unacceptable answer for the
expression N1*N2 for large values of N1 and N2, and when N1 and N2 are
used to specify array dimensions, I would expect SIZE() to result in a
similar kind of unacceptable answer for pretty much the same reason.

$.02 -Ron Shepard

[herrman...@gmail.com]

On Monday, March 28, 2016 at 6:19:00 PM UTC-7, James Van Buskirk wrote:

(snip, someone wrote)

> > Ever used y = sin (x) and got an unacceptable answer

> Yeah, y = sin(x) is quite controversial.

> x = 6.022d23
> y = sin(x)

> Now, some folks would like y = NaN because it effectively has negative
> significant figures. Others think that x should be reduced during
> execution of the SIN function modulo the exact value of 2*pi, no
> matter how long it takes.

Reminds me when I was in high school, and someone I knew had
a brand new HP-55 calculator. (You can figure out when that was.)

I tried sind(9.999999999e99) that is, sine in degrees, and the
result was not zero.

Extended with zeros, it is a multiple of 180, so should be zero.

The calculator owner was sure that HP would never make a
mistake, and it had to be right.

That was after I knew what the OS/360 Fortran library did with
large operands to SIN and DSIN, so I did understand the lack
of significance at the time.

[William Clodius]

Richard Maine <nos...@see.signature> wrote:

> steve kargl <s...@REMOVEtroutmask.apl.washington.edu> wrote:
>
> > Richard Maine wrote:
>
> > > Probably the main reason that man compilers keep default real and
> > > integer at 32 bits even for 64-bit address spaces is that there is an
> > > awful lot of code written with the non-portable assumption that the
> > > default sizes are 32 bits....
>
> > I suspect the main reason is the consequence of storage association.
> > A 64-bit REAL requires a 128-bit DOUBLE PRECSION. Common hardware
> > does not support a 128-bit hardware floating point; so, tall the DOUBLE
> > PRECISION with use a software implementation (which is slowwwwwww).
> > I suppose one could enforce the 128-bit storage association, and then
> > do computations in whatever extended precision is supported by the
> > hardware (e.g., Intel 80-bit extended double precision).
>
> Yes, it would be perfectly standard conforming for double to use 80-bit
> reals stored in 128-bits for storage association purposes. I've used
> compilers that did exactly that. I don't see that as probably having
> been a big driver against making default real 64 bits. Granted I could
> be wrong, but that would not be my guess as to a significant reason.

But then double precision would take up twice the memory footprint, with
the corresponding impact on cache performance for large double precision
arrays. There is also te desirable cabability to read binary files
produced by other "processors" on the same machine, and having it
inconvenient to read a single precision real out of a file would
discourage users.

[campbel...@gmail.com]

Richard,

Thanks for providing your reply. I agree with your doubt that this subthread will end up going anywhere useful.
My point is basically that we have lots of generic intrinsic functions and SIZE should be able to respond correctly to the properties of it's argument; the array "profile" and provide a suitable value accordingly.
As 64-bit usage becomes more common, I think this useability feature will be visited again.

Ron,
My use of N1*N2 was because not all compilers I use accept 2900000000 or 2900000000_8 as a value for ipos = 2900000000. A 23gb array is a valid use of 64-bit Fortran for the area where I work. Try using 64-bit integer constants in a readable way.

James,
I was only commenting on the generic property of function SIN, where it provides a result based on the KIND property of the argument. I don't think it is much of a leap to ask SIZE to consider the properties of it's argument in a similar way.

Based on the comments this thread has received, it appears there are few on this forum with a similar experience of Fortran usage that I have known and learnt from. Changing times !

[James Van Buskirk]

wrote in message
news:b88bfb4b-bbe1-47bd...@googlegroups.com...

> James,
> I was only commenting on the generic property of function SIN,
> where it provides a result based on the KIND property of the
> argument. I don't think it is much of a leap to ask SIZE to consider
> the properties of it's argument in a similar way.

> Based on the comments this thread has received, it appears there
> are few on this forum with a similar experience of Fortran usage
> that I have known and learnt from. Changing times !

Everybody had to make the change to 64 bits, and the Fortran
community was ahead of the rest by a few years. The issue with
SIZE(x) is that it was set up as a default integer function a
quarter century ago when many compilers didn't even have
bigger integers than default.

Its result kind can't be set according to the kind of its argument
because you can force the kind of its argument to be different
from any integer kind.

Its kind can't be adjusted according to the magnitude of its
result because Fortran doesn't allow variable kinds at runtime
so that it always knows what specific function to invoke
when the code contains something like y = f(size(x)).

But the situation is by no means confined to Fortran inquiry
functions. All the F77 style interfaces like LAPACK pass array
dimensions as arguments. Should a version for each possible
selection of these integer arguments be a part of an LAPACK
library?

[kargl]

Suppose in the following 'program foo' and 'subroutine bar' are in
different files (to prevent in-lining).

program foo
real,allocatable :: x(:)
integer(8) n
read(*,*) n ! n can exceed huge(1_4).
allocate(x(n))
call bar(size(x))
end program

subroutine bar(i)
integer i
print *, i
end subroutine bar

How is bar() suppose to deal with size() when it may or may not return
an INTEGER(8) instead of INTEGER(4). The Fortran standard can be
updated to impose a restriction that SIZE can't be used as an actual
argument. Should the standard also restrict UBOUND, LBOUND, etc.

Note, gfortran supports an INTEGER(16) (on some systems), should
size() automatically return an integer with the largest decimal
range?

--
steve

[Stefano Zaghi]

Il giorno martedì 29 marzo 2016 02:30:52 UTC+2, campbel...@gmail.com ha scritto:
> The standard is an ass !

I totally disagree. The standard is a great value of Fortran: I often do not like what the "standard says", but I like that there is a Standard says something. The Medioevo of Vendors Extensions is just a distant memory for me.

A standard conforming code is of greater quality than a code doing the same thing, but not being standard.

> Fortran must be useable and not just standard conforming. I hope there are some who understand this.

Your original code is unusable, it does not work as expected, moreover it is not standard conforming. I am sure there are many people who understand this.

My best regards.

[Clive Page]

On 29/03/2016 01:10, Richard Maine wrote:
> that supports 8-bit reals is pretty solid evidence that there is use for

Now that would be something to admire. But I suspect you meant 8-byte
reals.


--
Clive Page

[kargl]

IEEE754-2009 defines binary16 (sometimes referred to as a short float).
It has 11-bit precision and and min and max exponents of -14 and 15.
This is not quite Richard's 8-bit reals, but it is close.

--
steve

[Richard Maine]

Oops. :-)

I've actually seen 16-bit reals. The Space Shuttle backup system
telemetry stream had some. They were basically IBM mainframe 32-bit
reals truncated to 16 bits. Really horrid precision; just enough to give
you a general notion of about how big the value was.

[steve kargl]

James Van Buskirk wrote:

>
> CASE 2:
>
> x = atan(1.0d0)+(0,1)*nearest(1.5*log(2.0d0)+log(huge(1.0d0)),-1.0d0)
> y = sin(x)
>
> ifort gives y = (1.797693134862125E+308,1.797693134862125E+308),
> gfortran says y = ( Infinity, Infinity)
>

I get the ifort result with gfortran 4.9.4, 5.3.1, and 6.0 on x86_64-*-FreeBSD.
The RHS of x expression should be constant-folded with GMP/MPFR. Using
the -fdump-tree-original option shows this:

x = __complex__ (7.8539816339744827899949086713604629039764404296875e-1, 7.1082243366422380859148688614368438720703125e+2);
y = __builtin_csin (x);

So, either __builtin_csin() (or simply csin()) seems broke on your system.

--
steve

[Tim Prince]

Intel introduced hardware support for float16 in SandyBridge. The data
have to be unpacked to 32-bit single precision before performing
arithmetic and repacked for storage, so it would be used only for
vectorizable array operations. The only floating point instructions
with corresponding half precision (but nearly full single precision
exponent width) are the reciprocal and reciprocal sqrt approximations.

https://software.intel.com/en-us/blogs/2013/09/30/intel-half-precision-floating-point-format-conversion-instructions

hints that it was intended for computer imaging using C compiler simd
intrinsics.
The Fortran compiler support for half precision appears to consist of
the half-precision options -Qimf-precision

[Steve Lionel]

On 3/29/2016 11:49 AM, Tim Prince wrote:
> Intel introduced hardware support for float16 in SandyBridge. The data
> have to be unpacked to 32-bit single precision before performing
> arithmetic and repacked for storage, so it would be used only for
> vectorizable array operations. The only floating point instructions
> with corresponding half precision (but nearly full single precision
> exponent width) are the reciprocal and reciprocal sqrt approximations.
>
> https://software.intel.com/en-us/blogs/2013/09/30/intel-half-precision-floating-point-format-conversion-instructions
>
> hints that it was intended for computer imaging using C compiler simd
> intrinsics.
> The Fortran compiler support for half precision appears to consist of
> the half-precision options -Qimf-precision

The 16-bit float type was specifically for graphics processing. Intel
Fortran has no support for it. /Qimf-precision (-fimf-precision) has
nothing to do with 16-bit floating point - it selects an alternate math
library that has "fast but less accurate" versions of some intrinsics.

--
Steve Lionel
Developer Products Division
Intel Corporation
Merrimack, NH

For email address, replace "invalid" with "com"

User communities for Intel Software Development Products
http://software.intel.com/en-us/forums/
Intel Software Development Products Support
http://software.intel.com/sites/support/
My Fortran blog
http://www.intel.com/software/drfortran

Refer to http://software.intel.com/en-us/articles/optimization-notice
for more information regarding performance and optimization choices in
Intel software products.

[herrman...@gmail.com]

On Tuesday, March 29, 2016 at 2:04:37 AM UTC-7, Clive Page wrote:
> On 29/03/2016 01:10, Richard Maine wrote:
> > that supports 8-bit reals is pretty solid evidence that there is use for

> Now that would be something to admire.

VAX has six bit floating point for literal addressing mode.

Three bit fraction (binary point on the left), three bit unsigned exponent.

It allows for memory efficient loading of many common floating point values.

[FX]

> So, either __builtin_csin() (or simply csin()) seems broke on your
> system.

Windows' math library does not have the best quality complex routines.

--
FX

[Tim Prince]

You mentioned the built-in fsin of the x87 CPU. That gives up and
returns x for sin(x), |x| >= 2^64. A common workaround used to be to
write a function which reduced the argument recursively by multiplies of
the 2pi approximation for the case where fsin signals out of range, but
that doesn't solve the "negative significant figures" question; not as
good as any of the accepted methods for guarding range reduction. All it
accomplishes is to assure |y| <=1.
The internal Pi constant is an ordinary long double; it just happens
luckily that the next 2 bits are 0, so no special precautions are needed
for |x| < 8*pi or so (or considerably wider range if only IEEE double
accuracy is needed).
x87 fsin isn't considered fast enough on recent CPUs.

[herrman...@gmail.com]

On Wednesday, March 30, 2016 at 4:50:10 AM UTC-7, Tim Prince wrote:
> On 3/28/2016 9:18 PM, James Van Buskirk wrote:

(snip)

> > Yeah, y = sin(x) is quite controversial.

> > CASE 1:

> > x = 6.022d23
> > y = sin(x)

(snip)

> You mentioned the built-in fsin of the x87 CPU. That gives up and
> returns x for sin(x), |x| >= 2^64.

I believe that Intel expects the x87 instructions to be either used
inside called functions, or the tests to be generated inline.

Seems that many compilers don't do that.

> A common workaround used to be to
> write a function which reduced the argument recursively by multiplies of
> the 2pi approximation for the case where fsin signals out of range, but
> that doesn't solve the "negative significant figures" question; not as
> good as any of the accepted methods for guarding range reduction. All it
> accomplishes is to assure |y| <=1.

There might be some pure math problems that need accurate sin(x)
for very large x. Those are rare.

For engineering and experimental science, either x won't get that big,
or the program will do its own range reduction.

> The internal Pi constant is an ordinary long double; it just happens
> luckily that the next 2 bits are 0, so no special precautions are needed
> for |x| < 8*pi or so (or considerably wider range if only IEEE double
> accuracy is needed).
> x87 fsin isn't considered fast enough on recent CPUs.

Reminds of the VAX instructions POLY and INDEX, for computing
polynomials, and for array index calculations including bounds checking,
which turned out slower than doing either using ordinary instructions.
(Maybe fixed in later processors, though.)

Is FSIN too slow, as other instructions have been improved over
the years, but FSIN is still the same?

[James Van Buskirk]

"Tim Prince" wrote in message news:dm1sre...@mid.individual.net...

> The internal Pi constant is an ordinary long double; it just happens
> luckily that the next 2 bits are 0, so no special precautions are needed
> for |x| < 8*pi or so (or considerably wider range if only IEEE double
> accuracy is needed).
> x87 fsin isn't considered fast enough on recent CPUs.

Old Intel docs give a constant definitely with 66 significant bits.
Here is an example, gcc version 5.3.0, Core 17-4770:

D:\gfortran\clf\piacc>type piacc.f90
! Compare with 4*0.C90FDAA22168C234C
! Pentium Processor User's Manual
! Volume 3: Architecture and Programming Manual
! Intel Corporation 1994
program piacc
implicit none
integer, parameter :: ep = selected_real_kind(18,4931)
integer, parameter :: dp = kind(1.0d0)
integer, parameter :: ik8 = selected_int_kind(18)
real(ep) x,y,twopi
integer(ik8) i
write(*,'(*(g0))') 'Characteristics of x. KIND = ',kind(x),', &
PRECISION = ',precision(x),', RANGE = ',range(x)
write(*,'(*(g0))') 'Characteristics of i. KIND = ',kind(i),', RANGE =
',range
(i)
x = int(Z'C90FDAA2',ik8)
twopi = 8*atan(1.0_ep)
i = nint(x/twopi,ik8)
write(*,'((g0,z0.16))') 'Initial i = ',i
write(*,'(*(g0))') 'Initial x = ',x
y = sin(x)
x = cos(x)
x = atan2(y,x)/(2*i)
write(*,'(*(g0))') 'Final x = ',x
i = transfer(real(x,dp),i)
write(*,'(g0,Z16.16)') 'Low bits = ',ishft(i,1)
end program piacc

D:\gfortran\clf\piacc>gfortran piacc.f90 -opiacc

D:\gfortran\clf\piacc>piacc
Characteristics of x. KIND = 10, PRECISION = 18, RANGE = 4931
Characteristics of i. KIND = 8, RANGE = 18
Initial i = 0000000020000000
Initial x = 3373259426.00000000000
Final x = -.121542010126079319519E-0009
Low bits = 7BC168C234C00000

So I seem to be getting all those low bits documented for my old
Pentium Classic on my more current Haswell. Where is my analysis
wrong?

[Tim Prince]

Certainly, it's more of the latter, and the willingness to accept less
than the extended precision accuracy. I was on the crew attempting to
deal with performance of the original Pentium 4, which was a relative
disaster running gcc/g77 or CVF code using fsin compared with SSE2
library code which could use the full instruction issue rate and support
vector parallel math functions. Subsequent improvements in simd
parallel math libraries give at least another factor of 2 in performance
over x87 built-ins. Combining x87 with simd register code always
suffered from the extreme penalty of lacking any direct register to
register data transfer.

[James Van Buskirk]

"steve kargl" wrote in message news:nde772$vjn$1...@dont-email.me...

That doesn't seem unlikely. Microsoft seems more interested in
forcing or tricking users in installing the hated Windows 10 than
providing a usable system. But just to make sure we're on the
same page, here is an example that uses constant expressions,
ordinary expressions, and specification expressions:
gcc version 5.3.0 (x86_64-win32-seh-rev0, Built by MinGW-W64 project)

D:\gfortran\clf\piacc>type subs.f90
subroutine sub1(x,y)
implicit none

integer, parameter :: dp = kind(1.0d0)

complex(dp) x,y
y = sin(x)
end subroutine sub1

function fun2(x)
implicit none

integer, parameter :: dp = kind(1.0d0)

complex(dp) x
character(int(aimag(sin(x))/0.4e308_dp)) fun2
fun2 = 'Hello world'
end function fun2

D:\gfortran\clf\piacc>gfortran -c subs.f90

D:\gfortran\clf\piacc>type pirange.f90
program pirange
implicit none

integer, parameter :: dp = kind(1.0d0)

complex(dp), parameter :: x = atan(1.0_dp)+ &
(0,1)*nearest(1.5*log(2.0_dp)+log(huge(1.0_dp)),-1.0_dp)
complex(dp), parameter :: y = sin(x)
complex(dp) s,t
interface
function fun2(x)
implicit none

integer, parameter :: dp = kind(1.0d0)

complex(dp) x
character(int(aimag(sin(x))/0.4e308_dp)) fun2
end function fun2
end interface
s = x
t = y
write(*,*) 's=',s
write(*,*) 't=',t
call sub1(s,t)
write(*,*) 't=',t
write(*,*) len(fun2(s))
end program pirange

D:\gfortran\clf\piacc>gfortran pirange.f90 subs.o -opirange

D:\gfortran\clf\piacc>pirange
s= ( 0.78539816339744828 , 710.82243366422381 )
t= ( 1.7976931348621249E+308, 1.7976931348621251E+308)
t= ( Infinity, Infinity)
0

With ifort I get:

D:\gfortran\clf\piacc>pirange
s= (0.785398163397448,710.822433664224)
t= (1.797693134862125E+308,1.797693134862125E+308)
t= (1.797693134862125E+308,1.797693134862125E+308)
4

[Tim Prince]

gfortran on cygwin64 (either 5.3.0 or 6.0) gives me

s= ( 0.78539816339744828 , 710.82243366422381 )
t= ( 1.7976931348621249E+308, 1.7976931348621251E+308)

t= ( NaN, Infinity)
0

I suspect csin may use x87 code (even with ifort equivalent?) as a more
straightforward way to give full double precision exponent range. But I
don't think Microsoft would take any responsibility for what cygwin
gfortran is doing. I'm not coming up immediately with ideas on how to
debug. 'gfortran -O -march=native jb*.f90 -Bstatic -lm -Bdynamic' seems
to be a start on the way to linking some static math library without
failing on account of there not being full static library support.

[steve kargl]

James Van Buskirk wrote:

> "steve kargl" wrote in message news:nde772$vjn$1...@dont-email.me...
>
>> James Van Buskirk wrote:
>
>
>> > CASE 2:
>
>> > x = atan(1.0d0)+(0,1)*nearest(1.5*log(2.0d0)+log(huge(1.0d0)),-1.0d0)
>> > y = sin(x)
>
>> > ifort gives y = (1.797693134862125E+308,1.797693134862125E+308),
>> > gfortran says y = ( Infinity, Infinity)
>
>> I get the ifort result with gfortran 4.9.4, 5.3.1, and 6.0 on
>> x86_64-*-FreeBSD.
>> The RHS of x expression should be constant-folded with GMP/MPFR. Using
>> the -fdump-tree-original option shows this:
>
>> x = __complex__ (7.8539816339744827899949086713604629039764404296875e-1,
>> 7.1082243366422380859148688614368438720703125e+2);
>> y = __builtin_csin (x);
>
>> So, either __builtin_csin() (or simply csin()) seems broke on your
>> system.
>
> That doesn't seem unlikely. Microsoft seems more interested in
> forcing or tricking users in installing the hated Windows 10 than
> providing a usable system. But just to make sure we're on the
> same page, here is an example that uses constant expressions,
> ordinary expressions, and specification expressions:
> gcc version 5.3.0 (x86_64-win32-seh-rev0, Built by MinGW-W64 project)
>

(code removed for brevity)

> D:\gfortran\clf\piacc>gfortran pirange.f90 subs.o -opirange
>
> D:\gfortran\clf\piacc>pirange
> s= ( 0.78539816339744828 , 710.82243366422381 )
> t= ( 1.7976931348621249E+308, 1.7976931348621251E+308)
> t= ( Infinity, Infinity)
> 0
>
> With ifort I get:
>
> D:\gfortran\clf\piacc>pirange
> s= (0.785398163397448,710.822433664224)
> t= (1.797693134862125E+308,1.797693134862125E+308)
> t= (1.797693134862125E+308,1.797693134862125E+308)
> 4

We're on the same page to the extent that the 3 versiosn of
gfortran that I have on FreeBSD agrees with your ifort results.
It may be possible to get gfortran on Windows to use ifort's
math library, but I have no idea how.

--
steve

[steve kargl]

I would start with at least identifying what routine is called.

% nm z | grep csin
U csin@@FBSD_1.3

So,csin from libm.so.3 is being called on FreeBSD. Using
-static I get

% nm z | grep csin
000000000041c2d0 T csin
000000000041c050 T csinh

Openlibm supposedly compiles on Windows. Perhaps, it would be possible
to link with openlibm and intercept the calls to csin and csinh.

--
steve

[robin....@gmail.com]

On Tuesday, March 29, 2016 at 12:19:00 PM UTC+11, James Van Buskirk wrote:
> wrote in message
> news:19ad0ba4-f7bc-4188...@googlegroups.com...
>
> > Fortran must be useable and not just standard conforming.
> > I hope there are some who understand this.
>
> > Ever used y = sin (x) and got an unacceptable answer
>
> Yeah, y = sin(x) is quite controversial.
>
> CASE 1:
>
> x = 6.022d23
> y = sin(x)
>
> Now, some folks would like y = NaN because it effectively has negative
> significant figures. Others think that x should be reduced during
> execution of the SIN function modulo the exact value of 2*pi, no
> matter how long it takes. Others think that reduction via the machine
> representation of 2*pi in double precision should take place, and still
> others might count on reduction with a special 66-bit approximation to
> 2*pi like the x87 FPU does. In the last two instances probably the user
> hopes that there is a damping factor that wipes out these results.
>
> CASE 2:
>
> x = atan(1.0d0)+(0,1)*nearest(1.5*log(2.0d0)+log(huge(1.0d0)),-1.0d0)
> y = sin(x)
>
> ifort gives y = (1.797693134862125E+308,1.797693134862125E+308),
> gfortran says y = ( Infinity, Infinity)

Why is this important?

There are not 308 significant digits in the mantissa,
and on other hardware, such large values might not be achievable.

[Tim Prince]

Even if someone comes along with the Lahey version of gfortran, which
links against Visual Studio libraries, it's still unlikely that ifort
run-time could be used without a bunch of wrappers to implement
undocumented internal functionality.

[Tim Prince]

When I linked using the unusual command I quoted, I was able to go
slightly beyond this stage, and see the entry point in csin under both
gdb and objdump. No doubt, an expert in newlib could facilitate finding
the associated source code. Judging by file names, I guess there may be
some common origins shared between openlibm and newlib.
We're so far off the original topic of this thread that I can't see
whether anyone is interested in the question of why Windows Fortran has
less support than BSD for sin(ridiculously large argument). I only
jumped in to point out that the answer is not to be found in Microsoft's
corner.
openlibm is an interesting suggestion. There are so many places to
download it that I fear it will be as big a can of worms for the
beginner as mingw.