<html><body><div style="font-family: times new roman, new york, times, serif; font-size: 12pt; color: #000000"><div>The problem was difficult to reduce as reducing make things disappear... Luckily, I believe I got it (or at least, it looks "like" the one I "really" have...).<br></div><div><br></div><div>Seems that for square matrix, it works fine for csr and dense matrix. But, If I am not mistaken, it does not for dense rectangular matrix (still OK for csr).<br></div><div><br></div><div>matISCSRDenseSquare.cpp: 2 procs, global 3x3 matrix, each proc adds a 2x2 local matrix in the global matrix.<br></div><div>matISCSRDenseRect.cpp: 2 procs, global <strong>2</strong>x3 matrix, each proc adds a <strong>1</strong>x2 local <strong>vector</strong> in the global matrix.<br></div><div><br></div><div>reminder: running debian/testing with gcc-6.3 + petsc-3.7.6</div><div><br></div><div>Franck<br></div><div><br></div><div>>> mpirun -n 2 ./matISCSRDenseSquare.exe csr; mpirun -n 2 ./matISCSRDenseSquare.exe dense<br>csr<br>csr<br>Mat Object: 2 MPI processes<br> type: is<br> Mat Object: 1 MPI processes<br> type: seqaij<br>row 0: (0, 1.) (1, 0.) <br>row 1: (0, 0.) (1, 1.) <br> Mat Object: 1 MPI processes<br> type: seqaij<br>row 0: (0, 1.) (1, 0.) <br>row 1: (0, 0.) (1, 1.) <br>dense<br>dense<br>Mat Object: 2 MPI processes<br> type: is<br> Mat Object: 1 MPI processes<br> type: seqdense<br>1.0000000000000000e+00 0.0000000000000000e+00 <br>0.0000000000000000e+00 1.0000000000000000e+00 <br> Mat Object: 1 MPI processes<br> type: seqdense<br>1.0000000000000000e+00 0.0000000000000000e+00 <br>0.0000000000000000e+00 1.0000000000000000e+00 <br><br></div><div>>> mpirun -n 2 ./matISCSRDenseRect.exe csr; mpirun -n 2 ./matISCSRDenseRect.exe dense<br>csr<br>csr<br>Mat Object: 2 MPI processes<br> type: is<br> Mat Object: 1 MPI processes<br> type: seqaij<br>row 0: (0, 1.) (1, 0.) <br> Mat Object: 1 MPI processes<br> type: seqaij<br>row 0: (0, 1.) (1, 0.) <br>dense<br>dense<br>[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------<br>[0]PETSC ERROR: Argument out of range<br>[0]PETSC ERROR: [1]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------<br>[1]PETSC ERROR: Argument out of range<br>[1]PETSC ERROR: New nonzero at (0,1) caused a malloc<br>Use MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE) to turn off this check<br>[1]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting.<br>[1]PETSC ERROR: Petsc Release Version 3.7.6, Apr, 24, 2017 <br>[1]PETSC ERROR: ./matISCSRDenseRect.exe on a arch-linux2-c-debug named yoda by fghoussen Mon Jun 19 18:08:58 2017<br>New nonzero at (0,1) caused a malloc<br>Use MatSetOption(A, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE) to turn off this check<br>[0]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting.<br>[0]PETSC ERROR: Petsc Release Version 3.7.6, Apr, 24, 2017 <br>[1]PETSC ERROR: Configure options --prefix=/home/fghoussen/Documents/INRIA/petsc-3.7.6/local --with-mpi=1 --with-pthread=1 --download-f2cblaslapack=yes --download-mumps=yes --download-scalapack=yes --download-superlu=yes --download-suitesparse=yes<br>[1]PETSC ERROR: #1 MatSetValues_MPIAIJ() line 616 in /home/fghoussen/Documents/INRIA/petsc-3.7.6/src/mat/impls/aij/mpi/mpiaij.c<br>[1]PETSC ERROR: #2 MatSetValues() line 1190 in /home/fghoussen/Documents/INRIA/petsc-3.7.6/src/mat/interface/matrix.c<br>[1]PETSC ERROR: #3 MatSetValuesLocal() line 2053 in /home/fghoussen/Documents/INRIA/petsc-3.7.6/src/mat/interface/matrix.c<br>[1]PETSC ERROR: #4 MatISGetMPIXAIJ_IS() line 365 in /home/fghoussen/Documents/INRIA/petsc-3.7.6/src/mat/impls/is/matis.c<br>[1]PETSC ERROR: #5 MatISGetMPIXAIJ() line 437 in /home/fghoussen/Documents/INRIA/petsc-3.7.6/src/mat/impls/is/matis.c<br>[0]PETSC ERROR: ./matISCSRDenseRect.exe on a arch-linux2-c-debug named yoda by fghoussen Mon Jun 19 18:08:58 2017<br>[0]PETSC ERROR: Configure options --prefix=/home/fghoussen/Documents/INRIA/petsc-3.7.6/local --with-mpi=1 --with-pthread=1 --download-f2cblaslapack=yes --download-mumps=yes --download-scalapack=yes --download-superlu=yes --download-suitesparse=yes<br>[0]PETSC ERROR: #1 MatSetValues_MPIAIJ() line 582 in /home/fghoussen/Documents/INRIA/petsc-3.7.6/src/mat/impls/aij/mpi/mpiaij.c<br>[0]PETSC ERROR: #2 MatSetValues() line 1190 in /home/fghoussen/Documents/INRIA/petsc-3.7.6/src/mat/interface/matrix.c<br>[0]PETSC ERROR: #3 MatSetValuesLocal() line 2053 in /home/fghoussen/Documents/INRIA/petsc-3.7.6/src/mat/interface/matrix.c<br>[0]PETSC ERROR: #4 MatISGetMPIXAIJ_IS() line 365 in /home/fghoussen/Documents/INRIA/petsc-3.7.6/src/mat/impls/is/matis.c<br>[0]PETSC ERROR: #5 MatISGetMPIXAIJ() line 437 in /home/fghoussen/Documents/INRIA/petsc-3.7.6/src/mat/impls/is/matis.c<br>Mat Object: 2 MPI processes<br> type: is<br> Mat Object: 1 MPI processes<br> type: seqdense<br>1.0000000000000000e+00 0.0000000000000000e+00 <br> Mat Object: 1 MPI processes<br> type: seqdense<br>1.0000000000000000e+00 0.0000000000000000e+00 <br><br></div><div>>> diff matISCSRDenseSquare.cpp matISCSRDenseRect.cpp<br>3c3<br>< // ~> g++ -o matISCSRDenseSquare.exe matISCSRDenseSquare.cpp -lpetsc -lm; mpirun -n 2 matISCSRDenseSquare.exe<br>---<br>> // ~> g++ -o matISCSRDenseRect.exe matISCSRDenseRect.cpp -lpetsc -lm; mpirun -n 2 matISCSRDenseRect.exe<br>24c24<br>< ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD, 1, 2, localIdx, PETSC_COPY_VALUES, &rmap);<br>---<br>> ISLocalToGlobalMappingCreate(PETSC_COMM_WORLD, 1, 1, &rank, PETSC_COPY_VALUES, &rmap);<br>29c29<br>< MatCreateIS(PETSC_COMM_WORLD, 1, PETSC_DECIDE, PETSC_DECIDE, 3, 3, rmap, cmap, &A);<br>---<br>> MatCreateIS(PETSC_COMM_WORLD, 1, PETSC_DECIDE, PETSC_DECIDE, 2, 3, rmap, cmap, &A);<br>32,33c32,33<br>< if (matType == "csr") {cout << matType << endl; MatCreateSeqAIJ(PETSC_COMM_SELF, 2, 2, 2, NULL, &Aloc);}<br>< else {cout << matType << endl; MatCreateSeqDense(PETSC_COMM_SELF, 2, 2, NULL, &Aloc);}<br>---<br>> if (matType == "csr") {cout << matType << endl; MatCreateSeqAIJ(PETSC_COMM_SELF, 1, 2, 2, NULL, &Aloc);}<br>> else {cout << matType << endl; MatCreateSeqDense(PETSC_COMM_SELF, 1, 2, NULL, &Aloc);}<br>35,36c35,36<br>< PetscScalar localVal[4] = {1., 0., 0., 1.};<br>< MatSetValues(Aloc, 2, localIdx, 2, localIdx, localVal, ADD_VALUES); // Add local 2x2 matrix<br>---<br>> PetscScalar localVal[2] = {1., 0.}; PetscInt oneLocalRow = 0;<br>> MatSetValues(Aloc, 1, &oneLocalRow, 2, localIdx, localVal, ADD_VALUES); // Add local row<br><br></div><hr id="zwchr"><blockquote style="border-left:2px solid #1010FF;margin-left:5px;padding-left:5px;color:#000;font-weight:normal;font-style:normal;text-decoration:none;font-family:Helvetica,Arial,sans-serif;font-size:12pt;" data-mce-style="border-left: 2px solid #1010FF; margin-left: 5px; padding-left: 5px; color: #000; font-weight: normal; font-style: normal; text-decoration: none; font-family: Helvetica,Arial,sans-serif; font-size: 12pt;"><b>De: </b>"Stefano Zampini" <stefano.zampini@gmail.com><br><b>À: </b>"Franck Houssen" <franck.houssen@inria.fr><br><b>Cc: </b>"PETSc users list" <petsc-users@mcs.anl.gov><br><b>Envoyé: </b>Lundi 19 Juin 2017 15:25:35<br><b>Objet: </b>Re: [petsc-users] Building MatIS with dense local matrix ?<br><div><br></div><div dir="auto">Can you send a minimal working example so that I can fix the code?<div dir="auto"><br></div><div dir="auto">Thanks</div><div dir="auto">Stefano</div></div><div class="gmail_extra"><br><div class="gmail_quote">Il 19 Giu 2017 15:20, "Franck Houssen" <<a href="mailto:franck.houssen@inria.fr" target="_blank" data-mce-href="mailto:franck.houssen@inria.fr">franck.houssen@inria.fr</a>> ha scritto:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex" data-mce-style="margin: 0 0 0 .8ex; border-left: 1px #ccc solid; padding-left: 1ex;"><div><div style="font-family:times new roman,new york,times,serif;font-size:12pt;color:#000000" data-mce-style="font-family: times new roman,new york,times,serif; font-size: 12pt; color: #000000;"><div>Hi,<br></div><div><br></div><div>I try to call MatISGetMPIXAIJ on a MatIS (A) that has been feed locally by sequential (Aloc) dense matrix.</div><div>Seems this ends up with this error: [0]PETSC ERROR: New nonzero at (0,1) caused a malloc. Is this a known error / limitation ? (not supposed to work with dense matrix ?)<br></div><div><br></div><div>This (pseudo code) works fine:<br></div><div>MatCreateIS(..., A)</div><div>MatCreateSeqAIJ(..., Aloc)<br></div><div>MatISSetLocalMat(pcA, pcALoc) <br></div><div>MatISGetMPIXAIJ(A, ...) // OK !<br></div><div><br></div><div>When I try to replace MatCreateSeqAIJ(..., Aloc) with MatCreateSeqDense(..., Aloc), it does no more work. <br></div><div><br></div><div>Franck<br></div><div><br></div><div>PS: running debian/testing with gcc-6.3 + petsc-3.7.6<br></div></div></div></blockquote></div></div></blockquote><div><br></div></div></body></html>