使用scipy.sparse.block_diag时避免内存错误

2024-05-23 21:44:10 发布

您现在位置:Python中文网/ 问答频道 /正文

我有下面的代码来创建一个非常大的雅可比矩阵,我首先定义了一个关于三分之一相关向量的雅可比矩阵,然后将其转换为一个稀疏矩阵,该矩阵由块对角格式的三个子雅可比矩阵组成。这个“子”雅可比矩阵是密集创建的,以避免通过for循环填充时出现问题。for循环中的这一部分可能会被忽略/视为一个黑盒,但为了以防万一,我已经包含了完整的代码片段

然而,我遇到了一个问题,我遇到了一个错误“无法为sparse.block_diag分配6 gb”,即使我在集群上运行此代码,并且可以创建更大的密集矩阵。是否有一个简单的解决方案,我错过了,以避免这个问题

    J = np.zeros((nparts, 8*self.nnode))
    for elem in range(self.nelem):
        lnodes = self.elem2n[elem]
        lparts = self.elem2p[elem]
        lxc, lyc, lzc = self.xc[lnodes], self.yc[lnodes], self.zc[lnodes]
        lxpart, lypart, lzpart = xpart[lparts], ypart[lparts], zpart[lparts]
        for i in range(len(lnodes)):
            xl, xu = lxc[0], lxc[1]
            yl, yu = lyc[0], lyc[2]
            zl, zu = lzc[0], lzc[4]
            xpl = (lxpart - xl) / (xu - xl)
            ypl = (lypart - yl) / (yu - yl)
            zpl = (lzpart - zl) / (zu - zl)
            xcom = (xu - xl) / 1
            ycom = (yu - yl) / 1
            zcom = (zu - zl) / 1
            af = Amat.A[:,i]
            adfdx = Amat.A[:, 8 + i] * xcom
            adfdy = Amat.A[:, 16 + i] * ycom
            adfdz = Amat.A[:, 24 + i] * zcom
            ad2fdxdy = Amat.A[:, 32 + i] * xcom * ycom
            ad2fdxdz = Amat.A[:, 40 + i] * xcom * zcom
            ad2fdydz = Amat.A[:, 48 + i] * ycom * zcom
            ad3fdxdydz = Amat.A[:, 56 + i] * xcom * ycom * zcom
            J[lparts, 0*self.nnode + lnodes[i]] = cubic_mesh.__tricubic_eval(af, xpl, ypl, zpl)
            J[lparts, 1*self.nnode + lnodes[i]] = cubic_mesh.__tricubic_eval(adfdx, xpl, ypl, zpl)
            J[lparts, 2*self.nnode + lnodes[i]] = cubic_mesh.__tricubic_eval(adfdy, xpl, ypl, zpl)
            J[lparts, 3*self.nnode + lnodes[i]] = cubic_mesh.__tricubic_eval(adfdz, xpl, ypl, zpl)
            J[lparts, 4*self.nnode + lnodes[i]] = cubic_mesh.__tricubic_eval(ad2fdxdy, xpl, ypl, zpl)
            J[lparts, 5*self.nnode + lnodes[i]] = cubic_mesh.__tricubic_eval(ad2fdxdz, xpl, ypl, zpl)
            J[lparts, 6*self.nnode + lnodes[i]] = cubic_mesh.__tricubic_eval(ad2fdydz, xpl, ypl, zpl)
            J[lparts, 7*self.nnode + lnodes[i]] = cubic_mesh.__tricubic_eval(ad3fdxdydz, xpl, ypl, zpl)
    if lowmem:
        J = sparse.block_diag((J, J, J))

Tags: selfeval矩阵meshzplxplcubicxcom