Thank you for that fast and complete answer. sáb, 17 de jan de 2015 21:21, Barry Smith <[email protected]> escreveu:
> > As Dave points out the error is not detected because of a confluence of > options you happen to be using; normally a MatMult() in the Krylov methods > triggers an error but in your choices there is no MatMult(). The reason > PCSetUp() and PCSetUp_BJacobi() cannot detect the problem is because they > may not yet have available the vector that you pass into KSPSolve() hence > cannot detect the size difference. > > I have added an error check in PCApply() to make sure the sizes of the > objects passed in are correct; it is in the branch barry/fix-bjacobi-vector- > lengths and next (I learned my lesson and did not stick it immediately > into next and master :-(). I've also attached the patch file. > > Thanks > > Barry > > [0]PETSC ERROR: --------------------- Error Message > -------------------------------------------------------------- > [0]PETSC ERROR: Nonconforming object sizes > [0]PETSC ERROR: Preconditioner number of local rows 20 does not equal > resulting vector number of rows 30 > [0]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html > for trouble shooting. > [0]PETSC ERROR: Petsc Release Version 3.5.2, unknown > [0]PETSC ERROR: ./ex1 on a arch-maint named Barrys-MacBook-Pro-3.local by > barrysmith Sat Jan 17 17:15:12 2015 > [0]PETSC ERROR: Configure options --download-mpich --with-fc=0 > --with-cxx=0 --download-sowing=0 > [0]PETSC ERROR: #1 PCApply() line 439 in /Users/barrysmith/Src/petsc/sr > c/ksp/pc/interface/precon.c > [0]PETSC ERROR: #2 KSP_PCApply() line 230 in /Users/barrysmith/Src/petsc/ > include/petsc-private/kspimpl.h > [0]PETSC ERROR: #3 KSPSolve_PREONLY() line 26 in > /Users/barrysmith/Src/petsc/src/ksp/ksp/impls/preonly/preonly.c > [0]PETSC ERROR: #4 KSPSolve() line 460 in /Users/barrysmith/Src/petsc/sr > c/ksp/ksp/interface/itfunc.c > [0]PETSC ERROR: #5 main() line 59 in /Users/barrysmith/Src/petsc/te > st-dir/ex1.c > [0]PETSC ERROR: ----------------End of Error Message -------send entire > error message to [email protected] > > > On Jan 17, 2015, at 4:13 PM, Dave May <[email protected]> wrote: > > > > Many of the Mat-Vec operations will check that the local and global > sizes of the matrix and vector are compatible, for example MatMult. If you > changed your outer KSP type to anything other than PREONLY, an error about > the mismatched sizes would get thrown. > > > > It seems that none of the methods being called from within > PCSetUp_BJACOBI perform such a check (which is a bit of a surprise). I > cannot tell you exactly what the result of calling PCApply() with BJACOBI > will be in the case when the local/global sizes of x,y are larger than the > those for the operator. (I would guess that the "extra" entries might be > ignored by PCApply_BJACOBI, but without going through the function in > detail this really is a guess). > > > > The call to PCApply() does not check that the local/global sizes of the > input and output vectors are compatible with the sizes associated with the > preconditioned operator. It seems neither does PCSetUp(). > > > > Shouldn't PCSetUp() perform such a check to detect this error early on > and to not have the in-built assumption that a method called within > PCSetUp_XXX() will in fact catch the error?? > > > > > > Cheers, > > Dave > > > > > > > > On 17 January 2015 at 22:42, Michael Souza <[email protected]> > wrote: > > I'm not sure if this is a bug, but the PCJACOBI can be applied in a > vector with size different from its operators. > > In the code below, I define a PCJACOBI for a matrix A with size equal to > 20 and I apply it in a 30 length vector x. > > > > Is this just an unforeseen mistake (unexpected application)? If it > isn't, what exactly does PCJACOBI do in this situation? > > > > Cheers, > > Michael > > > > ------------------------------------------------------------ > ------------------------- > > int main(int argc, char **args){ > > PetscErrorCode ierr; > > PC pc; > > KSP ksp; > > MPI_Comm comm; > > Mat A; > > PetscInt i, na=20, nx=30; > > Vec x, y; > > > > ierr = PetscInitialize(&argc, &args, (char *) 0, help); > CHKERRQ(ierr); > > > > // matrix creation and setup > > ierr = MatCreate(PETSC_COMM_WORLD, &A); CHKERRQ(ierr); > > ierr = MatSetType(A, MATMPIAIJ); CHKERRQ(ierr); > > ierr = MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, na, na); > CHKERRQ(ierr); > > ierr = MatSetUp(A); CHKERRQ(ierr); > > for (i = 0; i < na; i++) { > > if(i > 0) > > ierr = MatSetValue(A, i, (i - 1), 1.0, INSERT_VALUES); > CHKERRQ(ierr); > > if(i < (na-1)) > > ierr = MatSetValue(A, i, i+1, 1.0, INSERT_VALUES); > CHKERRQ(ierr); > > ierr = MatSetValue(A, i, i, -2.0, INSERT_VALUES); CHKERRQ(ierr); > > } > > ierr = MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY); CHKERRQ(ierr); > > ierr = MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY); CHKERRQ(ierr); > > > > // vector creation and setup > > ierr = VecCreate(PETSC_COMM_WORLD, &x); CHKERRQ(ierr); > > ierr = VecSetSizes(x,PETSC_DECIDE,nx);CHKERRQ(ierr); > > ierr = VecSetFromOptions(x);CHKERRQ(ierr); > > ierr = VecSet(x,1.0);CHKERRQ(ierr); > > ierr = VecDuplicate(x,&y);CHKERRQ(ierr); > > > > // set KSP > > ierr = PetscObjectGetComm((PetscObject) A, &comm); CHKERRQ(ierr); > > ierr = KSPCreate(comm, &ksp); CHKERRQ(ierr); > > ierr = KSPSetType(ksp, KSPPREONLY); CHKERRQ(ierr); > > ierr = KSPSetOperators(ksp, A, A); CHKERRQ(ierr); > > ierr = KSPGetPC(ksp, &pc); CHKERRQ(ierr); > > ierr = PCSetType(pc, PCBJACOBI); CHKERRQ(ierr); > > ierr = PCSetUp(pc); CHKERRQ(ierr); > > { // setup the SubKSP > > PetscInt nlocal, dummy; > > KSP *subksp; > > PC subpc; > > ierr = PCBJacobiGetSubKSP(pc, &nlocal, &dummy, &subksp); > CHKERRQ(ierr); > > ierr = KSPSetType(*subksp, KSPPREONLY); CHKERRQ(ierr); > > ierr = KSPGetPC(*subksp, &subpc); CHKERRQ(ierr); > > ierr = PCSetType(subpc, PCILU); CHKERRQ(ierr); > > ierr = PCFactorSetLevels(subpc, 1); CHKERRQ(ierr); > > ierr = PCSetFromOptions(subpc); CHKERRQ(ierr); > > ierr = PCSetUp(subpc); CHKERRQ(ierr); > > } > > > > // solve > > ierr = VecView(x,PETSC_VIEWER_STDOUT_WORLD); > > ierr = KSPSolve(ksp,x,y); CHKERRQ(ierr); > > ierr = VecView(y,PETSC_VIEWER_STDOUT_WORLD); CHKERRQ(ierr); > > > > // destroy objects > > ierr = MatDestroy(&A); CHKERRQ(ierr); > > ierr = VecDestroy(&x); CHKERRQ(ierr); > > ierr = VecDestroy(&y); CHKERRQ(ierr); > > ierr = KSPDestroy(&ksp); CHKERRQ(ierr); > > > > return EXIT_SUCCESS; > > } > > > >
