Thanks for the answer.
Emmanuel:
> You can create a dense C with the required parallel layout without
> calling MatAssemblyBegin() and MatAssemblyEnd().
> Did you get error without calling these routines?
>
Yes, the output is (after create the C dense matrix and do not assembly it,
run1 - see attached file -):
[0]PETSC ERROR: --------------------- Error Message
--------------------------------------------------------------
[0]PETSC ERROR: Object is in wrong state
[0]PETSC ERROR: Not for unassembled matrix
[0]PETSC ERROR: See https://www.mcs.anl.gov/petsc/documentation/faq.html
for trouble shooting.
[0]PETSC ERROR: Petsc Release Version 3.12.4, Feb, 04, 2020
[0]PETSC ERROR: ./comp on a arch-linux-c-opt-O2 named ayala by ayala Thu
Feb 27 16:47:15 2020
[0]PETSC ERROR: Configure options --with-debugging=0 COPTFLAGS="-O2
-march=native -mtune=native" CXXOPTFLAGS="-O2 -march=native -mtune=native"
FOPTFLAGS="-O2 -march=native -mtune=native" --download-mpich=1
--download-fblaslapack=1 --with-cxx-dialect=C++11
[0]PETSC ERROR: #1 MatNorm() line 5123 in
/home/ayala/Documents/PETSc/petsc-3.12.4/src/mat/interface/matrix.c
We only updated the help manu, not internal implementation. In the next
> release, we'll introduce new set of API to consolidate the API of
> mat-mat-operations.
> Hong
>
I attach my test file, or maybe I'm doing something wrong. I tested this
file on my laptop ubuntu 18
Kind regards.
#include <petsc.h>
static char help[] = "MAT MAT MULT \n";
/*
This perform the matrix multiplication MatMatMult(A,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&C);
where A is sparse and B is dense, then C must be created with proper size,
then it's possible to use MAT_REUSE_MATRIX
*/
int main(int argc, char *argv[])
{
PetscErrorCode ierr;
PetscInitialize(&argc,&argv,PETSC_NULL,help);
Mat A, B, C;
PetscInt sz_rows = 1000, sz_cols = 1000, rows, cols, small = 12; // small is a small number of columns
PetscRandom rctx1, rctx2;
PetscScalar norma;
PetscBool create = PETSC_FALSE, assem = PETSC_FALSE;
PetscOptionsGetInt(NULL,NULL,"-rows",&sz_rows,NULL);
PetscOptionsGetInt(NULL,NULL,"-cols",&sz_cols,NULL);
PetscOptionsGetInt(NULL,NULL,"-small",&small,NULL);
PetscOptionsGetBool(NULL,NULL,"-create",&create,NULL);
PetscOptionsGetBool(NULL,NULL,"-assem",&assem,NULL);
PetscRandomCreate(PETSC_COMM_WORLD,&rctx1);
PetscRandomCreate(PETSC_COMM_WORLD,&rctx2);
ierr = MatCreate(PETSC_COMM_WORLD,&A); CHKERRQ(ierr);
ierr = MatSetType(A,MATAIJ);//by default use MATMPIAIJ when matirx is created
ierr = MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, sz_rows, sz_cols); CHKERRQ(ierr);
ierr = MatSetFromOptions(A); CHKERRQ(ierr);
ierr = MatSetUp(A); CHKERRQ(ierr);
MatSetRandom(A,rctx1);
ierr = MatCreate(PETSC_COMM_WORLD,&B); CHKERRQ(ierr);
ierr = MatSetType(B,MATDENSE);//by default use MATMPIAIJ when matirx is created
ierr = MatSetSizes(B, PETSC_DECIDE, PETSC_DECIDE, sz_cols, small); CHKERRQ(ierr);
ierr = MatSetFromOptions(B); CHKERRQ(ierr);
ierr = MatSetUp(B); CHKERRQ(ierr);
MatSetRandom(B,rctx2);
if (create)
{
ierr = MatCreate(PETSC_COMM_WORLD,&C); CHKERRQ(ierr);
ierr = MatSetType(C,MATDENSE);
ierr = MatSetSizes(C, PETSC_DECIDE, PETSC_DECIDE, sz_rows, small); CHKERRQ(ierr);
ierr = MatSetFromOptions(C); CHKERRQ(ierr);
ierr = MatSetUp(C); CHKERRQ(ierr);
MatZeroEntries(C);
}
MatGetSize(A,&rows,&cols);
PetscPrintf(PETSC_COMM_WORLD,"A: rows = %d, cols = %d\n",rows,cols);
MatGetSize(B,&rows,&cols);
PetscPrintf(PETSC_COMM_WORLD,"B: rows = %d, cols = %d\n",rows,cols);
if(create)
MatMatMult(A,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&C);
else
MatMatMult(A,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C);
if(assem)
{
ierr = MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY) ;CHKERRQ(ierr);
ierr = MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY) ;CHKERRQ(ierr);
}
MatGetSize(C,&rows,&cols);
PetscPrintf(PETSC_COMM_WORLD,"C: rows = %d, cols = %d\n",rows,cols);
MatNorm(A,NORM_INFINITY,&norma);
PetscPrintf(PETSC_COMM_WORLD,"A: norma = %f\n",norma);
MatNorm(B,NORM_INFINITY,&norma);
PetscPrintf(PETSC_COMM_WORLD,"B: norma = %f\n",norma);
MatNorm(C,NORM_INFINITY,&norma);
PetscPrintf(PETSC_COMM_WORLD,"C: norma = %f\n",norma);
ierr = PetscFinalize(); CHKERRQ(ierr);
return 0;
}
/*
run:
${MPIEXEC} -np 5 ./comp -rows 10000 -cols 10000 -small 12 -create 1 -assem 1
run1:
${MPIEXEC} -np 5 ./comp -rows 10000 -cols 10000 -small 12 -create 1 -assem 0
run2:
${MPIEXEC} -np 5 ./comp -rows 10000 -cols 10000 -small 12 -create 0 -assem 1
run3:
${MPIEXEC} -np 5 ./comp -rows 10000 -cols 10000 -small 12 -create 0 -assem 0
*/