2

我实现了一个有限差分算法来求解 PDE。

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

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

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

这很容易变得太大,并且出现“内存不足”错误。不幸的是sparse,它不适用于 N 维数组。

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

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

我是在打一场失败的战斗还是有更有效的方式来进行(其他类型的对象,矢量化......)?

谢谢你。

4

2 回答 2

4

您可以将数组存储为稀疏的线性形式;即长度等于维度乘积的列向量:

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
于 2022-02-10T16:04:41.127 回答
2

您可以创建一个稀疏矩阵的元胞数组来存储结果。但是,如果使用完整矩阵比使用稀疏矩阵更快并将完整矩阵转换为稀疏矩阵并将其放置在单元格中,则可以在完整矩阵上执行计算。

于 2022-02-10T17:34:34.463 回答