数值分析之牛顿法多项式求根

来源：互联网发布：现言推荐知乎编辑：程序博客网时间：2024/05/21 04:16

很朴素的方法，如果在区间[a,b]内有根，那么f(a)*f(b)<0

不难得到以下代码：

#include <iostream>#include <memory>using namespace std;double f(int m, double c [], double x){int i;double p = c[m];for (i = m; i > 0; i--)p = p*x + c[i - 1];return p;}int newton(double x0, double *r,double c [], double cp [], int n,double a, double b, double eps){int MAX_ITERATION = 1000;int i = 1;double x1, x2, fp, eps2 = eps / 10.0;x1 = x0;while (i < MAX_ITERATION) {x2 = f(n, c, x1);fp = f(n - 1, cp, x1);if ((fabs(fp)<0.000000001) && (fabs(x2)>1.0))return 0;x2 = x1 - x2 / fp;if (fabs(x1 - x2) < eps2) {if (x2<a || x2>b)return 0;*r = x2;return 1;}x1 = x2;i++;}return 0;}double Polynomial_Root(double c [], int n, double a, double b, double eps){double *cp;int i;double root;cp = (double *) calloc(n, sizeof(double) );for (i = n - 1; i >= 0; i--) {cp[i] = (i + 1)*c[i + 1];}if (a > b) {root = a; a = b; b = root;}if ((!newton(a, &root, c, cp, n, a, b, eps)) &&(!newton(b, &root, c, cp, n, a, b, eps)))newton((a + b)*0.5, &root, c, cp, n, a, b, eps);free(cp);if (fabs(root) < eps)return fabs(root);elsereturn root;}int main(){double c[2] = { 1,2 };int n = 1, a = -5, b = 5;double eps = 1e-8;double root = Polynomial_Root(c, n, a, b, eps);cout << root << endl;}