概念
LDL分解是经典Cholesky分解的一个变形：
- L 为下三角矩阵，且对角元素必须为1；
- D为对角矩阵。
转换
特征
- LDL变形如果得以有效运行，构造及使用时所需求的空间及计算的复杂性与经典Cholesky分解是相同的，但是可避免提取平方根。
- 某些不存在Cholesky分解的不定矩阵，也可以运行LDL分解，此时矩阵
中会出现负数元素。
因此人们更倾向于使用LDL分解。对于实数矩阵，该种分解的形式可被改写成
计算
减少乘法运算，引入辅助量,则公式为：
代码
#include <iostream>
#include <string.h>
#include <stdio.h>
#include <vector>
#include <math.h>
using namespace std;
const int N = 1005;
typedef double Type;
Type A[N][N], L[N][N];
/** 分解A得到A = L * L^T */
void Cholesky(Type A[][N], Type L[][N], int n)
{
for(int k = 0; k < n; k++)
{
Type sum = 0;
for(int i = 0; i < k; i++)
sum += L[k][i] * L[k][i];
sum = A[k][k] - sum;
L[k][k] = sqrt(sum > 0 ? sum : 0);
for(int i = k + 1; i < n; i++)
{
sum = 0;
for(int j = 0; j < k; j++)
sum += L[i][j] * L[k][j];
L[i][k] = (A[i][k] - sum) / L[k][k];
}
for(int j = 0; j < k; j++)
L[j][k] = 0;
}
}
/** 回带过程 */
vector<Type> Solve(Type L[][N], vector<Type> X, int n)
{
/** LY = B => Y */
for(int k = 0; k < n; k++)
{
for(int i = 0; i < k; i++)
X[k] -= X[i] * L[k][i];
X[k] /= L[k][k];
}
/** L^TX = Y => X */
for(int k = n - 1; k >= 0; k--)
{
for(int i = k + 1; i < n; i++)
X[k] -= X[i] * L[i][k];
X[k] /= L[k][k];
}
return X;
}
void Print(Type L[][N], const vector<Type> B, int n)
{
for(int i = 0; i < n; i++)
{
for(int j = 0; j < n; j++)
cout<<L[i][j]<<" ";
cout<<endl;
}
cout<<endl;
vector<Type> X = Solve(L, B, n);
vector<Type>::iterator it;
for(it = X.begin(); it != X.end(); it++)
cout<<*it<<" ";
cout<<endl;
}
int main()
{
int n;
cin>>n;
memset(L, 0, sizeof(L));
for(int i = 0; i < n; i++)
{
for(int j = 0; j < n; j++)
cin>>A[i][j];
}
vector<Type> B;
for(int i = 0; i < n; i++)
{
Type y;
cin>>y;
B.push_back(y);
}
Cholesky(A, L, n);
Print(L, B, n);
return 0;
}
/**data**
4
4 -2 4 2
-2 10 -2 -7
4 -2 8 4
2 -7 4 7
8 2 16 6
*/
