1dba47a55SKris Buschelman 24b9ad928SBarry Smith /* 33f457be1SHong Zhang This file defines a "solve the problem redundantly on each subgroup of processor" preconditioner. 44b9ad928SBarry Smith */ 5af0996ceSBarry Smith #include <petsc/private/pcimpl.h> 6c6db04a5SJed Brown #include <petscksp.h> /*I "petscksp.h" I*/ 74b9ad928SBarry Smith 84b9ad928SBarry Smith typedef struct { 93e065800SHong Zhang KSP ksp; 104b9ad928SBarry Smith PC pc; /* actual preconditioner used on each processor */ 11ce94432eSBarry Smith Vec xsub,ysub; /* vectors of a subcommunicator to hold parallel vectors of PetscObjectComm((PetscObject)pc) */ 123f457be1SHong Zhang Vec xdup,ydup; /* parallel vector that congregates xsub or ysub facilitating vector scattering */ 13b3804887SHong Zhang Mat pmats; /* matrix and optional preconditioner matrix belong to a subcommunicator */ 143f457be1SHong Zhang VecScatter scatterin,scatterout; /* scatter used to move all values to each processor group (subcommunicator) */ 15ace3abfcSBarry Smith PetscBool useparallelmat; 16c540e29cSHong Zhang PetscSubcomm psubcomm; 171fbd8f88SHong Zhang PetscInt nsubcomm; /* num of data structure PetscSubcomm */ 184b9ad928SBarry Smith } PC_Redundant; 194b9ad928SBarry Smith 204b9ad928SBarry Smith #undef __FUNCT__ 214b9ad928SBarry Smith #define __FUNCT__ "PCView_Redundant" 226849ba73SBarry Smith static PetscErrorCode PCView_Redundant(PC pc,PetscViewer viewer) 234b9ad928SBarry Smith { 244b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 25dfbe8321SBarry Smith PetscErrorCode ierr; 26ace3abfcSBarry Smith PetscBool iascii,isstring; 2703ccd0b4SBarry Smith PetscViewer subviewer; 284b9ad928SBarry Smith 294b9ad928SBarry Smith PetscFunctionBegin; 30251f4c67SDmitry Karpeev ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 31251f4c67SDmitry Karpeev ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);CHKERRQ(ierr); 3232077d6dSBarry Smith if (iascii) { 3303ccd0b4SBarry Smith if (!red->psubcomm) { 3403ccd0b4SBarry Smith ierr = PetscViewerASCIIPrintf(viewer," Redundant preconditioner: Not yet setup\n");CHKERRQ(ierr); 3503ccd0b4SBarry Smith } else { 363e065800SHong Zhang ierr = PetscViewerASCIIPrintf(viewer," Redundant preconditioner: First (color=0) of %D PCs follows\n",red->nsubcomm);CHKERRQ(ierr); 373f08860eSBarry Smith ierr = PetscViewerGetSubViewer(viewer,((PetscObject)red->pc)->comm,&subviewer);CHKERRQ(ierr); 38f5dd71faSBarry Smith if (!red->psubcomm->color) { /* only view first redundant pc */ 391575c14dSBarry Smith ierr = PetscViewerASCIIPushTab(subviewer);CHKERRQ(ierr); 403e065800SHong Zhang ierr = KSPView(red->ksp,subviewer);CHKERRQ(ierr); 411575c14dSBarry Smith ierr = PetscViewerASCIIPopTab(subviewer);CHKERRQ(ierr); 424b9ad928SBarry Smith } 433f08860eSBarry Smith ierr = PetscViewerRestoreSubViewer(viewer,((PetscObject)red->pc)->comm,&subviewer);CHKERRQ(ierr); 444b9ad928SBarry Smith } 4503ccd0b4SBarry Smith } else if (isstring) { 4603ccd0b4SBarry Smith ierr = PetscViewerStringSPrintf(viewer," Redundant solver preconditioner");CHKERRQ(ierr); 474b9ad928SBarry Smith } 484b9ad928SBarry Smith PetscFunctionReturn(0); 494b9ad928SBarry Smith } 504b9ad928SBarry Smith 514b9ad928SBarry Smith #undef __FUNCT__ 524b9ad928SBarry Smith #define __FUNCT__ "PCSetUp_Redundant" 536849ba73SBarry Smith static PetscErrorCode PCSetUp_Redundant(PC pc) 544b9ad928SBarry Smith { 554b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 56dfbe8321SBarry Smith PetscErrorCode ierr; 571b81debcSHong Zhang PetscInt mstart,mend,mlocal,M; 5813f74950SBarry Smith PetscMPIInt size; 59ce94432eSBarry Smith MPI_Comm comm,subcomm; 60ddc54837SHong Zhang Vec x; 613f457be1SHong Zhang 624b9ad928SBarry Smith PetscFunctionBegin; 63ce94432eSBarry Smith ierr = PetscObjectGetComm((PetscObject)pc,&comm);CHKERRQ(ierr); 64ddc54837SHong Zhang 65ddc54837SHong Zhang /* if pmatrix set by user is sequential then we do not need to gather the parallel matrix */ 66ddc54837SHong Zhang ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 67ddc54837SHong Zhang if (size == 1) red->useparallelmat = PETSC_FALSE; 681fbd8f88SHong Zhang 694b9ad928SBarry Smith if (!pc->setupcalled) { 701b81debcSHong Zhang PetscInt mloc_sub; 7175024027SHong Zhang if (!red->psubcomm) { /* create red->psubcomm, new ksp and pc over subcomm */ 7275024027SHong Zhang KSP ksp; 7375024027SHong Zhang ierr = PCRedundantGetKSP(pc,&ksp);CHKERRQ(ierr); 741b81debcSHong Zhang } 7575024027SHong Zhang subcomm = PetscSubcommChild(red->psubcomm); 761fbd8f88SHong Zhang 771b81debcSHong Zhang if (red->useparallelmat) { 781b81debcSHong Zhang /* grab the parallel matrix and put it into processors of a subcomminicator */ 7953cd1579SHong Zhang ierr = MatCreateRedundantMatrix(pc->pmat,red->psubcomm->n,subcomm,MAT_INITIAL_MATRIX,&red->pmats);CHKERRQ(ierr); 80b85f2e9bSHong Zhang 81b85f2e9bSHong Zhang ierr = MPI_Comm_size(subcomm,&size);CHKERRQ(ierr); 82b85f2e9bSHong Zhang if (size > 1) { 83b85f2e9bSHong Zhang PetscBool foundpack; 84b85f2e9bSHong Zhang ierr = MatGetFactorAvailable(red->pmats,NULL,MAT_FACTOR_LU,&foundpack);CHKERRQ(ierr); 85b85f2e9bSHong Zhang if (!foundpack) { /* reset default ksp and pc */ 86b85f2e9bSHong Zhang ierr = KSPSetType(red->ksp,KSPGMRES);CHKERRQ(ierr); 87b85f2e9bSHong Zhang ierr = PCSetType(red->pc,PCBJACOBI);CHKERRQ(ierr); 88*c1619fb6SBarry Smith } else { 89*c1619fb6SBarry Smith ierr = PCFactorSetMatSolverPackage(red->pc,NULL);CHKERRQ(ierr); 90b85f2e9bSHong Zhang } 91b85f2e9bSHong Zhang } 92b85f2e9bSHong Zhang 9323ee1639SBarry Smith ierr = KSPSetOperators(red->ksp,red->pmats,red->pmats);CHKERRQ(ierr); 944b9ad928SBarry Smith 951b81debcSHong Zhang /* get working vectors xsub and ysub */ 962a7a6963SBarry Smith ierr = MatCreateVecs(red->pmats,&red->xsub,&red->ysub);CHKERRQ(ierr); 972fa5cd67SKarl Rupp 988b96b0d2SHong Zhang /* create working vectors xdup and ydup. 998b96b0d2SHong Zhang xdup concatenates all xsub's contigously to form a mpi vector over dupcomm (see PetscSubcommCreate_interlaced()) 1008b96b0d2SHong Zhang ydup concatenates all ysub and has empty local arrays because ysub's arrays will be place into it. 101ce94432eSBarry Smith Note: we use communicator dupcomm, not PetscObjectComm((PetscObject)pc)! */ 1021b81debcSHong Zhang ierr = MatGetLocalSize(red->pmats,&mloc_sub,NULL);CHKERRQ(ierr); 10336be1a5eSBarry Smith ierr = VecCreateMPI(PetscSubcommContiguousParent(red->psubcomm),mloc_sub,PETSC_DECIDE,&red->xdup);CHKERRQ(ierr); 10436be1a5eSBarry Smith ierr = VecCreateMPIWithArray(PetscSubcommContiguousParent(red->psubcomm),1,mloc_sub,PETSC_DECIDE,NULL,&red->ydup);CHKERRQ(ierr); 1053f457be1SHong Zhang 106f68be91cSHong Zhang /* create vecscatters */ 107f68be91cSHong Zhang if (!red->scatterin) { /* efficiency of scatterin is independent from psubcomm_type! */ 1083f457be1SHong Zhang IS is1,is2; 1093f457be1SHong Zhang PetscInt *idx1,*idx2,i,j,k; 11045fc02eaSBarry Smith 1112a7a6963SBarry Smith ierr = MatCreateVecs(pc->pmat,&x,0);CHKERRQ(ierr); 1121b81debcSHong Zhang ierr = VecGetSize(x,&M);CHKERRQ(ierr); 1131b81debcSHong Zhang ierr = VecGetOwnershipRange(x,&mstart,&mend);CHKERRQ(ierr); 1141b81debcSHong Zhang mlocal = mend - mstart; 115dcca6d9dSJed Brown ierr = PetscMalloc2(red->psubcomm->n*mlocal,&idx1,red->psubcomm->n*mlocal,&idx2);CHKERRQ(ierr); 1163f457be1SHong Zhang j = 0; 1171fbd8f88SHong Zhang for (k=0; k<red->psubcomm->n; k++) { 1183f457be1SHong Zhang for (i=mstart; i<mend; i++) { 1193f457be1SHong Zhang idx1[j] = i; 120ddc54837SHong Zhang idx2[j++] = i + M*k; 1213f457be1SHong Zhang } 1223f457be1SHong Zhang } 12370b3c8c7SBarry Smith ierr = ISCreateGeneral(comm,red->psubcomm->n*mlocal,idx1,PETSC_COPY_VALUES,&is1);CHKERRQ(ierr); 12470b3c8c7SBarry Smith ierr = ISCreateGeneral(comm,red->psubcomm->n*mlocal,idx2,PETSC_COPY_VALUES,&is2);CHKERRQ(ierr); 125ddc54837SHong Zhang ierr = VecScatterCreate(x,is1,red->xdup,is2,&red->scatterin);CHKERRQ(ierr); 126fcfd50ebSBarry Smith ierr = ISDestroy(&is1);CHKERRQ(ierr); 127fcfd50ebSBarry Smith ierr = ISDestroy(&is2);CHKERRQ(ierr); 1283f457be1SHong Zhang 1296909748bSHong Zhang /* Impl below is good for PETSC_SUBCOMM_INTERLACED (no inter-process communication) and PETSC_SUBCOMM_CONTIGUOUS (communication within subcomm) */ 130ddc54837SHong Zhang ierr = ISCreateStride(comm,mlocal,mstart+ red->psubcomm->color*M,1,&is1);CHKERRQ(ierr); 1313f457be1SHong Zhang ierr = ISCreateStride(comm,mlocal,mstart,1,&is2);CHKERRQ(ierr); 132ddc54837SHong Zhang ierr = VecScatterCreate(red->xdup,is1,x,is2,&red->scatterout);CHKERRQ(ierr); 133fcfd50ebSBarry Smith ierr = ISDestroy(&is1);CHKERRQ(ierr); 134fcfd50ebSBarry Smith ierr = ISDestroy(&is2);CHKERRQ(ierr); 1351d79065fSBarry Smith ierr = PetscFree2(idx1,idx2);CHKERRQ(ierr); 136ddc54837SHong Zhang ierr = VecDestroy(&x);CHKERRQ(ierr); 1371b81debcSHong Zhang } 138ab661555SHong Zhang } else { /* !red->useparallelmat */ 13923ee1639SBarry Smith ierr = KSPSetOperators(red->ksp,pc->mat,pc->pmat);CHKERRQ(ierr); 1401b81debcSHong Zhang } 141ab661555SHong Zhang } else { /* pc->setupcalled */ 1424b9ad928SBarry Smith if (red->useparallelmat) { 143ab661555SHong Zhang MatReuse reuse; 1441b81debcSHong Zhang /* grab the parallel matrix and put it into processors of a subcomminicator */ 1451b81debcSHong Zhang /*--------------------------------------------------------------------------*/ 146ab661555SHong Zhang if (pc->flag == DIFFERENT_NONZERO_PATTERN) { 1474b9ad928SBarry Smith /* destroy old matrices */ 1486bf464f9SBarry Smith ierr = MatDestroy(&red->pmats);CHKERRQ(ierr); 149ab661555SHong Zhang reuse = MAT_INITIAL_MATRIX; 1504b9ad928SBarry Smith } else { 151ab661555SHong Zhang reuse = MAT_REUSE_MATRIX; 152ab661555SHong Zhang } 153306c2d5bSBarry Smith ierr = MatCreateRedundantMatrix(pc->pmat,red->psubcomm->n,PetscSubcommChild(red->psubcomm),reuse,&red->pmats);CHKERRQ(ierr); 15423ee1639SBarry Smith ierr = KSPSetOperators(red->ksp,red->pmats,red->pmats);CHKERRQ(ierr); 155ab661555SHong Zhang } else { /* !red->useparallelmat */ 15623ee1639SBarry Smith ierr = KSPSetOperators(red->ksp,pc->mat,pc->pmat);CHKERRQ(ierr); 1574b9ad928SBarry Smith } 158ab661555SHong Zhang } 1591b81debcSHong Zhang 1600c24e6a1SHong Zhang if (pc->setfromoptionscalled) { 1613e065800SHong Zhang ierr = KSPSetFromOptions(red->ksp);CHKERRQ(ierr); 1620c24e6a1SHong Zhang } 1633e065800SHong Zhang ierr = KSPSetUp(red->ksp);CHKERRQ(ierr); 1644b9ad928SBarry Smith PetscFunctionReturn(0); 1654b9ad928SBarry Smith } 1664b9ad928SBarry Smith 1674b9ad928SBarry Smith #undef __FUNCT__ 1684b9ad928SBarry Smith #define __FUNCT__ "PCApply_Redundant" 1696849ba73SBarry Smith static PetscErrorCode PCApply_Redundant(PC pc,Vec x,Vec y) 1704b9ad928SBarry Smith { 1714b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 172dfbe8321SBarry Smith PetscErrorCode ierr; 1733f457be1SHong Zhang PetscScalar *array; 1744b9ad928SBarry Smith 1754b9ad928SBarry Smith PetscFunctionBegin; 176ddc54837SHong Zhang if (!red->useparallelmat) { 177ddc54837SHong Zhang ierr = KSPSolve(red->ksp,x,y);CHKERRQ(ierr); 178ddc54837SHong Zhang PetscFunctionReturn(0); 179ddc54837SHong Zhang } 180ddc54837SHong Zhang 1813f457be1SHong Zhang /* scatter x to xdup */ 182ca9f406cSSatish Balay ierr = VecScatterBegin(red->scatterin,x,red->xdup,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 183ca9f406cSSatish Balay ierr = VecScatterEnd(red->scatterin,x,red->xdup,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1843f457be1SHong Zhang 1853f457be1SHong Zhang /* place xdup's local array into xsub */ 1863f457be1SHong Zhang ierr = VecGetArray(red->xdup,&array);CHKERRQ(ierr); 1873f457be1SHong Zhang ierr = VecPlaceArray(red->xsub,(const PetscScalar*)array);CHKERRQ(ierr); 1884b9ad928SBarry Smith 1894b9ad928SBarry Smith /* apply preconditioner on each processor */ 19083ab6a24SBarry Smith ierr = KSPSolve(red->ksp,red->xsub,red->ysub);CHKERRQ(ierr); 1913f457be1SHong Zhang ierr = VecResetArray(red->xsub);CHKERRQ(ierr); 1923f457be1SHong Zhang ierr = VecRestoreArray(red->xdup,&array);CHKERRQ(ierr); 1934b9ad928SBarry Smith 1943f457be1SHong Zhang /* place ysub's local array into ydup */ 1953f457be1SHong Zhang ierr = VecGetArray(red->ysub,&array);CHKERRQ(ierr); 1963f457be1SHong Zhang ierr = VecPlaceArray(red->ydup,(const PetscScalar*)array);CHKERRQ(ierr); 1973f457be1SHong Zhang 1983f457be1SHong Zhang /* scatter ydup to y */ 199ca9f406cSSatish Balay ierr = VecScatterBegin(red->scatterout,red->ydup,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 200ca9f406cSSatish Balay ierr = VecScatterEnd(red->scatterout,red->ydup,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 2013f457be1SHong Zhang ierr = VecResetArray(red->ydup);CHKERRQ(ierr); 2023f457be1SHong Zhang ierr = VecRestoreArray(red->ysub,&array);CHKERRQ(ierr); 2034b9ad928SBarry Smith PetscFunctionReturn(0); 2044b9ad928SBarry Smith } 2054b9ad928SBarry Smith 2064b9ad928SBarry Smith #undef __FUNCT__ 207d88bfacbSStefano Zampini #define __FUNCT__ "PCApplyTranspose_Redundant" 208d88bfacbSStefano Zampini static PetscErrorCode PCApplyTranspose_Redundant(PC pc,Vec x,Vec y) 209d88bfacbSStefano Zampini { 210d88bfacbSStefano Zampini PC_Redundant *red = (PC_Redundant*)pc->data; 211d88bfacbSStefano Zampini PetscErrorCode ierr; 212d88bfacbSStefano Zampini PetscScalar *array; 213d88bfacbSStefano Zampini 214d88bfacbSStefano Zampini PetscFunctionBegin; 215d88bfacbSStefano Zampini if (!red->useparallelmat) { 216d88bfacbSStefano Zampini ierr = KSPSolveTranspose(red->ksp,x,y);CHKERRQ(ierr); 217d88bfacbSStefano Zampini PetscFunctionReturn(0); 218d88bfacbSStefano Zampini } 219d88bfacbSStefano Zampini 220d88bfacbSStefano Zampini /* scatter x to xdup */ 221d88bfacbSStefano Zampini ierr = VecScatterBegin(red->scatterin,x,red->xdup,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 222d88bfacbSStefano Zampini ierr = VecScatterEnd(red->scatterin,x,red->xdup,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 223d88bfacbSStefano Zampini 224d88bfacbSStefano Zampini /* place xdup's local array into xsub */ 225d88bfacbSStefano Zampini ierr = VecGetArray(red->xdup,&array);CHKERRQ(ierr); 226d88bfacbSStefano Zampini ierr = VecPlaceArray(red->xsub,(const PetscScalar*)array);CHKERRQ(ierr); 227d88bfacbSStefano Zampini 228d88bfacbSStefano Zampini /* apply preconditioner on each processor */ 229d88bfacbSStefano Zampini ierr = KSPSolveTranspose(red->ksp,red->xsub,red->ysub);CHKERRQ(ierr); 230d88bfacbSStefano Zampini ierr = VecResetArray(red->xsub);CHKERRQ(ierr); 231d88bfacbSStefano Zampini ierr = VecRestoreArray(red->xdup,&array);CHKERRQ(ierr); 232d88bfacbSStefano Zampini 233d88bfacbSStefano Zampini /* place ysub's local array into ydup */ 234d88bfacbSStefano Zampini ierr = VecGetArray(red->ysub,&array);CHKERRQ(ierr); 235d88bfacbSStefano Zampini ierr = VecPlaceArray(red->ydup,(const PetscScalar*)array);CHKERRQ(ierr); 236d88bfacbSStefano Zampini 237d88bfacbSStefano Zampini /* scatter ydup to y */ 238d88bfacbSStefano Zampini ierr = VecScatterBegin(red->scatterout,red->ydup,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 239d88bfacbSStefano Zampini ierr = VecScatterEnd(red->scatterout,red->ydup,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 240d88bfacbSStefano Zampini ierr = VecResetArray(red->ydup);CHKERRQ(ierr); 241d88bfacbSStefano Zampini ierr = VecRestoreArray(red->ysub,&array);CHKERRQ(ierr); 242d88bfacbSStefano Zampini PetscFunctionReturn(0); 243d88bfacbSStefano Zampini } 244d88bfacbSStefano Zampini 245d88bfacbSStefano Zampini #undef __FUNCT__ 2461ea5a559SBarry Smith #define __FUNCT__ "PCReset_Redundant" 2471ea5a559SBarry Smith static PetscErrorCode PCReset_Redundant(PC pc) 2484b9ad928SBarry Smith { 2494b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 250dfbe8321SBarry Smith PetscErrorCode ierr; 2514b9ad928SBarry Smith 2524b9ad928SBarry Smith PetscFunctionBegin; 2531b81debcSHong Zhang if (red->useparallelmat) { 2546bf464f9SBarry Smith ierr = VecScatterDestroy(&red->scatterin);CHKERRQ(ierr); 2556bf464f9SBarry Smith ierr = VecScatterDestroy(&red->scatterout);CHKERRQ(ierr); 2566bf464f9SBarry Smith ierr = VecDestroy(&red->ysub);CHKERRQ(ierr); 2576bf464f9SBarry Smith ierr = VecDestroy(&red->xsub);CHKERRQ(ierr); 2586bf464f9SBarry Smith ierr = VecDestroy(&red->xdup);CHKERRQ(ierr); 2596bf464f9SBarry Smith ierr = VecDestroy(&red->ydup);CHKERRQ(ierr); 2601b81debcSHong Zhang } 2616bf464f9SBarry Smith ierr = MatDestroy(&red->pmats);CHKERRQ(ierr); 2621b81debcSHong Zhang ierr = KSPReset(red->ksp);CHKERRQ(ierr); 2631ea5a559SBarry Smith PetscFunctionReturn(0); 2641ea5a559SBarry Smith } 2651ea5a559SBarry Smith 2661ea5a559SBarry Smith #undef __FUNCT__ 2671ea5a559SBarry Smith #define __FUNCT__ "PCDestroy_Redundant" 2681ea5a559SBarry Smith static PetscErrorCode PCDestroy_Redundant(PC pc) 2691ea5a559SBarry Smith { 2701ea5a559SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 2711ea5a559SBarry Smith PetscErrorCode ierr; 2721ea5a559SBarry Smith 2731ea5a559SBarry Smith PetscFunctionBegin; 2741ea5a559SBarry Smith ierr = PCReset_Redundant(pc);CHKERRQ(ierr); 2756bf464f9SBarry Smith ierr = KSPDestroy(&red->ksp);CHKERRQ(ierr); 2766bf464f9SBarry Smith ierr = PetscSubcommDestroy(&red->psubcomm);CHKERRQ(ierr); 277c31cb41cSBarry Smith ierr = PetscFree(pc->data);CHKERRQ(ierr); 2784b9ad928SBarry Smith PetscFunctionReturn(0); 2794b9ad928SBarry Smith } 2804b9ad928SBarry Smith 2814b9ad928SBarry Smith #undef __FUNCT__ 2824b9ad928SBarry Smith #define __FUNCT__ "PCSetFromOptions_Redundant" 2834416b707SBarry Smith static PetscErrorCode PCSetFromOptions_Redundant(PetscOptionItems *PetscOptionsObject,PC pc) 2844b9ad928SBarry Smith { 285a98ce0f4SHong Zhang PetscErrorCode ierr; 286a98ce0f4SHong Zhang PC_Redundant *red = (PC_Redundant*)pc->data; 287a98ce0f4SHong Zhang 2884b9ad928SBarry Smith PetscFunctionBegin; 289e55864a3SBarry Smith ierr = PetscOptionsHead(PetscOptionsObject,"Redundant options");CHKERRQ(ierr); 29009a6bc64SHong Zhang ierr = PetscOptionsInt("-pc_redundant_number","Number of redundant pc","PCRedundantSetNumber",red->nsubcomm,&red->nsubcomm,0);CHKERRQ(ierr); 291a98ce0f4SHong Zhang ierr = PetscOptionsTail();CHKERRQ(ierr); 2924b9ad928SBarry Smith PetscFunctionReturn(0); 2934b9ad928SBarry Smith } 2944b9ad928SBarry Smith 2954b9ad928SBarry Smith #undef __FUNCT__ 29609a6bc64SHong Zhang #define __FUNCT__ "PCRedundantSetNumber_Redundant" 297f7a08781SBarry Smith static PetscErrorCode PCRedundantSetNumber_Redundant(PC pc,PetscInt nreds) 29809a6bc64SHong Zhang { 29909a6bc64SHong Zhang PC_Redundant *red = (PC_Redundant*)pc->data; 30009a6bc64SHong Zhang 30109a6bc64SHong Zhang PetscFunctionBegin; 30209a6bc64SHong Zhang red->nsubcomm = nreds; 30309a6bc64SHong Zhang PetscFunctionReturn(0); 30409a6bc64SHong Zhang } 30509a6bc64SHong Zhang 30609a6bc64SHong Zhang #undef __FUNCT__ 30709a6bc64SHong Zhang #define __FUNCT__ "PCRedundantSetNumber" 30809a6bc64SHong Zhang /*@ 30909a6bc64SHong Zhang PCRedundantSetNumber - Sets the number of redundant preconditioner contexts. 31009a6bc64SHong Zhang 3113f9fe445SBarry Smith Logically Collective on PC 31209a6bc64SHong Zhang 31309a6bc64SHong Zhang Input Parameters: 31409a6bc64SHong Zhang + pc - the preconditioner context 3159b21d695SBarry Smith - nredundant - number of redundant preconditioner contexts; for example if you are using 64 MPI processes and 3169b21d695SBarry Smith use an nredundant of 4 there will be 4 parallel solves each on 16 = 64/4 processes. 31709a6bc64SHong Zhang 31809a6bc64SHong Zhang Level: advanced 31909a6bc64SHong Zhang 32009a6bc64SHong Zhang .keywords: PC, redundant solve 32109a6bc64SHong Zhang @*/ 3227087cfbeSBarry Smith PetscErrorCode PCRedundantSetNumber(PC pc,PetscInt nredundant) 32309a6bc64SHong Zhang { 3244ac538c5SBarry Smith PetscErrorCode ierr; 32509a6bc64SHong Zhang 32609a6bc64SHong Zhang PetscFunctionBegin; 3270700a824SBarry Smith PetscValidHeaderSpecific(pc,PC_CLASSID,1); 328ce94432eSBarry Smith if (nredundant <= 0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONG, "num of redundant pc %D must be positive",nredundant); 3294ac538c5SBarry Smith ierr = PetscTryMethod(pc,"PCRedundantSetNumber_C",(PC,PetscInt),(pc,nredundant));CHKERRQ(ierr); 33009a6bc64SHong Zhang PetscFunctionReturn(0); 33109a6bc64SHong Zhang } 33209a6bc64SHong Zhang 33309a6bc64SHong Zhang #undef __FUNCT__ 3344b9ad928SBarry Smith #define __FUNCT__ "PCRedundantSetScatter_Redundant" 335f7a08781SBarry Smith static PetscErrorCode PCRedundantSetScatter_Redundant(PC pc,VecScatter in,VecScatter out) 3364b9ad928SBarry Smith { 3374b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 338dfbe8321SBarry Smith PetscErrorCode ierr; 3394b9ad928SBarry Smith 3404b9ad928SBarry Smith PetscFunctionBegin; 3414b9ad928SBarry Smith ierr = PetscObjectReference((PetscObject)in);CHKERRQ(ierr); 3426bf464f9SBarry Smith ierr = VecScatterDestroy(&red->scatterin);CHKERRQ(ierr); 3432fa5cd67SKarl Rupp 344c3122656SLisandro Dalcin red->scatterin = in; 3452fa5cd67SKarl Rupp 3464b9ad928SBarry Smith ierr = PetscObjectReference((PetscObject)out);CHKERRQ(ierr); 3476bf464f9SBarry Smith ierr = VecScatterDestroy(&red->scatterout);CHKERRQ(ierr); 348c3122656SLisandro Dalcin red->scatterout = out; 3494b9ad928SBarry Smith PetscFunctionReturn(0); 3504b9ad928SBarry Smith } 3514b9ad928SBarry Smith 3524b9ad928SBarry Smith #undef __FUNCT__ 3534b9ad928SBarry Smith #define __FUNCT__ "PCRedundantSetScatter" 3544b9ad928SBarry Smith /*@ 3554b9ad928SBarry Smith PCRedundantSetScatter - Sets the scatter used to copy values into the 3564b9ad928SBarry Smith redundant local solve and the scatter to move them back into the global 3574b9ad928SBarry Smith vector. 3584b9ad928SBarry Smith 3593f9fe445SBarry Smith Logically Collective on PC 3604b9ad928SBarry Smith 3614b9ad928SBarry Smith Input Parameters: 3624b9ad928SBarry Smith + pc - the preconditioner context 3634b9ad928SBarry Smith . in - the scatter to move the values in 3644b9ad928SBarry Smith - out - the scatter to move them out 3654b9ad928SBarry Smith 3664b9ad928SBarry Smith Level: advanced 3674b9ad928SBarry Smith 3684b9ad928SBarry Smith .keywords: PC, redundant solve 3694b9ad928SBarry Smith @*/ 3707087cfbeSBarry Smith PetscErrorCode PCRedundantSetScatter(PC pc,VecScatter in,VecScatter out) 3714b9ad928SBarry Smith { 3724ac538c5SBarry Smith PetscErrorCode ierr; 3734b9ad928SBarry Smith 3744b9ad928SBarry Smith PetscFunctionBegin; 3750700a824SBarry Smith PetscValidHeaderSpecific(pc,PC_CLASSID,1); 3760700a824SBarry Smith PetscValidHeaderSpecific(in,VEC_SCATTER_CLASSID,2); 3770700a824SBarry Smith PetscValidHeaderSpecific(out,VEC_SCATTER_CLASSID,3); 3784ac538c5SBarry Smith ierr = PetscTryMethod(pc,"PCRedundantSetScatter_C",(PC,VecScatter,VecScatter),(pc,in,out));CHKERRQ(ierr); 3794b9ad928SBarry Smith PetscFunctionReturn(0); 3804b9ad928SBarry Smith } 3814b9ad928SBarry Smith 3824b9ad928SBarry Smith #undef __FUNCT__ 38383ab6a24SBarry Smith #define __FUNCT__ "PCRedundantGetKSP_Redundant" 384f7a08781SBarry Smith static PetscErrorCode PCRedundantGetKSP_Redundant(PC pc,KSP *innerksp) 3854b9ad928SBarry Smith { 3865f06b7aaSBarry Smith PetscErrorCode ierr; 3874b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 38875024027SHong Zhang MPI_Comm comm,subcomm; 38975024027SHong Zhang const char *prefix; 3904b9ad928SBarry Smith 3914b9ad928SBarry Smith PetscFunctionBegin; 39275024027SHong Zhang if (!red->psubcomm) { 393e5acf8a4SHong Zhang ierr = PCGetOptionsPrefix(pc,&prefix);CHKERRQ(ierr); 394e5acf8a4SHong Zhang 39575024027SHong Zhang ierr = PetscObjectGetComm((PetscObject)pc,&comm);CHKERRQ(ierr); 39675024027SHong Zhang ierr = PetscSubcommCreate(comm,&red->psubcomm);CHKERRQ(ierr); 39775024027SHong Zhang ierr = PetscSubcommSetNumber(red->psubcomm,red->nsubcomm);CHKERRQ(ierr); 39875024027SHong Zhang ierr = PetscSubcommSetType(red->psubcomm,PETSC_SUBCOMM_CONTIGUOUS);CHKERRQ(ierr); 399e5acf8a4SHong Zhang 400e5acf8a4SHong Zhang ierr = PetscSubcommSetOptionsPrefix(red->psubcomm,prefix);CHKERRQ(ierr); 401e5acf8a4SHong Zhang ierr = PetscSubcommSetFromOptions(red->psubcomm);CHKERRQ(ierr); 40275024027SHong Zhang ierr = PetscLogObjectMemory((PetscObject)pc,sizeof(PetscSubcomm));CHKERRQ(ierr); 40375024027SHong Zhang 40475024027SHong Zhang /* create a new PC that processors in each subcomm have copy of */ 40575024027SHong Zhang subcomm = PetscSubcommChild(red->psubcomm); 40675024027SHong Zhang 40775024027SHong Zhang ierr = KSPCreate(subcomm,&red->ksp);CHKERRQ(ierr); 40875024027SHong Zhang ierr = KSPSetErrorIfNotConverged(red->ksp,pc->erroriffailure);CHKERRQ(ierr); 40975024027SHong Zhang ierr = PetscObjectIncrementTabLevel((PetscObject)red->ksp,(PetscObject)pc,1);CHKERRQ(ierr); 41075024027SHong Zhang ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)red->ksp);CHKERRQ(ierr); 41175024027SHong Zhang ierr = KSPSetType(red->ksp,KSPPREONLY);CHKERRQ(ierr); 41275024027SHong Zhang ierr = KSPGetPC(red->ksp,&red->pc);CHKERRQ(ierr); 41375024027SHong Zhang ierr = PCSetType(red->pc,PCLU);CHKERRQ(ierr); 41475024027SHong Zhang 41575024027SHong Zhang ierr = KSPSetOptionsPrefix(red->ksp,prefix);CHKERRQ(ierr); 41675024027SHong Zhang ierr = KSPAppendOptionsPrefix(red->ksp,"redundant_");CHKERRQ(ierr); 41775024027SHong Zhang } 41883ab6a24SBarry Smith *innerksp = red->ksp; 4194b9ad928SBarry Smith PetscFunctionReturn(0); 4204b9ad928SBarry Smith } 4214b9ad928SBarry Smith 4224b9ad928SBarry Smith #undef __FUNCT__ 42383ab6a24SBarry Smith #define __FUNCT__ "PCRedundantGetKSP" 4244b9ad928SBarry Smith /*@ 42583ab6a24SBarry Smith PCRedundantGetKSP - Gets the less parallel KSP created by the redundant PC. 4264b9ad928SBarry Smith 4274b9ad928SBarry Smith Not Collective 4284b9ad928SBarry Smith 4294b9ad928SBarry Smith Input Parameter: 4304b9ad928SBarry Smith . pc - the preconditioner context 4314b9ad928SBarry Smith 4324b9ad928SBarry Smith Output Parameter: 43383ab6a24SBarry Smith . innerksp - the KSP on the smaller set of processes 4344b9ad928SBarry Smith 4354b9ad928SBarry Smith Level: advanced 4364b9ad928SBarry Smith 4374b9ad928SBarry Smith .keywords: PC, redundant solve 4384b9ad928SBarry Smith @*/ 43983ab6a24SBarry Smith PetscErrorCode PCRedundantGetKSP(PC pc,KSP *innerksp) 4404b9ad928SBarry Smith { 4414ac538c5SBarry Smith PetscErrorCode ierr; 4424b9ad928SBarry Smith 4434b9ad928SBarry Smith PetscFunctionBegin; 4440700a824SBarry Smith PetscValidHeaderSpecific(pc,PC_CLASSID,1); 44583ab6a24SBarry Smith PetscValidPointer(innerksp,2); 44683ab6a24SBarry Smith ierr = PetscTryMethod(pc,"PCRedundantGetKSP_C",(PC,KSP*),(pc,innerksp));CHKERRQ(ierr); 4474b9ad928SBarry Smith PetscFunctionReturn(0); 4484b9ad928SBarry Smith } 4494b9ad928SBarry Smith 4504b9ad928SBarry Smith #undef __FUNCT__ 4514b9ad928SBarry Smith #define __FUNCT__ "PCRedundantGetOperators_Redundant" 452f7a08781SBarry Smith static PetscErrorCode PCRedundantGetOperators_Redundant(PC pc,Mat *mat,Mat *pmat) 4534b9ad928SBarry Smith { 4544b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 4554b9ad928SBarry Smith 4564b9ad928SBarry Smith PetscFunctionBegin; 457b3804887SHong Zhang if (mat) *mat = red->pmats; 458b3804887SHong Zhang if (pmat) *pmat = red->pmats; 4594b9ad928SBarry Smith PetscFunctionReturn(0); 4604b9ad928SBarry Smith } 4614b9ad928SBarry Smith 4624b9ad928SBarry Smith #undef __FUNCT__ 4634b9ad928SBarry Smith #define __FUNCT__ "PCRedundantGetOperators" 4644b9ad928SBarry Smith /*@ 4654b9ad928SBarry Smith PCRedundantGetOperators - gets the sequential matrix and preconditioner matrix 4664b9ad928SBarry Smith 4674b9ad928SBarry Smith Not Collective 4684b9ad928SBarry Smith 4694b9ad928SBarry Smith Input Parameter: 4704b9ad928SBarry Smith . pc - the preconditioner context 4714b9ad928SBarry Smith 4724b9ad928SBarry Smith Output Parameters: 4734b9ad928SBarry Smith + mat - the matrix 4744b9ad928SBarry Smith - pmat - the (possibly different) preconditioner matrix 4754b9ad928SBarry Smith 4764b9ad928SBarry Smith Level: advanced 4774b9ad928SBarry Smith 4784b9ad928SBarry Smith .keywords: PC, redundant solve 4794b9ad928SBarry Smith @*/ 4807087cfbeSBarry Smith PetscErrorCode PCRedundantGetOperators(PC pc,Mat *mat,Mat *pmat) 4814b9ad928SBarry Smith { 4824ac538c5SBarry Smith PetscErrorCode ierr; 4834b9ad928SBarry Smith 4844b9ad928SBarry Smith PetscFunctionBegin; 4850700a824SBarry Smith PetscValidHeaderSpecific(pc,PC_CLASSID,1); 4864482741eSBarry Smith if (mat) PetscValidPointer(mat,2); 4874482741eSBarry Smith if (pmat) PetscValidPointer(pmat,3); 4884ac538c5SBarry Smith ierr = PetscTryMethod(pc,"PCRedundantGetOperators_C",(PC,Mat*,Mat*),(pc,mat,pmat));CHKERRQ(ierr); 4894b9ad928SBarry Smith PetscFunctionReturn(0); 4904b9ad928SBarry Smith } 4914b9ad928SBarry Smith 4924b9ad928SBarry Smith /* -------------------------------------------------------------------------------------*/ 49337a17b4dSBarry Smith /*MC 49483ab6a24SBarry Smith PCREDUNDANT - Runs a KSP solver with preconditioner for the entire problem on subgroups of processors 49537a17b4dSBarry Smith 49683ab6a24SBarry Smith Options for the redundant preconditioners can be set with -redundant_pc_xxx for the redundant KSP with -redundant_ksp_xxx 49737a17b4dSBarry Smith 49809391456SBarry Smith Options Database: 4999b21d695SBarry Smith . -pc_redundant_number <n> - number of redundant solves, for example if you are using 64 MPI processes and 5009b21d695SBarry Smith use an n of 4 there will be 4 parallel solves each on 16 = 64/4 processes. 50109391456SBarry Smith 50237a17b4dSBarry Smith Level: intermediate 50337a17b4dSBarry Smith 50483ab6a24SBarry Smith Notes: The default KSP is preonly and the default PC is LU. 50583ab6a24SBarry Smith 50683ab6a24SBarry Smith Developer Notes: Note that PCSetInitialGuessNonzero() is not used by this class but likely should be. 5079cfaa89bSBarry Smith 50837a17b4dSBarry Smith .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PCRedundantSetScatter(), 50983ab6a24SBarry Smith PCRedundantGetKSP(), PCRedundantGetOperators(), PCRedundantSetNumber() 51037a17b4dSBarry Smith M*/ 51137a17b4dSBarry Smith 5124b9ad928SBarry Smith #undef __FUNCT__ 5134b9ad928SBarry Smith #define __FUNCT__ "PCCreate_Redundant" 5148cc058d9SJed Brown PETSC_EXTERN PetscErrorCode PCCreate_Redundant(PC pc) 5154b9ad928SBarry Smith { 516dfbe8321SBarry Smith PetscErrorCode ierr; 5174b9ad928SBarry Smith PC_Redundant *red; 51869db28dcSHong Zhang PetscMPIInt size; 5193f457be1SHong Zhang 5204b9ad928SBarry Smith PetscFunctionBegin; 521b00a9115SJed Brown ierr = PetscNewLog(pc,&red);CHKERRQ(ierr); 522ce94432eSBarry Smith ierr = MPI_Comm_size(PetscObjectComm((PetscObject)pc),&size);CHKERRQ(ierr); 5232fa5cd67SKarl Rupp 52469db28dcSHong Zhang red->nsubcomm = size; 5254b9ad928SBarry Smith red->useparallelmat = PETSC_TRUE; 5261fbd8f88SHong Zhang pc->data = (void*)red; 5274b9ad928SBarry Smith 5284b9ad928SBarry Smith pc->ops->apply = PCApply_Redundant; 529d88bfacbSStefano Zampini pc->ops->applytranspose = PCApplyTranspose_Redundant; 5304b9ad928SBarry Smith pc->ops->setup = PCSetUp_Redundant; 5314b9ad928SBarry Smith pc->ops->destroy = PCDestroy_Redundant; 5321ea5a559SBarry Smith pc->ops->reset = PCReset_Redundant; 5334b9ad928SBarry Smith pc->ops->setfromoptions = PCSetFromOptions_Redundant; 5344b9ad928SBarry Smith pc->ops->view = PCView_Redundant; 5352fa5cd67SKarl Rupp 536bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCRedundantSetScatter_C",PCRedundantSetScatter_Redundant);CHKERRQ(ierr); 537bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCRedundantSetNumber_C",PCRedundantSetNumber_Redundant);CHKERRQ(ierr); 538bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCRedundantGetKSP_C",PCRedundantGetKSP_Redundant);CHKERRQ(ierr); 539bdf89e91SBarry Smith ierr = PetscObjectComposeFunction((PetscObject)pc,"PCRedundantGetOperators_C",PCRedundantGetOperators_Redundant);CHKERRQ(ierr); 5404b9ad928SBarry Smith PetscFunctionReturn(0); 5414b9ad928SBarry Smith } 542b2573a8aSBarry Smith 543