Re: whole number constants

On Jul 21, 2:30=A0pm, Dick Hendrickson <dick.hendrick...@[EMAIL PROTECTED]
> wrote:
> dancerch...@[EMAIL PROTECTED]
 wrote:
> > This may be a dumb question but what is "best practice" when coding
> > whole numbers in real (kind specified) formulas. =A0I am looking for
> > accuracy, consistancy and prortability reagardless of what value dp
> > has below. =A0Many of my calcs have numerical cancellations or
accuracy
> > checks as per the error function below. =A0What are the pitfalls?
>
> > for example:
>
> > integer, parameter :: dp =3D selected_real_kind(15, 307)
>
> > real(dp) =A0:: x,xnew
>
> > real(dp), parameter :: one =3D 1.0_dp
> > real(dp), parameter :: onend =3D 1_dp
> > real(dp), parameter :: onenk =3D 1.0
> > real(dp), parameter :: onei =3D 1
>
> For real constants that are exact integers and are small enough
> to fit in one word, those 4 are (almost) exactly equivalent.
> No processor in the world will possibly mess up the conversion
> of 1 or 1.0 or 1.0_dp into some other kind. =A0My advice would
> be to do what looks natural, as long as you declare the result
> kind, via the "real(dp)".
>
> The interesting exception here is the "1_dp" form. =A0This
> declares an integer of kind dp. =A0It's not guaranteed that
> this will be ****table. =A0Many processors only sup****t one kind
> of integer and the selected_real_kind refereence you used
> to get the value of "dp" isn't likely to return a valid
> integer kind.
>
> As you almost certainly know, it's different for values that
> aren't exactly representable. =A03.14 and 3.14_dp are not the
> same thing (or at least aren't ****tably the same).
>
> > ! which should I use ?
>
> > 10 error =3D one-xnew/x
> > 20 error =3D onend-xnew/x
> > 30 error =3D onenk-xnew/x
> > 40 error =3D onei-xnew/x
>
> > ! or the more ledgable ?
>
> > 50 error =3D 1-xnew/x
>
> This would be my choice, because it's clearer about
> what you are doing.
>

I like label 50 as well but then what about cancellation, for example:

real(dp), parameter :: pi =3D 3.1415926535_dp
real(dp), parameter :: two =3D 2.0_dp

real(dp) :: x,y

!  this looks good/ledgable but.....

x =3D 2*pi
y =3D x/2

!  does the following result...

if (y =3D=3D pi) then
!   aha cancellation
else
!   integer value of two assigned different values of 'extra' digits
!   (beyond precision of dp) ergo not a constant
end if

Sorry to be so nitpicky.