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