方程的求解

一、非线性方程的求解

1、牛顿迭代法。牛顿迭代法的思想是根据函数f(x)的泰勒基数的前几项来寻找方程的根。牛顿迭代法具有平方收敛的速度。

f(x)=f(x₀)+f^‘(x₀)(x-x₀)+...(泰勒公式展开式）

根据泰勒公式展开式的前两项有：

f(x)=f(x₀)+f^‘(x₀)(x-x₀)=0

如果f^‘(x₀)不等于0则：

x₁=x₀-(f(x₀)/f^‘(x₀))

这样通过x₀得到了x₁,然后根据x₁进行下一次相似的计算，这样就得到牛顿法的一个迭代序列公式：

x_n+1=x_n-(f(x_n)/f^‘(x_n))

具体代码如下：

 1 /*
 2  * @author
 3  * 非线性方程的表达式
 4  */
 5 double function(double x) {
 6     return x * x * x * x - 3 * x * x * x + 1.5 * x * x - 4;
 7 }
 8
 9 /*
10  * @author
11  * 非线性方程的导函数
12  */
13 double dFunction(double x) {
14     return 4 * x * x * x - 9 * x * x + 3 * x;
15 }
16
17 /*
18  * @author
19  * 牛顿迭代算法
20  * precision为算法的精度，maxGen为最大迭代次数,x为一个地址变量,
21  * 算法刚开始赋值给x是一个初始值，并不是最终的结果，故利用地址变量，方便
22  * 将最后的结果赋值给x
23  */
24 int newtonMethod(double *x, double precision, int maxGen) {
25     double temp, t;
26     int i = 1;
27     temp = *x;
28     if (dFunction(temp) == 0.0) {
29         printf("Error! The derivative is zero\n");
30         return 0;
31     }
32     while (i <= maxGen) {
33         t = temp - (function(temp) / dFunction(temp));
34         if ((fabs(t - temp) < precision)
35                 || (fabs(function(temp)) < precision)) {
36             *x = t;
37             return 1;
38         } else {
39             temp = t;
40         }
41         i++;
42     }
43     printf("the number of iteration over maxGen!\n");
44     return 0;
45 }

2、二分法

二分法的思想：如果函数f(x)，在x=c时，f(c)=0,那么把x=c称为函数f(x)的零点。对于二分法，假定非线性方程f(x)在区间（x,y)上连续，如果存在两个实数a,b属于区间（x,y),使得满足如下式子：

f(a)*f(b)<0

也就是说f(a)，f(b)异号，这说明在区间（a,b)内一定有零点，可以根据f((a+b)/2)来判断方程解的位置，具体代码如下：

 1 double function(double x) {
 2     return 2 * x * x * x - 5 * x - 1;
 3 }
 4 /*
 5  * @author
 6  * 这里将negative 表示函数值小于零的自变量，positive表示函数值大于0的自变量
 7  */
 8 int dichotomy(double *x, double negative, double positive, double precision) {
 9     double mid, temp;
10     if (function(negative) > 0)
11         positive = negative;
12     if (function(positive) < 0)
13         negative = positive;
14     mid = (positive + negative) / 2;
15     temp = function(mid);
16     if (temp > 0) {
17         positive = mid;
18     } else {
19         negative = mid;
20     }
21     if (fabs(negative - positive) > precision) {
22         dichotomy(x, negative, positive, precision);
23     } else {
24         *x = temp;
25         printf("The result of x is %0.5f\n", mid);
26     }
27     return 0;
28 }
29 /*
30  * 另一种二分法的表示方法
31  */
32 double erfen(double a, double b, double precision) {
33     double c;
34     c = (a + b) / 2.0;
35     while ((fabs(function(c)) > precision) && (function(a - b) > precision)) {
36         if (function(c) * function(b) < 0) {
37             a = c;
38         }
39         if (function(c) * function(a) < 0) {
40             b = c;
41         }
42         c = (a + b) / 2.0;
43     }
44     return c;
45 }