xref: /petsc/src/ksp/pc/impls/redundant/redundant.c (revision fd7037dcf4346030654083f0b28040cd75512e8f)
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);
75d8a68f86SHong Zhang       ierr = PetscSubcommSetType(red->psubcomm,PETSC_SUBCOMM_INTERLACED);CHKERRQ(ierr);
76f68be91cSHong Zhang       /* enable runtime switch of psubcomm type, e.g., '-psubcomm_type contiguous */
77f68be91cSHong Zhang       ierr = PetscSubcommSetFromOptions(red->psubcomm);CHKERRQ(ierr);
78f68be91cSHong Zhang 
791fbd8f88SHong Zhang       ierr = PetscLogObjectMemory(pc,sizeof(PetscSubcomm));CHKERRQ(ierr);
801fbd8f88SHong Zhang 
811fbd8f88SHong Zhang       /* create a new PC that processors in each subcomm have copy of */
820d7810c8SBarry Smith       subcomm = red->psubcomm->comm;
832fa5cd67SKarl Rupp 
845f06b7aaSBarry Smith       ierr = KSPCreate(subcomm,&red->ksp);CHKERRQ(ierr);
855f06b7aaSBarry Smith       ierr = PetscObjectIncrementTabLevel((PetscObject)red->ksp,(PetscObject)pc,1);CHKERRQ(ierr);
865f06b7aaSBarry Smith       ierr = PetscLogObjectParent(pc,red->ksp);CHKERRQ(ierr);
875f06b7aaSBarry Smith       ierr = KSPSetType(red->ksp,KSPPREONLY);CHKERRQ(ierr);
885f06b7aaSBarry Smith       ierr = KSPGetPC(red->ksp,&red->pc);CHKERRQ(ierr);
89cf52b8b1SHong Zhang       ierr = PCSetType(red->pc,PCLU);CHKERRQ(ierr);
90cf52b8b1SHong Zhang 
911fbd8f88SHong Zhang       ierr = PCGetOptionsPrefix(pc,&prefix);CHKERRQ(ierr);
925f06b7aaSBarry Smith       ierr = KSPSetOptionsPrefix(red->ksp,prefix);CHKERRQ(ierr);
935f06b7aaSBarry Smith       ierr = KSPAppendOptionsPrefix(red->ksp,"redundant_");CHKERRQ(ierr);
941b81debcSHong Zhang     } else {
951b81debcSHong Zhang       subcomm = red->psubcomm->comm;
961b81debcSHong Zhang     }
971fbd8f88SHong Zhang 
981b81debcSHong Zhang     if (red->useparallelmat) {
991b81debcSHong Zhang       /* grab the parallel matrix and put it into processors of a subcomminicator */
100*fd7037dcSHong Zhang       ierr = MatGetRedundantMatrix(pc->pmat,red->psubcomm->n,MPI_COMM_NULL,red->psubcomm,MAT_INITIAL_MATRIX,&red->pmats);CHKERRQ(ierr);
1011b81debcSHong Zhang       ierr = KSPSetOperators(red->ksp,red->pmats,red->pmats,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
1024b9ad928SBarry Smith 
1031b81debcSHong Zhang       /* get working vectors xsub and ysub */
1041b81debcSHong Zhang       ierr = MatGetVecs(red->pmats,&red->xsub,&red->ysub);CHKERRQ(ierr);
1052fa5cd67SKarl Rupp 
1068b96b0d2SHong Zhang       /* create working vectors xdup and ydup.
1078b96b0d2SHong Zhang        xdup concatenates all xsub's contigously to form a mpi vector over dupcomm  (see PetscSubcommCreate_interlaced())
1088b96b0d2SHong Zhang        ydup concatenates all ysub and has empty local arrays because ysub's arrays will be place into it.
109ce94432eSBarry Smith        Note: we use communicator dupcomm, not PetscObjectComm((PetscObject)pc)! */
1101b81debcSHong Zhang       ierr = MatGetLocalSize(red->pmats,&mloc_sub,NULL);CHKERRQ(ierr);
1111fbd8f88SHong Zhang       ierr = VecCreateMPI(red->psubcomm->dupparent,mloc_sub,PETSC_DECIDE,&red->xdup);CHKERRQ(ierr);
1120298fd71SBarry Smith       ierr = VecCreateMPIWithArray(red->psubcomm->dupparent,1,mloc_sub,PETSC_DECIDE,NULL,&red->ydup);CHKERRQ(ierr);
1133f457be1SHong Zhang 
114f68be91cSHong Zhang       /* create vecscatters */
115f68be91cSHong Zhang       if (!red->scatterin) { /* efficiency of scatterin is independent from psubcomm_type! */
1163f457be1SHong Zhang         IS       is1,is2;
1173f457be1SHong Zhang         PetscInt *idx1,*idx2,i,j,k;
11845fc02eaSBarry Smith 
1191b81debcSHong Zhang         ierr = MatGetVecs(pc->pmat,&x,0);CHKERRQ(ierr);
1201b81debcSHong Zhang         ierr = VecGetSize(x,&M);CHKERRQ(ierr);
1211b81debcSHong Zhang         ierr = VecGetOwnershipRange(x,&mstart,&mend);CHKERRQ(ierr);
1221b81debcSHong Zhang         mlocal = mend - mstart;
1231714dc9eSHong Zhang         ierr = PetscMalloc2(red->psubcomm->n*mlocal,PetscInt,&idx1,red->psubcomm->n*mlocal,PetscInt,&idx2);CHKERRQ(ierr);
1243f457be1SHong Zhang         j    = 0;
1251fbd8f88SHong Zhang         for (k=0; k<red->psubcomm->n; k++) {
1263f457be1SHong Zhang           for (i=mstart; i<mend; i++) {
1273f457be1SHong Zhang             idx1[j]   = i;
128ddc54837SHong Zhang             idx2[j++] = i + M*k;
1293f457be1SHong Zhang           }
1303f457be1SHong Zhang         }
13170b3c8c7SBarry Smith         ierr = ISCreateGeneral(comm,red->psubcomm->n*mlocal,idx1,PETSC_COPY_VALUES,&is1);CHKERRQ(ierr);
13270b3c8c7SBarry Smith         ierr = ISCreateGeneral(comm,red->psubcomm->n*mlocal,idx2,PETSC_COPY_VALUES,&is2);CHKERRQ(ierr);
133ddc54837SHong Zhang         ierr = VecScatterCreate(x,is1,red->xdup,is2,&red->scatterin);CHKERRQ(ierr);
134fcfd50ebSBarry Smith         ierr = ISDestroy(&is1);CHKERRQ(ierr);
135fcfd50ebSBarry Smith         ierr = ISDestroy(&is2);CHKERRQ(ierr);
1363f457be1SHong Zhang 
137f68be91cSHong Zhang         /* efficiency of scatterout depends on psubcomm_type! Impl below is good for PETSC_SUBCOMM_INTERLACED */
138ddc54837SHong Zhang         ierr = ISCreateStride(comm,mlocal,mstart+ red->psubcomm->color*M,1,&is1);CHKERRQ(ierr);
1393f457be1SHong Zhang         ierr = ISCreateStride(comm,mlocal,mstart,1,&is2);CHKERRQ(ierr);
140ddc54837SHong Zhang         ierr = VecScatterCreate(red->xdup,is1,x,is2,&red->scatterout);CHKERRQ(ierr);
141fcfd50ebSBarry Smith         ierr = ISDestroy(&is1);CHKERRQ(ierr);
142fcfd50ebSBarry Smith         ierr = ISDestroy(&is2);CHKERRQ(ierr);
1431d79065fSBarry Smith         ierr = PetscFree2(idx1,idx2);CHKERRQ(ierr);
144ddc54837SHong Zhang         ierr = VecDestroy(&x);CHKERRQ(ierr);
1451b81debcSHong Zhang       }
146ab661555SHong Zhang     } else { /* !red->useparallelmat */
147ab661555SHong Zhang       ierr = KSPSetOperators(red->ksp,pc->mat,pc->pmat,pc->flag);CHKERRQ(ierr);
1481b81debcSHong Zhang     }
149ab661555SHong Zhang   } else { /* pc->setupcalled */
1504b9ad928SBarry Smith     if (red->useparallelmat) {
151ab661555SHong Zhang       MatReuse       reuse;
1521b81debcSHong Zhang       /* grab the parallel matrix and put it into processors of a subcomminicator */
1531b81debcSHong Zhang       /*--------------------------------------------------------------------------*/
154ab661555SHong Zhang       if (pc->flag == DIFFERENT_NONZERO_PATTERN) {
1554b9ad928SBarry Smith         /* destroy old matrices */
1566bf464f9SBarry Smith         ierr  = MatDestroy(&red->pmats);CHKERRQ(ierr);
157ab661555SHong Zhang         reuse = MAT_INITIAL_MATRIX;
1584b9ad928SBarry Smith       } else {
159ab661555SHong Zhang         reuse = MAT_REUSE_MATRIX;
160ab661555SHong Zhang       }
161*fd7037dcSHong Zhang       ierr = MatGetRedundantMatrix(pc->pmat,red->psubcomm->n,MPI_COMM_NULL,red->psubcomm,reuse,&red->pmats);CHKERRQ(ierr);
162ab661555SHong Zhang       ierr = KSPSetOperators(red->ksp,red->pmats,red->pmats,pc->flag);CHKERRQ(ierr);
163ab661555SHong Zhang     } else { /* !red->useparallelmat */
16490f1c854SHong Zhang       ierr = KSPSetOperators(red->ksp,pc->mat,pc->pmat,pc->flag);CHKERRQ(ierr);
1654b9ad928SBarry Smith     }
166ab661555SHong Zhang   }
1671b81debcSHong Zhang 
1680c24e6a1SHong Zhang   if (pc->setfromoptionscalled) {
1693e065800SHong Zhang     ierr = KSPSetFromOptions(red->ksp);CHKERRQ(ierr);
1700c24e6a1SHong Zhang   }
1713e065800SHong Zhang   ierr = KSPSetUp(red->ksp);CHKERRQ(ierr);
1724b9ad928SBarry Smith   PetscFunctionReturn(0);
1734b9ad928SBarry Smith }
1744b9ad928SBarry Smith 
1754b9ad928SBarry Smith #undef __FUNCT__
1764b9ad928SBarry Smith #define __FUNCT__ "PCApply_Redundant"
1776849ba73SBarry Smith static PetscErrorCode PCApply_Redundant(PC pc,Vec x,Vec y)
1784b9ad928SBarry Smith {
1794b9ad928SBarry Smith   PC_Redundant   *red = (PC_Redundant*)pc->data;
180dfbe8321SBarry Smith   PetscErrorCode ierr;
1813f457be1SHong Zhang   PetscScalar    *array;
1824b9ad928SBarry Smith 
1834b9ad928SBarry Smith   PetscFunctionBegin;
184ddc54837SHong Zhang   if (!red->useparallelmat) {
185ddc54837SHong Zhang     ierr = KSPSolve(red->ksp,x,y);CHKERRQ(ierr);
186ddc54837SHong Zhang     PetscFunctionReturn(0);
187ddc54837SHong Zhang   }
188ddc54837SHong Zhang 
1893f457be1SHong Zhang   /* scatter x to xdup */
190ca9f406cSSatish Balay   ierr = VecScatterBegin(red->scatterin,x,red->xdup,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
191ca9f406cSSatish Balay   ierr = VecScatterEnd(red->scatterin,x,red->xdup,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
1923f457be1SHong Zhang 
1933f457be1SHong Zhang   /* place xdup's local array into xsub */
1943f457be1SHong Zhang   ierr = VecGetArray(red->xdup,&array);CHKERRQ(ierr);
1953f457be1SHong Zhang   ierr = VecPlaceArray(red->xsub,(const PetscScalar*)array);CHKERRQ(ierr);
1964b9ad928SBarry Smith 
1974b9ad928SBarry Smith   /* apply preconditioner on each processor */
19883ab6a24SBarry Smith   ierr = KSPSolve(red->ksp,red->xsub,red->ysub);CHKERRQ(ierr);
1993f457be1SHong Zhang   ierr = VecResetArray(red->xsub);CHKERRQ(ierr);
2003f457be1SHong Zhang   ierr = VecRestoreArray(red->xdup,&array);CHKERRQ(ierr);
2014b9ad928SBarry Smith 
2023f457be1SHong Zhang   /* place ysub's local array into ydup */
2033f457be1SHong Zhang   ierr = VecGetArray(red->ysub,&array);CHKERRQ(ierr);
2043f457be1SHong Zhang   ierr = VecPlaceArray(red->ydup,(const PetscScalar*)array);CHKERRQ(ierr);
2053f457be1SHong Zhang 
2063f457be1SHong Zhang   /* scatter ydup to y */
207ca9f406cSSatish Balay   ierr = VecScatterBegin(red->scatterout,red->ydup,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
208ca9f406cSSatish Balay   ierr = VecScatterEnd(red->scatterout,red->ydup,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
2093f457be1SHong Zhang   ierr = VecResetArray(red->ydup);CHKERRQ(ierr);
2103f457be1SHong Zhang   ierr = VecRestoreArray(red->ysub,&array);CHKERRQ(ierr);
2114b9ad928SBarry Smith   PetscFunctionReturn(0);
2124b9ad928SBarry Smith }
2134b9ad928SBarry Smith 
2144b9ad928SBarry Smith #undef __FUNCT__
2151ea5a559SBarry Smith #define __FUNCT__ "PCReset_Redundant"
2161ea5a559SBarry Smith static PetscErrorCode PCReset_Redundant(PC pc)
2174b9ad928SBarry Smith {
2184b9ad928SBarry Smith   PC_Redundant   *red = (PC_Redundant*)pc->data;
219dfbe8321SBarry Smith   PetscErrorCode ierr;
2204b9ad928SBarry Smith 
2214b9ad928SBarry Smith   PetscFunctionBegin;
2221b81debcSHong Zhang   if (red->useparallelmat) {
2236bf464f9SBarry Smith     ierr = VecScatterDestroy(&red->scatterin);CHKERRQ(ierr);
2246bf464f9SBarry Smith     ierr = VecScatterDestroy(&red->scatterout);CHKERRQ(ierr);
2256bf464f9SBarry Smith     ierr = VecDestroy(&red->ysub);CHKERRQ(ierr);
2266bf464f9SBarry Smith     ierr = VecDestroy(&red->xsub);CHKERRQ(ierr);
2276bf464f9SBarry Smith     ierr = VecDestroy(&red->xdup);CHKERRQ(ierr);
2286bf464f9SBarry Smith     ierr = VecDestroy(&red->ydup);CHKERRQ(ierr);
2291b81debcSHong Zhang   }
2306bf464f9SBarry Smith   ierr = MatDestroy(&red->pmats);CHKERRQ(ierr);
2311b81debcSHong Zhang   ierr = KSPReset(red->ksp);CHKERRQ(ierr);
2321ea5a559SBarry Smith   PetscFunctionReturn(0);
2331ea5a559SBarry Smith }
2341ea5a559SBarry Smith 
2351ea5a559SBarry Smith #undef __FUNCT__
2361ea5a559SBarry Smith #define __FUNCT__ "PCDestroy_Redundant"
2371ea5a559SBarry Smith static PetscErrorCode PCDestroy_Redundant(PC pc)
2381ea5a559SBarry Smith {
2391ea5a559SBarry Smith   PC_Redundant   *red = (PC_Redundant*)pc->data;
2401ea5a559SBarry Smith   PetscErrorCode ierr;
2411ea5a559SBarry Smith 
2421ea5a559SBarry Smith   PetscFunctionBegin;
2431ea5a559SBarry Smith   ierr = PCReset_Redundant(pc);CHKERRQ(ierr);
2446bf464f9SBarry Smith   ierr = KSPDestroy(&red->ksp);CHKERRQ(ierr);
2456bf464f9SBarry Smith   ierr = PetscSubcommDestroy(&red->psubcomm);CHKERRQ(ierr);
246c31cb41cSBarry Smith   ierr = PetscFree(pc->data);CHKERRQ(ierr);
2474b9ad928SBarry Smith   PetscFunctionReturn(0);
2484b9ad928SBarry Smith }
2494b9ad928SBarry Smith 
2504b9ad928SBarry Smith #undef __FUNCT__
2514b9ad928SBarry Smith #define __FUNCT__ "PCSetFromOptions_Redundant"
2526849ba73SBarry Smith static PetscErrorCode PCSetFromOptions_Redundant(PC pc)
2534b9ad928SBarry Smith {
254a98ce0f4SHong Zhang   PetscErrorCode ierr;
255a98ce0f4SHong Zhang   PC_Redundant   *red = (PC_Redundant*)pc->data;
256a98ce0f4SHong Zhang 
2574b9ad928SBarry Smith   PetscFunctionBegin;
258a98ce0f4SHong Zhang   ierr = PetscOptionsHead("Redundant options");CHKERRQ(ierr);
25909a6bc64SHong Zhang   ierr = PetscOptionsInt("-pc_redundant_number","Number of redundant pc","PCRedundantSetNumber",red->nsubcomm,&red->nsubcomm,0);CHKERRQ(ierr);
260a98ce0f4SHong Zhang   ierr = PetscOptionsTail();CHKERRQ(ierr);
2614b9ad928SBarry Smith   PetscFunctionReturn(0);
2624b9ad928SBarry Smith }
2634b9ad928SBarry Smith 
2644b9ad928SBarry Smith #undef __FUNCT__
26509a6bc64SHong Zhang #define __FUNCT__ "PCRedundantSetNumber_Redundant"
266f7a08781SBarry Smith static PetscErrorCode  PCRedundantSetNumber_Redundant(PC pc,PetscInt nreds)
26709a6bc64SHong Zhang {
26809a6bc64SHong Zhang   PC_Redundant *red = (PC_Redundant*)pc->data;
26909a6bc64SHong Zhang 
27009a6bc64SHong Zhang   PetscFunctionBegin;
27109a6bc64SHong Zhang   red->nsubcomm = nreds;
27209a6bc64SHong Zhang   PetscFunctionReturn(0);
27309a6bc64SHong Zhang }
27409a6bc64SHong Zhang 
27509a6bc64SHong Zhang #undef __FUNCT__
27609a6bc64SHong Zhang #define __FUNCT__ "PCRedundantSetNumber"
27709a6bc64SHong Zhang /*@
27809a6bc64SHong Zhang    PCRedundantSetNumber - Sets the number of redundant preconditioner contexts.
27909a6bc64SHong Zhang 
2803f9fe445SBarry Smith    Logically Collective on PC
28109a6bc64SHong Zhang 
28209a6bc64SHong Zhang    Input Parameters:
28309a6bc64SHong Zhang +  pc - the preconditioner context
2849b21d695SBarry Smith -  nredundant - number of redundant preconditioner contexts; for example if you are using 64 MPI processes and
2859b21d695SBarry Smith                               use an nredundant of 4 there will be 4 parallel solves each on 16 = 64/4 processes.
28609a6bc64SHong Zhang 
28709a6bc64SHong Zhang    Level: advanced
28809a6bc64SHong Zhang 
28909a6bc64SHong Zhang .keywords: PC, redundant solve
29009a6bc64SHong Zhang @*/
2917087cfbeSBarry Smith PetscErrorCode  PCRedundantSetNumber(PC pc,PetscInt nredundant)
29209a6bc64SHong Zhang {
2934ac538c5SBarry Smith   PetscErrorCode ierr;
29409a6bc64SHong Zhang 
29509a6bc64SHong Zhang   PetscFunctionBegin;
2960700a824SBarry Smith   PetscValidHeaderSpecific(pc,PC_CLASSID,1);
297ce94432eSBarry Smith   if (nredundant <= 0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONG, "num of redundant pc %D must be positive",nredundant);
2984ac538c5SBarry Smith   ierr = PetscTryMethod(pc,"PCRedundantSetNumber_C",(PC,PetscInt),(pc,nredundant));CHKERRQ(ierr);
29909a6bc64SHong Zhang   PetscFunctionReturn(0);
30009a6bc64SHong Zhang }
30109a6bc64SHong Zhang 
30209a6bc64SHong Zhang #undef __FUNCT__
3034b9ad928SBarry Smith #define __FUNCT__ "PCRedundantSetScatter_Redundant"
304f7a08781SBarry Smith static PetscErrorCode  PCRedundantSetScatter_Redundant(PC pc,VecScatter in,VecScatter out)
3054b9ad928SBarry Smith {
3064b9ad928SBarry Smith   PC_Redundant   *red = (PC_Redundant*)pc->data;
307dfbe8321SBarry Smith   PetscErrorCode ierr;
3084b9ad928SBarry Smith 
3094b9ad928SBarry Smith   PetscFunctionBegin;
3104b9ad928SBarry Smith   ierr = PetscObjectReference((PetscObject)in);CHKERRQ(ierr);
3116bf464f9SBarry Smith   ierr = VecScatterDestroy(&red->scatterin);CHKERRQ(ierr);
3122fa5cd67SKarl Rupp 
313c3122656SLisandro Dalcin   red->scatterin  = in;
3142fa5cd67SKarl Rupp 
3154b9ad928SBarry Smith   ierr = PetscObjectReference((PetscObject)out);CHKERRQ(ierr);
3166bf464f9SBarry Smith   ierr = VecScatterDestroy(&red->scatterout);CHKERRQ(ierr);
317c3122656SLisandro Dalcin   red->scatterout = out;
3184b9ad928SBarry Smith   PetscFunctionReturn(0);
3194b9ad928SBarry Smith }
3204b9ad928SBarry Smith 
3214b9ad928SBarry Smith #undef __FUNCT__
3224b9ad928SBarry Smith #define __FUNCT__ "PCRedundantSetScatter"
3234b9ad928SBarry Smith /*@
3244b9ad928SBarry Smith    PCRedundantSetScatter - Sets the scatter used to copy values into the
3254b9ad928SBarry Smith      redundant local solve and the scatter to move them back into the global
3264b9ad928SBarry Smith      vector.
3274b9ad928SBarry Smith 
3283f9fe445SBarry Smith    Logically Collective on PC
3294b9ad928SBarry Smith 
3304b9ad928SBarry Smith    Input Parameters:
3314b9ad928SBarry Smith +  pc - the preconditioner context
3324b9ad928SBarry Smith .  in - the scatter to move the values in
3334b9ad928SBarry Smith -  out - the scatter to move them out
3344b9ad928SBarry Smith 
3354b9ad928SBarry Smith    Level: advanced
3364b9ad928SBarry Smith 
3374b9ad928SBarry Smith .keywords: PC, redundant solve
3384b9ad928SBarry Smith @*/
3397087cfbeSBarry Smith PetscErrorCode  PCRedundantSetScatter(PC pc,VecScatter in,VecScatter out)
3404b9ad928SBarry Smith {
3414ac538c5SBarry Smith   PetscErrorCode ierr;
3424b9ad928SBarry Smith 
3434b9ad928SBarry Smith   PetscFunctionBegin;
3440700a824SBarry Smith   PetscValidHeaderSpecific(pc,PC_CLASSID,1);
3450700a824SBarry Smith   PetscValidHeaderSpecific(in,VEC_SCATTER_CLASSID,2);
3460700a824SBarry Smith   PetscValidHeaderSpecific(out,VEC_SCATTER_CLASSID,3);
3474ac538c5SBarry Smith   ierr = PetscTryMethod(pc,"PCRedundantSetScatter_C",(PC,VecScatter,VecScatter),(pc,in,out));CHKERRQ(ierr);
3484b9ad928SBarry Smith   PetscFunctionReturn(0);
3494b9ad928SBarry Smith }
3504b9ad928SBarry Smith 
3514b9ad928SBarry Smith #undef __FUNCT__
35283ab6a24SBarry Smith #define __FUNCT__ "PCRedundantGetKSP_Redundant"
353f7a08781SBarry Smith static PetscErrorCode  PCRedundantGetKSP_Redundant(PC pc,KSP *innerksp)
3544b9ad928SBarry Smith {
3555f06b7aaSBarry Smith   PetscErrorCode ierr;
3564b9ad928SBarry Smith   PC_Redundant   *red = (PC_Redundant*)pc->data;
3575f06b7aaSBarry Smith   MPI_Comm       comm,subcomm;
3585f06b7aaSBarry Smith   const char     *prefix;
3594b9ad928SBarry Smith 
3604b9ad928SBarry Smith   PetscFunctionBegin;
3615f06b7aaSBarry Smith   if (!red->psubcomm) {
3625f06b7aaSBarry Smith     ierr = PetscObjectGetComm((PetscObject)pc,&comm);CHKERRQ(ierr);
363d8a68f86SHong Zhang     ierr = PetscSubcommCreate(comm,&red->psubcomm);CHKERRQ(ierr);
364d8a68f86SHong Zhang     ierr = PetscSubcommSetNumber(red->psubcomm,red->nsubcomm);CHKERRQ(ierr);
365d8a68f86SHong Zhang     ierr = PetscSubcommSetType(red->psubcomm,PETSC_SUBCOMM_INTERLACED);CHKERRQ(ierr);
3665f06b7aaSBarry Smith     ierr = PetscLogObjectMemory(pc,sizeof(PetscSubcomm));CHKERRQ(ierr);
3675f06b7aaSBarry Smith 
3685f06b7aaSBarry Smith     /* create a new PC that processors in each subcomm have copy of */
3695f06b7aaSBarry Smith     subcomm = red->psubcomm->comm;
3702fa5cd67SKarl Rupp 
3715f06b7aaSBarry Smith     ierr = KSPCreate(subcomm,&red->ksp);CHKERRQ(ierr);
3725f06b7aaSBarry Smith     ierr = PetscObjectIncrementTabLevel((PetscObject)red->ksp,(PetscObject)pc,1);CHKERRQ(ierr);
3735f06b7aaSBarry Smith     ierr = PetscLogObjectParent(pc,red->ksp);CHKERRQ(ierr);
3745f06b7aaSBarry Smith     ierr = KSPSetType(red->ksp,KSPPREONLY);CHKERRQ(ierr);
3755f06b7aaSBarry Smith     ierr = KSPGetPC(red->ksp,&red->pc);CHKERRQ(ierr);
3765f06b7aaSBarry Smith     ierr = PCSetType(red->pc,PCLU);CHKERRQ(ierr);
3775f06b7aaSBarry Smith 
3785f06b7aaSBarry Smith     ierr = PCGetOptionsPrefix(pc,&prefix);CHKERRQ(ierr);
3795f06b7aaSBarry Smith     ierr = KSPSetOptionsPrefix(red->ksp,prefix);CHKERRQ(ierr);
3805f06b7aaSBarry Smith     ierr = KSPAppendOptionsPrefix(red->ksp,"redundant_");CHKERRQ(ierr);
3815f06b7aaSBarry Smith   }
38283ab6a24SBarry Smith   *innerksp = red->ksp;
3834b9ad928SBarry Smith   PetscFunctionReturn(0);
3844b9ad928SBarry Smith }
3854b9ad928SBarry Smith 
3864b9ad928SBarry Smith #undef __FUNCT__
38783ab6a24SBarry Smith #define __FUNCT__ "PCRedundantGetKSP"
3884b9ad928SBarry Smith /*@
38983ab6a24SBarry Smith    PCRedundantGetKSP - Gets the less parallel KSP created by the redundant PC.
3904b9ad928SBarry Smith 
3914b9ad928SBarry Smith    Not Collective
3924b9ad928SBarry Smith 
3934b9ad928SBarry Smith    Input Parameter:
3944b9ad928SBarry Smith .  pc - the preconditioner context
3954b9ad928SBarry Smith 
3964b9ad928SBarry Smith    Output Parameter:
39783ab6a24SBarry Smith .  innerksp - the KSP on the smaller set of processes
3984b9ad928SBarry Smith 
3994b9ad928SBarry Smith    Level: advanced
4004b9ad928SBarry Smith 
4014b9ad928SBarry Smith .keywords: PC, redundant solve
4024b9ad928SBarry Smith @*/
40383ab6a24SBarry Smith PetscErrorCode  PCRedundantGetKSP(PC pc,KSP *innerksp)
4044b9ad928SBarry Smith {
4054ac538c5SBarry Smith   PetscErrorCode ierr;
4064b9ad928SBarry Smith 
4074b9ad928SBarry Smith   PetscFunctionBegin;
4080700a824SBarry Smith   PetscValidHeaderSpecific(pc,PC_CLASSID,1);
40983ab6a24SBarry Smith   PetscValidPointer(innerksp,2);
41083ab6a24SBarry Smith   ierr = PetscTryMethod(pc,"PCRedundantGetKSP_C",(PC,KSP*),(pc,innerksp));CHKERRQ(ierr);
4114b9ad928SBarry Smith   PetscFunctionReturn(0);
4124b9ad928SBarry Smith }
4134b9ad928SBarry Smith 
4144b9ad928SBarry Smith #undef __FUNCT__
4154b9ad928SBarry Smith #define __FUNCT__ "PCRedundantGetOperators_Redundant"
416f7a08781SBarry Smith static PetscErrorCode  PCRedundantGetOperators_Redundant(PC pc,Mat *mat,Mat *pmat)
4174b9ad928SBarry Smith {
4184b9ad928SBarry Smith   PC_Redundant *red = (PC_Redundant*)pc->data;
4194b9ad928SBarry Smith 
4204b9ad928SBarry Smith   PetscFunctionBegin;
421b3804887SHong Zhang   if (mat)  *mat  = red->pmats;
422b3804887SHong Zhang   if (pmat) *pmat = red->pmats;
4234b9ad928SBarry Smith   PetscFunctionReturn(0);
4244b9ad928SBarry Smith }
4254b9ad928SBarry Smith 
4264b9ad928SBarry Smith #undef __FUNCT__
4274b9ad928SBarry Smith #define __FUNCT__ "PCRedundantGetOperators"
4284b9ad928SBarry Smith /*@
4294b9ad928SBarry Smith    PCRedundantGetOperators - gets the sequential matrix and preconditioner matrix
4304b9ad928SBarry Smith 
4314b9ad928SBarry Smith    Not Collective
4324b9ad928SBarry Smith 
4334b9ad928SBarry Smith    Input Parameter:
4344b9ad928SBarry Smith .  pc - the preconditioner context
4354b9ad928SBarry Smith 
4364b9ad928SBarry Smith    Output Parameters:
4374b9ad928SBarry Smith +  mat - the matrix
4384b9ad928SBarry Smith -  pmat - the (possibly different) preconditioner matrix
4394b9ad928SBarry Smith 
4404b9ad928SBarry Smith    Level: advanced
4414b9ad928SBarry Smith 
4424b9ad928SBarry Smith .keywords: PC, redundant solve
4434b9ad928SBarry Smith @*/
4447087cfbeSBarry Smith PetscErrorCode  PCRedundantGetOperators(PC pc,Mat *mat,Mat *pmat)
4454b9ad928SBarry Smith {
4464ac538c5SBarry Smith   PetscErrorCode ierr;
4474b9ad928SBarry Smith 
4484b9ad928SBarry Smith   PetscFunctionBegin;
4490700a824SBarry Smith   PetscValidHeaderSpecific(pc,PC_CLASSID,1);
4504482741eSBarry Smith   if (mat)  PetscValidPointer(mat,2);
4514482741eSBarry Smith   if (pmat) PetscValidPointer(pmat,3);
4524ac538c5SBarry Smith   ierr = PetscTryMethod(pc,"PCRedundantGetOperators_C",(PC,Mat*,Mat*),(pc,mat,pmat));CHKERRQ(ierr);
4534b9ad928SBarry Smith   PetscFunctionReturn(0);
4544b9ad928SBarry Smith }
4554b9ad928SBarry Smith 
4564b9ad928SBarry Smith /* -------------------------------------------------------------------------------------*/
45737a17b4dSBarry Smith /*MC
45883ab6a24SBarry Smith      PCREDUNDANT - Runs a KSP solver with preconditioner for the entire problem on subgroups of processors
45937a17b4dSBarry Smith 
46083ab6a24SBarry Smith      Options for the redundant preconditioners can be set with -redundant_pc_xxx for the redundant KSP with -redundant_ksp_xxx
46137a17b4dSBarry Smith 
46209391456SBarry Smith   Options Database:
4639b21d695SBarry Smith .  -pc_redundant_number <n> - number of redundant solves, for example if you are using 64 MPI processes and
4649b21d695SBarry Smith                               use an n of 4 there will be 4 parallel solves each on 16 = 64/4 processes.
46509391456SBarry Smith 
46637a17b4dSBarry Smith    Level: intermediate
46737a17b4dSBarry Smith 
46883ab6a24SBarry Smith    Notes: The default KSP is preonly and the default PC is LU.
46983ab6a24SBarry Smith 
47083ab6a24SBarry Smith    Developer Notes: Note that PCSetInitialGuessNonzero()  is not used by this class but likely should be.
4719cfaa89bSBarry Smith 
47237a17b4dSBarry Smith .seealso:  PCCreate(), PCSetType(), PCType (for list of available types), PCRedundantSetScatter(),
47383ab6a24SBarry Smith            PCRedundantGetKSP(), PCRedundantGetOperators(), PCRedundantSetNumber()
47437a17b4dSBarry Smith M*/
47537a17b4dSBarry Smith 
4764b9ad928SBarry Smith #undef __FUNCT__
4774b9ad928SBarry Smith #define __FUNCT__ "PCCreate_Redundant"
4788cc058d9SJed Brown PETSC_EXTERN PetscErrorCode PCCreate_Redundant(PC pc)
4794b9ad928SBarry Smith {
480dfbe8321SBarry Smith   PetscErrorCode ierr;
4814b9ad928SBarry Smith   PC_Redundant   *red;
48269db28dcSHong Zhang   PetscMPIInt    size;
4833f457be1SHong Zhang 
4844b9ad928SBarry Smith   PetscFunctionBegin;
48538f2d2fdSLisandro Dalcin   ierr = PetscNewLog(pc,PC_Redundant,&red);CHKERRQ(ierr);
486ce94432eSBarry Smith   ierr = MPI_Comm_size(PetscObjectComm((PetscObject)pc),&size);CHKERRQ(ierr);
4872fa5cd67SKarl Rupp 
48869db28dcSHong Zhang   red->nsubcomm       = size;
4894b9ad928SBarry Smith   red->useparallelmat = PETSC_TRUE;
4901fbd8f88SHong Zhang   pc->data            = (void*)red;
4914b9ad928SBarry Smith 
4924b9ad928SBarry Smith   pc->ops->apply          = PCApply_Redundant;
4934b9ad928SBarry Smith   pc->ops->applytranspose = 0;
4944b9ad928SBarry Smith   pc->ops->setup          = PCSetUp_Redundant;
4954b9ad928SBarry Smith   pc->ops->destroy        = PCDestroy_Redundant;
4961ea5a559SBarry Smith   pc->ops->reset          = PCReset_Redundant;
4974b9ad928SBarry Smith   pc->ops->setfromoptions = PCSetFromOptions_Redundant;
4984b9ad928SBarry Smith   pc->ops->view           = PCView_Redundant;
4992fa5cd67SKarl Rupp 
500bdf89e91SBarry Smith   ierr = PetscObjectComposeFunction((PetscObject)pc,"PCRedundantSetScatter_C",PCRedundantSetScatter_Redundant);CHKERRQ(ierr);
501bdf89e91SBarry Smith   ierr = PetscObjectComposeFunction((PetscObject)pc,"PCRedundantSetNumber_C",PCRedundantSetNumber_Redundant);CHKERRQ(ierr);
502bdf89e91SBarry Smith   ierr = PetscObjectComposeFunction((PetscObject)pc,"PCRedundantGetKSP_C",PCRedundantGetKSP_Redundant);CHKERRQ(ierr);
503bdf89e91SBarry Smith   ierr = PetscObjectComposeFunction((PetscObject)pc,"PCRedundantGetOperators_C",PCRedundantGetOperators_Redundant);CHKERRQ(ierr);
5044b9ad928SBarry Smith   PetscFunctionReturn(0);
5054b9ad928SBarry Smith }
506b2573a8aSBarry Smith 
507