xref: /petsc/src/ksp/pc/impls/redundant/redundant.c (revision dcca6d9d80ebd869fe6029bd05a3aa9faafef49e)
1dba47a55SKris Buschelman 
24b9ad928SBarry Smith /*
33f457be1SHong Zhang   This file defines a "solve the problem redundantly on each subgroup of processor" preconditioner.
44b9ad928SBarry Smith */
507475bc1SBarry 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);
377adad957SLisandro Dalcin       ierr = PetscViewerGetSubcomm(viewer,((PetscObject)red->pc)->comm,&subviewer);CHKERRQ(ierr);
38f5dd71faSBarry Smith       if (!red->psubcomm->color) { /* only view first redundant pc */
394b9ad928SBarry Smith         ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
403e065800SHong Zhang         ierr = KSPView(red->ksp,subviewer);CHKERRQ(ierr);
414b9ad928SBarry Smith         ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
424b9ad928SBarry Smith       }
437adad957SLisandro Dalcin       ierr = PetscViewerRestoreSubcomm(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;
611fbd8f88SHong Zhang   const char     *prefix;
623f457be1SHong Zhang 
634b9ad928SBarry Smith   PetscFunctionBegin;
64ce94432eSBarry Smith   ierr = PetscObjectGetComm((PetscObject)pc,&comm);CHKERRQ(ierr);
65ddc54837SHong Zhang 
66ddc54837SHong Zhang   /* if pmatrix set by user is sequential then we do not need to gather the parallel matrix */
67ddc54837SHong Zhang   ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr);
68ddc54837SHong Zhang   if (size == 1) red->useparallelmat = PETSC_FALSE;
691fbd8f88SHong Zhang 
704b9ad928SBarry Smith   if (!pc->setupcalled) {
711b81debcSHong Zhang     PetscInt mloc_sub;
725f06b7aaSBarry Smith     if (!red->psubcomm) {
73d8a68f86SHong Zhang       ierr = PetscSubcommCreate(comm,&red->psubcomm);CHKERRQ(ierr);
74d8a68f86SHong Zhang       ierr = PetscSubcommSetNumber(red->psubcomm,red->nsubcomm);CHKERRQ(ierr);
75a23b1e67SHong Zhang       ierr = PetscSubcommSetType(red->psubcomm,PETSC_SUBCOMM_CONTIGUOUS);CHKERRQ(ierr);
76a23b1e67SHong Zhang       /* enable runtime switch of psubcomm type, e.g., '-psubcomm_type interlaced */
77a23b1e67SHong Zhang       ierr = PetscSubcommSetFromOptions(red->psubcomm);CHKERRQ(ierr);
783bb1ff40SBarry Smith       ierr = PetscLogObjectMemory((PetscObject)pc,sizeof(PetscSubcomm));CHKERRQ(ierr);
791fbd8f88SHong Zhang 
801fbd8f88SHong Zhang       /* create a new PC that processors in each subcomm have copy of */
810d7810c8SBarry Smith       subcomm = red->psubcomm->comm;
822fa5cd67SKarl Rupp 
835f06b7aaSBarry Smith       ierr = KSPCreate(subcomm,&red->ksp);CHKERRQ(ierr);
845f06b7aaSBarry Smith       ierr = PetscObjectIncrementTabLevel((PetscObject)red->ksp,(PetscObject)pc,1);CHKERRQ(ierr);
853bb1ff40SBarry Smith       ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)red->ksp);CHKERRQ(ierr);
865f06b7aaSBarry Smith       ierr = KSPSetType(red->ksp,KSPPREONLY);CHKERRQ(ierr);
875f06b7aaSBarry Smith       ierr = KSPGetPC(red->ksp,&red->pc);CHKERRQ(ierr);
88cf52b8b1SHong Zhang       ierr = PCSetType(red->pc,PCLU);CHKERRQ(ierr);
89cf52b8b1SHong Zhang 
901fbd8f88SHong Zhang       ierr = PCGetOptionsPrefix(pc,&prefix);CHKERRQ(ierr);
915f06b7aaSBarry Smith       ierr = KSPSetOptionsPrefix(red->ksp,prefix);CHKERRQ(ierr);
925f06b7aaSBarry Smith       ierr = KSPAppendOptionsPrefix(red->ksp,"redundant_");CHKERRQ(ierr);
931b81debcSHong Zhang     } else {
941b81debcSHong Zhang       subcomm = red->psubcomm->comm;
951b81debcSHong Zhang     }
961fbd8f88SHong Zhang 
971b81debcSHong Zhang     if (red->useparallelmat) {
981b81debcSHong Zhang       /* grab the parallel matrix and put it into processors of a subcomminicator */
99b2bf6370SHong Zhang       ierr = MatGetRedundantMatrix(pc->pmat,red->psubcomm->n,subcomm,MAT_INITIAL_MATRIX,&red->pmats);CHKERRQ(ierr);
1001b81debcSHong Zhang       ierr = KSPSetOperators(red->ksp,red->pmats,red->pmats,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
1014b9ad928SBarry Smith 
1021b81debcSHong Zhang       /* get working vectors xsub and ysub */
1031b81debcSHong Zhang       ierr = MatGetVecs(red->pmats,&red->xsub,&red->ysub);CHKERRQ(ierr);
1042fa5cd67SKarl Rupp 
1058b96b0d2SHong Zhang       /* create working vectors xdup and ydup.
1068b96b0d2SHong Zhang        xdup concatenates all xsub's contigously to form a mpi vector over dupcomm  (see PetscSubcommCreate_interlaced())
1078b96b0d2SHong Zhang        ydup concatenates all ysub and has empty local arrays because ysub's arrays will be place into it.
108ce94432eSBarry Smith        Note: we use communicator dupcomm, not PetscObjectComm((PetscObject)pc)! */
1091b81debcSHong Zhang       ierr = MatGetLocalSize(red->pmats,&mloc_sub,NULL);CHKERRQ(ierr);
1101fbd8f88SHong Zhang       ierr = VecCreateMPI(red->psubcomm->dupparent,mloc_sub,PETSC_DECIDE,&red->xdup);CHKERRQ(ierr);
1110298fd71SBarry Smith       ierr = VecCreateMPIWithArray(red->psubcomm->dupparent,1,mloc_sub,PETSC_DECIDE,NULL,&red->ydup);CHKERRQ(ierr);
1123f457be1SHong Zhang 
113f68be91cSHong Zhang       /* create vecscatters */
114f68be91cSHong Zhang       if (!red->scatterin) { /* efficiency of scatterin is independent from psubcomm_type! */
1153f457be1SHong Zhang         IS       is1,is2;
1163f457be1SHong Zhang         PetscInt *idx1,*idx2,i,j,k;
11745fc02eaSBarry Smith 
1181b81debcSHong Zhang         ierr = MatGetVecs(pc->pmat,&x,0);CHKERRQ(ierr);
1191b81debcSHong Zhang         ierr = VecGetSize(x,&M);CHKERRQ(ierr);
1201b81debcSHong Zhang         ierr = VecGetOwnershipRange(x,&mstart,&mend);CHKERRQ(ierr);
1211b81debcSHong Zhang         mlocal = mend - mstart;
122*dcca6d9dSJed Brown         ierr = PetscMalloc2(red->psubcomm->n*mlocal,&idx1,red->psubcomm->n*mlocal,&idx2);CHKERRQ(ierr);
1233f457be1SHong Zhang         j    = 0;
1241fbd8f88SHong Zhang         for (k=0; k<red->psubcomm->n; k++) {
1253f457be1SHong Zhang           for (i=mstart; i<mend; i++) {
1263f457be1SHong Zhang             idx1[j]   = i;
127ddc54837SHong Zhang             idx2[j++] = i + M*k;
1283f457be1SHong Zhang           }
1293f457be1SHong Zhang         }
13070b3c8c7SBarry Smith         ierr = ISCreateGeneral(comm,red->psubcomm->n*mlocal,idx1,PETSC_COPY_VALUES,&is1);CHKERRQ(ierr);
13170b3c8c7SBarry Smith         ierr = ISCreateGeneral(comm,red->psubcomm->n*mlocal,idx2,PETSC_COPY_VALUES,&is2);CHKERRQ(ierr);
132ddc54837SHong Zhang         ierr = VecScatterCreate(x,is1,red->xdup,is2,&red->scatterin);CHKERRQ(ierr);
133fcfd50ebSBarry Smith         ierr = ISDestroy(&is1);CHKERRQ(ierr);
134fcfd50ebSBarry Smith         ierr = ISDestroy(&is2);CHKERRQ(ierr);
1353f457be1SHong Zhang 
1366909748bSHong Zhang         /* Impl below is good for PETSC_SUBCOMM_INTERLACED (no inter-process communication) and PETSC_SUBCOMM_CONTIGUOUS (communication within subcomm) */
137ddc54837SHong Zhang         ierr = ISCreateStride(comm,mlocal,mstart+ red->psubcomm->color*M,1,&is1);CHKERRQ(ierr);
1383f457be1SHong Zhang         ierr = ISCreateStride(comm,mlocal,mstart,1,&is2);CHKERRQ(ierr);
139ddc54837SHong Zhang         ierr = VecScatterCreate(red->xdup,is1,x,is2,&red->scatterout);CHKERRQ(ierr);
140fcfd50ebSBarry Smith         ierr = ISDestroy(&is1);CHKERRQ(ierr);
141fcfd50ebSBarry Smith         ierr = ISDestroy(&is2);CHKERRQ(ierr);
1421d79065fSBarry Smith         ierr = PetscFree2(idx1,idx2);CHKERRQ(ierr);
143ddc54837SHong Zhang         ierr = VecDestroy(&x);CHKERRQ(ierr);
1441b81debcSHong Zhang       }
145ab661555SHong Zhang     } else { /* !red->useparallelmat */
146ab661555SHong Zhang       ierr = KSPSetOperators(red->ksp,pc->mat,pc->pmat,pc->flag);CHKERRQ(ierr);
1471b81debcSHong Zhang     }
148ab661555SHong Zhang   } else { /* pc->setupcalled */
1494b9ad928SBarry Smith     if (red->useparallelmat) {
150ab661555SHong Zhang       MatReuse       reuse;
1511b81debcSHong Zhang       /* grab the parallel matrix and put it into processors of a subcomminicator */
1521b81debcSHong Zhang       /*--------------------------------------------------------------------------*/
153ab661555SHong Zhang       if (pc->flag == DIFFERENT_NONZERO_PATTERN) {
1544b9ad928SBarry Smith         /* destroy old matrices */
1556bf464f9SBarry Smith         ierr  = MatDestroy(&red->pmats);CHKERRQ(ierr);
156ab661555SHong Zhang         reuse = MAT_INITIAL_MATRIX;
1574b9ad928SBarry Smith       } else {
158ab661555SHong Zhang         reuse = MAT_REUSE_MATRIX;
159ab661555SHong Zhang       }
160b2bf6370SHong Zhang       ierr = MatGetRedundantMatrix(pc->pmat,red->psubcomm->n,red->psubcomm->comm,reuse,&red->pmats);CHKERRQ(ierr);
161ab661555SHong Zhang       ierr = KSPSetOperators(red->ksp,red->pmats,red->pmats,pc->flag);CHKERRQ(ierr);
162ab661555SHong Zhang     } else { /* !red->useparallelmat */
16390f1c854SHong Zhang       ierr = KSPSetOperators(red->ksp,pc->mat,pc->pmat,pc->flag);CHKERRQ(ierr);
1644b9ad928SBarry Smith     }
165ab661555SHong Zhang   }
1661b81debcSHong Zhang 
1670c24e6a1SHong Zhang   if (pc->setfromoptionscalled) {
1683e065800SHong Zhang     ierr = KSPSetFromOptions(red->ksp);CHKERRQ(ierr);
1690c24e6a1SHong Zhang   }
1703e065800SHong Zhang   ierr = KSPSetUp(red->ksp);CHKERRQ(ierr);
1714b9ad928SBarry Smith   PetscFunctionReturn(0);
1724b9ad928SBarry Smith }
1734b9ad928SBarry Smith 
1744b9ad928SBarry Smith #undef __FUNCT__
1754b9ad928SBarry Smith #define __FUNCT__ "PCApply_Redundant"
1766849ba73SBarry Smith static PetscErrorCode PCApply_Redundant(PC pc,Vec x,Vec y)
1774b9ad928SBarry Smith {
1784b9ad928SBarry Smith   PC_Redundant   *red = (PC_Redundant*)pc->data;
179dfbe8321SBarry Smith   PetscErrorCode ierr;
1803f457be1SHong Zhang   PetscScalar    *array;
1814b9ad928SBarry Smith 
1824b9ad928SBarry Smith   PetscFunctionBegin;
183ddc54837SHong Zhang   if (!red->useparallelmat) {
184ddc54837SHong Zhang     ierr = KSPSolve(red->ksp,x,y);CHKERRQ(ierr);
185ddc54837SHong Zhang     PetscFunctionReturn(0);
186ddc54837SHong Zhang   }
187ddc54837SHong Zhang 
1883f457be1SHong Zhang   /* scatter x to xdup */
189ca9f406cSSatish Balay   ierr = VecScatterBegin(red->scatterin,x,red->xdup,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
190ca9f406cSSatish Balay   ierr = VecScatterEnd(red->scatterin,x,red->xdup,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
1913f457be1SHong Zhang 
1923f457be1SHong Zhang   /* place xdup's local array into xsub */
1933f457be1SHong Zhang   ierr = VecGetArray(red->xdup,&array);CHKERRQ(ierr);
1943f457be1SHong Zhang   ierr = VecPlaceArray(red->xsub,(const PetscScalar*)array);CHKERRQ(ierr);
1954b9ad928SBarry Smith 
1964b9ad928SBarry Smith   /* apply preconditioner on each processor */
19783ab6a24SBarry Smith   ierr = KSPSolve(red->ksp,red->xsub,red->ysub);CHKERRQ(ierr);
1983f457be1SHong Zhang   ierr = VecResetArray(red->xsub);CHKERRQ(ierr);
1993f457be1SHong Zhang   ierr = VecRestoreArray(red->xdup,&array);CHKERRQ(ierr);
2004b9ad928SBarry Smith 
2013f457be1SHong Zhang   /* place ysub's local array into ydup */
2023f457be1SHong Zhang   ierr = VecGetArray(red->ysub,&array);CHKERRQ(ierr);
2033f457be1SHong Zhang   ierr = VecPlaceArray(red->ydup,(const PetscScalar*)array);CHKERRQ(ierr);
2043f457be1SHong Zhang 
2053f457be1SHong Zhang   /* scatter ydup to y */
206ca9f406cSSatish Balay   ierr = VecScatterBegin(red->scatterout,red->ydup,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
207ca9f406cSSatish Balay   ierr = VecScatterEnd(red->scatterout,red->ydup,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
2083f457be1SHong Zhang   ierr = VecResetArray(red->ydup);CHKERRQ(ierr);
2093f457be1SHong Zhang   ierr = VecRestoreArray(red->ysub,&array);CHKERRQ(ierr);
2104b9ad928SBarry Smith   PetscFunctionReturn(0);
2114b9ad928SBarry Smith }
2124b9ad928SBarry Smith 
2134b9ad928SBarry Smith #undef __FUNCT__
2141ea5a559SBarry Smith #define __FUNCT__ "PCReset_Redundant"
2151ea5a559SBarry Smith static PetscErrorCode PCReset_Redundant(PC pc)
2164b9ad928SBarry Smith {
2174b9ad928SBarry Smith   PC_Redundant   *red = (PC_Redundant*)pc->data;
218dfbe8321SBarry Smith   PetscErrorCode ierr;
2194b9ad928SBarry Smith 
2204b9ad928SBarry Smith   PetscFunctionBegin;
2211b81debcSHong Zhang   if (red->useparallelmat) {
2226bf464f9SBarry Smith     ierr = VecScatterDestroy(&red->scatterin);CHKERRQ(ierr);
2236bf464f9SBarry Smith     ierr = VecScatterDestroy(&red->scatterout);CHKERRQ(ierr);
2246bf464f9SBarry Smith     ierr = VecDestroy(&red->ysub);CHKERRQ(ierr);
2256bf464f9SBarry Smith     ierr = VecDestroy(&red->xsub);CHKERRQ(ierr);
2266bf464f9SBarry Smith     ierr = VecDestroy(&red->xdup);CHKERRQ(ierr);
2276bf464f9SBarry Smith     ierr = VecDestroy(&red->ydup);CHKERRQ(ierr);
2281b81debcSHong Zhang   }
2296bf464f9SBarry Smith   ierr = MatDestroy(&red->pmats);CHKERRQ(ierr);
2301b81debcSHong Zhang   ierr = KSPReset(red->ksp);CHKERRQ(ierr);
2311ea5a559SBarry Smith   PetscFunctionReturn(0);
2321ea5a559SBarry Smith }
2331ea5a559SBarry Smith 
2341ea5a559SBarry Smith #undef __FUNCT__
2351ea5a559SBarry Smith #define __FUNCT__ "PCDestroy_Redundant"
2361ea5a559SBarry Smith static PetscErrorCode PCDestroy_Redundant(PC pc)
2371ea5a559SBarry Smith {
2381ea5a559SBarry Smith   PC_Redundant   *red = (PC_Redundant*)pc->data;
2391ea5a559SBarry Smith   PetscErrorCode ierr;
2401ea5a559SBarry Smith 
2411ea5a559SBarry Smith   PetscFunctionBegin;
2421ea5a559SBarry Smith   ierr = PCReset_Redundant(pc);CHKERRQ(ierr);
2436bf464f9SBarry Smith   ierr = KSPDestroy(&red->ksp);CHKERRQ(ierr);
2446bf464f9SBarry Smith   ierr = PetscSubcommDestroy(&red->psubcomm);CHKERRQ(ierr);
245c31cb41cSBarry Smith   ierr = PetscFree(pc->data);CHKERRQ(ierr);
2464b9ad928SBarry Smith   PetscFunctionReturn(0);
2474b9ad928SBarry Smith }
2484b9ad928SBarry Smith 
2494b9ad928SBarry Smith #undef __FUNCT__
2504b9ad928SBarry Smith #define __FUNCT__ "PCSetFromOptions_Redundant"
2516849ba73SBarry Smith static PetscErrorCode PCSetFromOptions_Redundant(PC pc)
2524b9ad928SBarry Smith {
253a98ce0f4SHong Zhang   PetscErrorCode ierr;
254a98ce0f4SHong Zhang   PC_Redundant   *red = (PC_Redundant*)pc->data;
255a98ce0f4SHong Zhang 
2564b9ad928SBarry Smith   PetscFunctionBegin;
257a98ce0f4SHong Zhang   ierr = PetscOptionsHead("Redundant options");CHKERRQ(ierr);
25809a6bc64SHong Zhang   ierr = PetscOptionsInt("-pc_redundant_number","Number of redundant pc","PCRedundantSetNumber",red->nsubcomm,&red->nsubcomm,0);CHKERRQ(ierr);
259a98ce0f4SHong Zhang   ierr = PetscOptionsTail();CHKERRQ(ierr);
2604b9ad928SBarry Smith   PetscFunctionReturn(0);
2614b9ad928SBarry Smith }
2624b9ad928SBarry Smith 
2634b9ad928SBarry Smith #undef __FUNCT__
26409a6bc64SHong Zhang #define __FUNCT__ "PCRedundantSetNumber_Redundant"
265f7a08781SBarry Smith static PetscErrorCode  PCRedundantSetNumber_Redundant(PC pc,PetscInt nreds)
26609a6bc64SHong Zhang {
26709a6bc64SHong Zhang   PC_Redundant *red = (PC_Redundant*)pc->data;
26809a6bc64SHong Zhang 
26909a6bc64SHong Zhang   PetscFunctionBegin;
27009a6bc64SHong Zhang   red->nsubcomm = nreds;
27109a6bc64SHong Zhang   PetscFunctionReturn(0);
27209a6bc64SHong Zhang }
27309a6bc64SHong Zhang 
27409a6bc64SHong Zhang #undef __FUNCT__
27509a6bc64SHong Zhang #define __FUNCT__ "PCRedundantSetNumber"
27609a6bc64SHong Zhang /*@
27709a6bc64SHong Zhang    PCRedundantSetNumber - Sets the number of redundant preconditioner contexts.
27809a6bc64SHong Zhang 
2793f9fe445SBarry Smith    Logically Collective on PC
28009a6bc64SHong Zhang 
28109a6bc64SHong Zhang    Input Parameters:
28209a6bc64SHong Zhang +  pc - the preconditioner context
2839b21d695SBarry Smith -  nredundant - number of redundant preconditioner contexts; for example if you are using 64 MPI processes and
2849b21d695SBarry Smith                               use an nredundant of 4 there will be 4 parallel solves each on 16 = 64/4 processes.
28509a6bc64SHong Zhang 
28609a6bc64SHong Zhang    Level: advanced
28709a6bc64SHong Zhang 
28809a6bc64SHong Zhang .keywords: PC, redundant solve
28909a6bc64SHong Zhang @*/
2907087cfbeSBarry Smith PetscErrorCode  PCRedundantSetNumber(PC pc,PetscInt nredundant)
29109a6bc64SHong Zhang {
2924ac538c5SBarry Smith   PetscErrorCode ierr;
29309a6bc64SHong Zhang 
29409a6bc64SHong Zhang   PetscFunctionBegin;
2950700a824SBarry Smith   PetscValidHeaderSpecific(pc,PC_CLASSID,1);
296ce94432eSBarry Smith   if (nredundant <= 0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONG, "num of redundant pc %D must be positive",nredundant);
2974ac538c5SBarry Smith   ierr = PetscTryMethod(pc,"PCRedundantSetNumber_C",(PC,PetscInt),(pc,nredundant));CHKERRQ(ierr);
29809a6bc64SHong Zhang   PetscFunctionReturn(0);
29909a6bc64SHong Zhang }
30009a6bc64SHong Zhang 
30109a6bc64SHong Zhang #undef __FUNCT__
3024b9ad928SBarry Smith #define __FUNCT__ "PCRedundantSetScatter_Redundant"
303f7a08781SBarry Smith static PetscErrorCode  PCRedundantSetScatter_Redundant(PC pc,VecScatter in,VecScatter out)
3044b9ad928SBarry Smith {
3054b9ad928SBarry Smith   PC_Redundant   *red = (PC_Redundant*)pc->data;
306dfbe8321SBarry Smith   PetscErrorCode ierr;
3074b9ad928SBarry Smith 
3084b9ad928SBarry Smith   PetscFunctionBegin;
3094b9ad928SBarry Smith   ierr = PetscObjectReference((PetscObject)in);CHKERRQ(ierr);
3106bf464f9SBarry Smith   ierr = VecScatterDestroy(&red->scatterin);CHKERRQ(ierr);
3112fa5cd67SKarl Rupp 
312c3122656SLisandro Dalcin   red->scatterin  = in;
3132fa5cd67SKarl Rupp 
3144b9ad928SBarry Smith   ierr = PetscObjectReference((PetscObject)out);CHKERRQ(ierr);
3156bf464f9SBarry Smith   ierr = VecScatterDestroy(&red->scatterout);CHKERRQ(ierr);
316c3122656SLisandro Dalcin   red->scatterout = out;
3174b9ad928SBarry Smith   PetscFunctionReturn(0);
3184b9ad928SBarry Smith }
3194b9ad928SBarry Smith 
3204b9ad928SBarry Smith #undef __FUNCT__
3214b9ad928SBarry Smith #define __FUNCT__ "PCRedundantSetScatter"
3224b9ad928SBarry Smith /*@
3234b9ad928SBarry Smith    PCRedundantSetScatter - Sets the scatter used to copy values into the
3244b9ad928SBarry Smith      redundant local solve and the scatter to move them back into the global
3254b9ad928SBarry Smith      vector.
3264b9ad928SBarry Smith 
3273f9fe445SBarry Smith    Logically Collective on PC
3284b9ad928SBarry Smith 
3294b9ad928SBarry Smith    Input Parameters:
3304b9ad928SBarry Smith +  pc - the preconditioner context
3314b9ad928SBarry Smith .  in - the scatter to move the values in
3324b9ad928SBarry Smith -  out - the scatter to move them out
3334b9ad928SBarry Smith 
3344b9ad928SBarry Smith    Level: advanced
3354b9ad928SBarry Smith 
3364b9ad928SBarry Smith .keywords: PC, redundant solve
3374b9ad928SBarry Smith @*/
3387087cfbeSBarry Smith PetscErrorCode  PCRedundantSetScatter(PC pc,VecScatter in,VecScatter out)
3394b9ad928SBarry Smith {
3404ac538c5SBarry Smith   PetscErrorCode ierr;
3414b9ad928SBarry Smith 
3424b9ad928SBarry Smith   PetscFunctionBegin;
3430700a824SBarry Smith   PetscValidHeaderSpecific(pc,PC_CLASSID,1);
3440700a824SBarry Smith   PetscValidHeaderSpecific(in,VEC_SCATTER_CLASSID,2);
3450700a824SBarry Smith   PetscValidHeaderSpecific(out,VEC_SCATTER_CLASSID,3);
3464ac538c5SBarry Smith   ierr = PetscTryMethod(pc,"PCRedundantSetScatter_C",(PC,VecScatter,VecScatter),(pc,in,out));CHKERRQ(ierr);
3474b9ad928SBarry Smith   PetscFunctionReturn(0);
3484b9ad928SBarry Smith }
3494b9ad928SBarry Smith 
3504b9ad928SBarry Smith #undef __FUNCT__
35183ab6a24SBarry Smith #define __FUNCT__ "PCRedundantGetKSP_Redundant"
352f7a08781SBarry Smith static PetscErrorCode  PCRedundantGetKSP_Redundant(PC pc,KSP *innerksp)
3534b9ad928SBarry Smith {
3545f06b7aaSBarry Smith   PetscErrorCode ierr;
3554b9ad928SBarry Smith   PC_Redundant   *red = (PC_Redundant*)pc->data;
3565f06b7aaSBarry Smith   MPI_Comm       comm,subcomm;
3575f06b7aaSBarry Smith   const char     *prefix;
3584b9ad928SBarry Smith 
3594b9ad928SBarry Smith   PetscFunctionBegin;
3605f06b7aaSBarry Smith   if (!red->psubcomm) {
3615f06b7aaSBarry Smith     ierr = PetscObjectGetComm((PetscObject)pc,&comm);CHKERRQ(ierr);
362d8a68f86SHong Zhang     ierr = PetscSubcommCreate(comm,&red->psubcomm);CHKERRQ(ierr);
363d8a68f86SHong Zhang     ierr = PetscSubcommSetNumber(red->psubcomm,red->nsubcomm);CHKERRQ(ierr);
364d8a68f86SHong Zhang     ierr = PetscSubcommSetType(red->psubcomm,PETSC_SUBCOMM_INTERLACED);CHKERRQ(ierr);
3653bb1ff40SBarry Smith     ierr = PetscLogObjectMemory((PetscObject)pc,sizeof(PetscSubcomm));CHKERRQ(ierr);
3665f06b7aaSBarry Smith 
3675f06b7aaSBarry Smith     /* create a new PC that processors in each subcomm have copy of */
3685f06b7aaSBarry Smith     subcomm = red->psubcomm->comm;
3692fa5cd67SKarl Rupp 
3705f06b7aaSBarry Smith     ierr = KSPCreate(subcomm,&red->ksp);CHKERRQ(ierr);
3715f06b7aaSBarry Smith     ierr = PetscObjectIncrementTabLevel((PetscObject)red->ksp,(PetscObject)pc,1);CHKERRQ(ierr);
3723bb1ff40SBarry Smith     ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)red->ksp);CHKERRQ(ierr);
3735f06b7aaSBarry Smith     ierr = KSPSetType(red->ksp,KSPPREONLY);CHKERRQ(ierr);
3745f06b7aaSBarry Smith     ierr = KSPGetPC(red->ksp,&red->pc);CHKERRQ(ierr);
3755f06b7aaSBarry Smith     ierr = PCSetType(red->pc,PCLU);CHKERRQ(ierr);
3765f06b7aaSBarry Smith 
3775f06b7aaSBarry Smith     ierr = PCGetOptionsPrefix(pc,&prefix);CHKERRQ(ierr);
3785f06b7aaSBarry Smith     ierr = KSPSetOptionsPrefix(red->ksp,prefix);CHKERRQ(ierr);
3795f06b7aaSBarry Smith     ierr = KSPAppendOptionsPrefix(red->ksp,"redundant_");CHKERRQ(ierr);
3805f06b7aaSBarry Smith   }
38183ab6a24SBarry Smith   *innerksp = red->ksp;
3824b9ad928SBarry Smith   PetscFunctionReturn(0);
3834b9ad928SBarry Smith }
3844b9ad928SBarry Smith 
3854b9ad928SBarry Smith #undef __FUNCT__
38683ab6a24SBarry Smith #define __FUNCT__ "PCRedundantGetKSP"
3874b9ad928SBarry Smith /*@
38883ab6a24SBarry Smith    PCRedundantGetKSP - Gets the less parallel KSP created by the redundant PC.
3894b9ad928SBarry Smith 
3904b9ad928SBarry Smith    Not Collective
3914b9ad928SBarry Smith 
3924b9ad928SBarry Smith    Input Parameter:
3934b9ad928SBarry Smith .  pc - the preconditioner context
3944b9ad928SBarry Smith 
3954b9ad928SBarry Smith    Output Parameter:
39683ab6a24SBarry Smith .  innerksp - the KSP on the smaller set of processes
3974b9ad928SBarry Smith 
3984b9ad928SBarry Smith    Level: advanced
3994b9ad928SBarry Smith 
4004b9ad928SBarry Smith .keywords: PC, redundant solve
4014b9ad928SBarry Smith @*/
40283ab6a24SBarry Smith PetscErrorCode  PCRedundantGetKSP(PC pc,KSP *innerksp)
4034b9ad928SBarry Smith {
4044ac538c5SBarry Smith   PetscErrorCode ierr;
4054b9ad928SBarry Smith 
4064b9ad928SBarry Smith   PetscFunctionBegin;
4070700a824SBarry Smith   PetscValidHeaderSpecific(pc,PC_CLASSID,1);
40883ab6a24SBarry Smith   PetscValidPointer(innerksp,2);
40983ab6a24SBarry Smith   ierr = PetscTryMethod(pc,"PCRedundantGetKSP_C",(PC,KSP*),(pc,innerksp));CHKERRQ(ierr);
4104b9ad928SBarry Smith   PetscFunctionReturn(0);
4114b9ad928SBarry Smith }
4124b9ad928SBarry Smith 
4134b9ad928SBarry Smith #undef __FUNCT__
4144b9ad928SBarry Smith #define __FUNCT__ "PCRedundantGetOperators_Redundant"
415f7a08781SBarry Smith static PetscErrorCode  PCRedundantGetOperators_Redundant(PC pc,Mat *mat,Mat *pmat)
4164b9ad928SBarry Smith {
4174b9ad928SBarry Smith   PC_Redundant *red = (PC_Redundant*)pc->data;
4184b9ad928SBarry Smith 
4194b9ad928SBarry Smith   PetscFunctionBegin;
420b3804887SHong Zhang   if (mat)  *mat  = red->pmats;
421b3804887SHong Zhang   if (pmat) *pmat = red->pmats;
4224b9ad928SBarry Smith   PetscFunctionReturn(0);
4234b9ad928SBarry Smith }
4244b9ad928SBarry Smith 
4254b9ad928SBarry Smith #undef __FUNCT__
4264b9ad928SBarry Smith #define __FUNCT__ "PCRedundantGetOperators"
4274b9ad928SBarry Smith /*@
4284b9ad928SBarry Smith    PCRedundantGetOperators - gets the sequential matrix and preconditioner matrix
4294b9ad928SBarry Smith 
4304b9ad928SBarry Smith    Not Collective
4314b9ad928SBarry Smith 
4324b9ad928SBarry Smith    Input Parameter:
4334b9ad928SBarry Smith .  pc - the preconditioner context
4344b9ad928SBarry Smith 
4354b9ad928SBarry Smith    Output Parameters:
4364b9ad928SBarry Smith +  mat - the matrix
4374b9ad928SBarry Smith -  pmat - the (possibly different) preconditioner matrix
4384b9ad928SBarry Smith 
4394b9ad928SBarry Smith    Level: advanced
4404b9ad928SBarry Smith 
4414b9ad928SBarry Smith .keywords: PC, redundant solve
4424b9ad928SBarry Smith @*/
4437087cfbeSBarry Smith PetscErrorCode  PCRedundantGetOperators(PC pc,Mat *mat,Mat *pmat)
4444b9ad928SBarry Smith {
4454ac538c5SBarry Smith   PetscErrorCode ierr;
4464b9ad928SBarry Smith 
4474b9ad928SBarry Smith   PetscFunctionBegin;
4480700a824SBarry Smith   PetscValidHeaderSpecific(pc,PC_CLASSID,1);
4494482741eSBarry Smith   if (mat)  PetscValidPointer(mat,2);
4504482741eSBarry Smith   if (pmat) PetscValidPointer(pmat,3);
4514ac538c5SBarry Smith   ierr = PetscTryMethod(pc,"PCRedundantGetOperators_C",(PC,Mat*,Mat*),(pc,mat,pmat));CHKERRQ(ierr);
4524b9ad928SBarry Smith   PetscFunctionReturn(0);
4534b9ad928SBarry Smith }
4544b9ad928SBarry Smith 
4554b9ad928SBarry Smith /* -------------------------------------------------------------------------------------*/
45637a17b4dSBarry Smith /*MC
45783ab6a24SBarry Smith      PCREDUNDANT - Runs a KSP solver with preconditioner for the entire problem on subgroups of processors
45837a17b4dSBarry Smith 
45983ab6a24SBarry Smith      Options for the redundant preconditioners can be set with -redundant_pc_xxx for the redundant KSP with -redundant_ksp_xxx
46037a17b4dSBarry Smith 
46109391456SBarry Smith   Options Database:
4629b21d695SBarry Smith .  -pc_redundant_number <n> - number of redundant solves, for example if you are using 64 MPI processes and
4639b21d695SBarry Smith                               use an n of 4 there will be 4 parallel solves each on 16 = 64/4 processes.
46409391456SBarry Smith 
46537a17b4dSBarry Smith    Level: intermediate
46637a17b4dSBarry Smith 
46783ab6a24SBarry Smith    Notes: The default KSP is preonly and the default PC is LU.
46883ab6a24SBarry Smith 
46983ab6a24SBarry Smith    Developer Notes: Note that PCSetInitialGuessNonzero()  is not used by this class but likely should be.
4709cfaa89bSBarry Smith 
47137a17b4dSBarry Smith .seealso:  PCCreate(), PCSetType(), PCType (for list of available types), PCRedundantSetScatter(),
47283ab6a24SBarry Smith            PCRedundantGetKSP(), PCRedundantGetOperators(), PCRedundantSetNumber()
47337a17b4dSBarry Smith M*/
47437a17b4dSBarry Smith 
4754b9ad928SBarry Smith #undef __FUNCT__
4764b9ad928SBarry Smith #define __FUNCT__ "PCCreate_Redundant"
4778cc058d9SJed Brown PETSC_EXTERN PetscErrorCode PCCreate_Redundant(PC pc)
4784b9ad928SBarry Smith {
479dfbe8321SBarry Smith   PetscErrorCode ierr;
4804b9ad928SBarry Smith   PC_Redundant   *red;
48169db28dcSHong Zhang   PetscMPIInt    size;
4823f457be1SHong Zhang 
4834b9ad928SBarry Smith   PetscFunctionBegin;
48438f2d2fdSLisandro Dalcin   ierr = PetscNewLog(pc,PC_Redundant,&red);CHKERRQ(ierr);
485ce94432eSBarry Smith   ierr = MPI_Comm_size(PetscObjectComm((PetscObject)pc),&size);CHKERRQ(ierr);
4862fa5cd67SKarl Rupp 
48769db28dcSHong Zhang   red->nsubcomm       = size;
4884b9ad928SBarry Smith   red->useparallelmat = PETSC_TRUE;
4891fbd8f88SHong Zhang   pc->data            = (void*)red;
4904b9ad928SBarry Smith 
4914b9ad928SBarry Smith   pc->ops->apply          = PCApply_Redundant;
4924b9ad928SBarry Smith   pc->ops->applytranspose = 0;
4934b9ad928SBarry Smith   pc->ops->setup          = PCSetUp_Redundant;
4944b9ad928SBarry Smith   pc->ops->destroy        = PCDestroy_Redundant;
4951ea5a559SBarry Smith   pc->ops->reset          = PCReset_Redundant;
4964b9ad928SBarry Smith   pc->ops->setfromoptions = PCSetFromOptions_Redundant;
4974b9ad928SBarry Smith   pc->ops->view           = PCView_Redundant;
4982fa5cd67SKarl Rupp 
499bdf89e91SBarry Smith   ierr = PetscObjectComposeFunction((PetscObject)pc,"PCRedundantSetScatter_C",PCRedundantSetScatter_Redundant);CHKERRQ(ierr);
500bdf89e91SBarry Smith   ierr = PetscObjectComposeFunction((PetscObject)pc,"PCRedundantSetNumber_C",PCRedundantSetNumber_Redundant);CHKERRQ(ierr);
501bdf89e91SBarry Smith   ierr = PetscObjectComposeFunction((PetscObject)pc,"PCRedundantGetKSP_C",PCRedundantGetKSP_Redundant);CHKERRQ(ierr);
502bdf89e91SBarry Smith   ierr = PetscObjectComposeFunction((PetscObject)pc,"PCRedundantGetOperators_C",PCRedundantGetOperators_Redundant);CHKERRQ(ierr);
5034b9ad928SBarry Smith   PetscFunctionReturn(0);
5044b9ad928SBarry Smith }
505b2573a8aSBarry Smith 
506