共轭梯度

View on Github 线性共轭梯度算法库（C++ Library of Linear Conjugate Gradient，LIBLCG） 张壹（zhangyiss@icloud.com） 浙江大学地球科学学院·地球物理研究所 简介 liblcg是一个简单的C++线性共轭梯度算法库，包含了一般形式的共轭梯度与预优共轭梯度算法。可用于求解如下形式的线性方程组： Ax = B （可用共轭梯度求解）或者 PAx = PB（可用预优共轭梯度求解） 其中，A是一个N阶的对称矩阵、x为N*1的待求解的模型向量，B为N*1需拟合的目标向量。P则为预优算法中的预优矩阵，一般是一个N阶的对角阵。共轭梯度法广泛应用于无约束的线性最优化问题，拥有优良收敛与计算效率。 安装 GCC（MacOS） mkdir build cd build cmake .. make make install Windows 请自行拷贝代码新建项目并编译…嘿嘿😁！ 使用说明 自定义Ax计算函数 通常我们在使用共轭梯度法求解线性方程组Ax=B时A的维度可能会很大，直接储存A将消耗大量的内存空间，因此一般并不直接计算并储存A而是在需要的时候计算Ax的乘积。因此用户在使用liblcg时需要定义Ax的计算函数，同时提供初始解x与共轭梯度的B项（即拟合的对象）。如果使用预优方法还需要提供预优矩阵P项。Ax计算函数的形式必须满足算法库定义的一般形式： typedef void (*lcg_axfunc_ptr)(void* instance, const lcg_float* x, lcg_float* prod_Ax, const int n_size); 函数接受4个参数，分别为： void *instance 传入的实例对象； const lcg_float *input_array Ax计算中的x数组的指针； lcg_float *output_array Ax的乘积； const int n_size 矩阵的大小。 自定义进程监控函数 用户可以以下面的模版创建函数来显示共轭梯度迭代中的参数，并可以在适当的情况下停止迭代的进程。 typedef int (*lcg_progress_ptr)(void* instance, const lcg_float* m, const lcg_float converge, const lcg_para* param, const int n_size, const int k); 函数接收6个参数，分别为： ...

非线性共轭梯度法-DY法 算法可用于求解形如$min(f(x))$的非线性最优化问题，解的迭代形式为$x_{k+1}=x_k+a_kd_k$。其中$a_k$为迭代步长，$d_k$为迭代方向。（不同的$d_k$计算方法即为不同的非线性共轭梯度算法，它们对于$a_k$计算的要求也各有区别。） 要求目标函数$f(x)$连续可微且梯度$g(x)$已知。 DY法具有全局收敛性，但有一定几率陷入局部极小值。 DY法具有二次终止性，迭代终止于零梯度值位置或者$f(x) \leq threshold$处。 算法描述 step 1 -> 给出解的初值$f(x_0)$并计算其梯度$g(x_0)$，如果$\left\lVert g(x_0) \right\rVert = 0$则算法终止，说明$x_0$是$f(x)$的全局或局部极小解。否则设置$d_0 = -g_0$作为初始迭代方向（即解的最速下降方向）。 step 2 -> 执行满足如下所示的标准Wolfe条件的非精确线性搜索，计算得到$a_k$。 $f(x_k+a_k d_k)-f(x_k) \leq c_1 a_k {g_k}^T d_k$， ${g(x_k+a_k d_k)}^T d_K > c_2 {g_k}^T d_k$ 其中$0< c_1 < c_2 < 1$，一般取$c_1=0.1$，$c_2 \in (0.6,0.8)$。对于重磁非线形反演一般取$c_1=1e-4$，$c_2=0.9$。 step 3 -> 更新$x_{k+1}=x_k+a_kd_k $，若$\left\lVert g(x_{k+1}) \right\rVert = 0 $则算法终止，说明$x_{k+1}$是$f(x)$的全局或局部极小解。否则继续step 4。 step 4 -> 更新$d_{k+1}=-g_{k+1}+\beta_k d_k$，其中$\beta_k = \left\lVert g_{k+1} \right\rVert / {d_k}^T (g_{k+1}-g_k)$。$k:=k+1$返回step 2。 C++代码 func.h ...