On 31/03/2011 20:05, Barry Smith wrote: > On Mar 31, 2011, at 12:38 PM, gouarin wrote: > >> On 31/03/2011 19:06, Matthew Knepley wrote: >>> On Thu, Mar 31, 2011 at 10:10 AM, gouarin<loic.gouarin at math.u-psud.fr> >>> wrote: >>> Thanks for the reply. >>> >>> >>> On 31/03/2011 16:48, Hong Zhang wrote: >>> Loic, >>> >>> I want to do the same for MATMPI. This is what I do: >>> >>> MatGetLocalMat(A,MAT_INITIAL_MATRIX,&Aloc); >>> MatGetRowIJ(Aloc, 1, PETSC_FALSE, PETSC_FALSE,&(Acsr.n),&(Acsr.ia), >>> &(Acsr.ja),&done); >>> MatGetArray(Aloc,&(Acsr.a)); >>> >>> Now, I need to know on what process the non local points are. >>> Each process get its own local submatrix. See >>> http://www.mcs.anl.gov/petsc/petsc-as/snapshots/petsc-current/docs/manualpages/Mat/MatCreateMPIAIJ.html >>> on data structure of MPIAIJ format. Applying your statements above to >>> the example listed on the page (a 8x8 matrix), Proc0 will get >>> Aloc=A(0:2,:), ..., Porc2 gets Aloc=A(6:7,:). >>> Is there a petsc function that gives this partition ? >>> MatGetOwnershipRange() >>> I use this function to know the local points index. In the same way, with >>> this function I know the non local index. >>> >>> But this function doesn't give the process. I give you an example. I have >>> this matrix: >>> >>> 1 0 | 0 1 >>> 1 2 | 1 0 >>> ----------- >>> 0 0 | 1 0 >>> 0 2 | 0 1 >>> >>> I want this kind of format >>> >>> for proc 0 >>> a = {1,1,1,2,1} >>> ia = {0,2,5} >>> ja = {0,3,0,1,2} >>> rank = {0,0,1,1} >>> >>> for proc 1 >>> a = {1,2,1} >>> ia = {0,1,3} >>> ja = {0,2,1} >>> rank = {1,1,0} >>> >>> rank give me the process for each non zero column. Local rows should have >>> numbers 1,...N and offdiagonal entries with j>N. >>> >>> The ownership range is the row partition. However, MPIAIJ matrices are not >>> stored as you indicate above. >>> >>> Different storage formats make the process you describe fragile. WHy do you >>> need the CSR form of the matrix? If you really >>> do, use MatGetSubMatrix() to get a SEQAIJ matrix on each process and pull >>> out the values. >>> >> I want to use this format because I try to create a new multigrid PC with >> this method: >> >> http://homepages.ulb.ac.be/~ynotay/AGMG/userguide/index.html >> >> I don't understand why it's not possible to have the process partition. When >> you create a MPIAIJ matrix and you use a KSP to solve the system, you know >> what values to exchange. No? Where is this information ? > The point is not that the information is not available, it is just that > it is not trivially available since we don't use the data in that format. > > Here is what you have to do. Include "src/mat/impls/aij/mpi/mpiaij.h" > directly in your interface files and access the data directly, do not try to > find some PETSc function that provides exactly what you want directly, that > cannot exist because we cannot know what format each other package wants. In > particular inside Mat_MPIAJ are the two submatrices called A and B which are > Mat_SeqAIJ matrices. The data in > #if defined (PETSC_USE_CTABLE) > PetscTable colmap; > #else > PetscInt *colmap; /* local col number of off-diag col > */ > #endif > PetscInt *garray; /* global index of all off-processor > columns */ > provides the information about the "off process" columns for the matrices > stored in B. > > But you will have to understand the data structures somewhat in order to > get the data in the form you want. > Thanks Matt and Barry. I looked at the interface for hypre to see what you have done and I saw this include file.
I'll try to do that in a simple way with MatGetOwnershipRanges for my first PC version. There are a lot of functions in PETSc and it's sometime difficult to find the good one. Thanks again. Loic > Barry > > > > > >> I can also use MatGetOwnershipRange() and MPI_Allgather to have this >> information. >> >> Loic >> >> >> >> >>> Matt >>> >>> I suppose that I have no information about the structure of MPIAIJ. >>> >>> Loic >>> >>> >>> >>> >>> -- >>> 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
