Prev: Call for participants
Next: (paypal payment) wholesale Lacoste shirts from china,and free shipping
From: Kai-Uwe Bux on 8 May 2010 19:25
Chad wrote:

> On May 8, 1:28 pm, Kai-Uwe Bux <jkherci...(a)gmx.net> wrote:
>> Chad wrote:
>> > The question is about the harmonic() function in the following
>> > code....
>>
>> > #include <stdio.h>
>>
>> > double harmonic(unsigned long n)
>> > {
>> > double acc = 0;
>> > int i;
>>
>> > for (i = 1; i <= n; i++) {
>> > acc += 1.0/i;
>> > }
>>
>> > return acc;
>> > }
>>
>> Hint: summing up the terms in reverse order (starting with the smallest)
>> handles rounding errors much better and is numerically preferable.
>>
>
> Why would summing up the terms in reverse order handle rounding errors
> much better?

This is due to the floating point format. It uses a mantissa (think
significant digits) and an exponent (think position of the decimal point).
For the sake of exposition, let us do the computation in base 10 and have
just two significant digits (so that the effects become visible early). Let
us compute harmonic(10000) two ways. If you sum from the front

1.0
+ .50
+ .33
+ .25
+ .20
+ .17
+ .14
+ .13
+ .11
+ .10
+ .091
....
+ .0001

When you add 1/21 = 0.048, the current value of the sum is roughly 3.x (with
only two significant digits!). But, 3.x + 0.048 =3.x when rounded to two
significant digits. That means that from i=21, the current value in the
summation does not change anymore. The value is stuck.

If you do it the other way around, the small values are actually accounted
for.

The precise numerical analysis is, of course, more subtle; but you should
get the idea. Generally, large values overpower small values in addition. To
minimize rounding errors, it pays to add numbers of roughly the same order
of magnitude: this way, the significant digits of the summands are lined up
properly and only few of them are ignored.


Best

Kai-Uwe Bux
From: Pascal J. Bourguignon on 8 May 2010 19:34
Chad <cdalten(a)gmail.com> writes:

> On May 8, 1:28�pm, Kai-Uwe Bux <jkherci...(a)gmx.net> wrote:
>> Chad wrote:
>> > The question is about the harmonic() function in the following
>> > code....
>>
>> > #include <stdio.h>
>>
>> > double harmonic(unsigned long n)
>> > {
>> > � double acc = 0;
>> > � int i;
>>
>> > � for (i = 1; i <= n; i++) {
>> > � � acc += 1.0/i;
>> > � }
>>
>> > � return acc;
>> > }
>>
>> Hint: summing up the terms in reverse order (starting with the smallest)
>> handles rounding errors much better and is numerically preferable.
>>
>
> Why would summing up the terms in reverse order handle rounding errors
> much better?

How much is 1.0e0 + 1.0e-8 ?
How much is 1.0e0 + 1.0e-7 ?

What is the result of:
{
float f=1.0e0;
int i;
for(i=0;i<10;i++){
f+=1.0e-8;
}
return(f);
}
vs.
{
float f=0.0e0;
int i;
for(i=0;i<10;i++){
f+=1.0e-8;
}
return(1.0e0+f);
}
?

(For double, s/8/16/ and s/7/15/).

--
__Pascal Bourguignon__
From: Chad on 8 May 2010 20:15
On May 8, 4:25 pm, Kai-Uwe Bux <jkherci...(a)gmx.net> wrote:
> Chad wrote:
> > On May 8, 1:28 pm, Kai-Uwe Bux <jkherci...(a)gmx.net> wrote:
> >> Chad wrote:
> >> > The question is about the harmonic() function in the following
> >> > code....
>
> >> > #include <stdio.h>
>
> >> > double harmonic(unsigned long n)
> >> > {
> >> > double acc = 0;
> >> > int i;
>
> >> > for (i = 1; i <= n; i++) {
> >> > acc += 1.0/i;
> >> > }
>
> >> > return acc;
> >> > }
>
> >> Hint: summing up the terms in reverse order (starting with the smallest)
> >> handles rounding errors much better and is numerically preferable.
>
> > Why would summing up the terms in reverse order handle rounding errors
> > much better?
>
> This is due to the floating point format. It uses a mantissa (think
> significant digits) and an exponent (think position of the decimal point)..
> For the sake of exposition, let us do the computation in base 10 and have
> just two significant digits (so that the effects become visible early). Let
> us compute harmonic(10000) two ways. If you sum from the front
>
>   1.0
> +  .50
> +  .33
> +  .25
> +  .20
> +  .17
> +  .14
> +  .13
> +  .11
> +  .10
> +  .091
> ...
> +  .0001
>
> When you add 1/21 = 0.048, the current value of the sum is roughly 3.x (with
> only two significant digits!). But, 3.x + 0.048 =3.x when rounded to two
> significant digits. That means that from i=21, the current value in the
> summation does not change anymore. The value is stuck.
>
> If you do it the other way around, the small values are actually accounted
> for.
>
> The precise numerical analysis is, of course, more subtle; but you should
> get the idea. Generally, large values overpower small values in addition. To
> minimize rounding errors, it pays to add numbers of roughly the same order
> of magnitude: this way, the significant digits of the summands are lined up
> properly and only few of them are ignored.
>

Maybe I'm missing the broader concept here, I tried this in python. I
get the same results.

Here is what I get when I go one way...

[cdalten(a)localhost oakland]$ python
Python 2.6.2 (r262:71600, May 3 2009, 17:04:44)
[GCC 4.1.1 20061011 (Red Hat 4.1.1-30)] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> acc = 1.0
>>> for x in range(1,50):
.... acc = acc + 1.0/x
.... print acc
....
2.0
2.5
2.83333333333
3.08333333333
3.28333333333
3.45
3.59285714286
3.71785714286
3.82896825397
3.92896825397
4.01987734488
4.10321067821
4.18013375513
4.25156232656
4.31822899323
4.38072899323
4.43955252264
4.4951080782
4.54773965714
4.59773965714
4.64535870476
4.69081325022
4.73429151109
4.77595817775
4.81595817775
4.85441971622
4.89145675325
4.92717103897
4.96165379759
4.99498713092
5.02724519544
5.05849519544
5.08879822574
5.11820999045
5.14678141902
5.17455919679
5.20158622382
5.2279020133
5.25354303894
5.27854303894
5.30293328284
5.32674280665
5.3499986206
5.37272589333
5.39494811555
5.41668724599
5.43796384173
5.45879717506
5.47920533833
>>>

And here is what I get when I reverse it...

[cdalten(a)localhost oakland]$ python
Python 2.6.2 (r262:71600, May 3 2009, 17:04:44)
[GCC 4.1.1 20061011 (Red Hat 4.1.1-30)] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> acc = 1.0
>>> for x in reversed(range(1,50)):
.... acc = acc + 1.0/x
.... print acc
....
1.02040816327
1.0412414966
1.06251809234
1.08425722278
1.106479445
1.12920671773
1.15246253168
1.17627205549
1.20066229939
1.22566229939
1.25130332503
1.27761911451
1.30464614153
1.33242391931
1.36099534788
1.39040711259
1.42071014289
1.45196014289
1.48421820741
1.51755154074
1.55203429936
1.58774858508
1.62478562211
1.66324716058
1.70324716058
1.74491382724
1.78839208811
1.83384663357
1.88146568119
1.93146568119
1.98409726013
2.03965281569
2.0984763451
2.1609763451
2.22764301177
2.2990715832
2.37599466012
2.45932799345
2.55023708436
2.65023708436
2.76134819547
2.88634819547
3.02920533833
3.195872005
3.395872005
3.645872005
3.97920533833
4.47920533833
5.47920533833
>>>


From: bart.c on 8 May 2010 21:03
Chad wrote:
> On May 8, 4:25 pm, Kai-Uwe Bux <jkherci...(a)gmx.net> wrote:
>> Chad wrote:
>>> On May 8, 1:28 pm, Kai-Uwe Bux <jkherci...(a)gmx.net> wrote:
>>>> Chad wrote:

>>>>> for (i = 1; i <= n; i++) {
>>>>> acc += 1.0/i;
>>>>> }

>>>> Hint: summing up the terms in reverse order (starting with the
>>>> smallest) handles rounding errors much better and is numerically
>>>> preferable.
>>
>>> Why would summing up the terms in reverse order handle rounding
>>> errors much better?

> Maybe I'm missing the broader concept here, I tried this in python. I
> get the same results.
>
> Here is what I get when I go one way...

>>>> acc = 1.0
>>>> for x in range(1,50):
> ... acc = acc + 1.0/x
> ... print acc

> 5.47920533833


> And here is what I get when I reverse it...

>>>> acc = 1.0
>>>> for x in reversed(range(1,50)):
> ... acc = acc + 1.0/x
> ... print acc

> 5.47920533833

I think you need to print more decimal places. Where you have ...833...
above, I get differences after ...832942...

Also you're starting with acc=1.0, which itself will hide some accuracy (so
when adding 0.02, you're losing another couple of decimals). Perhaps start
with acc=0.0 and add the 1.0 a bit later.

--
Bartc

From: Chad on 8 May 2010 23:51
On May 8, 6:03 pm, "bart.c" <ba...(a)freeuk.com> wrote:
> Chad wrote:
> > On May 8, 4:25 pm, Kai-Uwe Bux <jkherci...(a)gmx.net> wrote:
> >> Chad wrote:
> >>> On May 8, 1:28 pm, Kai-Uwe Bux <jkherci...(a)gmx.net> wrote:
> >>>> Chad wrote:
> >>>>> for (i = 1; i <= n; i++) {
> >>>>> acc += 1.0/i;
> >>>>> }
> >>>> Hint: summing up the terms in reverse order (starting with the
> >>>> smallest) handles rounding errors much better and is numerically
> >>>> preferable.
>
> >>> Why would summing up the terms in reverse order handle rounding
> >>> errors much better?
> > Maybe I'm missing the broader concept here, I tried this in python. I
> > get the same results.
>
> > Here is what I get when I go one way...
> >>>> acc = 1.0
> >>>> for x in range(1,50):
> > ...   acc = acc + 1.0/x
> > ...   print acc
> > 5.47920533833
> > And here is what I get when I reverse it...
> >>>> acc = 1.0
> >>>> for x in reversed(range(1,50)):
> > ...   acc = acc + 1.0/x
> > ...   print acc
> > 5.47920533833
>
> I think you need to print more decimal places. Where you have ...833...
> above, I get differences after ...832942...
>
> Also you're starting with acc=1.0, which itself will hide some accuracy (so
> when adding 0.02, you're losing another couple of decimals). Perhaps start
> with acc=0.0 and add the 1.0 a bit later.
>

I saw the difference when I took the output to 40 decimal places..

[cdalten(a)localhost oakland]$ python
Python 2.6.2 (r262:71600, May 3 2009, 17:04:44)
[GCC 4.1.1 20061011 (Red Hat 4.1.1-30)] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> acc = 1.0
>>> for x in range(1,50):
.... acc = acc + 1.0/x
....
>>> print "(%.40f)" % acc
(5.4792053383294225810118405206594616174698)
>>>
[cdalten(a)localhost oakland]$ python
Python 2.6.2 (r262:71600, May 3 2009, 17:04:44)
[GCC 4.1.1 20061011 (Red Hat 4.1.1-30)] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> acc = 1.0
>>> for x in reversed(range(1,50)):
.... acc = acc + 1.0/x
....
>>> print "(%.40f)" % acc
(5.4792053383294252455470996210351586341858)
>>>

I'm still not sure how I would start off with acc = 0.0 and then add
1.0 later on. Who knows, maybe I will have it figured out before
someone tells me how to do it.
First  |  Prev  |  Next  |  Last
Pages: 1 2 3
Prev: Call for participants
Next: (paypal payment) wholesale Lacoste shirts from china,and free shipping