二分法求一个非线性方程的多个解

一、问题引人
工作上遇到一个算法问题，需要从一张图片的中心点开始，生成一条螺旋线并沿着螺旋线进行特征提取。参考一些资料，我选择了一条比较常用的螺旋线方程，如下所示：

上式中x和y表示坐标点，t表示时间，可见x和y都是关于t的函数。当我对t按照步长为1来递增时，效果如下：

t以步长为1进行增长

在生成螺旋线时，发现相邻两个点之间的距离越来越大，而我需要寻找的是螺旋线上等间距分布的点集。

二、数学建模
为了解决上述的问题，要设计了一个算法，在螺旋线方程的基础上引入一个限定条件，求出t的“合理值”使得相邻两个点之间距离相等，算法思路如下：
步骤一：给定t一个初始值，比如t等于4，带入螺旋线方程求出p(x,y)；
步骤二：以p为圆心，r为半径生成一个小圆，求出圆与螺旋线之间的交点；
步骤三：选择距离上一个圆心最近的交点作为新的圆心，再次生成一个圆，圆的半径等于r，继续求圆与螺旋线之间的交点，以此类推，直到生成的点集个数达到一定阈值。

如图所示，红色点表示圆心，蓝色代表小圆，为了更好的示意，我将圆缩小了，实际上每个红色的点恰好分布在前一个圆的边上，所有红色的点集合就是要求的螺旋线上等间距分布的点。假设圆心的坐标用（a,b)表示，半径为r，每次求圆与螺旋线的交点都是求一个方程组的解。

乍一看觉得这个方程很好解，只要把x和y带入最后一个式子，只剩一个自变量，然而经过漫长的求解最终幡然醒悟，这个方程没有解析解（即无法用解析式可以表达的解），只能用数值逼近的方法求近似解。世界上有很多方程都没有解析解，但是解确实存在，这点确实也让人感到好奇和感叹，并非所有的问题都能依靠公式来表达，反而复杂的问题却可以用简单的二分法来找答案。

三、二分法求解方程
将x、y带入最后一个式子，化简得到关于t的方程：
f(t) = t * t - 2 * t * (a * cos(t) + b * sin(t)) + a * a + b * b - r * r;

导函数为
f’(t) = 2 * (t - (a + b * t) * cos(t) - (b - a * t) * sin(t));

观察函数图像和导函数图像，发现
函数的性质
1.存在多个根
2.随着x增大，y震荡递增且振幅越来越大
3.极值点均匀分布，周期为π

算法步骤
步骤一：从t = 0开始，以固定步长的形式搜索第一个极值点t1和第二个极值点t2（导函数为0作为判断条件）；
步骤二：T = t2 - t1；
步骤三：第三个极值点t3 = t2 + T，第四个极值点t4 = t3 + T,以此类推，求出前k个极值点，直到连续3个极值点大于都大于0为止；
步骤四：将前k个极值点位置分成k - 1个区间，[t1,t2)、[t2,t3) …
步骤五：对每个区间应用二分查找法来找根。

二分法参考链接：https://blog.csdn.net/u013232740/article/details/43342165/

四、代码

/*************************************************************//**
 * @file        Bisection.cpp
 * @brief       this class provides interfaces for bisection algorithm
 * @author      Qingquan Zhou (zhouqingquan118@163.com)
 * @date        10/20/2022 (M/D/Y)
 *
 ***************************************************************/

#include <iostream>
#include <vector>
#include <cmath>
#include "Bisection.hpp"


double Bisection::function(double x, const double a, const double b, const double r