Floating point error: Overflow (abnormal program termination)

[easynhappy]

CAN SOMEONE BE KIND ENOUGH TO TELL ME WHY I HAVE A

FLOATING POINT ERROR: OVERFLOW
 ABNORMAL PROGRAM TERMINATION

IN THE CODE BELOW?

#include<stdio.h>
#include <math.h>
#include<conio.h>
#define A 0
#define B 1
float F(float x0, float y0)
{
return (0.2*x0)+(y0*y0);
       }
float Rung4(float x0, float y0, float h)
    {
	float k1 = F(x0,y0);
	float k2 = F(x0+h/2,y0+k1/2);
	float k3 = F(x0+h/2,y0+k2/2);
	float k4 = F(x0+h/2,y0+k2/2);
	float y1 = y0 + ( k1 + 2*k2 + 2*k3 + k4)/6;
    /*	printf("\n\n  k1 = %.4lf  ",k1);
	printf("\n\n  k2 = %.4lf ",k2);
	printf("\n\n  k3 = %.4lf ",k3);
	printf("\n\n  k4 = %.4lf ",k4);
	printf("\n\n  y(%.4lf) = %.3lf ",x0+h,y1);*/

 return y1;
    }

float Maks(float *y, float *y2, float n)
{
   int i;
   float Maks=fabs(y[0]-y2[0]);
	for(i=1; i<n; i++)
	 if (fabs(y[i]-y2[2*i])>Maks)
      Maks= fabs(y[i]- y2[2*i]);
return Maks;
}

void Rk4()
{
  int n=10, i;
  float h=(float)(B-A)/n, *x, *y, *x1, *y2, esp=0.001;
    do
      {
	 x=new float [n+1];
	 y=new float [n+1];
	 x1=new float [n*2+1];
	 y2=new float [n*2+1];
	      for (i=0; i<n; i++)
	       {
		   y[0]=0.1;
		   x[n]=B;
		   x[i]=A+i*h;
		   y[i+1]=Rung4(x[i], y[i], h);
		   }
	 n=2*n;
	 h=(float)(B-A)/n;
	 for (i=0; i<n; i++)
	  {
	     y2[0]=0.1;
	     x1[n]=B;
	     x1[i]=A+i*h;
	     y2[i+1]=Rung4(x1[i], y2[i], h);
	     }
      }
  while ((Maks(y, y2, n/2)/15)>=esp);
  puts("Runge- Kutta 4th Order");
    for (i=0; i<n; i++)
printf("x[%d]=%.2f       y[%d]=%6f\n", i, x[i], i, y[i]);

}
int main ()
{
// clrscr ();

  Rk4();
  getch ();
  return 0;
}

[VictorN]

Well, reading about Floating point would be a good start!

[easynhappy]

Well, reading about Floating point would be a good start!

I already did that. But am still kind of lost.

[Paul McKenzie]

CAN SOMEONE BE KIND ENOUGH TO TELL ME WHY I HAVE A

FLOATING POINT ERROR: OVERFLOW
 ABNORMAL PROGRAM TERMINATION

IN THE CODE BELOW?

Also, why are you not debugging your own code?

1) Your program leaks memory while in that loop. You are calling new[] without any call to delete[], and to make it worse, you're doing this while that do loop is running. Each time you loop, your leak gets bigger and bigger.

Either deallocate that memory, or use a container instead of new/delete, such as std::vector or std::list.

2)

for (i=0; i<n; i++)

The i is an integer, the n is a computed floating point value. At the start, n is 0.1. So if you loop only once, this loop at the end of the do-loop doesn't execute. If n is greater than 1, then there is no guarantee how many times you will loop, since floating point calculations are not exact. What if n is 9.99999 when it should be 10? So you loop only 9 times instead of 10.

Floating point calculations are not exact. You should not be using floating point types when specifying the for() loop limits. Always use integer, so that the loop runs consistently.

These are the glaring errors, all without running your program. Given all of this, it isn't surprising that you may be dividing by zero or a number close to 0 somewhere.

But again, you are supposed to debug your own code by using your compiler's debugger. You should be single-stepping through your program and seeing what the values are. You wrote the program, therefore you should be debugging your own code.

Edit:

I made a mistake in reading your code. I see that n is an integer, so your use of it in the loop and in computing the number of elements seems to be consistent. However, the memory leak still exists.

Regards,

Paul McKenzie

[Paul McKenzie]

To add to point 2 I mentioned:

What you should be doing is always use integers, and only use floating point for calculations. If this means dividing the integer by some value before you use it for the calculation, then so be it.

For example:

float n = 1.0;
for (float i = 0; i < n; i += 0.1)
{ 
   float someComputedValue = i / 20.0F;  // whatever
}

How many times will that loop execute? You think it will be 10 times? It sure looks like it, but there is no guarantee it will be 10 times, since adding 0.1 to a float is not exact.

The proper way to write this loop is:

int n = 10;
for (int i = 0 ; i < n; ++i)
{
    float myValue = (float)i / 10;
    float someComputedValue  = myValue / 20.0F;  // whatever
}

The loop will always execute 10 times. Note that I delay any usage of float until I need to make the calculation by dividing i by 10.

The bottom line is floating point computation using a binary machine such as a computer is not like high-school algebra, where you have a formula and on paper, things work.

I know that you're using the Runge-Kutta algorithm, but it is only an algorithm -- it's your job to convert that algorithm into something that a binary computing machine will handle correctly, given all the flaws that floating point calculations will give you. Courses such as Numerical Analysis teaches how to write such programs.

The other thing with your code is that you make no check to see if you're accessing elements beyond the boundaries of the array. If you're accessing invalid elements, there is a chance that the bogus value from the array is NAN, a number close to 0 but not 0, or some other "bad" floating point value.

So your program has at least one flaw and two possible flaws:

1) Memory leak (this is definitely a flaw)
2) Accessing element beyond the boundaries of an array
3) Possible incorrect usage of floating point values (i.e. assuming that they are exact computations).

Regards,

Paul McKenzie

[TheGreatCthulhu]

The proper way to write this loop is:

int n = 10;
for (int i = 0 ; i < n; ++i)
{
    float myValue = (float)i / 10;
    float someComputedValue  = myValue / 20.0F;  // whatever
}

To add to that - if performance is important, you can also avoid division operations on each iteration, by replacing them with less costly multiplication. All you need to do is to precalculate some values, like this:

float inv10 = 1.0F / 10.0F;
float inv20 = 1.0F / 20.0F;

int n = 10;
for (int i = 0 ; i < n; ++i)
{
    float myValue = i * inv10;
    float someComputedValue  = myValue * inv20;  // whatever
}

[/QUOTE]

Or even doing it all in one step:

float invDenominator = 1.0F / 200.0F;        // 1/(10*20)

int n = 10;
for (int i = 0 ; i < n; ++i)
{
    float someComputedValue  = i * invDenominator;  // whatever
}

[/QUOTE]

[easynhappy]

Am getting even more lost.

[Paul McKenzie]

Am getting even more lost.

Did you write this program? If so, then how could you even begin to write it if you couldn't debug it if something goes wrong? I don't get it.

When you wrote the program, you had something in mind, else it would not exist. So now that it is written, where does the program break down? Or is the problem you don't understand what you're doing with this algorithm at all?

I'll restate my first post to you:

But again, you are supposed to debug your own code by using your compiler's debugger. You should be single-stepping through your program and seeing what the values are. You wrote the program, therefore you should be debugging your own code.

Regards,

Paul McKenzie