二分根查找分段错误
bisection root finding segmentation fault
我使用以下代码通过简单的二分算法找到函数的根
#include <stdio.h>
#include <stdlib.h>
typedef float (*continous_function)(float);
static const float epsilon = 0.0000001;
#if __STDC__
static __inline float fabsf(float x)
{
return x >=0 ? x : -1*x;
}
#else
#include <math.h>
#endif
float poly_1(float x)
{
return x*x*x+2*x*x+3;
}
float poly_2(float x)
{
return 2*x + 3;
}
float poly_3(float x)
{
return 5*x*x*x*x*x+3*x*x*x+2*x+7;
}
static __inline void swap_in_place(float *a, float *b)
{
float tmp = *a;
*a = *b;
*b = tmp;
}
float bisection_root(continous_function f, float a, float b)
{
float neg = f(a);
float pos = f(b);
float c = 0.5 * (a + b);
float mid;
if (neg * pos > 0) {
/* neg and pos should have different sizes */
abort();
}
if (neg > 0) {
/* Ensure f(a)=neg is negative */
swap_in_place(&neg, &pos);
swap_in_place(&a, &b);
}
mid = f(c);
if (fabsf(mid) < epsilon) {
return c;
} else {
if (mid > 0)
return bisection_root(f, a, c);
else
return bisection_root(f, c, b);
}
}
int main()
{
{
float a = -5;
float b = 2;
printf("Root of x^3+2*x^2+3 is %f\n", bisection_root(poly_1,a,b));
}
{
float a = -2;
float b = 0;
printf("Root of 2*x+3 is %f\n", bisection_root(poly_2,a,b));
}
{
float a = -1;
float b = 0;
printf("Root of 5*x^5+3*x^3+2*x+7 is %f\n", bisection_root(poly_3,a,b));
}
return 0;
}
此程序在 windows xp(32 位)上使用 mingw gcc 编译时会导致段错误。
当epsilon
的小数位数减少时,可以避免分段错误。因此,我断定分段错误与溢出有关。
我想知道此错误发生的原因和发生方式,以便我可以找到一种可靠的方法来设置 epsilon 或修复可能导致该问题的其他错误。
static const float epsilon = 0.0000001;
多项式的值与 0 的距离为 epsilon
时,不需要浮点数。特别是如果根是 x
的较大值,则多项式可能从小于 -epsilon
跳到大于 +epsilon
在一个浮点数和它的后继者之间。
在这种情况下,以及在许多其他情况下,您的算法将循环,因为您使用递归实现它并且 C 编译器通常不保证尾调用优化,这种无限循环可能导致分段错误(当堆栈变满)。
无论使用递归还是简单的 while
循环,解决方案都是限制迭代次数。
注意:您应该阅读有关计算多项式的 Horner's scheme。
我使用以下代码通过简单的二分算法找到函数的根
#include <stdio.h>
#include <stdlib.h>
typedef float (*continous_function)(float);
static const float epsilon = 0.0000001;
#if __STDC__
static __inline float fabsf(float x)
{
return x >=0 ? x : -1*x;
}
#else
#include <math.h>
#endif
float poly_1(float x)
{
return x*x*x+2*x*x+3;
}
float poly_2(float x)
{
return 2*x + 3;
}
float poly_3(float x)
{
return 5*x*x*x*x*x+3*x*x*x+2*x+7;
}
static __inline void swap_in_place(float *a, float *b)
{
float tmp = *a;
*a = *b;
*b = tmp;
}
float bisection_root(continous_function f, float a, float b)
{
float neg = f(a);
float pos = f(b);
float c = 0.5 * (a + b);
float mid;
if (neg * pos > 0) {
/* neg and pos should have different sizes */
abort();
}
if (neg > 0) {
/* Ensure f(a)=neg is negative */
swap_in_place(&neg, &pos);
swap_in_place(&a, &b);
}
mid = f(c);
if (fabsf(mid) < epsilon) {
return c;
} else {
if (mid > 0)
return bisection_root(f, a, c);
else
return bisection_root(f, c, b);
}
}
int main()
{
{
float a = -5;
float b = 2;
printf("Root of x^3+2*x^2+3 is %f\n", bisection_root(poly_1,a,b));
}
{
float a = -2;
float b = 0;
printf("Root of 2*x+3 is %f\n", bisection_root(poly_2,a,b));
}
{
float a = -1;
float b = 0;
printf("Root of 5*x^5+3*x^3+2*x+7 is %f\n", bisection_root(poly_3,a,b));
}
return 0;
}
此程序在 windows xp(32 位)上使用 mingw gcc 编译时会导致段错误。
当epsilon
的小数位数减少时,可以避免分段错误。因此,我断定分段错误与溢出有关。
我想知道此错误发生的原因和发生方式,以便我可以找到一种可靠的方法来设置 epsilon 或修复可能导致该问题的其他错误。
static const float epsilon = 0.0000001;
多项式的值与 0 的距离为 epsilon
时,不需要浮点数。特别是如果根是 x
的较大值,则多项式可能从小于 -epsilon
跳到大于 +epsilon
在一个浮点数和它的后继者之间。
在这种情况下,以及在许多其他情况下,您的算法将循环,因为您使用递归实现它并且 C 编译器通常不保证尾调用优化,这种无限循环可能导致分段错误(当堆栈变满)。
无论使用递归还是简单的 while
循环,解决方案都是限制迭代次数。
注意:您应该阅读有关计算多项式的 Horner's scheme。