Rounding off double precision

I'm trying some calculations in double precision, and I'd like the
output rounded off so that errors in calculation are not shown. For
instance, the result of COS(90 * 3.141592653589793238d0 / 180.) is
-4.37113886E-06. I'd like it to be zero. Thanks for any help.

On Mar 29, 11:43 am, Bamm <bamms...@gmail.com> wrote:
> I'm trying some calculations in double precision, and I'd like the
> output rounded off so that errors in calculation are not shown. For
> instance, the result of COS(90 * 3.141592653589793238d0 / 180.) is
> -4.37113886E-06. I'd like it to be zero. Thanks for any help.

For that matter, I also wonder why the result of Cosine is real and
not double precision. Thanks for any clarification.

Bamm <bammster@gmail.com> wrote:

> On Mar 29, 11:43 am, Bamm <bamms...@gmail.com> wrote: 
>
> For that matter, I also wonder why the result of Cosine is real and
> not double precision. Thanks for any clarification.

Well...

1. You can't get "output rounded so that errors in calculation are not
shown." That is equivalent to asking for there to be no errors in
calculations. Ain't gonna happen. If you can establish what the order of
magnitude of the errors will be, then you can round the output
appropriately if you use an F edit descriptor. You aren't ever going to
get exactly zero out of an E edit descriptor (or list-directed output)
unless the number is indeed exactly zero.

2. The result of that particular COS *IS* double precision. If you are
seeing something else, then either a) you are misinterpreting the cause
of what you are seeing, or b) it is a compiler bug. I'd place my money
on (a). Quite a lot of it, and at good odds. :-) You need to show
exactly what is leading you to this conclusion. That means more than
just an isolated expression and an alleged result. I'd need to see a
complete program, most definitely including all declarations and
whatever output statement you used to get this result.  For example,
when I run the complete program

write (*,*) COS(90 * 3.141592653589793238d0 / 180.)
end

using g95 on this Mac, I get

6.123256244561421E-17

which is about the expected accuracy for double precision.

As an example of using an F edit descriptor to round, if I change the
program to

write (*,'(f15.12)') COS(90 * 3.141592653589793238d0 / 180.)
end

I then get

0.000000000000

--
Richard Maine                    | Good judgement comes from experience;
email: last name at domain . net | experience comes from bad judgement.
domain: summertriangle           |  -- Mark Twain

In article
<3277fc93-726f-428f-82cd-aae971d0711f@d4g2000prg.googlegroups.com>,
Bamm <bammster@gmail.com> wrote:

> I'm trying some calculations in double precision, and I'd like the
> output rounded off so that errors in calculation are not shown. For
> instance, the result of COS(90 * 3.141592653589793238d0 / 180.) is
> -4.37113886E-06. I'd like it to be zero. Thanks for any help.

The closest number is "never" going to be zero because there are
simply too many nonzero numbers near zero compared to floating point
numbers near the exact value of pi/2.  However, I get
6.12323399573676604E-017 for the result of your expression, which is
a lot closer to zero than your value.

$.02 -Ron Shepard

> when I run the complete program
>
>     write (*,*) COS(90 * 3.141592653589793238d0 / 180.)
>     end
>
> using g95 on this Mac, I get
>
>     6.123256244561421E-17
>
> which is about the expected accuracy for double precision.

Here's the answer I get for the same program above using g77 on Linux:

6.12303177E-17

Am I doing anything wrong?

In article <c97c106c-801f-48fa-8751-2e4e91cb028f@e6g2000prf.googlegroups.com
>,
Bamm  <bammster@gmail.com> wrote:

>Am I doing anything wrong?

Yes, you failed to show us the source code of the original program
which showed a much larger value.

-- greg

On Mar 29, 12:48=A0am, Bamm <bamms...@gmail.com> wrote: 
> 
> 
> 
> 
>
> Here's the answer I get for the same program above using g77 on Linux:
>
> =A0 6.12303177E-17
>
> Am I doing anything wrong?

Yes! You are expecting calculations done with floating point numbers
to be EXACT. Such calculations are actually approximate. Each
operation has a small quantity of error built in. That's the nature of
floating point.

To start with, you can not represent the EXACT value of PI with a
finite number of bits. It only gets worse from there. :-).

- e

Bamm wrote:
> On Mar 29, 11:43 am, Bamm <bamms...@gmail.com> wrote: 

You may want to investigate using fixed-point arithmetic?

http://en.wikipedia.org/wiki/Fixed-point_arithmetic

> For that matter, I also wonder why the result of Cosine is real and
> not double precision. Thanks for any clarification.

Don't know; here's what I'm seeing:

program precision_trial
integer, parameter :: dp = selected_real_kind(15,307)
integer, parameter :: sp = selected_real_kind(6,37)
write(*,*) cos( 90 * 3.141592653589793238_dp / 180. )
write(*,*) cos( 90 * 3.141592653589793238_sp / 180. )
end program precision_trial

$ g95 precision_trial.f90 && ./a.out
6.123256244561421E-17
-4.371139E-8

$ gfortran precision_trial.f90 && ./a.out
6.123233995736766E-017
-4.3711388E-08

$ g95 -v
gcc version 4.0.3 (g95 0.90!) Mar  7 2008

$ gfortran --version
GNU Fortran (GCC) 4.3.0 20070713 (experimental)

Regards,
--
http://twitter.com/bil_kleb

Bamm wrote: 
>
> Here's the answer I get for the same program above using g77 on Linux:
>
>   6.12303177E-17
>
> Am I doing anything wrong?
Probably not, different compilers print out a different number
of digits for list directed output.  The standard doesn't
specify how many digits to print out for an * format.  If
you want to experiment, try an explicit format.

double precision T
real R
T = COS(90 * 3.141592653589793238d0 / 180.)
R = T
write (*, '(E20.15)' )  T, R
write (*, '(F10.5)'  )  T, R
end

Try different format widths (the 10 or 20) and number of digits
(the 5 or 15) until you get something you like.

Also, it's not obvious, but for output, double precision values
are printed with an E exponent when you use the E edit descriptor.
You could try using D20.15.  Some processors print this with a D
exponent, others use an E.


Dick Hendrickson

> Probably not, different compilers print out a different number
> of digits for list directed output.  The standard doesn't
> specify how many digits to print out for an * format.  If
> you want to experiment, try an explicit format.
>
>        double precision T
>        real R
>        T = COS(90 * 3.141592653589793238d0 / 180.)
>        R = T
>        write (*, '(E20.15)' )  T, R
>        write (*, '(F10.5)'  )  T, R
>        end
>
> Try different format widths (the 10 or 20) and number of digits
> (the 5 or 15) until you get something you like.
>
> Also, it's not obvious, but for output, double precision values
> are printed with an E exponent when you use the E edit descriptor.
> You could try using D20.15.  Some processors print this with a D
> exponent, others use an E.
>
> Dick Hendrickson

Yey! That explains why I have an E exponent instead of D, and why I
have fewer digits displayed. You are right, it's not a real error, but
that the compiler just didn't show the entire value when using the *
format.