2

我正在尝试从已经存在的 csc-matrix 中创建一个 petsc-matrix。考虑到一点,我创建了以下示例代码:

import numpy as np
import scipy.sparse as sp
import math as math
from petsc4py import PETSc
n=100

A = sp.csc_matrix((n,n),dtype=np.complex128)
print A.shape
A[1:5,:]=1+1j*5*math.pi

p1=A.indptr
p2=A.indices
p3=A.data
petsc_mat = PETSc.Mat().createAIJ(size=A.shape,csr=(p1,p2,p3))

只要矩阵只包含实数值,这种方法就可以很好地工作。当矩阵很复杂时,运行这段代码会产生 TypeError: Cannot cast array data from dtype('complex128') to dtype('float64') according to the rule 'safe'. 我试图找出错误发生的确切位置,但无法理解回溯:

petsc_mat = PETSc.Mat().createAIJ(size=A.shape,csr=(p1,p2,p3))  File "Mat.pyx", line 265, in petsc4py.PETSc.Mat.createAIJ (src/petsc4py.PETSc.c:98970)
File "petscmat.pxi", line 662, in petsc4py.PETSc.Mat_AllocAIJ (src/petsc4py.PETSc.c:24264)
File "petscmat.pxi", line 633, in petsc4py.PETSc.Mat_AllocAIJ_CSR (src/petsc4py.PETSc.c:23858)
File "arraynpy.pxi", line 136, in petsc4py.PETSc.iarray_s (src/petsc4py.PETSc.c:8048)
File "arraynpy.pxi", line 117, in petsc4py.PETSc.iarray (src/petsc4py.PETSc.c:7771)

有没有一种有效的方法从复杂的 scipy csc 矩阵中创建一个 petsc 矩阵(我想稍后检索一些特征对)?

如果你们能帮我找到我(希望不是太明显)的错误,我会非常高兴。

4

1 回答 1

2

我在让 PETSc 工作时遇到了麻烦,所以我不止一次配置了它,在最后一次运行中我显然忘记了选项--with-scalar-type=complex

这是我应该做的:

  • 检查日志文件$PETSC_DIR/arch-linux2-c-opt/conf/configure.log

  • 或者看看reconfigure-arch-linux2-c-opt.py

在那里您可以找到用于配置 PETSc 的所有选项。如果您也使用 SLEPc,您还需要重新编译它。现在,由于我将选项 ( --with-scalar-type=complex) 添加到重新配置脚本并运行它,所以一切正常。

于 2015-04-09T18:48:55.970 回答