主要内容

ichol

不完全Cholesky分解

折叠所有页面

语法

L = ichol(A)
L = ichol(A，选项)

描述

L = ichol(一个）的不完全Cholesky分解一个补零。

L = ichol(一个，选项）的不完全Cholesky分解一个使用结构指定的选项选项．

默认情况下,ichol的下三角形一个得到下三角因子。

输入参数

一个

稀疏矩阵

选项

结构，最多包含五个字段:

字段名 总结 描述
类型 分解的类型

指示要执行的不完整Cholesky的哪种风味。该字段的有效值为“nofill”而且“信息与通信技术”．的“nofill”变体执行零填充不完全Cholesky（集成电路(0))．的“信息与通信技术”变体执行阈值下降的不完全Cholesky(ICT)。默认值为“nofill”．
droptol 类型为时的落差公差“信息与通信技术”

执行ICT时使用的非负标量作为落差公差。从结果因子中删除大小小于局部容错的元素，但对角线元素永远不会删除。步进时的局部落差公差j因式分解是规范((结束,j), 1) * droptol．“droptol”被忽略，如果“类型”是“nofill”．默认值为0．
michol 是否执行修改的不完全Cholesky

指示是否改良的不完全Cholesky(MIC)执行。这个领域可能是“上”或“关闭”．在执行MIC时，对角线会对删除的元素进行补偿，以加强这种关系A*e = L*L'*e在哪里e = ones(size(A,2)，1)．默认值为“关闭”．
diagcomp 用指定的系数进行补偿不完全Cholesky

实非负标量用作全局对角线移位α形成不完全Cholesky因子。也就是说，不执行不完全的Cholesky一个的因式分解A + alpha*diag(diag(A))就形成了。默认值为0．
形状 确定引用和返回哪个三角形

有效值为“上”而且“低”．如果“上”时，只指定上三角形一个是引用的R是这样构造的一个近似为‘* R．如果“低”时，只指定了一个是引用的l是这样构造的一个近似为L * L '．默认值为“低”．

例子

全部折叠

打开实时脚本

这个例子生成了一个不完整的Cholesky分解。

从一个对称正定矩阵开始，一个：

N = 100;A = delsq(numgrid(“年代”, N));

一个是带有狄利克雷边界条件的100 × 100平方网格上的二维五点离散负拉普拉斯算子。的大小一个是98*98 = 9604(不是10000，因为网格的边界是用来施加狄利克雷条件的)。

无填充不完全Cholesky分解是在相同位置上只包含非零的一种分解一个包含非零。这种因式分解的计算成本非常低。尽管产品L * L '是非常不同的一个，乘积L * L '将匹配一个直到四舍五入为止。

L = ichol(A);规范(L * L ',“摇来摇去”)。/规范(,“摇来摇去”）
Ans = 0.0916
规范(A - (L * L ')。* spones (A),“摇来摇去”)。/规范(,“摇来摇去”）
Ans = 4.9606e-17

ichol也可用于生成阈值下降的不完全Cholesky分解。随着耐滴度的降低，因子趋于致密，产物趋于致密L * L '是一个更好的近似一个．下面的图显示了不完全因子分解的相对误差与下降容忍度的关系，以及不完全因子的密度与完全Cholesky因子的密度之比。

n = size(A,1);nols = 20;Droptol = logspace(-8,0, nols);NRM = 0 (1, nols);Nz = 0 (1, nols);nzComplete = nnz(chol(A，“低”));为k = 1: L = ichol(A,struct(“类型”，“信息与通信技术”，“droptol”droptol (k)));nz(k) = nnz(L);nrm(k) = norm(A-L*L'，“摇来摇去”)。/规范(,“摇来摇去”）;结束图重对数(droptol全国抵抗运动,“线宽”2)标题(放的宽容与规范(L * L”、“向后”)。/规范(A,“来回”)）

图中包含一个轴对象。标题为Drop tolerance vs norm(A-L*L'，'fro')./norm(A，'fro')的axis对象包含一个line类型的对象。

图semilogx (droptol nz. / nzComplete,“线宽”2)标题(“跌落容忍度vs填充比ichol/chol”）

图中包含一个轴对象。标题为Drop tolerance vs fill ratio ichol/chol的axes对象包含一个类型为line的对象。

相对误差通常与落差容差处于同一量级，尽管不能保证这一点。

打开实时脚本

这个例子展示了如何使用不完全Cholesky分解作为改善收敛性的预处理。

创建一个对称正定矩阵，一个．

N = 100;A = delsq(numgrid(“年代”, N));

创建一个不完整的Cholesky分解作为预处理pcg．右边用一个常数向量。作为基线，执行pcg不需要预处理。

b = ones(size(A,1)，1);Tol = 1e-6;Maxit = 100;[x0,fl0,rr0,it0,rv0] = pcg(A,b,tol,maxit);

请注意,Fl0 = 1这表明pcg在允许的最大迭代中，没有将相对残差驱动到所要求的公差。尝试无填充不完全Cholesky分解作为预处理条件。

L1 = ichol(A);[x1,fl1,rr1,it1,rv1] = pcg(A,b,tol,maxit,L1,L1');

Fl1 = 0，表明pcg收敛到所请求的公差，并在59次迭代中完成(的值it1)．由于这个矩阵是一个离散的拉普拉斯矩阵，使用改进的不完全Cholesky可以得到一个更好的预条件。改进的不完全Cholesky分解构造了一个近似分解，它保留了算子对常数向量的作用。也就是说,规范(* e-L * (L ' * e))近似为零e = ones(size(A,2)，1)尽管规范(L * L”、“向后”)/规范(A,“来回”)不接近于零。没有必要为此语法指定类型，因为nofill是默认值，但这是很好的实践。

选项。类型=“nofill”；选项。michol =“上”；L2 = ichol(A,options);e = ones(size(A,2)，1);规范(* e-L2 * (L2‘* e))
Ans = 3.7983e-14
[x2,fl2,rr2,it2,rv2] = pcg(A,b,tol,maxit,L2,L2');

pcg收敛(Fl2 = 0)但只进行了38次迭代。绘制所有三个收敛历史可以显示收敛。

semilogy(0:麦克斯特,rv0. /规范(b),“b”。）;持有在semilogy (0: it1 rv1. /规范(b),“r”。）;semilogy (0: it2 rv2. /规范(b),“k”。）;传奇(“没有预调节器”，“集成电路(0)”，“麦克风(0)”）;

图中包含一个轴对象。axis对象包含3个line类型的对象。这些对象表示No preconditioners, IC(0)， MIC(0)。

图中显示，修改后的不完全Cholesky预条件创造了更快的收敛速度。

您还可以尝试使用阈值降低的不完全Cholesky分解。下面的图显示的收敛pcg预调节装置具有不同的下降公差。

L3 = ichol(A, struct(“类型”，“信息与通信技术”，“droptol”1 e 1));[x3, fl3 rr3、it3 rv3] = pcg (A, b,托尔,麦克斯特,L3, L3”);L4 = ichol(A, struct(“类型”，“信息与通信技术”，“droptol”, 1依照));[x4, fl4 rr4、it4 rv4] = pcg (A, b,托尔,麦克斯特、L4、L4的);L5 = ichol(A, struct(“类型”，“信息与通信技术”，“droptol”1 e - 3));[x5, fl5 rr5、it5 rv5] = pcg (A, b,托尔,麦克斯特、L5、L5”);图semilogy(0:麦克斯特,rv0. /规范(b),“b -”，“线宽”2);持有在semilogy (0: it3 rv3. /规范(b),b -。，“线宽”2);semilogy (0: it4 rv4. /规范(b),“b——”，“线宽”2);semilogy (0: it5 rv5. /规范(b),”乙:“，“线宽”2);传奇(“没有预调节器”，“ICT (1 e 1)”，“ICT(1依照)”，.．.“ICT (1 e - 3)”，“位置”，“东南”）;

图中包含一个轴对象。axis对象包含4个line类型的对象。这些对象代表No preconditioners, ICT(1e-1)， ICT(1e-2)， ICT(1e-3)。

请注意，不完整的Cholesky预调理器构造了降容1)依照表示为ICT(1飞行)．

与零填充不完全Cholesky一样，阈值下降分解可以受益于修改(即。选项。michol =“上”)因为矩阵是由椭圆偏微分方程产生的。与麦克风(0)，基于降低不完全Cholesky的修正阈值将保留预处理条件对常数向量的作用，即规范(* e-L * (L ' * e))近似为零。

打开实时脚本

的用法diagcomp选择ichol．

正定矩阵的不完全Cholesky分解并不总是存在的。下面的代码构造了一个随机对称正定矩阵，并尝试用pcg．

S = rng(“默认”）;A = sprandsym(1000,1e-2,1e-4,1);rng(年代);b = full(sum(A,2));[x0,fl0,rr0,it0,rv0] = pcg(A,b,1e-6,100);

由于收敛性没有达到，尝试构造一个不完全的Cholesky预处理条件。

L = ichol(A,struct(“类型”，“信息与通信技术”，“droptol”1 e - 3));
遇到非正枢轴。

如果ichol如上所述，您可以使用diagcomp构造移位不完全Cholesky分解的选项。也就是说，不是构造l这样L * L '接近一个，ichol对角补偿结构l这样L * L '接近M = A + alpha*diag(diag(A))没有明确地形成米．由于对角占优矩阵总是存在不完全因式分解，α可以找到做什么米对角占优。

alpha = max(sum(abs(A)，2)./diag(A)))-2
Alpha = 62.9341
L1 = ichol(A, struct(“类型”，“信息与通信技术”，“droptol”1 e - 3,“diagcomp”、α));[x1,fl1,rr1,it1,rv1] = pcg(A,b,1e-6,100,L1,L1');

在这里,pcg在期望的迭代次数内仍然不能收敛到期望的容差，但是如下图所示，收敛对pcg有这个预处理的比没有预处理的要多。选择小一点的α也许会有帮助。通过一些实验，我们可以确定一个合适的值α．

Alpha = .1;L2 = ichol(A, struct(“类型”，“信息与通信技术”，“droptol”1 e - 3,“diagcomp”、α));[x2,fl2,rr2,it2,rv2] = pcg(A,b,1e-6,100,L2,L2');

现在,pcg收敛和一个图可以显示每个前置条件的收敛历史。

semilogy (0:100 rv0. /规范(b),“b”。）;持有在；semilogy (0:100 rv1. /规范(b),“r”。）;semilogy (0: it2 rv2. /规范(b),“k”。）;传奇(“没有预调节器”，“\alpha \约63”，'\alpha = .1'）;包含(的迭代次数）;ylabel (的相对剩余的）;

图中包含一个轴对象。axis对象包含3个line类型的对象。这些对象表示No preconditioning， \alpha \约63，\alpha = .1。

提示

  • 这个例程给出的因子可以作为线性方程组用迭代方法求解的预条件，例如pcg或minres．

  • ichol仅适用于稀疏方阵

参考文献

Saad, Yousef。“预处理技术。”稀疏线性系统的迭代方法．PWS出版公司，1996年。

[2] Manteuffel T.A. <正定线性系统的不完全因式分解技术>。数学。第一版。34, 473-497, 1980。

扩展功能

另请参阅

|||