type
status
date
slug
summary
tags
category
icon
password
不动点迭代
定义:当 ,实数 是函数 的不动点
从不同视角来看,当输入和输出相等时,不动点 就是方程 的解
所有的方程求解问题 都可以转换为形如 的不动点迭代问题
收敛条件:,导数越接近0,收敛速度越快
重根问题
定义:函数 有根 ,且它在该点出现不止一次,则称 是重根
如果 有 重根,则满足
例如:
图像意义
- 单根:图像穿过 x 轴
- 偶数重根:图像在根处与 x 轴相切,但不穿过 x 轴
- 奇数重根:图像在根处变的平滑或s形,穿过 x 轴
重根会导致数值求解变难
牛顿法
牛顿法的线性收敛
这意味着每次迭代,下一步缩小为原来的一半、1/3等
重根会导致牛顿法收敛为线性,通过改进可修正为二次收敛
牛顿法的二次收敛
这意味着每次迭代,下一步的误差是 ,有效数字翻倍增长,远快于线性收敛
方程组直接求解
高斯消元法:将矩阵通过行变换变为好解的上三角形式,再从最后一行向上回代求解
LU分解:将一个方阵 分解成下三角矩阵 和上三角矩阵 的乘积
其中
- L:下三角矩阵,主对角线为1,表示消元的系数
- U: 上三角矩阵,表示高斯消元后的结果
LU分解将 转化为 ,则求解变为:
- 先解
- 再解
有些矩阵不能直接做LU分解,如主对角线元素为0时,可以引入行交换来避免
其中, 是置换矩阵,用于表示行交换
分解将 转化为 ,则求解变为:
- 先解
- 再解
误差
范数:一种衡量向量长度或大小的函数,对于矩阵衡量的是这个矩阵对输入向量的放大能力
常见范数类型
名称 | 向量 | 矩阵 |
1-范数(曼哈顿范数) | 所有分量绝对值之和 | 每列元素绝对值求和中的最大值 |
2-范数(欧几里范数) | 所有分量平方和开根 | 最大奇异值 |
∞-范数 | 所有分量最大绝对值 | 每行元素绝对值求和中的最大值 |
条件数:输入发生微小变化时,输出会被放大的倍数
换句话说,对于线性系统 求解, 是输入, 是输出, 的微小变化对应的误差放大因子即是 的条件数
其中, 是矩阵范数(如2-范数、∞-范数)
条件数越大,解的相对误差越可能被放大
- ,非常稳定,几乎无误差放大
- ,稳定,常见系统
- ,病态,对误差敏感
- ,极度病态,结果几乎不可用
- ,矩阵奇异,无解
方程组迭代求解
雅可比法:上一轮的所有旧值计算下一轮的新值
将矩阵 分解为
其中
- D:仅包含A主对角线元素的对角矩阵
- L:主对角线元素为0的严格下三角矩阵
- U:主对角线元素为0的严格上三角矩阵
则迭代公式为
高斯-赛德尔法:当前轮已更新的值和剩下的旧值计算下一轮的值
采用同样的分解,迭代公式为
连续过松弛法:在高斯-赛德尔基础上增加松弛参数 控制更新幅度
迭代公式为
- :变成了高斯-塞德尔法
- :欠松弛,稳定收敛但慢
- :过松弛,加速收敛
适用于以上三种方法的收敛条件为
- 严格对角占优矩阵:对角元素绝对值大于该行其他元素绝对值之和
共轭梯度法:专门为对称正定矩阵线性求解设计的高效迭代算法
条件1-对称
条件2-正定,对任意非零向量 有
共轭梯度法基于最速下降法改进,使每一步都尽可能朝正确方向前进
迭代步骤
- 根据初始猜测,计算残差
- 计算步长, 表示正确前进方向
- 更新残差
- 计算共轭方向更新系数
- 更新方向
重复以上步骤,直到满足收敛条件
使用预条件技术可以降低问题的条件数,加快共轭梯度法迭代的收敛速度
构造形式为
其中, 是预处理矩阵,也预条件子
一种特别简单的方式是雅可比预处理,即令
迭代步骤
- 计算步长
- 更新解向量
- 更新残差
- 解预处理方程
- 计算共轭方向更新系数
- 更新共轭方向
问题1:既然线性方程组可以用高斯消元法直接求解,为什么还要用迭代法?
比较维度 | 直接法(高斯消元) | 迭代法 |
计算复杂度 | O(n^3) | 每轮低,适合大规模系统 |
存储需求 | 高(特别是稠密矩阵) | 低(适合稀疏矩阵) |
精度 | 精确(理论上) | 可控误差 |
并行计算 | 不好并行 | 雅可比易并行,高斯-赛德尔差一些 |
适用规模 | 小到中等规模系统 | 大规模系统(上百万个未知数) |
鲁棒性 | 对病态矩阵敏感 | 有时更稳定(加松弛因子) |
实际应用 | 一般用于小系统、基准测试 | 工程计算/科学仿真常用 |
牛顿拉夫逊法:求解多元非线性方程组
雅可比矩阵
雅可比矩阵是一个多变量函数的偏导数组成的矩阵,描述局部变化率,常用于非线性系统的牛顿法求解
给定一个多元向量函数:
它的导数对应的雅可比矩阵表示如下:
牛顿法
函数 的根是其与 轴的交点,迭代求解根的步骤如下
- 给定一个初始值
- 画出 在 处的切线
- 切线与 轴的交点
- 画出 在 处的切线
- …迭代逼近精确解
点 的切线斜率方程是
令
因此单变量中,牛顿法迭代公式表示为
在多元中,用泰勒展开的一阶近似
泰勒展开公式
核心思想:用多项式函数逼近任意光滑函数
设 在 处 阶可导,则
实际使用中,通常选前 项作为近似,如牛顿法取前两项(切线)近似
误差由余项 表示
- 拉格朗日余项:
- 佩亚诺余项:
令
残差
残差(Residual)表示近似解 与真实解 的误差
对于函数 ,残差是 ,即
对于方程组 ,残差是一个向量
其中
- 是第 次迭代的近似解
- 是第 次迭代的残差向量
采用牛顿法时,收敛判据常用残差范数小于某个容差 (工程上常用 )
举例一:在管网中,直接采用经典牛顿法,根据能量守恒和质量守恒,构建出非线性方程组
构造雅可比矩阵
举例二:求解二维非线性方程组
- 定义
- 雅可比矩阵
- 每次迭代做
- 更新
⚠️ 注意事项
- 初值必须“足够接近”真正解,否则可能发散
- 雅可比矩阵可能稀疏、奇异,需要数值稳定的解法(LU 分解、SVD、迭代法)
- 可使用数值微分或自动微分近似雅可比
函数 的根是其与 轴的交点,迭代求解根的步骤如下
- 给定一个初始值
- 画出 在 处的切线
- 切线与 轴的交点
- 画出 在 处的切线
- …迭代逼近精确解
点 的切线斜率方程是
令
因此单变量中,牛顿法迭代公式表示为
在多元中,用泰勒展开的一阶近似
泰勒展开公式
核心思想:用多项式函数逼近任意光滑函数
设 在 处 阶可导,则
实际使用中,通常选前 项作为近似,如牛顿法取前两项(切线)近似
误差由余项 表示
- 拉格朗日余项:
- 佩亚诺余项:
令
残差
残差(Residual)表示近似解 与真实解 的误差
对于函数 ,残差是 ,即
对于方程组 ,残差是一个向量
其中
- 是第 次迭代的近似解
- 是第 次迭代的残差向量
采用牛顿法时,收敛判据常用残差范数小于某个容差 (工程上常用 )
雅可比矩阵
雅可比矩阵是一个多变量函数的偏导数组成的矩阵,描述局部变化率,常用于非线性系统的牛顿法求解
给定一个多元向量函数:
它的导数对应的雅可比矩阵表示如下:
举例:求解二维非线性方程组
- 定义
- 雅可比矩阵
- 每次迭代做
- 更新
⚠️ 注意事项
- 初值必须“足够接近”真正解,否则可能发散
- 雅可比矩阵可能稀疏、奇异,需要数值稳定的解法(LU 分解、SVD、迭代法)
- 可使用数值微分或自动微分近似雅可比
插值
插值是指已知一组离散点 构造一个多项式
拉格朗日插值
构造一组基函数 使每个函数在对应点上为1,其他点上为0,精确复现
设有n+1个已知点,插值多项式表示为
其中
牛顿插值
基于差商构造插值多项式,具备递推性,适合程序实现与动态添加点
差商的定义
- 一阶差商
- 二阶差商
- 更高阶以此类推
插值多项式表示为
误差
插值多项式只能在插值点精确重建,在区间外和区间内的其他点是逼近原函数,误差由以下公式描述:
若原函数 是 次可导函数,插值多项式为 ,则对任意 ,存在某个 ,使得
这个公式表明
- 误差大小与 相关,原函数越平滑,误差越小
- 与插值点的分布有关,通过合理选择插值点可以减小误差
龙格(Runge)现象:当使用高阶多项式在等距点上进行插值时,插值多项式在区间端点附近会产生剧烈震荡和大误差,即使原函数非常光滑,而且插值点越多反而区间边界产生振荡与误差更大。
切比雪夫(Chebyshev)节点:选择最优的插值点,显著减少最大误差抑制边界振荡
对于区间 上的 个插值点,切比雪夫节点为
切比雪夫插值
切比雪夫多项式是一个正交多项式序列,记作 ,定义在区间 上:
前几个切比雪夫多项式为:
它们在数值上非常稳定,适合做高阶逼近。
三次样条插值
在每两个数据点之间分别构造一个三次多项式,通过一系列三次多项式来光滑连接数据点的插值方法
给定 个已知数据点 其中升序
我们希望在每个 子区间构造一个三次多项式
最终构成的这个分段函数满足以下条件
- 插值条件
- 一阶导数连续
- 二阶导数连续
边界条件:三次样条插值总共有4n个未知系数,还需补充2个额外的条件才能求解
常见边界条件类型有
- 自然边界
- 固定导数边界
- 周期性边界
最小二乘
给定一组数据点,但方程数量多于未知数,因而解不能精确满足所有方程,则寻找一个最优近似函数逼近它们,使得满足近似解和精确解间误差的平方和最小
目标由精确解转为衡量与精确解的接近程度,即误差最小
方程数多于未知数称为超定线性方程组 ,定义误差向量
我们用欧几里距离衡量,使问题变为求解
正规方程法
对x求导后令导数为0,得到正规方程
最小二乘解为
对于每个数据点 要找到一条直线尽量靠近所有点,直线模型为
表示为矩阵形式
问题就转为了最小二乘标准形式
拟合曲线同理,多项式模型
把 的每列扩展成 即可
周期数据采用傅里叶级数拟合,模型如下
其中
- ,T是数据的周期
- 是拟合系数,通过最小二乘法求解
- 是阶数,表示用了多少对正/余弦波来组成曲线
数据点数量 | 推荐阶数范围 |
少于 20 个 | 1 ~ 2 阶 |
20 ~ 100 个 | 2 ~ 5 阶 |
100 ~ 500 个 | 3 ~ 10 阶 |
很大数据量 | 可视化调参 + 验证法选择阶数 |
QR 分解法
把 分解成两个矩阵
其中
- 是正交矩阵,它的列是单位正交向量
- 是上三角矩阵
施密特(Schmidt)正交变换
施密特正交化,可以把n个线性无关的向量(如A的列)变成n个彼此垂直的单位向量
假设取 中的 作为基底,想让 变成垂直于 的新方向
- 把 归一化,即除以它的长度获取单位向量,得到
- 减去其在 方向上的投影,得到垂直方向的分量
- 把 归一化,得到
就是两个正交单位向量,即 矩阵的列
表示向量 在正交向量 上的投影
豪斯霍尔德(Householder)反射子变换
在 QR 分解中,Householder 反射常作为 Gram-Schmidt 正交化的替代方法,因为它在数值上更高效稳定
构造一个镜子,把任意向量 反射到跟 (第一坐标轴)对齐的方向

图中把向量 通过构造的镜子(虚线),反射为与其等长的向量 ,其与坐标轴重合
这个镜子就是 Householder 矩阵,记为
通过一系列反射,逐步消去每列下面的元素,形成一个上三角矩阵
高斯-牛顿法
当有一组数据点,想对一个非线性模型函数拟合,如
作为非线性最小二乘问题,目标是最小化残差平方和函数
模仿牛顿法的思路,最小化函数就是极值问题,得出迭代公式为
其中
- 是梯度,即向量导数
- 是海森矩阵,即二阶导数
对残差平方和函数求梯度
构造海森矩阵
代入牛顿迭代公式
示例:
设有一组实验数据,表示不同雷诺数
Re
下的摩擦阻力系数 λ
,怀疑他们之间满足如下非线性模型通过非线性最小二乘,拟合出参数 a,b,c
数值微分
在很多实际问题中,我们知道一个函数的值(比如测量、仿真结果),但无法获得其解析表达式。这时候,如果想知道函数的导数(比如速度、加速度、应力梯度等),我们需要通过数值方法来近似求导
有限差分
有限差分是用函数在离散点上的值,来近似导数的常用方法,是数值计算、科学计算和仿真中的基础工具
令 为步长, 为待求导函数,则
一阶导数近似的中心差分公式
二阶导数近似的中心差分公式
数值积分
数值积分的目的是在无法解析求积分时,用一个插值多项式去逼近函数,然后对这个多项式积分
梯形法
用直线逼近曲线,等价于将面积近似为梯形面积
复合梯形法将区间[a, b]分成n个小区间,每个区间用梯形法则,然后累加
其中,,
辛普森法
用二次多项式逼近函数,计算面积
复合辛普森法将区间分成偶数个小区间,对每对相邻子区间使用辛普森法,然后累加
常微分方程
常微分方程用于理解以及预测随时间变化的系统,是一种含有未知函数及其导数的方程,只有一个自变量,通常是时间 ,形如
初值问题
给定微分方程和某个时间点上的函数初始值,求函数在一段时间上的解
欧拉法
欧拉法类似牛顿法,也是用切线近似曲线,预测下一个时间点的函数值
根据一阶泰勒展开
欧拉法公式表示为
其中
- 是当前点近似值
- 是时间步长
- 是当前导数值(切线斜率)
假设现在有一个微分方程,如何知道它的解是存在且唯一的?
利普西茨连续:存在一个常数 使得对所有 都有
只要函数 对 是 Lipschitz 连续,这个微分方程就一定有唯一解
显式梯形法
显式梯形法是对欧拉法采用梯形积分公式思想的改进,先用欧拉预测出下一点,再估算这条梯形斜边的平均斜率,用这个斜率进行更新,此方法是二阶精度收敛更快
中点方法
中点方法通过在每一个时间步上,先用欧拉法预测一个中间点的位置,再用该中点的斜率来更新主变量,此方法是二阶精度收敛更快
自适应步长
普通求解方法每一步都有固定的步长 ,RKF45 是当前最好的自适应步长的 Runge-Kutta 方法,它同时计算四阶和五阶两个近似解,利用两者差值估计误差,误差太大减小步长
Fehlberg 专门构造的一组系数,使得四阶和五阶解共享斜率
五阶近似为
四阶近似为
局部截断误差估计
如果 ,减小 重算
隐式方法
普通的欧拉法是一种显式方法,每一步都用上一步的已知量直接计算新值
刚性问题:一类解中有快速变化成分的ODE,例如
它的解析解是
如果用欧拉法显式求解,步长必须 ,否则解会震荡甚至爆炸
隐式欧拉法求解刚性问题非常稳定,形式如下
下一步的值出现在等式的右侧,必须解方程(一般是非线性方程)才能得到
边值问题
初值问题只给了左边界的值,边值问题则给定两端的值
射线法
将边值问题转化为初值问题,通过尝试不同的初始导数值 ,用初值问题求解 ,看得到的 是否等于目标值 ,如果不等于就调整 直到满足
有限差分法
有限差分法的思路将微分算子用差分公式近似替代,把原来的微分方程转化为一个代数方程组,从而在离散点上求解函数的值
将区间[a, b]分成 N 个等分,记步长为
构造网格点
使用中心差分公式近似二阶导数
替换微分方程为
根据 对 线性相关或非线性相关求解方程组
有限元法
将复杂问题分而治之,把连续域划分成许多小的子区域(称为元素),在每个子区域上用简单函数近似解,然后拼接这些局部近似以构建整个域上的全局近似解
偏微分方程
偏微分方程是含有多个变量的偏导数的方程,形如
其中
- 为待求解的函数
- 为自变量,如时间、空间等
- 涉及一阶、二阶等偏导数
例 1
一维热传导方程
表示在位置 和 时间 时的温度, 是热扩散系数
例 2
二维拉普拉斯方程
例 3
一维波动方程
表示弦在位置 和时间 的位移, 是波速
二阶PDE的分类
二阶线性偏微分方程的一般形式为
主导阶项不同,解的性质不同,对具有两个独立变量的二阶PDE通过判别式 来分类
- 椭圆型 (如拉普拉斯方程)
描述平衡态,解在整个区域内都很"光滑",边界条件影响整个求解域
- 抛物型 (如热传导方程)
描述扩散过程,具有时间方向性,信息向所有方向传播
- 双曲型 (如波动方程)
描述波动现象,信息以有限速度传播
椭圆型方程求解
椭圆型偏微分方程一般不含时间项,描述的是稳态现象,因此非常适合用于稳态建模
有限差分法
使用有限差分代替导数,把连续的偏微分方程变成离散的代数方程系统,在网格节点上求解函数近似值
在前面的边值问题中我们使用有限差分方法构造一维网格求解全微分方程,现在通过构建二维网格并使用有限差分求解偏微分方程
对于矩形区域
- Author:风之旅人
- URL:https://www.hrmi.fun//article/numerical-analysis
- Copyright:All articles in this blog, except for special statements, adopt BY-NC-SA agreement. Please indicate the source!