0.0000000567 becomes 0.000000000, why?

0.0000000567 becomes 0.000000000, why?

Hi, there is a parameter sigma which should be equal to 5.67*10**(-8) in my programming.  When the code is like this

          sigma= 5.67*10**(-8)

         L=0.1 !m

         a=1.0

         e=1.0

         ua=2000.0 !K

         dx=L/21.0

         t=sigma*dx*(a*ua**4-e*298.0**4)

        print*, t

the results is t=0.00000000

However, when the code is like this

         sigma= 0.0000000567

         L=0.1 !m

         a=1.0

         e=1.0

         ua=2000.0 !K

         dx=L/21.0

         t=sigma*dx*(a*ua**4-e*298.0**4)

        print*, t

I get the exact result that t=4317.871

I don't know why I got two different results because of the sigma function. Could you tell me why? 

                 

31 posts / 0 new
Last post
For more complete information about compiler optimizations, see our Optimization Notice.
Best Reply

You don't tell us what compiler you are using, or compiler version, operating system with version, or hardware.  Therefore, I can only guess as to what is going on.  Right now I don't have the time to run your source code through a compiler to see what, if anything, changes.

It also appears that you have little understanding of how Fortran represents numbers or how Fortran evaluates arithmetic expressions.

Here is my guess.  It appears that the root cause of your problem is roundoff in operations on floating-point numbers.  With numbers that can include a fractional component, the computer cannot represent most real numbers exactly, so they are represented by the closest approximation.  When you perform arithmetic operations on such numbers, the computer also needs to use the closest approximation for the result.  Hence, with every operation, there is a small amount of error between the exact result and the closest approximation that the computer can generate.

All Fortran compilers must provide at least two methods of representing floating-point, also known as real, numbers.  These are known as single precision and double precision.  Many Fortran compilers provide more than two.  For example, the GFortran compiler on Windows provides four different precisions.  Most, but not all, Fortran compilers use the IEEE definitions for single and double precision.  With these definitions, single precision is 7 decimal digits and an exponent range of 38, double precision is 15 decimal digits and an exponent range of 308.

By rule, in Fortran, default real numbers are single precision.  If you want any other precision, you need to specify it using a kind type parameter when you define a variable or you use a numeric constant.  In your example, all of the numeric constants do not have kind type parameters after them, so they are single precision.  Most likely, they are accurate to only 7 decimal digits.

In the first version, the program is computing sigma by a series of arithmetic operations.  The compiler may be computing a value that is too small to be represented, so it may be assigning a value of zero to sigma.  Again, this is only a guess.

I need to point out that in Fortran, and many other languages, using notation like 10**(-8) to specify exponential notation is not correct.  So instead of writing:

sigma = 5.67*10**(-8)

you should write:

sigma = 5.67E-08

The "E" in the notation stands for 10 raised to the power of the exponent after the "E".

I strongly recommend that programmers use double precision for calculating results whenever possible because you get far less roundoff with each arithmetic operation.  I have a module of commonly used named constants that I use in all programs and procedures.  Among many other things I define the named constants for all of the precisions that I need for real number.  E.g., I use these definitions:

Integer, Parameter :: SP = Selected_Real_Kind (P=6,R=35)
Integer, Parameter :: DP = Selected_Real_Kind (P=15,R=300)

Using these definitions you would define sigma as double precision by writing:

sigma = 5.67E-08_DP

Hope this helps.

Quote:

Baosheng B. wrote:

          0.0000000567 becomes 0.000000000, why?

Try this very short Fortran program, and consult Fortran textbooks and manuals to understand why the result printed is correct, even if it is not what you expect.

program real

	real sigma

	sigma=10**(-1)

	write(*,*)sigma

	end

Thank you Craig Dedo, it works. The fortran I'm using is fortran 90, intel compiler.

Baosheng,

Mecej4 indicated the problem with your coding, perhaps the following will more clearly show the problem.

    real*8 sigma, s

	    sigma = 5.8*10**(-8)

	    s     = 5.8*10.**(-8)

	    write (*,*) sigma,s

	    end

The following perhaps better shows what is required to improve the precison due to the implied type and kind of each operation.

    real*8 sigma, s, s8,ss8

	    sigma = 5.8   * 10   **(-8)

	    s     = 5.8   * 10.  **(-8)

	    s8    = 5.8d0 * 10.  **(-8)

	    ss8   = 5.8d0 * 10.d0**(-8)

	    write (*,*) sigma

	    write (*,*) s

	    write (*,*) s8

	    write (*,*) ss8

	    write (*,*) 5.8d-8

	    end

Craig,

I can understand the use of kind qualifiers, as in:
    Integer, Parameter :: I8 = Selected_Int_Kind  (R=10)
    write (*,*) 10_i8**10, 10**10

But, Why "5.67E-08_DP" !!

John

Hi Mecej4 I tried what you said, but it didn't work,  sigma is still written as 0.00000 on the screen. Thank  you all the same.

Hi John, I tried what you said. they worked except this one:

03

        sigma = 5.8   * 10   **(-8)

Without a "." or "d0" after the values, it still shows 0.0000. Thank you very much. 

Baosheng,

I wonder whether this last case is a matter of integer arithmetic.  You are raising an integer to an integer power, and so the result is an integer.  In this case, 1E-8 as an integer is zero.

David

The problem is integer arithmetic, due to the precidence of **

5.8*10**(-8)  becomes   5.8 * ( 10**(-8) ) and 10**(-8) becomes integer = 0

The problem with ( 10.**(-8) ) is that this becomes real(4) so 1.e-8 with only 7 figures of accuracy.

To further explain the alternatives I provided, to identify both type and kind problems, when comparing 10.d0 ** (-8) to 10. ** (-8), even though 10._4 converts to the same value as 10._8 ( the extra ~32 bits are all zero), the resulting value becomes real(4) and 1.e-8 differs from 1.d-8, as the extra ~32 bits of 1.d-8 are not zero. I'd expect that the same 64-bit calculation would be performed, but for 10._4 the result returned would be converted/truncated to 4 bytes. ( I think there are compiler options to retain the register values which may retain the precision in calculating the s8 example I provided)

Similarly s = 5.8_4 and s = 5.8_8 are different, while s = 58._4 and s = 58._8 would result in the same value for s, as the value 58. is a "rational" binary number.

John

Hi David,

Thank you for your reminding.  it's a matter of integer arithmetic. I tried other examples, let sigma1=30**(-3), sigma2=30.**(-3), sigma3=30**(-3.0), sigma4=30.**(-3.0), respectively. Only does sigma1 equal to zero, the rest of other three sigmas get the same right value.

Hi John,

I see. If I want to get a real number, not a integer or zero, there should be at least a real number(following a dot "." or something) before or after "**".

Baosheng

Tips: Whenever you need to represent any "exact" number (say sigma = 5.67E-08, tau=0.99991200) in fortran, just type sigma = 5.67D-8 and tau=0.999912d0 for short and you get the accurate recognition by compilers. 

David and Baosheng:

Yes, it is a matter of integer arithmetic.  I completely overlooked this very important issue in my rather lengthy answer yesterday.  My bad.  Still, everything I said in my post of yesterday is correct, as far as it goes.

In the assignment statement, sigma = 5.67*10**(-8), exponentiation has a higher precedence than multiplication or division.  Therefore, the exponentiation is performed first.  In exponentiation between two integers, the Fortran standard requires the use of integer arithmetic.  Therefore, 10**(-8) is evaluated as 1 / 100000000, using integer division, not floating-point division.  Thus the result is zero, not 1.0E-08.  This is specifically required by the rules in sections 7.1.5.2.1 and 7.1.5.2.2 of the Fortran 2008 standard.  If you want the quotation of these rules, I will be happy to provide them in a separate message.

Thus, if you want to use the form of 10**exponent instead of the standard exponential notation, you will need to include a decimal point in the base and possibly also in the exponent.  E.g., you would write sigma = 5.67*10.0**(-8) or sigma = 5.67*10.0**(-8.0).  Of course, if you did not specify any kind type parameter, you would by default get single precision, which on most Fortran compilers is 7 decimal digits and an exponent range of 38.  It is still very much preferred to use the standard method of exponential notation using E+nn or E-nn.

By the way, in case you are wondering, in an operation of raising a integer base to a real or complex exponent, the base is first converted to real or complex.  Raising a negative-valued real base to a real power is prohibited.

Hope this helps.

Quote:

John Campbell wrote:

Craig,

I can understand the use of kind qualifiers, as in:
    Integer, Parameter :: I8 = Selected_Int_Kind  (R=10)
    write (*,*) 10_i8**10, 10**10

But, Why "5.67E-08_DP" !!

John

Using constant_kind, where constant is a numeric, character, or logical constant and kind is the kind type parameter, is the standard method of specifying kind type parameters for constants since Fortran 90 was published in June 1991.  Thus, 5.67E-08_DP means exactly the same thing as 5.67D-08.  Without a trailing kind parameter, if you use E for exponent notation, you get single precision and if you use D for exponent notation, you get double precision.  In order to specify kinds other than single or double precision, you have to use the kind parameter at the end.  I find that using kind parameters for everything is neater and more consistent that sticking with the older practice from Fortran 77.

We've travelled a long way on this post, where the answer was probably just: use  "sigma = 5.67E-8"

Quote:

John Campbell wrote:

We've travelled a long way on this post, where the answer was probably just: use  "sigma = 5.67E-8"

Why do you compound the confusion now by dropping the "_DP" part; Craig Dedo's post was a highly valuable lesson in its usefulness.

FortranFan,

I suggested sigma = 5.67e-8 as there was no initial indication of the kind of sigma, although t was reported to 7 figures.

While 5.67E-08_DP may provide a valuable lesson in the definition of precision, I regard this as an "overloaded" constant, where 5.67d-8 is more typical. I don't think we should adopt 5.67E-08_8 ? Can I describe this as the more typical way of defining a real(8) constant ?

Craig Dedo's post is a valuable lesson as 5.67E-08_16 is a necessary way of defining a real*16 constant, which should probably be provided as 5.67E-08_QP.

John

Quote:

Craig Dedo wrote:

...

Among many other things I define the named constants for all of the precisions that I need for real number.  E.g., I use these definitions:


...
   Integer, Parameter :: DP = Selected_Real_Kind (P=15,R=300)

...

Note the above-mentioned SELECTED_REAL_KIND function standardized with Fortran 95 has at least two parts: the precision part with the P argument and the range part with the R argument.

Also, note most FORTRAN 77 coders are used to "DOUBLE PRECISION" attribute to their data types for real numbers requiring more than single precision.  Now, note on most systems (IA32 processors, etc.), types declared with "DOUBLE PRECISION" attribute have a precision of 15 digits and range of 307 (roughly speaking, the type can represent numbers between -1E+307 and -1E-307 and +1E-307 and +1E+307).

So if one wants to be completely consistent with code written with DOUBLE PRECISION attribute in FORTRAN 77 (e.g., procedure calls involving older libraries/DLLs) , one should use KIND and RANGE functions on a given system to determine the precision and range of such types.  One may find a fully consistent named constant is

   INTEGER, PARAMETER :: DP = SELECTED_REAL_KIND(P=15, R=307)

However, one can also take a shortcut as

   INTEGER, PARAMETER :: DP = KIND(1D0)

which is a style I'm not fond of.

Quote:

John Campbell wrote:

FortranFan,

I suggested sigma = 5.67e-8 as there was no initial indication of the kind of sigma, although t was reported to 7 figures.

While 5.67E-08_DP may provide a valuable lesson in the definition of precision, I regard this as an "overloaded" constant, where 5.67d-8 is more typical. I don't think we should adopt 5.67E-08_8 ? Can I describe this as the more typical way of defining a real(8) constant ?

Craig Dedo's post is a valuable lesson as 5.67E-08_16 is a necessary way of defining a real*16 constant, which should probably be provided as 5.67E-08_QP.

John

A sidebar: as you know, asterisk (*) declarations such as INTEGER*2, REAL*4, REAL*8, CHARACTER* are marked for obsolescence in the newer standards even though most compilers still support them.  I prefer not to mention them in communications at all since it's a convention I don't like but which is too sticky.

Re: "5.67d-8 is more typical. I don't think we should adopt 5.67E-08_8", well most of our code adopts the convention explained by Craig above and use either "_DP/_dp".  Use of "D" (0D0, etc.) is generally avoided.  For generic code intended to work with a variety of real kinds, we also use "_WP/_wp" convention denoting working precision.  You'll also find it used widely in a recent book, "Numerical Computing with Modern Fortran" by Hanson and Hopkins that was reviewed by Steve in his Dr Fortran blog.

Fortranfan,

I am not a fan of SELECTED_REAL_KIND, as I find what it implies to be very misleading.

The concept that you can request a precision or exponent range that is required for the numerical calculation to be carried out is far removed from reality, as all you have is real*4, real*8 and either real*10 or real*16. Whenever I use SELECTED_REAL_KIND, I have to go and check what are the values for P= or R= that are consistent with the real*x I am going to use.

Even then; real*4 is never good enough, real*8 should work and if you require real*16, then you should rewrite your algorithm. How do you check a real*16 calculation ?

Have you ever developed a numerical computation by deciding I need 13 significant figures so set P=13 ? You don't consider the precision of the inputs, but the expected precision of the result and estimate what *4, *8 or better might provide. Given that consideration, you select the KIND required and then have to check which P will get it right. My background is numerical methods computation, so I don't have an interest in a verbose definition that does not define the accuracy of the result.

Enough of my runblings, I'll let you get back to Computing with Modern Fortran.

John

John, you misunderstand what SELECTED_REAL_KIND is for. The problem it solves is that there are no guarantees what the precision and range of various real kinds are, across implementations. Also, syntax such as real*4, etc., is non-standard. You use SELECTED_REAL_KIND to specify a minimum acceptable precision and range and let the implementation give you the kind that best meets your needs.

Steve - Intel Developer Support

Steve,

My apologies for continuing this discussion, as the standard and it's direction is not going to change, but to me SELECTED_REAL_KIND misunderstands what Fortran is for. There are few number precisions available on the existing hardware platforms in use, and the consideration of what precision to use is more based on what number formats are available and not the P= or R=, as byte size is often the primary consideration.
Having experienced most hardware types over the last 40 years, and received and converted code from many sources, the specification of real*8 is a much preferable syntax that any of the other offerings. Certainly, when this is not used in a code, it flags that there can be many problems ahead in using this code. While use of "DP = Selected_Real_Kind (P=15,R=300)" may be elegant, often it can be difficult to find with code listings, as it can be located in modules removed from the main code. Look at all the KIND definitions associated with ifort.
That there was no guidance for the value of kind in the F90 standard was a shortcoming of the standard, as with an intrinsic function to return the byte count, although I suppose the existence of votes from hardware representatives not based on 8-bit bytes didn't help.
I start from the position that Fortran is used for numerical computation and the standard should support that. I think that the recent standards have moved away from this and what Fortran is used for has either been forgotten or has changed significantly.

John

Quote:

John Campbell wrote:

The concept that you can request a precision or exponent range that is required for the numerical calculation to be carried out is far removed from reality, as all you have is real*4, real*8 and either real*10 or real*16. Whenever I use SELECTED_REAL_KIND, I have to go and check what are the values for P= or R= that are consistent with the real*x I am going to use.

John,

You are assuming that Fortran compilers support only three kinds:  4-byte, 8-byte, and either 10-byte or 16-byte.  This is a false assumption.

It is possible for  a Fortran compiler to support more than three (3) kinds of real.  This point is not a hypothetical.  GFortran for Windows supports four (4) kinds of real:  4-byte, 8-byte, 10-byte and 16-byte.  If you disbelieve this, then please look at the documentation and/or obtain the compiler and run some demo programs.  I have actually worked with both the 10-byte and 16-byte reals in the same program.  Except for a few isolated defects, they work as designed.

I should also point out that, by rule, the Fortran standard does not limit the number of kinds of a particular data type, with only two exceptions: (1) a Fortran compiler must provide at least two kinds of real; and (2) a Fortran compiler must provide a kind of complex for each kind of real.

 

Quote:

John Campbell wrote:

My apologies for continuing this discussion, as the standard and it's direction is not going to change, but to me SELECTED_REAL_KIND misunderstands what Fortran is for. There are few number precisions available on the existing hardware platforms in use, and the consideration of what precision to use is more based on what number formats are available and not the P= or R=, as byte size is often the primary consideration.
 

Byte size may often be the primary consideration, but it is not always the primary consideration.  There are real Fortran compilers that implement two different kinds of real in the same number of bytes.  In the past I worked with one of these.  It was the former DEC Fortran for VMS, now HP Fortran for OpenVMS.

This Fortran compiler had two different kinds of real that were both denoted by REAL*8.  One kind had 16 digits of decimal precision, an exponent range of 38, and was known as D_FLOAT.  The other had 15 decimal digits of precision, an exponent range of 308, and was known as G_FLOAT.  Which kind you got if you used REAL*8 notation depended solely on which compiler option you used, /D_FLOAT or /G_FLOAT.  However, at the operating system level, the data type designations were different, so the OS considered them (rightly) to be different data types.

 

Quote:

John Campbell wrote:

Having experienced most hardware types over the last 40 years, and received and converted code from many sources, the specification of real*8 is a much preferable syntax that any of the other offerings. Certainly, when this is not used in a code, it flags that there can be many problems ahead in using this code. While use of "DP = Selected_Real_Kind (P=15,R=300)" may be elegant, often it can be difficult to find with code listings, as it can be located in modules removed from the main code. Look at all the KIND definitions associated with ifort.
 

I strongly disagree with both of the points in this paragraph.

1.  I find the use of kind type parameters, preferably by using named constants defined using the Selected_*_Kind() intrinsic functions, to be far better than using the old, non-standard *n notation.  The use of kind type parameters that are named constants is very flexible and adapts well to situations in which there are more than common number of kind parameters for a particular data type.  This point is not a hypothetical.  There are real Fortran compilers, of commercial significance, in the market today, that support more than the common number of kinds of a particular data type.  E.g., GFortran supports 4 kinds of real:" 4, 8, 10, and 16 bytes.  GFortran supports 2 kinds of character:  ASCII (Kind=1) and ISO 10646 UCS-4 (Kind=4).  NAG supports 3 kinds of character:  ASCII (Kind=1), ISO 10646 UCS-2 (Kind=2), and ISO 10646 UCS-4 (Kind=4).

2.  It is not difficult to find named constants for kind type parameters in source code if the software developer has a minimal amount of discipline in how he/she writes the source code.  Almost always all of the kind type parameters are defined in one module and used everywhere else.  I define mine in a module I name Common_Definitions_M.  It's quite obvious what that module is about, isn't it?  Besides defining kind type parameters, I also define a few other useful items in this module.  There are named constants for pi, Golden Ratio, Silver Ratio, Apery Constant, Euler-Mascheroni Constant, and Yes and No as synonyms for .True. and .False. .  I also define generic procedures and operators for a few common operations such as Greatest Common Divisor (GCD), Least Common Multiple (LCM), Remainder, Modulo, and Adjusted Modulo.

 

Quote:

John Campbell wrote:

That there was no guidance for the value of kind in the F90 standard was a shortcoming of the standard, as with an intrinsic function to return the byte count, although I suppose the existence of votes from hardware representatives not based on 8-bit bytes didn't help.

The development of the kind type parameter for all intrinsic data types was part of the development of Fortran 90, so it occurred before my time on the Fortran Standard Technical Committee.  I joined the committee in November 1993.  I had to leave after November 2002 because I could no longer afford to attend the meetings.

The kind type parameter feature was actually very well thought through.  It originated in the Japanese comment on the Fortran 8X draft, which was the first public draft of the Fortran standard that was to become Fortran 90.  The Japanese suggested a kind type parameter for the Character data type because they wanted to support more than one character set, e.g., Hiragana, Katakana, and Kanji, in addition to ASCII.  The committee decided to accept this idea and extend it to all other data types.

The committee, very deliberately, decided to leave the number of kinds of each data type and the kind numbers as processor dependent.  As part of its standard operating procedures, the committee conscientiously tries to avoid any and all hardware dependencies.  This core value comes up constantly in committee deliberations.  Thus, there is no tie-in to byte counts.  What if some hardware is not addressable by 8-bit octets?  This has happened in the past and may happen again in the future.  Although right now memory is byte-addressable and floating-point arithmetic is according to the IEEE standard, this situation may not continue indefinitely into the future.

BTW, there are some Fortran compilers that do not use the byte count as the kind number.  The Salford, now Silverfrost, Fortran compiler uses kind numbers 1, 2, and 3 for 32-bit, 64-bit, and 80-bit reals.

 

Quote:

John Campbell wrote:

I start from the position that Fortran is used for numerical computation and the standard should support that. I think that the recent standards have moved away from this and what Fortran is used for has either been forgotten or has changed significantly.

Although primarily used for numerical computation, Fortran is not used exclusively for numerical computation.  Although Fortran is not used nearly to the extent that it could be in other application areas, it is a well-designed language for many application areas outside of the field of numeric computation.  I have used it in many areas that do not use much numeric computation.  I like Fortran because it has an easy-to-use syntactic structure.  The (non-numeric) CHARACTER data type is especially robust and easy to use, much more so than the character facilities of the C-like languages.

The Fortran Standards Technical Committee has not forgotten what the core market of Fortran is nor has it moved Fortran away from it.  Quite the contrary.  The major feature of the Fortran 2008 standard was Co-Arrays, a form of parallel processing.  This was done precisely to facilitate high-powered numeric computation by doing operations in parallel across a large number of CPUs.

 

@Craig. I can't disagree with the "correctness" of your extensive comments however I am still inclined to sympathise with John's point of view. When I see Integer, Parameter :: DP = Selected_Real_Kind (P=15,R=300) I have to go and look it up. This is something that might occur once in an application hidden away with some constants in some code we generally don't need to look at...

This is an Intel Fortran forum for Windows. The practical choices are DP=4 or DP=8. That might not work on a Salford compiler for the  "Deep thought" platform but if I want to port to that then changing DP=8 to DP=42 in one location is the least of my problems.... 

The  P=15,R=300 implies that there is a whole range of choices available, but that is misleading there aren't. I think in this case it is each to their own preference.

My advice is to write portable code whenever possible. That means using the SELECTED_xxx_KIND intrinsics and named constants for your kinds. I know this looks strange to those who have programmed Fortran for many years, but it works, is portable and also doesn't get you in trouble when the user starts throwing switches such as /r8.

Steve - Intel Developer Support

For any of us that happen to have a PDP-10, it was a 36-bit computer, which by convention held seven 7-bit bytes/characters per 36-bit word. Using (4) and (8) would be quite different on this machine. I think early Cray computers were 60-bit.

FWIW I have too much code using real(4) and real(8) to make a cosmetic change to use SELECTED_REAL_KIIND(...).

I wish to compliment Craig on his contributions to this thread.

Jim Dempsey

www.quickthreadprogramming.com

Quote:

Steve Lionel (Intel) wrote:
My advice is to write portable code whenever possible. That means using the SELECTED_xxx_KIND intrinsics and named constants for your kinds. I know this looks strange to those who have programmed Fortran for many years, but it works, is portable and also doesn't get you in trouble when the user starts throwing switches such as /r8.

I would agree entirely with striving for portability, and  using global named constants for KIND archives this. The single occurrence you would  of SELECTED_xxx_KIND doesn't really bring much to the party. If the real representations available on one platform are radically different you are going to have to apply some thinking, evaluation, testing and possibly rewriting of algorithms. SELECTED_xxx_KIND isn't a global panacea that is going so save you from grief unless the application is fairly trivial in the first place.

Quote:

app4619 wrote:

... The single occurrence you would  of SELECTED_xxx_KIND doesn't really bring much to the party. ... SELECTED_xxx_KIND isn't a global panacea ...

I think this is missing Craig's and Steve's points entirely.

Regardless, I've had to do a fair share of porting several non-trivial applications, mostly for process simulation involving partial differential equations, that were developed for IBM mainframes during 1960s and 1970s over to Windows and Linux platforms:

  • My experience was once the kinds of programs variables and constants were established correctly, I'd a solid foundation that made the rest of the process much easier; the intricacies related to any hardware representation, compiler options, linker optimizations, etc. were mostly taken care of.
  • I could then focus strictly on separating business logic (calculation algorithms and like) from user input/output i.e., pre- and post-processing aspects; also on my test suite and file management which can vary from platform to platform.
  • For any issues, I could back to the underlying math fairly independent of any Fortran language or compiler option considerations to figure out the path forward.

Leave a Comment

Please sign in to add a comment. Not a member? Join today