On 6 April 2016 at 16:20, FRANCAVILLA MATTEO ALESSANDRO <d019...@polito.it> wrote:
> Fantastic, it definitely solved my problem! Thank you! > > Should I always create the vectors with MatCreateVecs when I use them for > matrix-vector multiplications? I used to create vectors with VecCreate, and > it has always worked fine with "standard" matrices (matrices not created > with MATNEST structure). > Yes > > In the example at > http://www.mcs.anl.gov/petsc/petsc-current/src/ksp/ksp/examples/tutorials/ex2.c.html > I read the following: > > "PETSc automatically generates appropriately partitioned matrices and > vectors when MatCreate() and VecCreate() are used with the same > communicator." > This comment also assumes that you use PETSC_DECIDE for the local sizes. So if you create a matrix with global size M x M and a vector of global size M and you (i) use PETSC_DECIDE for the local size for the row/col with MatCreate() and PETSC_DECIDE for the local row for VecCreate(); and (ii) you use the same communicator - then the row partitioning adopted by the matrix will be consistent with the row partitioning for the vector. > It is not clear to me when I can simply use VecCreate (as I've always done > so far), and when to use MatCreateVecs... > You should always use MatCreateVecs(). The implementation of Mat will give you back a consistent vector. Thanks, Dave > > > > On Wed, 6 Apr 2016 14:38:12 +0200 > > Dave May <dave.mayhe...@gmail.com> wrote: > >> On 6 April 2016 at 14:12, FRANCAVILLA MATTEO ALESSANDRO < >> d019...@polito.it> >> wrote: >> >> Thanks Dave, that's exactly what I was looking for, and it (partially) >>> solved my problem. >>> >>> The MatMult operation now works fine, unfortunately not for every number >>> of MPI processes (while the MatMult involving each submatrix works fine >>> for >>> any number of processes): >>> >>> .... >>> matArray(1) = L_discretized >>> matArray(2) = K_discretized >>> matArray(3) = L_MR >>> matArray(4) = K_MR >>> CALL >>> >>> MatCreateNest(PETSC_COMM_WORLD,2,PETSC_NULL_OBJECT,2,PETSC_NULL_OBJECT,matArray,TheMatrix,IERR) >>> CALL MatSetFromOptions(TheMatrix,IERR) >>> ! >>> CALL VecCreate(PETSC_COMM_WORLD,x,IERR) >>> CALL VecSetSizes(x,PETSC_DECIDE,2*NumberOfUnknowns,IERR) >>> CALL VecSetFromOptions(x,IERR) >>> CALL VecSet(x,0.+0.*PETSC_i,IERR) >>> CALL VecSetValues(x,1,0,1.+0.*PETSC_i,INSERT_VALUES,IERR) >>> CALL VecAssemblyBegin(x,IERR) >>> CALL VecAssemblyEnd(x,IERR) >>> CALL VecCreate(PETSC_COMM_WORLD,y,IERR) >>> CALL VecSetSizes(y,PETSC_DECIDE,NumberOfUnknowns+NumberOfMeasures,IERR) >>> CALL VecSetFromOptions(y,IERR) >>> CALL MatMult(TheMatrix,x,y,IERR) >>> .... >>> >>> >>> I guess there is something wrong with the distribution of vectors between >>> the MPI processes. With the specific sizes of my problem everything works >>> fine with 1-2-3-6 processes, otherwise I get "Nonconforming object sizes" >>> error, e.g. with 8 processes: >>> >>> [2]PETSC ERROR: --------------------- Error Message >>> -------------------------------------------------------------- >>> [2]PETSC ERROR: Nonconforming object sizes >>> [2]PETSC ERROR: Mat mat,Vec y: local dim 383 382 >>> [2]PETSC ERROR: Petsc Release Version 3.6.3, Dec, 03, 2015 >>> [2]PETSC ERROR: Nonconforming object sizes >>> ./ISS on a arch-linux2-c-debug named alessandro-HP-ProDesk-490-G2-MT by >>> alessandro Wed Apr 6 13:31:45 2016 >>> [2]PETSC ERROR: Configure options >>> --with-blas-lapack-dir=/opt/intel/composer_xe_2015.2.164/mkl >>> --with-mpi-cc=/home/alessandro/.openmpi-1.8.3/bin/mpicc >>> --with-mpi-f90=/home/alessandro/.openmpi-1.8.3/bin/mpif90 >>> --with-mpiexec=/home/alessandro/.openmpi-1.8.3/bin/mpiexec >>> --with-scalar-type=complex >>> [2]PETSC ERROR: #1 MatMult() line 2216 in >>> /home/alessandro/CODE/PETSc/petsc-3.6.3/src/mat/interface/matrix.c >>> ... >>> [6]PETSC ERROR: Mat mat,Vec y: local dim 381 382 >>> ... >>> [7]PETSC ERROR: Mat mat,Vec y: local dim 381 382 >>> ... >>> [1]PETSC ERROR: Mat mat,Vec y: local dim 383 382 >>> ... >>> >>> >>> I noticed that the code works (runs and provides correct results) when >>> the >>> number of processes divides the length of y2 (number of rows of A21 and >>> A22); however when I apply MatMult separately on A21 and A22 everything >>> works fine with any number of processes >>> >>> Any idea where I'm messing up here? >>> >>> >> Sounds like the large vector you are using to multiple with the matnest >> object wasn't created via MatCreateVecs(). If you use this method you will >> obtain a vector with the correct local and global sizes. >> >> >> >> >> >>> >>> Thanks, >>> Alessandro >>> >>> >>> >>> >>> On Wed, 6 Apr 2016 11:16:43 +0200 >>> Dave May <dave.mayhe...@gmail.com> wrote: >>> >>> On 6 April 2016 at 11:08, FRANCAVILLA MATTEO ALESSANDRO < >>>> d019...@polito.it> >>>> wrote: >>>> >>>> Hi, >>>> >>>>> >>>>> I'm trying to set a linear problem with a 2x2 block matrix (say A=[A11 >>>>> A12; A21 A22] in matlab notation). Obviously I could generate a single >>>>> matrix, but I would like to keep it as a 2x2 block matrix with the 4 >>>>> blocks >>>>> generated and used independently (e.g., because at some point I need to >>>>> change some of those blocks to matrix-free). My plan was then to >>>>> generate a >>>>> MatShell and set a MATOP_MULT shell operation to explicitly compute >>>>> y=A*x >>>>> by splitting x=[x1;x2] and y=[y1;y2] and carrying out the single >>>>> matrix-vector multiplications, and reassembling the final result. >>>>> >>>>> >>>> >>>> Take a look at MatNest - it does exactly what you want. >>>> >>>> >>>> >>>> http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateNest.html >>>> >>>> >>>> >>>> http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MATNEST.html#MATNEST >>>> >>>> Thanks, >>>> Dave >>>> >>>> >>>> My problem is that my efforts in extracting x1 and x2 from the original >>>> x >>>> >>>>> vector have failed. The following Fortran code is what seemed to me the >>>>> most reasonable thing to do to set the matrix-vector multiplication, >>>>> but >>>>> it >>>>> does not work >>>>> >>>>> >>>>> !====================================================== >>>>> subroutine MyMult(A,x,y,IERR) >>>>> >>>>> implicit none >>>>> >>>>> #include <petsc/finclude/petscsys.h> >>>>> #include <petsc/finclude/petscvec.h> >>>>> #include <petsc/finclude/petscvec.h90> >>>>> #include <petsc/finclude/petscmat.h> >>>>> >>>>> Mat A >>>>> Vec x,y, x1, x2, y1, y2 >>>>> IS :: ix1, ix2, iy1, iy2 >>>>> PetscErrorCode ierr >>>>> PetscInt M, N, II, NumberOfUnknowns, NumberOfMeasures >>>>> PetscInt, ALLOCATABLE :: INDICES(:) >>>>> PetscMPIInt :: RANK >>>>> >>>>> CALL VecGetSize(x,N,IERR) >>>>> NumberOfUnknowns = N / 2 >>>>> CALL VecGetSize(y,N,IERR) >>>>> NumberOfMeasures = N - NumberOfUnknowns >>>>> >>>>> N = NumberOfUnknowns + MAX(NumberOfUnknowns,NumberOfMeasures) >>>>> ALLOCATE(INDICES(N)) >>>>> INDICES = (/ (ii, ii=0, N-1) /) >>>>> >>>>> CALL >>>>> >>>>> >>>>> ISCreateGeneral(MPI_COMM_SELF,NumberOfUnknowns,INDICES(1:NumberOfUnknowns),PETSC_COPY_VALUES,ix1,IERR) >>>>> CALL VecGetSubVector(x,ix1,x1,IERR) >>>>> CALL >>>>> >>>>> >>>>> ISCreateGeneral(MPI_COMM_SELF,NumberOfUnknowns,INDICES(NumberOfUnknowns+1:2*NumberOfUnknowns),PETSC_COPY_VALUES,ix2,IERR) >>>>> CALL VecGetSubVector(x,ix2,x2,IERR) >>>>> CALL >>>>> >>>>> >>>>> ISCreateGeneral(MPI_COMM_SELF,NumberOfUnknowns,INDICES(1:NumberOfUnknowns),PETSC_COPY_VALUES,iy1,IERR) >>>>> CALL VecGetSubVector(y,iy1,y1,IERR) >>>>> CALL >>>>> >>>>> >>>>> ISCreateGeneral(MPI_COMM_SELF,NumberOfUnknowns,INDICES(NumberOfUnknowns+1:NumberOfUnknowns+NumberOfMeasures),PETSC_COPY_VALUES,iy2,IERR) >>>>> CALL VecGetSubVector(y,iy2,y2,IERR) >>>>> >>>>> CALL MatMult(L_discretized,x1,y1,IERR) >>>>> CALL MatMultAdd(K_discretized,x2,y1,y1,IERR) >>>>> CALL MatMult(L_MR,x1,y2,IERR) >>>>> CALL MatMultAdd(K_MR,x2,y2,y2,IERR) >>>>> >>>>> CALL VecRestoreSubVector(y,iy1,y1,IERR) >>>>> CALL VecRestoreSubVector(y,iy2,y2,IERR) >>>>> >>>>> return >>>>> end subroutine MyMult >>>>> !====================================================== >>>>> >>>>> >>>>> Obviously the sequence of calls to ISCreateGeneral and VecGetSubVector >>>>> does not do what I expect, as the errors I'm getting are in the >>>>> following >>>>> MatMult multiplications (the global sizes of x1 and x2 SHOULD BE both >>>>> 771, >>>>> while y1 and y2 should be 771 and 2286): >>>>> >>>>> 1) executed with 1 MPI process: >>>>> >>>>> [0]PETSC ERROR: --------------------- Error Message >>>>> -------------------------------------------------------------- >>>>> [0]PETSC ERROR: Nonconforming object sizes >>>>> [0]PETSC ERROR: Mat mat,Vec y: global dim 2286 771 >>>>> [0]PETSC ERROR: See >>>>> http://www.mcs.anl.gov/petsc/documentation/faq.html >>>>> for trouble shooting. >>>>> [0]PETSC ERROR: Petsc Release Version 3.6.3, Dec, 03, 2015 >>>>> [0]PETSC ERROR: ./ISS on a arch-linux2-c-debug named >>>>> alessandro-HP-ProDesk-490-G2-MT by alessandro Wed Apr 6 10:54:03 2016 >>>>> [0]PETSC ERROR: Configure options >>>>> --with-blas-lapack-dir=/opt/intel/composer_xe_2015.2.164/mkl >>>>> --with-mpi-cc=/home/alessandro/.openmpi-1.8.3/bin/mpicc >>>>> --with-mpi-f90=/home/alessandro/.openmpi-1.8.3/bin/mpif90 >>>>> --with-mpiexec=/home/alessandro/.openmpi-1.8.3/bin/mpiexec >>>>> --with-scalar-type=complex >>>>> [0]PETSC ERROR: #1 MatMult() line 2215 in >>>>> /home/alessandro/CODE/PETSc/petsc-3.6.3/src/mat/interface/matrix.c >>>>> [0]PETSC ERROR: --------------------- Error Message >>>>> -------------------------------------------------------------- >>>>> [0]PETSC ERROR: Nonconforming object sizes >>>>> [0]PETSC ERROR: Mat mat,Vec v3: local dim 2286 771 >>>>> [0]PETSC ERROR: See >>>>> http://www.mcs.anl.gov/petsc/documentation/faq.html >>>>> for trouble shooting. >>>>> [0]PETSC ERROR: Petsc Release Version 3.6.3, Dec, 03, 2015 >>>>> [0]PETSC ERROR: ./ISS on a arch-linux2-c-debug named >>>>> alessandro-HP-ProDesk-490-G2-MT by alessandro Wed Apr 6 10:54:03 2016 >>>>> [0]PETSC ERROR: Configure options >>>>> --with-blas-lapack-dir=/opt/intel/composer_xe_2015.2.164/mkl >>>>> --with-mpi-cc=/home/alessandro/.openmpi-1.8.3/bin/mpicc >>>>> --with-mpi-f90=/home/alessandro/.openmpi-1.8.3/bin/mpif90 >>>>> --with-mpiexec=/home/alessandro/.openmpi-1.8.3/bin/mpiexec >>>>> --with-scalar-type=complex >>>>> [0]PETSC ERROR: #2 MatMultAdd() line 2396 in >>>>> /home/alessandro/CODE/PETSc/petsc-3.6.3/src/mat/interface/matrix.c >>>>> >>>>> >>>>> >>>>> 2) executed with 4 MPI processes (I'm copying the output of process 0 >>>>> only, they're all the same): >>>>> >>>>> .... >>>>> [0]PETSC ERROR: --------------------- Error Message >>>>> -------------------------------------------------------------- >>>>> [0]PETSC ERROR: Nonconforming object sizes >>>>> [0]PETSC ERROR: Mat mat,Vec x: global dim 771 3084 >>>>> [0]PETSC ERROR: See >>>>> http://www.mcs.anl.gov/petsc/documentation/faq.html >>>>> for trouble shooting. >>>>> [0]PETSC ERROR: Petsc Release Version 3.6.3, Dec, 03, 2015 >>>>> [0]PETSC ERROR: ./ISS on a arch-linux2-c-debug named >>>>> alessandro-HP-ProDesk-490-G2-MT by alessandro Wed Apr 6 10:49:13 2016 >>>>> [0]PETSC ERROR: Configure options >>>>> --with-blas-lapack-dir=/opt/intel/composer_xe_2015.2.164/mkl >>>>> --with-mpi-cc=/home/alessandro/.openmpi-1.8.3/bin/mpicc >>>>> --with-mpi-f90=/home/alessandro/.openmpi-1.8.3/bin/mpif90 >>>>> --with-mpiexec=/home/alessandro/.openmpi-1.8.3/bin/mpiexec >>>>> --with-scalar-type=complex >>>>> [0]PETSC ERROR: #1 MatMult() line 2214 in >>>>> /home/alessandro/CODE/PETSc/petsc-3.6.3/src/mat/interface/matrix.c >>>>> [0]PETSC ERROR: --------------------- Error Message >>>>> -------------------------------------------------------------- >>>>> [0]PETSC ERROR: Nonconforming object sizes >>>>> [0]PETSC ERROR: Mat mat,Vec v1: global dim 771 3084 >>>>> [0]PETSC ERROR: See >>>>> http://www.mcs.anl.gov/petsc/documentation/faq.html >>>>> for trouble shooting. >>>>> [0]PETSC ERROR: Petsc Release Version 3.6.3, Dec, 03, 2015 >>>>> [0]PETSC ERROR: ./ISS on a arch-linux2-c-debug named >>>>> alessandro-HP-ProDesk-490-G2-MT by alessandro Wed Apr 6 10:49:13 2016 >>>>> [0]PETSC ERROR: Configure options >>>>> --with-blas-lapack-dir=/opt/intel/composer_xe_2015.2.164/mkl >>>>> --with-mpi-cc=/home/alessandro/.openmpi-1.8.3/bin/mpicc >>>>> --with-mpi-f90=/home/alessandro/.openmpi-1.8.3/bin/mpif90 >>>>> --with-mpiexec=/home/alessandro/.openmpi-1.8.3/bin/mpiexec >>>>> --with-scalar-type=complex >>>>> [0]PETSC ERROR: #2 MatMultAdd() line 2393 in >>>>> /home/alessandro/CODE/PETSc/petsc-3.6.3/src/mat/interface/matrix.c >>>>> [0]PETSC ERROR: --------------------- Error Message >>>>> -------------------------------------------------------------- >>>>> [0]PETSC ERROR: Nonconforming object sizes >>>>> [0]PETSC ERROR: Mat mat,Vec x: global dim 771 3084 >>>>> [0]PETSC ERROR: See >>>>> http://www.mcs.anl.gov/petsc/documentation/faq.html >>>>> for trouble shooting. >>>>> [0]PETSC ERROR: Petsc Release Version 3.6.3, Dec, 03, 2015 >>>>> [0]PETSC ERROR: ./ISS on a arch-linux2-c-debug named >>>>> alessandro-HP-ProDesk-490-G2-MT by alessandro Wed Apr 6 10:49:13 2016 >>>>> [0]PETSC ERROR: Configure options >>>>> --with-blas-lapack-dir=/opt/intel/composer_xe_2015.2.164/mkl >>>>> --with-mpi-cc=/home/alessandro/.openmpi-1.8.3/bin/mpicc >>>>> --with-mpi-f90=/home/alessandro/.openmpi-1.8.3/bin/mpif90 >>>>> --with-mpiexec=/home/alessandro/.openmpi-1.8.3/bin/mpiexec >>>>> --with-scalar-type=complex >>>>> [0]PETSC ERROR: #3 MatMult() line 2214 in >>>>> /home/alessandro/CODE/PETSc/petsc-3.6.3/src/mat/interface/matrix.c >>>>> [0]PETSC ERROR: --------------------- Error Message >>>>> -------------------------------------------------------------- >>>>> [0]PETSC ERROR: Nonconforming object sizes >>>>> [0]PETSC ERROR: Mat mat,Vec v1: global dim 771 3084 >>>>> [0]PETSC ERROR: See >>>>> http://www.mcs.anl.gov/petsc/documentation/faq.html >>>>> for trouble shooting. >>>>> [0]PETSC ERROR: Petsc Release Version 3.6.3, Dec, 03, 2015 >>>>> [0]PETSC ERROR: ./ISS on a arch-linux2-c-debug named >>>>> alessandro-HP-ProDesk-490-G2-MT by alessandro Wed Apr 6 10:49:13 2016 >>>>> [0]PETSC ERROR: Configure options >>>>> --with-blas-lapack-dir=/opt/intel/composer_xe_2015.2.164/mkl >>>>> --with-mpi-cc=/home/alessandro/.openmpi-1.8.3/bin/mpicc >>>>> --with-mpi-f90=/home/alessandro/.openmpi-1.8.3/bin/mpif90 >>>>> --with-mpiexec=/home/alessandro/.openmpi-1.8.3/bin/mpiexec >>>>> --with-scalar-type=complex >>>>> [0]PETSC ERROR: #4 MatMultAdd() line 2393 in >>>>> /home/alessandro/CODE/PETSc/petsc-3.6.3/src/mat/interface/matrix.c >>>>> .... >>>>> >>>>> >>>>> Any idea? I would like to keep the flexibility of having a sequence of >>>>> mat-vec multiplications as in >>>>> >>>>> CALL MatMult(L_discretized,x1,y1,IERR) >>>>> CALL MatMultAdd(K_discretized,x2,y1,y1,IERR) >>>>> CALL MatMult(L_MR,x1,y2,IERR) >>>>> CALL MatMultAdd(K_MR,x2,y2,y2,IERR) >>>>> >>>>> but for doing that I need to be able to extract x1 and x2, and to >>>>> reassemble the partial results into y. >>>>> >>>>> Thanks in advance, >>>>> Alessandro >>>>> >>>>> >>>>> ___________________________________ >>> Matteo Alessandro Francavilla, PhD >>> Antenna and EMC Lab (LACE) >>> Istituto Superiore Mario Boella (ISMB) >>> Torino, Italy >>> phone +39 011 2276 711 >>> >>> > ___________________________________ > Matteo Alessandro Francavilla, PhD > Antenna and EMC Lab (LACE) > Istituto Superiore Mario Boella (ISMB) > Torino, Italy > phone +39 011 2276 711 >