C语言利用函数指针进行非线性方程求根

下面给出了两种方法求解非线性方程零点的程序，即求解：
$$f(x)=0$$形式的方程的根$x^{\ast }$。

黄金分割法

计算原理为在区间范围$[a,b]$内试探，每次并非去中点而是取黄金分割点。也是一种区间优化方法，只需要提供真实解$x^{\ast }$区间范围。
迭代公式为$$m_k=a_k+0.618\left(b_k-a_k\right)$$若$f(a_k)$与$f(m_k)$同号，则$a_{k+1}=(赋值)m_k$；否则$b_{k+1}=(赋值)m_k$. 直至$b_k-a_k<ϵ$.

牛顿迭代法

计算原理为切线近似非线性函数，然后求根。是一种局部优化和求根方法，需要接近真实解$x^{\ast }$的初始猜测。
迭代公式为$$x_{k+1}=x_k-\frac{f(x_k)}{f^{\prime }(x_k)}$$

需要提供函数的导数。终止条件为达到迭代上限$k>k_{max}$或精度足够高：$x_{k+1}-x_k<ϵ$。

函数指针

函数指针用于指向一个函数，函数名是函数体的入口地址。下面是几个例子：

// 最简单的函数及其对应的函数指针：
void f();
void (*f_ptr)();
// 复杂点的函数和他的函数指针
int f(double b, int i);
int (*f_ptr)(double b, int i);

为了避免指针的使用，我的习惯是用typedef重新给它改个名字：
typedef int(*FUNC_P)(int, int)。以下都是FUNC_P的实例:

int func(int a, int b) {
    cout << "1997年写的func" << endl;
    return 0;
}
int func2(int a, int b) {
    cout << "1999年写的func2" << endl;
}
//2018年想添加一个新的子业务
int func2021(int a, int b) {
    cout << "2021年写的func2021" << endl;
}

显然，由于上面几个函数的输出、形参表、参数类型都完全一样，所以可以用同一个函数指针类型。使用时用以下声明和定义

    FUNC_P f1 = NULL, f2 = NULL;// 声明
    f1 = func;		// 赋值
    f2 = func2021;
    f1(10, 20);    // 等价 func(10,20);

此种定义类似MATLAB中的函数句柄@fun（handle），在多个源文件中调用时非常方便。

源文件

// functionHandleTest.cpp : 此文件包含 "main" 函数。程序执行将在此处开始并结束。
//

#include <iostream>
#include <cmath>

typedef double(*FUNC_P)(double);//定义函数句柄类型

double f1(double x) {
    return x * x - 2.0;// 求解根号2
}
double f1p(double x) {
    return 2.0* x ;
}

double f2(double x) {
    return cos(x)-x;// 求解cos(x)-x==0
}
double f2p(double x) {
    return -sin(x)-1;
}

static int iterCount = 0;

double goldenSearch(FUNC_P fun,double left, double right, double precision) {
    double mid = 0.618*right+0.382*left;
    iterCount++;
    printf("%d\t%1.16f\n", iterCount, mid);
    double result = fun(mid);
    if (right - left < precision) {
        return mid;
    }
    else if (result * fun(left) < 0) {
        return goldenSearch(fun , left, mid, precision);
    }
    else if (result * fun(right) < 0) {
        return goldenSearch(fun , mid, right, precision);
    }
    else {
        return goldenSearch(fun , mid, right, precision);
    }
    return 0;
}

double newtonIter(FUNC_P func, FUNC_P diffF,double x0, double MAX_ITER) {
    double preX;
    double x;

    preX = x0;
    x = preX - func(preX) / diffF(preX);
    double precision = 1.e-10;
    while (iterCount < MAX_ITER &&abs(preX-x)>precision) {
        preX = x;
        x = preX - func(preX) / diffF(preX);
        iterCount++;
        printf("%d\t%1.16f\n", iterCount,x);
    }
    return x;
}

int main()
{  
    FUNC_P f = NULL,fp=NULL;
    f = f2; fp = f2p;
    //f3 = fun2;
    std::cout << "Iteration \t Function Value" << std::endl;
    std::cout << goldenSearch(f,0,2,1.0e-10);

    iterCount = 0;
    std::cout << std::endl;

    std::cout << "Newton \t Function Value" << std::endl;
    std::cout << newtonIter(f,fp, 1.5,20);
}

// 运行程序: Ctrl + F5 或调试 >“开始执行(不调试)”菜单
// 调试程序: F5 或调试 >“开始调试”菜单

求解结果如图