5

有一个用于Ax=b类型问题的大矩阵。 A是对称的。有什么算法可以让我们只保存一半的矩阵并对其进行操作x=A\b吗?

4

3 回答 3

6

您只会节省一半的内存,但您可以通过创建矩阵的平面版本来做到这一点,保存它,然后将其展开。请注意,所需的额外时间可能并不值得节省:

% pretend this is symettric...
A = rand(10, 10);

% store it as a flat list
flatA = [];
for k = 1:size(A,1)
    flatA = [flatA; A([1:k],k)];
end

save('flatA.mat', 'flatA');

% undo
N = (-1 + (1+4*2*length(flatA))^0.5)/2;
newA = NaN(N,N); 
cur = 1;
for k = 1:N
    len = k;
    newA(1:len, k) = flatA(cur:cur+len-1);
    cur = cur + len;
end
% real A cannot have NaNs or this trick fails
newA(isnan(newA)) = newA(isnan(newA'))';
于 2011-06-07T09:25:32.370 回答
3

这是一个想法,但我还没有测试过。如果您的对称矩阵是正定的,请对对称矩阵 A 进行 Cholesky 分解,得到 A = U*U'。如果您使用 MATLAB 的内置稀疏矩阵将 U 存储为稀疏矩阵,则您拥有所需的一切,并且使用了大约一半的内存。由于它使用 MATLAB 的稀疏矩阵类型,因此您可以使用标准 MATLAB 函数对其进行操作,只要您记住 A = U*U'

For example, to compute A*x = b, use x = U'\U\b. Unlike the other proposed solutions, MATLAB will never be actually using a full matrix internally, and will even use accelerated solvers that will take advantage of the fact that you're only solving with triangular. The cost is that to solve a single system, you've actually running the backslash operator twice (see above). However, that's the price you pay for never instantiating the full matrix.

于 2011-06-07T17:38:07.697 回答
1

如果您提取上三角部分并转换为稀疏矩阵,它应该可以节省内存。这种技术相当快。

% Create a symmetric matrix
A = rand(1000);
A = A + A';

% Extract the upper triangular part
B = sparse(triu(A))              % This is the important line, the rest is just checking.

% Undo
BB = full(B);
C = BB + BB' - diag(diag(BB));

% Check correct
assert(all(A(:) == C(:)));

% Save matrices
save('A.mat', 'A');
save('B.mat', 'B');

% Get file sizes
infoA = dir('A.mat'); infoA.bytes
infoB = dir('B.mat'); infoB.bytes

编辑以澄清木片的事情

上述解决方案演示了一种使用较少文件空间保存矩阵的方法。该矩阵B也比原始矩阵占用更少的内存A。如果你想对 进行线性代数运算B,那效果很好。比较

b = rand(1000);
fullTriUA = triu(A);
sparseTriUA = sparse(fullTriUA);  %same as B above
fullResult = fullTriUA\b;
sparseResult = sparseTriUA\b;
assert(all(fullResult(:) == sparseResult(:)))
于 2011-06-07T10:35:32.987 回答