[petsc-users] PCJacobi applied to incompatible vector
Michael Souza
souza.michael at gmail.com
Sat Jan 17 15:42:28 CST 2015
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;
}
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20150117/3783e72d/attachment-0001.html>
More information about the petsc-users
mailing list