【原】Python实现所有算法-雅可比方法(Jacobian)

云深无际 2022-07-09 发布于内蒙古

展开全文

Python实现所有算法-二分法

Python实现所有算法-力系统是否静态平衡

Python实现所有算法-力系统是否静态平衡(补篇)

Python实现所有算法-高斯消除法

Python实现所有算法-牛顿-拉夫逊（拉弗森）方法

断断续续的写了五篇了，夸我！

可恶，优不见了

有一说一，矩阵的数值算法不是那么简单的写，我这里会推荐一些学习的资源假如你愿意学的话。

首先是最出名（国内的嗷）的一个教材

其中一个封面是这个

还有一个是这样的

当然了，还有一个黄皮的数值计算的书，和数学建模的封面是一样的。其实我收集了很多相关的图书，只是硬盘不在手边，有机会可以分享一下。

这里推荐一个课程，中科大的，当然没有什么视频资料，不过PPT做的很好。

科大的数值计算

PPT一观

课程的特点

课程研究的东西

关于代码也推荐一本C数值代码，就平平无奇的名字

对于视频的话，这里推荐苏州大学：

妈妈，我在B大学上课

反正万国造

你以后可以这样说，看的清华的教程，翻的中科大的PPT，研究的苏州大学的视频，抄的云深无际的公众号的代码。

不管怎么说，学会就行，一个确定的知识和结论不会随着载体的不一样而便便，因为你学的是本质。

我相信很多人的学校都寡，但是又没有剥夺你学习的权力，所以你不学才是蛋疼的地方。

开始正文了：

对于矩阵求解，我们大体分为，A稠密不稠密：

那么就会演化出来两个解决的办法

感谢CSDN一位作者的总结，绘制了一观漂亮的思维导图

雅克比迭代法的优点明显，计算公式简单，每迭代一次只需计算一次矩阵和向量的乘法，且计算过程中原始矩阵 A （A在上面说明了）始终不变，比较容易并行计算。然而这种迭代方式收敛速度较慢，而且占据的存储空间较大，所以工程中一般不直接用雅克比迭代法，而用其改进方法。

这个迭代法又称为辗转法，是用计算机解决问题的一种基本方法，为一种不断用变量的旧值递推新值的过程，与直接法相对应，一次性解决问题。迭代法分为精确迭代和近似迭代，“二分法”和“牛顿迭代法”属于近似迭代法。迭代法利用计算机运算速度快、适合做重复性操作的特点，让计算机对一组指令（或一定步骤）进行重复执行，在每次执行这组指令（或这些步骤）时，都从变量的原值推出它的一个新值。

再说矩阵的求解：

考虑线性方程组Ax = b时，一般当A为低阶稠密矩阵时，用主元消去法解此方程组是有效方法。但是，对于由工程技术中产生的大型稀疏矩阵方程组（A的阶数很高，但零元素较多，例如求某些偏微分方程数值解所产生的线性方程组），利用迭代法求解此方程组就是合适的，在计算机内存和运算两方面，迭代法通常都可利用A中有大量零元素的特点。雅克比迭代法就是众多迭代法中比较早且较简单的一种，其命名也是为纪念普鲁士著名数学家雅可比。

这里总结一句话，雅可比算法是：用于确定严格对角占优的线性方程组的解。

在数学中，如果对于矩阵的每一行，一行中对角线条目的大小大于或等于所有其他（非对角线）的大小之和，则称方阵为对角占优该行中的条目。更准确地说，矩阵A是对角占优的，如果：

定义给出来了

多说无疑，你可以参考这个学习对角占优矩阵

所以这里的A是指非奇异的大规模稀疏矩阵。

什么是稀疏矩阵？？？毕竟一开始就写了。

概念：在实际问题中，特别是微分方程数值解法中，出现的线性代数方程组的系数矩阵往往系数很高，但其非零元素所占的比例很小，我们常把这类矩阵成为大型稀疏矩阵。

理解：零元素很多的多阶矩阵。

注意：求解此类系数矩阵若使用Gauss消元法常常会破坏矩阵稀疏性，另分解过程中出现大量非零元素。

再插一个：

什么是非奇异阵呢？非奇异矩阵是行列式不为 0 的矩阵，也就是可逆矩阵。意思是n 阶方阵 A 是非奇异方阵的充要条件是 A 为可逆矩阵，也即A的行列式不为零。即矩阵（方阵）A可逆与矩阵A非奇异是等价的概念。

理论的东西先说上面那么多，都是概念，说计算的时候的样子。

首先将方程组中的系数矩阵A分解成三部分，即：A = L+D+U，如下图所示，其中D为对角阵，L为下三角矩阵，U为上三角矩阵。

之后确定迭代格式，X^(k+1) = B*X^(k) +f ，（这里^表示的是上标，括号内数字即迭代次数），如下图所示，其中B称为迭代矩阵，雅克比迭代法中一般记为J。（k = 0,1，......）

再选取初始迭代向量X^(0)，开始逐次迭代。

理论是这样的

也就是说，对角线分量D，下三角部分L，和上三角U

伪算法是这样的

其实这个算法除了矩阵符合要求以外，最重要的迭代格式的转换。

我们的线性方程组是这样的

把A可以分解成三个矩阵的加和

最终的迭代格式是这样的

对具体元素的迭代公式是这样

在代码实现的部分，频繁的会看到__future__这个库，SO？是什么。

是这样的

关于annotation其实是一个PEP的提案。

这也就是它的作用

我们对函数的设计是需要两个多维的数组，对应在我们的线性方程组，迭代的初始值也要设计，不然从哪里开始迭代呢？以及想迭代的次数。

对于此函数的一个输入来说是这样的

coefficient = np.array([[4, 1, 1], [1, 5, 2], [1, 2, 4]])constant = np.array([[2], [-6], [-4]])init_val = [0.5, -0.5, -0.5]iterations = 3
jacobi_iteration_method(coefficient, constant, init_val, iterations)[0.909375, -1.14375, -0.7484375]