|Re: Just reasoning...|
Message #6 Posted by Dave Shaffer on 10 Aug 2003, 7:18 p.m.,
in response to message #4 by Vieira, Luiz C. (Brazil)
Here is my fortran program (with proper line breaks, I hope!):
C CALCULATES SUM OF POWERS OF RECIPROCAL INTEGERS
IMPLICIT REAL*8 (A-H,O-Z)
C DO A LOT OF NUMBERS TO GET TIMING
C TRY IT TWO WAYS - MULTIPLY & EPONENTIATION
SUM = 0D0
DO 100 N = 1,10000000
X = 1D0/DFLOAT(N)
C RUN IT ONE WAY OR THE OTHER, BY MOVING THE COMMENT
SUM = SUM + X*X*X
C SUM = SUM + X^3
C IF(MOD(N,50000) .EQ. 0) WRITE (*,2000) N,SUM
write (*,2000) n,sum
Nothing fancy. The line (commented out here, because I was doing the timing tests) just before the CONTINUE statement gives me the value after every 50,000 terms. I was NOT rounding off; just summing the results to the precision of the compiler/computer.
I just looked at those values to make my statement that it was between 50k and 100k terms that the sum settled down in the tenth decimal place. (I was too lazy to try to find exactly how many terms were needed!) It was looking at the printout that also indicated that somewhere around 200k-250k terms, the tenth decimal place should be rounded up (IF you were interested in only the "correct" properly rounded-off value to ten decimal places).
Double precision fortran should be good to 15 or 16 significant figures - plenty enough for this test. If amybody cares, there is quadruple precision, too, I think.
This doesn't offer much insight into HP calculator precision or round-off. Conversely, if you can define the algorithm that the calculator uses, it would be easy to check it this way.
In any case, Hugh's solution with alternating terms is far more elegant!