# minres
求解线性方程-最小残差法
函数库: TyMath
# 语法
(x,stats)=minres(A,b::AbstrctVector{T};atol::T=√eps(T)/100, rtol::T=√eps(T)/100)
# 说明
(x,stats)=minres(A,b::AbstrctVector{T};atol::T=√eps(T)/100, rtol::T=√eps(T)/100) 解决移位线性最小二乘问题，minimize‖b - (A + λI)x‖₂² 或平移线性系统 (A + λI) x = b，使用 MINRES 方法，其中 λ ≥ 0 是移位参数，其中 A 是正方形且对称的。
当 A 为正定时，MINRES 在形式上等同于将 CR 应用于 Ax=b，但通常更稳定。
也适用于这种情况，其中 A 是不定的，MINRES 产生单调残差 ‖r‖₂ 和最优残差 ‖Aᵀr‖₂。可以以线性算子的形式提供预处理器 M 并假定它是对称的和正定的。示例
# 示例
正定矩阵 - 预条件共轭梯度法求解
创建线性对称系统，利用预条件共轭梯度法求解。
using TyMath
cg_tol=1e-6
A,b=[[1 2 3;2 8 10;3.0 10 5],[1e-9,5.5 ,9.998]]##确定系数矩阵和b
(x,stats)=minres(A,b)
([-2.1877499981250006, 1.9372499993750012, -0.562249999874999], SimpleStats
niter: 3
solved: true
inconsistent: true
residuals: []
Aresiduals: []
κ₂(A): []
timer: 98.50μs
status: found approximate minimum least-squares solution
)
对角矩阵 - 预条件共轭梯度法求解
创建线性对称系统，利用预条件共轭梯度法求解。
using TyMath
cg_tol=1e-6
A, b = [1.0 0 0;0 5 0;0 0 8],[5.0,5.0,5.0]
(x,stats)=minres(A,b)
([4.999999999999998, 1.0000000000000009, 0.6250000000000004], SimpleStats
niter: 3
solved: true
inconsistent: false
residuals: []
Aresiduals: []
κ₂(A): []
timer: 11.80μs
status: found approximate zero-residual solution
)
线性方程组的迭代解
使用采用默认设置的 minres 求解系数矩阵为方阵的线性方程组，然后在求解过程中调整使用的容差和迭代次数。
创建一个稀疏三对角矩阵 A 作为系数矩阵。使用 A 的行总和作为 Ax=b 右侧的向量 b，以便 x 的预期解是由 1 组成的向量。
using TyMath
n = 400
on = ones(n)
A = spdiagm(-1 => -2*on[1:n-1], 0 => 4*on, 1 => -2*on[1:n-1])
b = sum(A, dims=2)
b = vec(b);
使用 ty_minres 函数求解线性方程组 Ax = b。
x, = ty_minres(A, b);
x
400-element Vector{Float64}:
0.9302325581395354
0.8607671398369076
0.7919057686499551
0.7239504681365149
0.6572032618544249
0.591966173361522
0.5285412262156448
0.4672304439746301
⋮
0.4672304439746301
0.5285412262156448
0.591966173361522
0.6572032618544249
0.7239504681365149
0.7919057686499551
0.8607671398369076
0.9302325581395354
使用指定了预条件子的 minres
检查使用指定了预条件子矩阵的 minres 求解线性方程组的效果。
创建一个对称正定带状系数矩阵。
using TyMath
A = delsq(numgrid("S", 102));
定义 b 以使实际解是全为 1 的向量。
b = sum(A,dims = 2);
设置容差和最大迭代次数。
tol = 1e-12;
maxit = 100;
使用 minres 根据请求的容差和迭代次数求解。指定六个输出以返回有关求解过程的信息：
x0 是计算 A*x0 = b 所得的解；
fl0 是指示算法是否收敛的标志；
rr0 是计算的解 x0 的残差；
it0 是计算出 x0 时的迭代次数；
rv0 是 ‖Ax−b‖ 的残差历史记录组成的向量；
rvcg0 是 ‖AᵀAx-Aᵀb‖ 的共轭梯度残差历史记录组成的向量。
x0, fl0, rr0, it0, rv0, rvcg0 = ty_minres(A, b, tol, maxit)
fl0
minres 未在请求的 100 次迭代内收敛至请求的容差 1e-12，因此 fl0 为 1。
使用函数句柄代替数值矩阵
定义 afun 函数。
using TyMath
function afun(x::AbstractVector{T}) where T
weights = vcat((10:-1:0), (1:10)) # 确保 weights 是一个 21 元素的列向量
y = [0.0; x[1:20]] .+
weights .* x .+
[x[2:21]; 0.0]
return y
end
创建矩阵 A。
A = TestArrays.wilk( 21)
21×21 Matrix{Float64}:
10.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 … 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
1.0 9.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
0.0 1.0 8.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
0.0 0.0 1.0 7.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
0.0 0.0 0.0 1.0 6.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 1.0 5.0 1.0 0.0 … 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 0.0 1.0 4.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 0.0 0.0 1.0 3.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
⋮ ⋮ ⋱ ⋮ ⋮
0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 3.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 4.0 1.0 0.0 0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 … 0.0 1.0 5.0 1.0 0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 6.0 1.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 7.0 1.0 0.0 0.0
0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 8.0 1.0 0.0
0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 9.0 1.0
0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 … 0.0 0.0 0.0 0.0 0.0 0.0 1.0 10.0
创建向量 b。
b = ones(21);
设置 tol 和 maxit 参数。
tol = 1e-12;
maxit = 50;
求解线性系统 Ax = b。
x1, stats = ty_minres(afun, b, tol, maxit)
minres converged at iteration 11 to a solution with relative residual 2.448730691563479e-15.
([0.09101018481309177, 0.08989815186908864, 0.09990644836511511, 0.11085026120999213, 0.12414172316494039, 0.14429939980036555, 0.15436127783323175, 0.23825548886670725, 0.13087225556664642, 0.49999999999999983 … 0.5, 0.13087225556664636, 0.23825548886670725, 0.15436127783323172, 0.14429939980036555, 0.12414172316494036, 0.11085026120999214, 0.09990644836511513, 0.08989815186908862, 0.09101018481309177], 0, 2.448730691563479e-15, 11, [4.58257569495584, 1.7237035261766585, 0.7618748106770391, 0.42338436765117127, 0.35684160390923325, 0.354688037699217, 0.277144483796505, 0.1800779584395271, 0.10833378356337581, 0.05280991429796305, 0.01604965982657888, 1.1221493750651204e-14], [4.58257569495584, 3.9430113443329655, 0.0, 4.469902172583435, 4.552096394379473, 4.5888598150794255, 5.574660329063579, 4.578951403107975, 4.578082374773965, 4.580740577447539, 4.582061929322622, 4.58252234365156, 4.58257569495584])
提供初始估计值
检查向 minres 提供解的初始估计值的效果。
创建一个三对角稀疏矩阵。使用每行的总和作为 Ax=b 右侧的向量，以便 x 的预期解是由 1 组成的向量。
using TyMath
n = 900
e = ones(n)
A = spdiagm(-1 => e[1:n-1], 0 => 2*e, 1 => e[1:n-1])
b = sum(A, dims=2);
使用 minres 求解 Ax=b 两次：一次是使用默认的初始估计值，一次是使用解的良好初始估计值。对这两个解都使用 200 次迭代，并将初始估计值指定为所有元素均等于 0.99 的向量。
maxit = 200
x1, = ty_minres(A, b, nothing ,maxit)
x1
minres converged at iteration 27 to a solution with relative residual NaN.
900-element Vector{Float64}:
0.9998580901231674
1.0002817226507053
⋮
1.0002817226507053
0.9998580901231674
# 输入参数
A - 输入矩阵对称正定矩阵
minres 方法要求线性系统为对称方形系统，A 为系数矩阵。
数据类型： Int64 | Int32 | Int16 | Int128 | Float16 | Float32 | Float64
b - 输入矩阵矩阵
常数向量，一般用来表示 Ax=b 中 b 向量。
数据类型： Int64 | Int32 | Int16 | Int128 | Float16 | Float32 | Float64
atol - 绝对容差正标量
绝对误差（Absolute error）=测量值-真值，是测量值（单一测量值或多次测量值的均值）与真值之差。若测量结果大于真值，误差为正，反之为负。
数据类型：Float16 | Float32 | Float64
rtol - 相对容差正标量
相对误差（Relative error）=绝对误差÷真值，为绝对误差与真值的比值（可以用百分比（%）、千分比（ppt）、百万分比（ppm）表示，但常以百分比表示）。一般来说，相对误差更能反映测量的可信程度。
数据类型：Float16 | Float32 | Float64
# 算法
minres 和 symmlq 方法是支撑共轭梯度法 PCG 的 Lanczos 方法的变体。像 PCG 一样，系数矩阵仍需是对称矩阵，但 minres 和 symmlq 允许它是不定矩阵（不要求所有特征值都必须为正值）。这是通过避免 Lanczos方法中通常存在的隐式 LU 分解来实现的。当不定矩阵遇到主元为零的情况时，该方法容易出现故障。minres 最小化 2 - 范数残差，而 symmlq 使用 LQ 分解求解投影方程组，并使残差与所有先前的残差正交。GMRES 方法旨在将 minres 推广到非对称问题。
# 参考文献
[1] C. C. Paige and M. A. Saunders, Solution of Sparse Indefinite Systems of Linear Equations, SIAM Journal on Numerical Analysis, 12(4), pp. 617–629, 1975.