对于一个非线性方程f(x)=0求改方程的根,我们的思路可以这么想:
1.根的存在性。若该方程没有根,何必徒劳想法设法去求它的解呢?对于一个方程,我们怎么去找他的根,有连续函数零点定理可知:若有f(a)f(b)<0,则在(a, b)区间有解,究竟是一个还是多少个,还是要看具体的方程。
2.根的分布。这个方程的根分布在哪个区间,我们在程序实现的时候就可以一一搜索,用什么方法呢?我们可以采用一个不怎么高效的方法,等步长扫描法,在扫描的期间如果出现f(x1)(fy1)<0,则说明(x1, y1)区间有根。
等步长扫描法:
设定h>0为给定的步长,基础区间为(a, b),取x0=a,x1=x0+h,若f(x0)(x1)<0,则扫描成功,有根区间锁定在(x0, x1),否则,有x0=x1, x1=x0+h,然后再进行扫描,直到x1>b为止,但这不等于该方程没有根,因为你的步长如果很大,误差就大,很容易错过了有根的区间,所以当然建议采用尽量小的步长扫描。
#include <iostream> #include <list> using namespace std; /* Value类: 用来存储一个区间的左边界值和右边界值 */ class Value { private: double leftBound; double rightBound; public: double getLeftBound() { return leftBound; } void setLeftBound(double leftBound) { this->leftBound = leftBound; } double getRightBound() { return rightBound; } void setRightBound(double rightBound) { this->rightBound = rightBound; } }; /* Array类: 利用list类库作为模式定义一个存储容器来存储结果 */ typedef list<Value> Array; /* f(x)=0这个函数的逻辑实现 */ double f(double x) { return x*x-5*x+6; } /* 等步长扫描法实现: SameStepScan方法: 参数: left:左边界 right:右边界 stepLength:步长值 array: 结果集 */ void SameStepScan(double left, double right, double stepLength, Array *array) { double x0 = left; double x1 = x0 + stepLength; while(x1 <= right) { if(f(x0)*f(x1)<0) { Value value; value.setLeftBound(x0); value.setRightBound (x1); (*array).push_back(value); } x0 = x1; x1 = x0 + stepLength; } } /* main方法测试 */ int main() { Array *array = new Array(); SameStepScan(1, 10, 0.3, array); for(Array::iterator it = (*array).begin(); it != (*array).end(); it++) { cout<<"("<<(*it).getLeftBound()<<", "<<(*it).getRightBound()<<")"<<endl; } return 0; }
3.根的精确化。我们可以通过上面的方法得到一个区间,相当于已知了一个根的近似值,最后我们最迫切的就是让这个近似值尽量可能的靠近真值。我们通过等步长扫描法找到区间,然后通过二分法在区间中找到近似值。
二分法:二分法的原理也是基于连续函数的零点定理,设定f(x)=0在(a, b)区间有唯一的实根,令a1=a, b1=b, x=(a1+b1)/2, 如果f(a)f(x)<0,则新的区间为[a1, x],否则为[x, b1],然后再重复上面的步骤,最后a1-b1足够小的时候,x就是为近似值了。
/* 二分法实现: BinaryDivide: left:左边界 right:右边界 accuracy:精度控制量 value:近似值 */ void BinaryDivide(double left, double right, double accuracy, double &value) { double a1 = left; double b1 = right; double x = (a1+b1)/2; while(b1-a1 > accuracy) { if(f(a1)*f(x) < 0) { b1 = x; x = (a1+b1)/2; } else { a1 = x; x = (a1+b1)/2; } } value = x; } /* main方法测试 */ int main() { double value; Array *array = new Array(); SameStepScan(1, 10, 0.3, array); for(Array::iterator it = (*array).begin(); it != (*array).end(); it++) { cout<<"("<<(*it).getLeftBound()<<", "<<(*it).getRightBound()<<")"<<endl; BinaryDivide((*it).getLeftBound(), (*it).getRightBound(), 0.00001, value); cout<<"value : "<<value<<endl; } return 0; }