Commits

Graham Markall committed 935313a

Add getLocalSubMatrix, restoreLocalSubMatrix, and createNest to Mat

Comments (0)

Files changed (2)

src/PETSc/Mat.pyx

             Mat_AllocDense_DEFAULT(self.mat)
         return self
 
+    def createNest(self, isrow not None, iscol not None, submats not None, comm=None):
+        cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT)
+        cdef int nr = isrow.getSize(), nc = iscol.getSize()
+        cdef PetscMat newmat = NULL
+        matcount = len(submats)
+        cdef PetscMat *mats = <PetscMat *> malloc(matcount*sizeof(PetscMat))
+        cdef PetscMat mat
+        cdef Mat m
+
+        
+        for i in range(matcount):
+            m = submats[i]
+            mat = m.mat
+            mats[i] = mat
+
+        #cdef PetscIS cis
+        #cdef IS ind
+        #cdef PetscIS *cisrows = <PetscIS *> malloc(nr*sizeof(PetscIS))
+        #cdef PetscIS *ciscols = <PetscIS *> malloc(nc*sizeof(PetscIS))
+        #for i in range(nr):
+        #    ind = isrow[i]
+        #    cis = ind.iset
+        #    cisrows[i] = cis
+        #for i in range(nc):
+        #    ind = iscol[i]
+        #    cis = ind.iset
+        #    ciscols[i] = cis
+
+        #CHKERR( MatCreateNest(ccomm, nr, cisrows, nc, ciscols, mats, &newmat) )
+        CHKERR( MatCreateNest(ccomm, nr, NULL, nc, NULL, mats, &newmat) )
+        self.mat = newmat
+        return self
+
     def setPreallocationNNZ(self, nnz, bsize=None):
         cdef PetscInt bs = PETSC_DECIDE
         CHKERR( Mat_BlockSize(bsize, &bs) )
                                 reuse, &submat.mat) )
         return submat
 
+    def getLocalSubMatrix(self, IS isrow not None, IS iscol not None,
+                          Mat submat=None):
+        cdef PetscIS ciscol = NULL
+        if submat is None: submat = Mat()
+        CHKERR( MatGetLocalSubMatrix(self.mat, isrow.iset, iscol.iset, 
+                                     &submat.mat) )
+        return submat
+
+    def restoreLocalSubMatrix(self, IS isrow not None, IS iscol not None, 
+                              Mat submat not None):
+        CHKERR( MatRestoreLocalSubMatrix(self.mat, isrow.iset, iscol.iset, 
+                                         &submat.mat) )
+
     def increaseOverlap(self, IS iset not None, overlap=1):
         cdef PetscInt ival = asInt(overlap)
         CHKERR( MatIncreaseOverlap(self.mat, 1, &iset.iset, ival) )

src/PETSc/petscmat.pxi

     int MatCreateLRC(PetscMat,PetscMat,PetscMat,PetscMat*)
     int MatCreateSubMatrix(PetscMat,PetscIS,PetscIS,PetscMat*)
     int MatCreateShell(MPI_Comm,PetscInt,PetscInt,PetscInt,PetscInt,void*,PetscMat*)
+    int MatCreateNest(MPI_Comm,PetscInt,PetscIS[],PetscInt,PetscIS[],PetscMat[],PetscMat*)
 
     int MatCreateSeqAIJWithArrays(MPI_Comm,PetscInt,PetscInt,PetscInt[],PetscInt[],PetscScalar[],PetscMat*)
     int MatCreateSeqBAIJWithArrays(MPI_Comm,PetscInt,PetscInt,PetscInt,PetscInt[],PetscInt[],PetscScalar[],PetscMat*)
     int MatMerge(MPI_Comm,PetscMat,PetscInt,PetscMatReuse,PetscMat*)
     int MatGetSubMatrix(PetscMat,PetscIS,PetscIS,PetscMatReuse,PetscMat*)
     int MatGetSubMatrices(PetscMat,PetscInt,PetscIS[],PetscIS[],PetscMatReuse,PetscMat*[])
+    int MatGetLocalSubMatrix(PetscMat,PetscIS,PetscIS,PetscMat*)
+    int MatRestoreLocalSubMatrix(PetscMat,PetscIS,PetscIS,PetscMat*)
     int MatIncreaseOverlap(PetscMat,PetscInt,PetscIS[],PetscInt)
     int MatGetDiagonalBlock(PetscMat,PetscMat*)