本征和乔尔斯基分解：不稳定问题？

发布于 2024-11-26 15:42:08 字数 1580 浏览 0 评论 0原文

我正在研究 3-D 表面网格，并尝试通过将二次多项式拟合到每个顶点邻域来获取高斯曲率值。为了获得多项式系数，我使用 Eigen 包 (C++) 中包含的标准乔列斯基分解 (LL^T) 来求解我的系统 Ax = b，其中 A 是对称正定矩阵。

当求解该系统时，我只得到其中一个解的最奇怪的值（大约 2000 个顶点！）。这是我的对称正定矩阵的逆矩阵的示例：

Vertex 85: 
288.413   6.45563   3.95006  -48.7131  -17.8428  -30.1609
6.45563  0.199845  0.103744  -1.09986 -0.409136 -0.673235
3.95006  0.103744 0.0979785 -0.666895 -0.241934  -0.41527
-48.7131  -1.09986 -0.666895   8.24633   3.03219   5.09083
-17.8428 -0.409136 -0.241934   3.03219   1.14272    1.8648
-30.1609 -0.673235  -0.41527   5.09083    1.8648   3.16057

Vertex 86: 
496411  2072.94 -8090.48 -48680.5  23267.4 -88125.4
2072.94  8.69566 -33.8001 -203.285  97.1605 -367.998
-8090.48 -33.8001  131.916  793.391 -379.213  1436.27
-48680.5 -203.285  793.391  4773.85 -2281.71  8641.99
23267.4  97.1605 -379.213 -2281.71   1090.6 -4130.56
-88125.4 -367.998  1436.27  8641.99 -4130.56  15644.5

Vertex 87: 
523.131  -1.91535    -7.666  -49.6631   17.1292   -96.191
-1.91535 0.0430769 0.0135821  0.181351 -0.061675  0.352979
-7.666 0.0135821  0.174082   0.72835 -0.249803   1.41043
-49.6631  0.181351   0.72835   4.72212  -1.62466   9.12386
17.1292 -0.061675 -0.249803  -1.62466   0.58246  -3.16276
-96.191  0.352979   1.41043   9.12386  -3.16276   17.7163

只需比较每个矩阵的第一个元素，您就可以看到差异：288（顶点 85），496411（顶点 86） 并返回顶点 87 上的正常值：523.131...

这是我正在使用的代码的一部分（B 是我的对称正定矩阵）：

MatrixXd B = A.adjoint()*A;
MatrixXd x = B.llt().solve(A.adjoint()*b)

我不知道如果我遗漏了一些东西......这可能是一个不稳定问题，也许是舍入错误？如何摆脱这些意想不到的结果？谢谢，米格尔

原文

I'm working on a 3-D surface mesh and trying to obtain gaussian curvature values by fitting a quadric polynomial to each vertex neighborhood. In order to obtain the polynomial coefficients I'm using the standard cholesky decomposition (LL^T) included on Eigen package (C++), for solving my system Ax = b, where A is a symmetric, positive definite matrix.

When solving the system I'm getting the strangest values only for one of the solution (among about 2000 vertices!). Here's a sample of the what would be the inverse of my symmetric, positive definite matrix:

Vertex 85: 
288.413   6.45563   3.95006  -48.7131  -17.8428  -30.1609
6.45563  0.199845  0.103744  -1.09986 -0.409136 -0.673235
3.95006  0.103744 0.0979785 -0.666895 -0.241934  -0.41527
-48.7131  -1.09986 -0.666895   8.24633   3.03219   5.09083
-17.8428 -0.409136 -0.241934   3.03219   1.14272    1.8648
-30.1609 -0.673235  -0.41527   5.09083    1.8648   3.16057

Vertex 86: 
496411  2072.94 -8090.48 -48680.5  23267.4 -88125.4
2072.94  8.69566 -33.8001 -203.285  97.1605 -367.998
-8090.48 -33.8001  131.916  793.391 -379.213  1436.27
-48680.5 -203.285  793.391  4773.85 -2281.71  8641.99
23267.4  97.1605 -379.213 -2281.71   1090.6 -4130.56
-88125.4 -367.998  1436.27  8641.99 -4130.56  15644.5

Vertex 87: 
523.131  -1.91535    -7.666  -49.6631   17.1292   -96.191
-1.91535 0.0430769 0.0135821  0.181351 -0.061675  0.352979
-7.666 0.0135821  0.174082   0.72835 -0.249803   1.41043
-49.6631  0.181351   0.72835   4.72212  -1.62466   9.12386
17.1292 -0.061675 -0.249803  -1.62466   0.58246  -3.16276
-96.191  0.352979   1.41043   9.12386  -3.16276   17.7163

Just by comparing the first element of each matrix you can see the differences: 288 (vertex 85), 496411 (vertex 86) and back to normal values on vertex 87: 523.131...

Here's part of the code I'm using (B is my symmetric, positive definite matrix):

MatrixXd B = A.adjoint()*A;
MatrixXd x = B.llt().solve(A.adjoint()*b)

I don't know if I'm missing something... It might be an instability issue, maybe round-off errors? How to get rid of these unexpected results?
Thanks,
Miguel

分享到QQ

分享到微博