在 Matlab 中有效存储大部分为零的 N 维数组

Question

我实现了有限差分算法来求解偏微分方程。

网格是大小为 [Nx, Nz] 的结构化二维域，已求解 Nt 次。

我预先分配包含所有解决方案的对象：

sol = zeros(Nx, Nz, Nt, 'single') ;

这很容易变得太大，并且我收到“内存不足”错误。不幸的sparse是不适用于 N 维数组。

对于这个问题来说，知道这些值并不重要，不言而喻，随着网格间距的减小和模拟时间的增加，RAM 使用量呈指数增长。

我知道我不需要为了改进解决方案而存储每个时刻。只需存储前两个时间步就足够了。然而，出于后处理的原因，我需要在所有时间步长（或至少是总数的约数）访问解决方案。这可能有助于指定，即使在解决方案之后，网格仍然主要由零填充。

我是在打一场失败的战斗还是有更有效的方法来继续（其他类型的对象，矢量化......）？

谢谢。

Answer 1

您可以以稀疏的线性形式存储数组；即，长度等于维度乘积的列向量：

sol = sparse([], [], [], Nx*Nz*Nt, 1); % sparse column vector containing zeros

然后，而不是正常索引，

sol(x, z, t),

您需要将索引x, z,t转换为相应的线性索引：

• 对于您使用的标量索引

  sol(x + Nx*(z-1) + Nx*Nz*(t-1))
  

  为了方便起见，您可以定义一个辅助函数：

  ind = @(sol, x, y, t) sol(x + Nx*(z-1) + Nx*Nz*(t-1))
  

  因此索引变得更具可读性：

  ind(sol, x, z, t)
  

• 对于一般（数组）索引，您需要reshape沿不同维度的索引，以便隐式扩展产生适当的线性索引：

  sol(reshape(x,[],1,1) + Nx*(reshape(z,1,[],1)-1) + Nx*Nz*(reshape(t,1,1,[])-1))
  

  当然也可以封装成一个函数。

检查线性索引的转换是否有效（一般情况，使用非稀疏数组与普通索引进行比较）：

Nx = 15; Nz = 18; Nt = 11;
sol = randi(9, Nx, Nz, Nt);
x = [5 6; 7 8]; z = 7; t = [4 9 1];
isequal(sol(x, z, t), ...
    sol(reshape(x,[],1,1) + Nx*(reshape(z,1,[],1)-1) + Nx*Nz*(reshape(t,1,1,[])-1)))

给出

ans =
  logical
   1