Hi:
I have a simple code myself to learn the basic usage of nest matrix. It
assembles a 3x3 nest matrix and solve with that nest matrix. I can run it
without a preconditioner. Running the code with -pc_type fieldsplit just
leads to an error message. I have no clue of what goes wrong. The message
is in below. I also have my code attached.
Thanks,
M
[0]PETSC ERROR: PetscTrFreeDefault() called from VecRestoreArray_Nest()
line 678 in /home/lib/petsc-3.9.1-src/src/vec/vec/impls/nest/vecnest.c
[0]PETSC ERROR: Block at address 0x2784ec8 is corrupted; cannot free;
may be block not allocated with PetscMalloc()
[0]PETSC ERROR: --------------------- Error Message
--------------------------------------------------------------
[0]PETSC ERROR: Memory corruption:
http://www.mcs.anl.gov/petsc/documentation/installation.html#valgrind
[0]PETSC ERROR: Bad location or corrupted memory
[0]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html for
trouble shooting.
[0]PETSC ERROR: Petsc Release Version 3.9.1, Apr, 29, 2018
[0]PETSC ERROR: ./matnest on a named bacon by mike Tue Aug 14 10:03:56 2018
[0]PETSC ERROR: Configure options --with-x=0 -with-pic --with-make-np=12
--with-mpi-compilers=1 --with-mpi-dir=/home/lib/mpich-3.2/
--with-scalar-type=real --with-shared-libraries=1 --with-precision=double
--with-chaco=1 --download-chaco --with-hypre=1 --download-hypre
--with-plapack=1 --download-plapack --with-spai=1 --download-spai
--with-sundials=1 --download-sundials --with-mumps=1 --download-mumps
--with-scalapack=1 --download-scalapack --with-blacs=1 --download-blacs
--with-spooles=1 --download-spooles --with-suitesparse=1
--download-suitesparse --with-superlu_dist=1 --download-superlu_dist
--with-superlu=1 --download-superlu --with-zoltan=1 --download-zoltan=1
--with-debugging=yes --download-fblaslapack --with-ml=1 --download-ml
--with-eigen=1 --download-eigen --prefix=~/lib/petsc-3.9.1-debug
[0]PETSC ERROR: #1 PetscTrFreeDefault() line 269 in
/home/lib/petsc-3.9.1-src/src/sys/memory/mtr.c
[0]PETSC ERROR: #2 VecRestoreArray_Nest() line 678 in
/home/lib/petsc-3.9.1-src/src/vec/vec/impls/nest/vecnest.c
[0]PETSC ERROR: #3 VecRestoreArrayRead() line 1835 in
/home/lib/petsc-3.9.1-src/src/vec/vec/interface/rvector.c
[0]PETSC ERROR: #4 VecRestoreArrayPair() line 511 in
/home/lib/petsc-3.9.1-src/include/petscvec.h
[0]PETSC ERROR: #5 VecScatterBegin_SSToSS() line 671 in
/home/lib/petsc-3.9.1-src/src/vec/vscat/impls/vscat.c
[0]PETSC ERROR: #6 VecScatterBegin() line 1779 in
/home/lib/petsc-3.9.1-src/src/vec/vscat/impls/vscat.c
[0]PETSC ERROR: #7 PCApply_FieldSplit() line 1040 in
/home/lib/petsc-3.9.1-src/src/ksp/pc/impls/fieldsplit/fieldsplit.c
[0]PETSC ERROR: #8 PCApply() line 457 in
/home/lib/petsc-3.9.1-src/src/ksp/pc/interface/precon.c
[0]PETSC ERROR: #9 KSP_PCApply() line 276 in
/home/lib/petsc-3.9.1-src/include/petsc/private/kspimpl.h
[0]PETSC ERROR: #10 KSPInitialResidual() line 67 in
/home/lib/petsc-3.9.1-src/src/ksp/ksp/interface/itres.c
[0]PETSC ERROR: #11 KSPSolve_GMRES() line 233 in
/home/lib/petsc-3.9.1-src/src/ksp/ksp/impls/gmres/gmres.c
[0]PETSC ERROR: #12 KSPSolve() line 669 in
/home/lib/petsc-3.9.1-src/src/ksp/ksp/interface/itfunc.c
// This is a study code for PETSc NestMat data structure
// Aug 5 2018
#include "PETSc_Tools.hpp"
#include "petscksp.h"
int main(int argc, char **argv)
{
PetscInitialize(&argc,&argv,(char*)0, PETSC_NULL);
Mat A; // block matrix
Mat subA[9]; // sub block matrices
IS isg[3]; // index sets
PetscInt row, col, mstart, mend;
PetscScalar val;
// subA 0
MatCreateAIJ( PETSC_COMM_WORLD, 5, 5, PETSC_DETERMINE, PETSC_DETERMINE,
3, NULL, 0, NULL, &subA[0] );
MatGetOwnershipRange(subA[0], &mstart, &mend);
for(row = mstart; row < mend; ++row)
{
val = 1.0; col = row;
MatSetValues(subA[0], 1, &row, 1, &col, &val, INSERT_VALUES);
}
MatAssemblyBegin(subA[0], MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(subA[0], MAT_FINAL_ASSEMBLY);
// subA 1
MatCreateAIJ( PETSC_COMM_WORLD, 5, 3, PETSC_DETERMINE, PETSC_DETERMINE,
1, NULL, 1, NULL, &subA[1] );
MatGetOwnershipRange(subA[1], &mstart, &mend);
for(row = mstart; row < mend; ++row)
{
col = 1;
val = 0.0;
MatSetValues(subA[1], 1, &row, 1, &col, &val, INSERT_VALUES);
}
MatAssemblyBegin(subA[1], MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(subA[1], MAT_FINAL_ASSEMBLY);
// subA 2
MatCreateAIJ( PETSC_COMM_WORLD, 5, 4, PETSC_DETERMINE, PETSC_DETERMINE,
1, NULL, 1, NULL, &subA[2] );
MatGetOwnershipRange(subA[2], &mstart, &mend);
for(row = mstart; row < mend; ++row)
{
col = 1;
val = 0.0;
MatSetValues(subA[2], 1, &row, 1, &col, &val, INSERT_VALUES);
}
MatAssemblyBegin(subA[2], MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(subA[2], MAT_FINAL_ASSEMBLY);
// subA 3
MatCreateAIJ( PETSC_COMM_WORLD, 3, 5, PETSC_DETERMINE, PETSC_DETERMINE,
1, NULL, 1, NULL, &subA[3] );
MatGetOwnershipRange(subA[3], &mstart, &mend);
for(row = mstart; row < mend; ++row)
{
col = row;
val = 0.0;
MatSetValues(subA[3], 1, &row, 1, &col, &val, INSERT_VALUES);
}
MatAssemblyBegin(subA[3], MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(subA[3], MAT_FINAL_ASSEMBLY);
// subA 4
MatCreateAIJ( PETSC_COMM_WORLD, 3, 3, PETSC_DETERMINE, PETSC_DETERMINE,
1, NULL, 1, NULL, &subA[4] );
MatGetOwnershipRange(subA[4], &mstart, &mend);
for(row = mstart; row < mend; ++row)
{
col = row;
val = 4.0;
MatSetValues(subA[4], 1, &row, 1, &col, &val, INSERT_VALUES);
}
MatAssemblyBegin(subA[4], MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(subA[4], MAT_FINAL_ASSEMBLY);
// subA 5
MatCreateAIJ( PETSC_COMM_WORLD, 3, 4, PETSC_DETERMINE, PETSC_DETERMINE,
1, NULL, 1, NULL, &subA[5] );
MatGetOwnershipRange(subA[5], &mstart, &mend);
for(row = mstart; row < mend; ++row)
{
col = row;
val = 0.0;
MatSetValues(subA[5], 1, &row, 1, &col, &val, INSERT_VALUES);
}
MatAssemblyBegin(subA[5], MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(subA[5], MAT_FINAL_ASSEMBLY);
// subA 6
MatCreateAIJ( PETSC_COMM_WORLD, 4, 5, PETSC_DETERMINE, PETSC_DETERMINE,
1, NULL, 1, NULL, &subA[6] );
MatGetOwnershipRange(subA[6], &mstart, &mend);
for(row = mstart; row < mend; ++row)
{
col = 2;
val = 0.0;
MatSetValues(subA[6], 1, &row, 1, &col, &val, INSERT_VALUES);
}
MatAssemblyBegin(subA[6], MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(subA[6], MAT_FINAL_ASSEMBLY);
// subA 7
MatCreateAIJ( PETSC_COMM_WORLD, 4, 3, PETSC_DETERMINE, PETSC_DETERMINE,
1, NULL, 1, NULL, &subA[7] );
MatGetOwnershipRange(subA[7], &mstart, &mend);
for(row = mstart; row < mend; ++row)
{
col = 1;
val = 0.0;
MatSetValues(subA[7], 1, &row, 1, &col, &val, INSERT_VALUES);
}
MatAssemblyBegin(subA[7], MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(subA[7], MAT_FINAL_ASSEMBLY);
// subA 8
MatCreateAIJ( PETSC_COMM_WORLD, 4, 4, PETSC_DETERMINE, PETSC_DETERMINE,
1, NULL, 1, NULL, &subA[8] );
MatGetOwnershipRange(subA[8], &mstart, &mend);
for(row = mstart; row < mend; ++row)
{
col = row;
val = 8.0;
MatSetValues(subA[8], 1, &row, 1, &col, &val, INSERT_VALUES);
}
MatAssemblyBegin(subA[8], MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(subA[8], MAT_FINAL_ASSEMBLY);
// Create A the big nest matrix
MatCreateNest(PETSC_COMM_WORLD, 3, NULL, 3, NULL, subA, &A);
// Generate the index sets
MatNestGetISs(A, isg, NULL);
//ISView(isg[0], PETSC_VIEWER_STDOUT_WORLD);
//ISView(isg[1], PETSC_VIEWER_STDOUT_WORLD);
//ISView(isg[2], PETSC_VIEWER_STDOUT_WORLD);
// Create sub vectors
Vec subb[3];
Vec b;
VecCreate(PETSC_COMM_WORLD, &subb[0]);
VecSetSizes(subb[0], 5, PETSC_DECIDE);
VecSetFromOptions(subb[0]);
VecSet(subb[0], 1.0);
VecCreate(PETSC_COMM_WORLD, &subb[1]);
VecSetSizes(subb[1], 3, PETSC_DECIDE);
VecSetFromOptions(subb[1]);
VecSet(subb[1], 2.0);
VecCreate(PETSC_COMM_WORLD, &subb[2]);
VecSetSizes(subb[2], 4, PETSC_DECIDE);
VecSetFromOptions(subb[2]);
VecSet(subb[2], 3.0);
VecCreateNest(PETSC_COMM_WORLD, 3, NULL, subb, &b);
// duplicate
Vec r;
VecDuplicate(b, &r);
VecCopy(b, r);
// mat mult
MatMult(A, b, r);
PETSc_T::print(b);
VecSet(b, 0.0);
PETSc_T::print(r);
KSP ksp;
PC pc;
KSPCreate(PETSC_COMM_WORLD,&ksp);
KSPSetOperators(ksp,A,A);
KSPSetFromOptions(ksp);
KSPGetPC(ksp,&pc);
PCFieldSplitSetIS(pc,"a",isg[0]);
PCFieldSplitSetIS(pc,"b",isg[1]);
PCFieldSplitSetIS(pc,"c",isg[2]);
KSPSolve(ksp,r,b);
KSPView(ksp,PETSC_VIEWER_STDOUT_WORLD);
PETSc_T::print(b);
//VecAXPY(r, -1.0, b);
//MatSetOption(subA[0], MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE);
//MatSetOption(subA[0], MAT_NEW_NONZERO_LOCATIONS, PETSC_FALSE);
//row = 0; col = 4; val = -100.0;
//MatSetValues(subA[0], 1, &row, 1, &col, &val, INSERT_VALUES);
//MatAssemblyBegin(subA[0], MAT_FINAL_ASSEMBLY);
//MatAssemblyEnd(subA[0], MAT_FINAL_ASSEMBLY);
MatDestroy(&subA[0]);
MatDestroy(&subA[1]);
MatDestroy(&subA[2]);
MatDestroy(&subA[3]);
MatDestroy(&subA[4]);
MatDestroy(&subA[5]);
MatDestroy(&subA[6]);
MatDestroy(&subA[7]);
MatDestroy(&subA[8]);
MatDestroy(&A);
VecDestroy(&subb[0]);
VecDestroy(&subb[1]);
VecDestroy(&subb[2]);
VecDestroy(&b);
VecDestroy(&r);
KSPDestroy(&ksp);
PetscFinalize();
return 0;
}
// EOF
