Kahan summation....

BlitzMax Forums/BlitzMax Programming/Kahan summation....

BinaryBurst

(Posted 2015) [#1]

I'm trying to run a code implementing Kahan summation and.... it doesn't want to work. It's very frustrating because when try to manually compile the code with bmk from the command line, it works(sometimes).
How do I solve this iritating problem? :)

Graphics 800,800
h#=1
n#=1
While n<1e9
	x#=0    
	y#=1
	h=1.0/n
	For i#=1 To n
		dx#=h+cx#
		new_x#=x+dx
		cx=(x-new_x)+dx
		x=new_x
		
		dy#=h*(-y)+cy#
		new_y#=y+dy
		cy=(y-new_y)+dy 
		y=new_y
	Next
	Print "KAH "+x+" "+y+" "+Exp(-1)+" "+(y-Exp(-1))
	SetColor(0,255,0)
	Plot(3e1*Log(n),-3e1*Log((y-Exp(-1)))) ;Flip(0)
	
	x#=0
	y#=1
	For i#=1 To n 
		x=x+h
		y=y+h*(-y)
	Next
	Print "NOR "+x+" "+y+" "+Exp(-1)+" "+(y-Exp(-1))
	SetColor(255,0,0)
	Plot(3e1*Log(n),-3e1*Log((y-Exp(-1)))) ;Flip(0)
	
	n=n*1.1
Wend

WaitKey()

Floyd

(Posted 2015) [#2]

Imagine a floating point system which has three decimal digits of accuracy, and an exponent. 999 would be 9.99e2. If we add 1 we get 1000, which is 1.00e3. Note this is really only using one digit for the mantissa.

The snag comes when we try to add 1 to 1000. The result would ideally be 1.001e3. But that needs four digit accuracy and only three digits are available. Hence 1 cannot be added correctly to 1000. Attempting to add 1 to 1000 gives a result of 1000, no change. The smallest amount that could be added is 10 since 1000 + 10 = 1010 = 1.01e3 which can be represented with three digit mantissa.

Real floating point does essentially the same thing in binary. Single precision has 24 bits of accuracy. Adding 1 to 2^24 fails because the result needs 25 bits to be exact. You can see it here.

i# = 2^24 - 5

For n = 1 To 10
	Print i
	i = i + 1
Next

Your code attempts a loop with i# as the counter and an upper limit much bigger than 2^24. You need a data type which can correctly increment such a large number. Integer, Long and Double can all do it.

EDIT: I forgot to mention that sometimes BlitzMax will not properly enforce the fact that i# is only single precision. It has no choice when dealing with numbers in memory. But if the variable is held in a register it will really be in a higher precision. So generated code occasionally does calculations in higher precision than you specified, depending on how registers are allocated by the compiler. Usually this is of no consequence, but might affect you since you are trying to compare relatively tiny differences in results.