Re: [petsc-dev] Using PETSc MatIS, how to get local matrix (= one domain) before and after assembly ?

Wed, 24 May 2017 04:42:51 -0700


On May 24, 2017, at 11:46 AM, Franck Houssen <[email protected]> wrote:

The code I sent compile and run at my side with petsc-3.7.6 (on debian/testing with gcc-6.3). The code you sent does not compile at my side. Anyway, no big deal.


MatGetSubMatrix/MatGetSubMatrices have been renamed  to MatCreateSubMatrix/MatCreateSubMatrices in petsc-dev. I thought you were using the master branch and not the latest release. Sorry for the confusion. 

To compile the code I have sent, just rename MatCreateSubMatrices with MatGetSubMatrices and it should work. 


The modification you propose as far as I understand is to replace "ISCreateGeneral(PETSC_COMM_WORLD" with "ISCreateGeneral(PETSC_COMM_SELF" : still not working at my side (empty dirichlet local matrix).
I will try to get that with a MPI matrix (that would contain same data that MatIS : that's what I tried to avoid as this doubles allocations - anyway, no big deal).


In the code, you are already extracting submatrices from MPIAIJ format, not from MATIS. Attached a code that compiles and runs with petsc-3.7.6 


Franck


De: "Stefano Zampini" <[email protected]>
À: "Franck Houssen" <[email protected]>
Cc: "petsc-dev" <[email protected]>, "PETSc users list" <[email protected]>, "petsc-maint" <[email protected]>
Envoyé: Mardi 23 Mai 2017 20:23:49
Objet: Re: [petsc-dev] Using PETSc MatIS, how to get local matrix (= one domain) before and after assembly ?


On May 23, 2017, at 6:34 PM, Franck Houssen <[email protected]> wrote:

OK. I am supposed to destroy the matrix returned by MatISGetMPIXAIJ ?
Yes
Also, my example still not get the final assembled local matrix (the MatCreateSubMatrix returns an empty matrix) but as far as I understand my (global) index set is OK: what did I miss ?

I really doubt you can use the example you have sent. It doesn’t compile, as MatCreateSubMatrix needs an extra argument.
Attached a modified version that does what I guess is what you are looking for (sequential Dirichlet problems on the subdomains).


Franck


De: "Stefano Zampini" <[email protected]>
À: "Franck Houssen" <[email protected]>
Cc: "petsc-dev" <[email protected]>, "PETSc users list" <[email protected]>, "petsc-maint" <[email protected]>
Envoyé: Mardi 23 Mai 2017 13:16:18
Objet: Re: [petsc-dev] Using PETSc MatIS, how to get local matrix (= one domain) before and after assembly ?

MatISGetMPIXAIJ is collective, as it assembles the global operator. To get the matrices you are looking for, you should call MatCreateSubMatrix on the assembled global operator, with the global indices representing the subdomain problem. Each process needs to call both functions

Stefano

Il 23 Mag 2017 11:41, "Franck Houssen" <[email protected]> ha scritto:
I have a 3x3 global matrix made of two overlapping 2x2 local matrix (= diagonal with 1.). Each local matrix correspond to one domain (each domain is delegated to one MPI proc, so, I have 2 MPI procs because I have 2 domains).
This is the simplest possible example: I have two 2x2 (local) diag matrix that overlap so that the global matrix built from them is 1, 2, 1 on the diagonal (local contributions add up in the middle).

Now, I need for each MPI proc to get the assembled local matrix (sometimes called the dirichlet matrix) : this is a local matrix (sequential - not distributed with MPI) that accounts for contribution of neighboring domains (MPI proc).

How to get the local assembled matrix ? MatGetLocalSubMatrix does not work (throw error - see example attached). MatGetSubMatrix returns a MPI distributed matrix, not a local (sequential) one.
  1. My understanding is that MatISGetMPIXAIJ should return a local matrix (sequential AIJ matrix) : the MPI in the name recall that you get the assembled matrix (with contributions from the shared border) from the other MPI processus. Correct ? In my simple example, I replaced MatGetLocalSubMatrix with MatISGetMPIXAIJ : I get a deadlock which was surprising to me... Is MatISGetMPIXAIJ a collective call ?
  2. Supposing this is a collective call (and that point 1 is not correct), I ride up  MatISGetMPIXAIJ before the "if (rank > 0)" : I don't deadlock now, but it seems I get a global matrix which is not the assembled local matrix I am looking for.
  3. I am supposed to destroy the matrix returned by MatISGetMPIXAIJ ? (I believe yes - not sure as AFAIU wording should associate Destroy methods to Create methods)
Franck

The git diff illustrate modifications I tried to add to the initial file attached to this thread:
--- a/matISLocalMat.cpp
+++ b/matISLocalMat.cpp
@@ -31,6 +31,8 @@ int main(int argc,char **argv) {
   MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY); MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
   MatView(A, PETSC_VIEWER_STDOUT_WORLD); PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD); // Diag: 1, 2, 1
 
+  Mat assembledLocalMat;
+  MatISGetMPIXAIJ(A, MAT_INITIAL_MATRIX, &assembledLocalMat);
   if (rank > 0) { // Do not pollute stdout: print only 1 proc
     std::cout << std::endl << "non assembled local matrix:" << std::endl << std::endl;
     Mat nonAssembledLocalMat;
@@ -38,11 +40,10 @@ int main(int argc,char **argv) {
     MatView(nonAssembledLocalMat, PETSC_VIEWER_STDOUT_SELF); // Diag: 1, 1
 
     std::cout << std::endl << "assembled local matrix:" << std::endl << std::endl;
-    Mat assembledLocalMat;
-    IS is; ISCreateGeneral(PETSC_COMM_SELF, localSize, localIdx, PETSC_COPY_VALUES, &is);
-    MatGetLocalSubMatrix(A, is, is, &assembledLocalMat); // KO ?!...
-    MatView(assembledLocalMat, PETSC_VIEWER_STDOUT_SELF); // Would like to get => Diag: 2, 1
+    //IS is; ISCreateGeneral(PETSC_COMM_SELF, localSize, localIdx, PETSC_COPY_VALUES, &is);
+    //MatGetLocalSubMatrix(A, is, is, &assembledLocalMat); // KO ?!...
   }
+  MatView(assembledLocalMat, PETSC_VIEWER_STDOUT_WORLD); // Would like to get => Diag: 2, 1


De: "Stefano Zampini" <[email protected]>
À: "petsc-maint" <[email protected]>
Cc: "petsc-dev" <[email protected]>, "PETSc users list" <[email protected]>, "Franck Houssen" <[email protected]>
Envoyé: Dimanche 21 Mai 2017 22:51:34
Objet: Re: [petsc-dev] Using PETSc MatIS, how to get local matrix (= one domain) before and after assembly ?

To assemble the operator in aij format, use 
MatISGetMPIXAIJ

Il 21 Mag 2017 18:43, "Matthew Knepley" <[email protected]> ha scritto:
On Sun, May 21, 2017 at 11:23 AM, Franck Houssen <[email protected]> wrote:
I have a 3x3 global matrix is built (diag: 1, 2, 1): it's made of 2 overlapping 2x2 local matrix (diag: 1, 1).
Getting non assembled local matrix is OK with MatISGetLocalMat.
How to get assembled local matrix (initial local matrix + neigbhor contributions on the borders) ? (expected result is diag: 2, 1)

You can always use


to get copies, but if you just want to build things, you can use


  Thanks,

     Matt
 
Franck



--
What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.
-- Norbert Wiener

http://www.caam.rice.edu/~mk51/


<matISLocalMat.cpp>


Attachment: out.log
Description: Binary data

Attachment: matISLocalMat.cpp
Description: Binary data

Reply via email to