1dba47a55SKris Buschelman #define PETSCKSP_DLL 2dba47a55SKris Buschelman 34b9ad928SBarry Smith /* 43f457be1SHong Zhang This file defines a "solve the problem redundantly on each subgroup of processor" preconditioner. 54b9ad928SBarry Smith */ 66356e834SBarry Smith #include "private/pcimpl.h" /*I "petscpc.h" I*/ 74b9ad928SBarry Smith #include "petscksp.h" 84b9ad928SBarry Smith 94b9ad928SBarry Smith typedef struct { 101fbd8f88SHong Zhang MPI_Comm parent; /* parent communicator */ 111fbd8f88SHong Zhang MPI_Comm dupparent; /* duplicate parent communicator, under which the processors of this subcomm have contiguous rank */ 121fbd8f88SHong Zhang MPI_Comm comm; /* this communicator */ 131fbd8f88SHong Zhang PetscInt n; /* num of subcommunicators under the parent communicator */ 141fbd8f88SHong Zhang PetscInt color; /* color of processors belong to this communicator */ 151fbd8f88SHong Zhang } PetscSubcomm; 161fbd8f88SHong Zhang 171fbd8f88SHong Zhang #undef __FUNCT__ 181fbd8f88SHong Zhang #define __FUNCT__ "PetscSubcommDestroy" 191fbd8f88SHong Zhang PetscErrorCode PetscSubcommDestroy(PetscSubcomm *psubcomm) 201fbd8f88SHong Zhang { 211fbd8f88SHong Zhang PetscErrorCode ierr; 221fbd8f88SHong Zhang 231fbd8f88SHong Zhang PetscFunctionBegin; 241fbd8f88SHong Zhang ierr = PetscFree(psubcomm);CHKERRQ(ierr); 251fbd8f88SHong Zhang PetscFunctionReturn(0); 261fbd8f88SHong Zhang } 271fbd8f88SHong Zhang 281fbd8f88SHong Zhang /*-------------------------------------------------------------------------------------------------- 291fbd8f88SHong Zhang To avoid data scattering from subcomm back to original comm, we create subcomm by iteratively taking a 301fbd8f88SHong Zhang processe into a subcomm. 311fbd8f88SHong Zhang An example: size=4, nsubcomm=3 321fbd8f88SHong Zhang pc->comm: 331fbd8f88SHong Zhang rank: [0] [1] [2] [3] 341fbd8f88SHong Zhang color: 0 1 2 0 351fbd8f88SHong Zhang 361fbd8f88SHong Zhang subcomm: 371fbd8f88SHong Zhang subrank: [0] [0] [0] [1] 381fbd8f88SHong Zhang 391fbd8f88SHong Zhang dupcomm: 401fbd8f88SHong Zhang duprank: [0] [2] [3] [1] 411fbd8f88SHong Zhang 421fbd8f88SHong Zhang Here, subcomm[color = 0] has subsize=2, owns process [0] and [3] 431fbd8f88SHong Zhang subcomm[color = 1] has subsize=1, owns process [1] 441fbd8f88SHong Zhang subcomm[color = 2] has subsize=1, owns process [2] 451fbd8f88SHong Zhang dupcomm has same number of processes as pc->comm, and its duprank maps 461fbd8f88SHong Zhang processes in subcomm contiguously into a 1d array: 471fbd8f88SHong Zhang duprank: [0] [1] [2] [3] 481fbd8f88SHong Zhang rank: [0] [3] [1] [2] 491fbd8f88SHong Zhang subcomm[0] subcomm[1] subcomm[2] 501fbd8f88SHong Zhang ----------------------------------------------------------------------------------------*/ 511fbd8f88SHong Zhang #undef __FUNCT__ 521fbd8f88SHong Zhang #define __FUNCT__ "PetscSubcommCreate" 531fbd8f88SHong Zhang PetscErrorCode PetscSubcommCreate(MPI_Comm comm,PetscInt nsubcomm,PetscSubcomm **psubcomm) 541fbd8f88SHong Zhang { 551fbd8f88SHong Zhang PetscErrorCode ierr; 561fbd8f88SHong Zhang PetscMPIInt rank,size,*subsize,duprank,subrank; 571fbd8f88SHong Zhang PetscInt np_subcomm,nleftover,i,j,color; 581fbd8f88SHong Zhang MPI_Comm subcomm=0,dupcomm=0; 591fbd8f88SHong Zhang const char *prefix; 601fbd8f88SHong Zhang PetscSubcomm *psubcomm_tmp; 611fbd8f88SHong Zhang 621fbd8f88SHong Zhang PetscFunctionBegin; 631fbd8f88SHong Zhang ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 641fbd8f88SHong Zhang ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 651fbd8f88SHong Zhang if (nsubcomm < 1 || nsubcomm > size) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE, "Num of subcommunicators %D cannot be < 1 or > input comm size %D",nsubcomm,size); 661fbd8f88SHong Zhang 671fbd8f88SHong Zhang /* get size of each subcommunicator */ 681fbd8f88SHong Zhang ierr = PetscMalloc((1+nsubcomm)*sizeof(PetscMPIInt),&subsize);CHKERRQ(ierr); 691fbd8f88SHong Zhang np_subcomm = size/nsubcomm; 701fbd8f88SHong Zhang nleftover = size - nsubcomm*np_subcomm; 711fbd8f88SHong Zhang for (i=0; i<nsubcomm; i++){ 721fbd8f88SHong Zhang subsize[i] = np_subcomm; 731fbd8f88SHong Zhang if (i<nleftover) subsize[i]++; 741fbd8f88SHong Zhang } 751fbd8f88SHong Zhang 761fbd8f88SHong Zhang /* find color for this proc */ 771fbd8f88SHong Zhang color = rank%nsubcomm; 781fbd8f88SHong Zhang subrank = rank/nsubcomm; 791fbd8f88SHong Zhang 801fbd8f88SHong Zhang ierr = MPI_Comm_split(comm,color,subrank,&subcomm);CHKERRQ(ierr); 811fbd8f88SHong Zhang 821fbd8f88SHong Zhang j = 0; duprank = 0; 831fbd8f88SHong Zhang for (i=0; i<nsubcomm; i++){ 841fbd8f88SHong Zhang if (j == color){ 851fbd8f88SHong Zhang duprank += subrank; 861fbd8f88SHong Zhang break; 871fbd8f88SHong Zhang } 881fbd8f88SHong Zhang duprank += subsize[i]; j++; 891fbd8f88SHong Zhang } 901fbd8f88SHong Zhang 911fbd8f88SHong Zhang /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */ 921fbd8f88SHong Zhang ierr = MPI_Comm_split(comm,0,duprank,&dupcomm);CHKERRQ(ierr); 931fbd8f88SHong Zhang ierr = PetscFree(subsize);CHKERRQ(ierr); 941fbd8f88SHong Zhang 951fbd8f88SHong Zhang ierr = PetscNew(PetscSubcomm,&psubcomm_tmp);CHKERRQ(ierr); 961fbd8f88SHong Zhang psubcomm_tmp->parent = comm; 971fbd8f88SHong Zhang psubcomm_tmp->dupparent = dupcomm; 981fbd8f88SHong Zhang psubcomm_tmp->comm = subcomm; 991fbd8f88SHong Zhang psubcomm_tmp->n = nsubcomm; 1001fbd8f88SHong Zhang psubcomm_tmp->color = color; 1011fbd8f88SHong Zhang *psubcomm = psubcomm_tmp; 1021fbd8f88SHong Zhang PetscFunctionReturn(0); 1031fbd8f88SHong Zhang } 1041fbd8f88SHong Zhang 1051fbd8f88SHong Zhang typedef struct { 1064b9ad928SBarry Smith PC pc; /* actual preconditioner used on each processor */ 1073f457be1SHong Zhang Vec xsub,ysub; /* vectors of a subcommunicator to hold parallel vectors of pc->comm */ 1083f457be1SHong Zhang Vec xdup,ydup; /* parallel vector that congregates xsub or ysub facilitating vector scattering */ 109b3804887SHong Zhang Mat pmats; /* matrix and optional preconditioner matrix belong to a subcommunicator */ 1103f457be1SHong Zhang VecScatter scatterin,scatterout; /* scatter used to move all values to each processor group (subcommunicator) */ 1114b9ad928SBarry Smith PetscTruth useparallelmat; 1121fbd8f88SHong Zhang PetscSubcomm *psubcomm; 1131fbd8f88SHong Zhang PetscInt nsubcomm; /* num of data structure PetscSubcomm */ 1144b9ad928SBarry Smith } PC_Redundant; 1154b9ad928SBarry Smith 1164b9ad928SBarry Smith #undef __FUNCT__ 1174b9ad928SBarry Smith #define __FUNCT__ "PCView_Redundant" 1186849ba73SBarry Smith static PetscErrorCode PCView_Redundant(PC pc,PetscViewer viewer) 1194b9ad928SBarry Smith { 1204b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 121dfbe8321SBarry Smith PetscErrorCode ierr; 12213f74950SBarry Smith PetscMPIInt rank; 12332077d6dSBarry Smith PetscTruth iascii,isstring; 124a47c9f9aSHong Zhang PetscViewer sviewer,subviewer; 1251fbd8f88SHong Zhang PetscInt color = red->psubcomm->color; 1264b9ad928SBarry Smith 1274b9ad928SBarry Smith PetscFunctionBegin; 1284b9ad928SBarry Smith ierr = MPI_Comm_rank(pc->comm,&rank);CHKERRQ(ierr); 12932077d6dSBarry Smith ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr); 1304b9ad928SBarry Smith ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_STRING,&isstring);CHKERRQ(ierr); 13132077d6dSBarry Smith if (iascii) { 132a98ce0f4SHong Zhang ierr = PetscViewerASCIIPrintf(viewer," Redundant solver preconditioner: First PC (color=0) follows\n");CHKERRQ(ierr); 133a98ce0f4SHong Zhang ierr = PetscViewerGetSubcomm(viewer,red->pc->comm,&subviewer);CHKERRQ(ierr); 134a98ce0f4SHong Zhang if (!color) { /* only view first redundant pc */ 1354b9ad928SBarry Smith ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 136a47c9f9aSHong Zhang ierr = PCView(red->pc,subviewer);CHKERRQ(ierr); 1374b9ad928SBarry Smith ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 1384b9ad928SBarry Smith } 139a98ce0f4SHong Zhang ierr = PetscViewerRestoreSubcomm(viewer,red->pc->comm,&subviewer);CHKERRQ(ierr); 140a98ce0f4SHong Zhang } else if (isstring) { /* not test it yet! */ 1414b9ad928SBarry Smith ierr = PetscViewerStringSPrintf(viewer," Redundant solver preconditioner");CHKERRQ(ierr); 1424b9ad928SBarry Smith ierr = PetscViewerGetSingleton(viewer,&sviewer);CHKERRQ(ierr); 1434b9ad928SBarry Smith if (!rank) { 1444b9ad928SBarry Smith ierr = PCView(red->pc,sviewer);CHKERRQ(ierr); 1454b9ad928SBarry Smith } 1464b9ad928SBarry Smith ierr = PetscViewerRestoreSingleton(viewer,&sviewer);CHKERRQ(ierr); 1474b9ad928SBarry Smith } else { 14879a5c55eSBarry Smith SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PC redundant",((PetscObject)viewer)->type_name); 1494b9ad928SBarry Smith } 1504b9ad928SBarry Smith PetscFunctionReturn(0); 1514b9ad928SBarry Smith } 1524b9ad928SBarry Smith 153b9147fbbSdalcinl #include "include/private/matimpl.h" /*I "petscmat.h" I*/ 1544b9ad928SBarry Smith #undef __FUNCT__ 1554b9ad928SBarry Smith #define __FUNCT__ "PCSetUp_Redundant" 1566849ba73SBarry Smith static PetscErrorCode PCSetUp_Redundant(PC pc) 1574b9ad928SBarry Smith { 1584b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 159dfbe8321SBarry Smith PetscErrorCode ierr; 1603f457be1SHong Zhang PetscInt mstart,mend,mlocal,m; 16113f74950SBarry Smith PetscMPIInt size; 1624b9ad928SBarry Smith MatReuse reuse = MAT_INITIAL_MATRIX; 1634b9ad928SBarry Smith MatStructure str = DIFFERENT_NONZERO_PATTERN; 1641fbd8f88SHong Zhang MPI_Comm comm = pc->comm; 16523ce1328SBarry Smith Vec vec; 1663f457be1SHong Zhang PetscInt mlocal_sub; 1673f457be1SHong Zhang PetscMPIInt subsize,subrank; 1681fbd8f88SHong Zhang PetscInt rstart_sub,rend_sub,mloc_sub,nsubcomm; 1691fbd8f88SHong Zhang const char *prefix; 1703f457be1SHong Zhang 1714b9ad928SBarry Smith PetscFunctionBegin; 17223ce1328SBarry Smith ierr = MatGetVecs(pc->pmat,&vec,0);CHKERRQ(ierr); 17323ce1328SBarry Smith ierr = VecGetSize(vec,&m);CHKERRQ(ierr); 1741fbd8f88SHong Zhang 1754b9ad928SBarry Smith if (!pc->setupcalled) { 1761fbd8f88SHong Zhang ierr = PetscSubcommCreate(comm,red->nsubcomm,&red->psubcomm);CHKERRQ(ierr); 1771fbd8f88SHong Zhang ierr = PetscLogObjectMemory(pc,sizeof(PetscSubcomm));CHKERRQ(ierr); 1781fbd8f88SHong Zhang 1791fbd8f88SHong Zhang /* create a new PC that processors in each subcomm have copy of */ 1801fbd8f88SHong Zhang MPI_Comm subcomm = red->psubcomm->comm; 1811fbd8f88SHong Zhang ierr = PCCreate(subcomm,&red->pc);CHKERRQ(ierr); 1821fbd8f88SHong Zhang ierr = PCSetType(red->pc,PCLU);CHKERRQ(ierr); 1831fbd8f88SHong Zhang ierr = PCGetOptionsPrefix(pc,&prefix);CHKERRQ(ierr); 1841fbd8f88SHong Zhang ierr = PCSetOptionsPrefix(red->pc,prefix);CHKERRQ(ierr); 1851fbd8f88SHong Zhang ierr = PCAppendOptionsPrefix(red->pc,"redundant_");CHKERRQ(ierr); 1861fbd8f88SHong Zhang ierr = PCSetFromOptions(red->pc);CHKERRQ(ierr); 1871fbd8f88SHong Zhang 1883f457be1SHong Zhang /* create working vectors xsub/ysub and xdup/ydup */ 18923ce1328SBarry Smith ierr = VecGetLocalSize(vec,&mlocal);CHKERRQ(ierr); 1903f457be1SHong Zhang ierr = VecGetOwnershipRange(vec,&mstart,&mend);CHKERRQ(ierr); 1914b9ad928SBarry Smith 1923f457be1SHong Zhang /* get local size of xsub/ysub */ 1931fbd8f88SHong Zhang ierr = MPI_Comm_size(subcomm,&subsize);CHKERRQ(ierr); 1941fbd8f88SHong Zhang ierr = MPI_Comm_rank(subcomm,&subrank);CHKERRQ(ierr); 1951fbd8f88SHong Zhang rstart_sub = pc->pmat->rmap.range[red->psubcomm->n*subrank]; /* rstart in xsub/ysub */ 1963f457be1SHong Zhang if (subrank+1 < subsize){ 1971fbd8f88SHong Zhang rend_sub = pc->pmat->rmap.range[red->psubcomm->n*(subrank+1)]; 1983f457be1SHong Zhang } else { 1993f457be1SHong Zhang rend_sub = m; 2003f457be1SHong Zhang } 2013f457be1SHong Zhang mloc_sub = rend_sub - rstart_sub; 2021fbd8f88SHong Zhang ierr = VecCreateMPI(subcomm,mloc_sub,PETSC_DECIDE,&red->ysub);CHKERRQ(ierr); 2033f457be1SHong Zhang /* create xsub with empty local arrays, because xdup's arrays will be placed into it */ 2041fbd8f88SHong Zhang ierr = VecCreateMPIWithArray(subcomm,mloc_sub,PETSC_DECIDE,PETSC_NULL,&red->xsub);CHKERRQ(ierr); 2053f457be1SHong Zhang 2063f457be1SHong Zhang /* create xdup and ydup. ydup has empty local arrays because ysub's arrays will be place into it. 2073f457be1SHong Zhang Note: we use communicator dupcomm, not pc->comm! */ 2081fbd8f88SHong Zhang ierr = VecCreateMPI(red->psubcomm->dupparent,mloc_sub,PETSC_DECIDE,&red->xdup);CHKERRQ(ierr); 2091fbd8f88SHong Zhang ierr = VecCreateMPIWithArray(red->psubcomm->dupparent,mloc_sub,PETSC_DECIDE,PETSC_NULL,&red->ydup);CHKERRQ(ierr); 2103f457be1SHong Zhang 2113f457be1SHong Zhang /* create vec scatters */ 2123f457be1SHong Zhang if (!red->scatterin){ 2133f457be1SHong Zhang IS is1,is2; 2143f457be1SHong Zhang PetscInt *idx1,*idx2,i,j,k; 2151fbd8f88SHong Zhang ierr = PetscMalloc(2*red->psubcomm->n*mlocal*sizeof(PetscInt),&idx1);CHKERRQ(ierr); 2161fbd8f88SHong Zhang idx2 = idx1 + red->psubcomm->n*mlocal; 2173f457be1SHong Zhang j = 0; 2181fbd8f88SHong Zhang for (k=0; k<red->psubcomm->n; k++){ 2193f457be1SHong Zhang for (i=mstart; i<mend; i++){ 2203f457be1SHong Zhang idx1[j] = i; 2213f457be1SHong Zhang idx2[j++] = i + m*k; 2223f457be1SHong Zhang } 2233f457be1SHong Zhang } 2241fbd8f88SHong Zhang ierr = ISCreateGeneral(comm,red->psubcomm->n*mlocal,idx1,&is1);CHKERRQ(ierr); 2251fbd8f88SHong Zhang ierr = ISCreateGeneral(comm,red->psubcomm->n*mlocal,idx2,&is2);CHKERRQ(ierr); 2263f457be1SHong Zhang ierr = VecScatterCreate(vec,is1,red->xdup,is2,&red->scatterin);CHKERRQ(ierr); 2273f457be1SHong Zhang ierr = ISDestroy(is1);CHKERRQ(ierr); 2283f457be1SHong Zhang ierr = ISDestroy(is2);CHKERRQ(ierr); 2293f457be1SHong Zhang 2301fbd8f88SHong Zhang ierr = ISCreateStride(comm,mlocal,mstart+ red->psubcomm->color*m,1,&is1);CHKERRQ(ierr); 2313f457be1SHong Zhang ierr = ISCreateStride(comm,mlocal,mstart,1,&is2);CHKERRQ(ierr); 2323f457be1SHong Zhang ierr = VecScatterCreate(red->xdup,is1,vec,is2,&red->scatterout);CHKERRQ(ierr); 2333f457be1SHong Zhang ierr = ISDestroy(is1);CHKERRQ(ierr); 2343f457be1SHong Zhang ierr = ISDestroy(is2);CHKERRQ(ierr); 2353f457be1SHong Zhang ierr = PetscFree(idx1);CHKERRQ(ierr); 2364b9ad928SBarry Smith } 2374b9ad928SBarry Smith } 23823ce1328SBarry Smith ierr = VecDestroy(vec);CHKERRQ(ierr); 2394b9ad928SBarry Smith 2404b9ad928SBarry Smith /* if pmatrix set by user is sequential then we do not need to gather the parallel matrix */ 2413f457be1SHong Zhang ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 2424b9ad928SBarry Smith if (size == 1) { 2434b9ad928SBarry Smith red->useparallelmat = PETSC_FALSE; 2444b9ad928SBarry Smith } 2454b9ad928SBarry Smith 2464b9ad928SBarry Smith if (red->useparallelmat) { 2474b9ad928SBarry Smith if (pc->setupcalled == 1 && pc->flag == DIFFERENT_NONZERO_PATTERN) { 2484b9ad928SBarry Smith /* destroy old matrices */ 2494b9ad928SBarry Smith if (red->pmats) { 250b3804887SHong Zhang ierr = MatDestroy(red->pmats);CHKERRQ(ierr); 2514b9ad928SBarry Smith } 2524b9ad928SBarry Smith } else if (pc->setupcalled == 1) { 2534b9ad928SBarry Smith reuse = MAT_REUSE_MATRIX; 2544b9ad928SBarry Smith str = SAME_NONZERO_PATTERN; 2554b9ad928SBarry Smith } 2564b9ad928SBarry Smith 2573f457be1SHong Zhang /* grab the parallel matrix and put it into processors of a subcomminicator */ 258f664ae05SHong Zhang /*--------------------------------------------------------------------------*/ 259f664ae05SHong Zhang ierr = VecGetLocalSize(red->ysub,&mlocal_sub);CHKERRQ(ierr); 260*69db28dcSHong Zhang ierr = MatGetRedundantMatrix(pc->pmat,red->psubcomm->n,red->psubcomm->comm,mlocal_sub,reuse,&red->pmats);CHKERRQ(ierr); 261*69db28dcSHong Zhang 2623f457be1SHong Zhang /* tell PC of the subcommunicator its operators */ 263b3804887SHong Zhang ierr = PCSetOperators(red->pc,red->pmats,red->pmats,str);CHKERRQ(ierr); 2644b9ad928SBarry Smith } else { 2654b9ad928SBarry Smith ierr = PCSetOperators(red->pc,pc->mat,pc->pmat,pc->flag);CHKERRQ(ierr); 2664b9ad928SBarry Smith } 2674b9ad928SBarry Smith ierr = PCSetFromOptions(red->pc);CHKERRQ(ierr); 2684b9ad928SBarry Smith ierr = PCSetUp(red->pc);CHKERRQ(ierr); 2694b9ad928SBarry Smith PetscFunctionReturn(0); 2704b9ad928SBarry Smith } 2714b9ad928SBarry Smith 2724b9ad928SBarry Smith #undef __FUNCT__ 2734b9ad928SBarry Smith #define __FUNCT__ "PCApply_Redundant" 2746849ba73SBarry Smith static PetscErrorCode PCApply_Redundant(PC pc,Vec x,Vec y) 2754b9ad928SBarry Smith { 2764b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 277dfbe8321SBarry Smith PetscErrorCode ierr; 2783f457be1SHong Zhang PetscScalar *array; 2794b9ad928SBarry Smith 2804b9ad928SBarry Smith PetscFunctionBegin; 2813f457be1SHong Zhang /* scatter x to xdup */ 2823f457be1SHong Zhang ierr = VecScatterBegin(x,red->xdup,INSERT_VALUES,SCATTER_FORWARD,red->scatterin);CHKERRQ(ierr); 2833f457be1SHong Zhang ierr = VecScatterEnd(x,red->xdup,INSERT_VALUES,SCATTER_FORWARD,red->scatterin);CHKERRQ(ierr); 2843f457be1SHong Zhang 2853f457be1SHong Zhang /* place xdup's local array into xsub */ 2863f457be1SHong Zhang ierr = VecGetArray(red->xdup,&array);CHKERRQ(ierr); 2873f457be1SHong Zhang ierr = VecPlaceArray(red->xsub,(const PetscScalar*)array);CHKERRQ(ierr); 2884b9ad928SBarry Smith 2894b9ad928SBarry Smith /* apply preconditioner on each processor */ 2903f457be1SHong Zhang ierr = PCApply(red->pc,red->xsub,red->ysub);CHKERRQ(ierr); 2913f457be1SHong Zhang ierr = VecResetArray(red->xsub);CHKERRQ(ierr); 2923f457be1SHong Zhang ierr = VecRestoreArray(red->xdup,&array);CHKERRQ(ierr); 2934b9ad928SBarry Smith 2943f457be1SHong Zhang /* place ysub's local array into ydup */ 2953f457be1SHong Zhang ierr = VecGetArray(red->ysub,&array);CHKERRQ(ierr); 2963f457be1SHong Zhang ierr = VecPlaceArray(red->ydup,(const PetscScalar*)array);CHKERRQ(ierr); 2973f457be1SHong Zhang 2983f457be1SHong Zhang /* scatter ydup to y */ 2993f457be1SHong Zhang ierr = VecScatterBegin(red->ydup,y,INSERT_VALUES,SCATTER_FORWARD,red->scatterout);CHKERRQ(ierr); 3003f457be1SHong Zhang ierr = VecScatterEnd(red->ydup,y,INSERT_VALUES,SCATTER_FORWARD,red->scatterout);CHKERRQ(ierr); 3013f457be1SHong Zhang ierr = VecResetArray(red->ydup);CHKERRQ(ierr); 3023f457be1SHong Zhang ierr = VecRestoreArray(red->ysub,&array);CHKERRQ(ierr); 3034b9ad928SBarry Smith PetscFunctionReturn(0); 3044b9ad928SBarry Smith } 3054b9ad928SBarry Smith 3064b9ad928SBarry Smith #undef __FUNCT__ 3074b9ad928SBarry Smith #define __FUNCT__ "PCDestroy_Redundant" 3086849ba73SBarry Smith static PetscErrorCode PCDestroy_Redundant(PC pc) 3094b9ad928SBarry Smith { 3104b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 311dfbe8321SBarry Smith PetscErrorCode ierr; 3124b9ad928SBarry Smith 3134b9ad928SBarry Smith PetscFunctionBegin; 3144b9ad928SBarry Smith if (red->scatterin) {ierr = VecScatterDestroy(red->scatterin);CHKERRQ(ierr);} 3154b9ad928SBarry Smith if (red->scatterout) {ierr = VecScatterDestroy(red->scatterout);CHKERRQ(ierr);} 3163f457be1SHong Zhang if (red->ysub) {ierr = VecDestroy(red->ysub);CHKERRQ(ierr);} 3173f457be1SHong Zhang if (red->xsub) {ierr = VecDestroy(red->xsub);CHKERRQ(ierr);} 3183f457be1SHong Zhang if (red->xdup) {ierr = VecDestroy(red->xdup);CHKERRQ(ierr);} 3193f457be1SHong Zhang if (red->ydup) {ierr = VecDestroy(red->ydup);CHKERRQ(ierr);} 320b3804887SHong Zhang if (red->pmats) { 321b3804887SHong Zhang ierr = MatDestroy(red->pmats);CHKERRQ(ierr); 3223f457be1SHong Zhang } 3231fbd8f88SHong Zhang ierr = PetscSubcommDestroy(red->psubcomm);CHKERRQ(ierr); 3244b9ad928SBarry Smith ierr = PCDestroy(red->pc);CHKERRQ(ierr); 3254b9ad928SBarry Smith ierr = PetscFree(red);CHKERRQ(ierr); 3264b9ad928SBarry Smith PetscFunctionReturn(0); 3274b9ad928SBarry Smith } 3284b9ad928SBarry Smith 3294b9ad928SBarry Smith #undef __FUNCT__ 3304b9ad928SBarry Smith #define __FUNCT__ "PCSetFromOptions_Redundant" 3316849ba73SBarry Smith static PetscErrorCode PCSetFromOptions_Redundant(PC pc) 3324b9ad928SBarry Smith { 333a98ce0f4SHong Zhang PetscErrorCode ierr; 334a98ce0f4SHong Zhang PC_Redundant *red = (PC_Redundant*)pc->data; 3351fbd8f88SHong Zhang PetscMPIInt size; 336a98ce0f4SHong Zhang 3374b9ad928SBarry Smith PetscFunctionBegin; 338a98ce0f4SHong Zhang ierr = PetscOptionsHead("Redundant options");CHKERRQ(ierr); 339*69db28dcSHong Zhang ierr = PetscOptionsInt("-pc_redundant_number_comm","Number of subcommunicators","PCRedundantSetNumComm",red->nsubcomm,&red->nsubcomm,0);CHKERRQ(ierr); 340a98ce0f4SHong Zhang ierr = PetscOptionsTail();CHKERRQ(ierr); 3414b9ad928SBarry Smith PetscFunctionReturn(0); 3424b9ad928SBarry Smith } 3434b9ad928SBarry Smith 3444b9ad928SBarry Smith EXTERN_C_BEGIN 3454b9ad928SBarry Smith #undef __FUNCT__ 3464b9ad928SBarry Smith #define __FUNCT__ "PCRedundantSetScatter_Redundant" 347dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCRedundantSetScatter_Redundant(PC pc,VecScatter in,VecScatter out) 3484b9ad928SBarry Smith { 3494b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 350dfbe8321SBarry Smith PetscErrorCode ierr; 3514b9ad928SBarry Smith 3524b9ad928SBarry Smith PetscFunctionBegin; 3534b9ad928SBarry Smith red->scatterin = in; 3544b9ad928SBarry Smith red->scatterout = out; 3554b9ad928SBarry Smith ierr = PetscObjectReference((PetscObject)in);CHKERRQ(ierr); 3564b9ad928SBarry Smith ierr = PetscObjectReference((PetscObject)out);CHKERRQ(ierr); 3574b9ad928SBarry Smith PetscFunctionReturn(0); 3584b9ad928SBarry Smith } 3594b9ad928SBarry Smith EXTERN_C_END 3604b9ad928SBarry Smith 3614b9ad928SBarry Smith #undef __FUNCT__ 3624b9ad928SBarry Smith #define __FUNCT__ "PCRedundantSetScatter" 3634b9ad928SBarry Smith /*@ 3644b9ad928SBarry Smith PCRedundantSetScatter - Sets the scatter used to copy values into the 3654b9ad928SBarry Smith redundant local solve and the scatter to move them back into the global 3664b9ad928SBarry Smith vector. 3674b9ad928SBarry Smith 3684b9ad928SBarry Smith Collective on PC 3694b9ad928SBarry Smith 3704b9ad928SBarry Smith Input Parameters: 3714b9ad928SBarry Smith + pc - the preconditioner context 3724b9ad928SBarry Smith . in - the scatter to move the values in 3734b9ad928SBarry Smith - out - the scatter to move them out 3744b9ad928SBarry Smith 3754b9ad928SBarry Smith Level: advanced 3764b9ad928SBarry Smith 3774b9ad928SBarry Smith .keywords: PC, redundant solve 3784b9ad928SBarry Smith @*/ 379dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCRedundantSetScatter(PC pc,VecScatter in,VecScatter out) 3804b9ad928SBarry Smith { 381dfbe8321SBarry Smith PetscErrorCode ierr,(*f)(PC,VecScatter,VecScatter); 3824b9ad928SBarry Smith 3834b9ad928SBarry Smith PetscFunctionBegin; 3844482741eSBarry Smith PetscValidHeaderSpecific(pc,PC_COOKIE,1); 3854482741eSBarry Smith PetscValidHeaderSpecific(in,VEC_SCATTER_COOKIE,2); 3864482741eSBarry Smith PetscValidHeaderSpecific(out,VEC_SCATTER_COOKIE,3); 3874b9ad928SBarry Smith ierr = PetscObjectQueryFunction((PetscObject)pc,"PCRedundantSetScatter_C",(void (**)(void))&f);CHKERRQ(ierr); 3884b9ad928SBarry Smith if (f) { 3894b9ad928SBarry Smith ierr = (*f)(pc,in,out);CHKERRQ(ierr); 3904b9ad928SBarry Smith } 3914b9ad928SBarry Smith PetscFunctionReturn(0); 3924b9ad928SBarry Smith } 3934b9ad928SBarry Smith 3944b9ad928SBarry Smith EXTERN_C_BEGIN 3954b9ad928SBarry Smith #undef __FUNCT__ 3964b9ad928SBarry Smith #define __FUNCT__ "PCRedundantGetPC_Redundant" 397dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCRedundantGetPC_Redundant(PC pc,PC *innerpc) 3984b9ad928SBarry Smith { 3994b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 4004b9ad928SBarry Smith 4014b9ad928SBarry Smith PetscFunctionBegin; 4024b9ad928SBarry Smith *innerpc = red->pc; 4034b9ad928SBarry Smith PetscFunctionReturn(0); 4044b9ad928SBarry Smith } 4054b9ad928SBarry Smith EXTERN_C_END 4064b9ad928SBarry Smith 4074b9ad928SBarry Smith #undef __FUNCT__ 4084b9ad928SBarry Smith #define __FUNCT__ "PCRedundantGetPC" 4094b9ad928SBarry Smith /*@ 4104b9ad928SBarry Smith PCRedundantGetPC - Gets the sequential PC created by the redundant PC. 4114b9ad928SBarry Smith 4124b9ad928SBarry Smith Not Collective 4134b9ad928SBarry Smith 4144b9ad928SBarry Smith Input Parameter: 4154b9ad928SBarry Smith . pc - the preconditioner context 4164b9ad928SBarry Smith 4174b9ad928SBarry Smith Output Parameter: 4184b9ad928SBarry Smith . innerpc - the sequential PC 4194b9ad928SBarry Smith 4204b9ad928SBarry Smith Level: advanced 4214b9ad928SBarry Smith 4224b9ad928SBarry Smith .keywords: PC, redundant solve 4234b9ad928SBarry Smith @*/ 424dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCRedundantGetPC(PC pc,PC *innerpc) 4254b9ad928SBarry Smith { 426dfbe8321SBarry Smith PetscErrorCode ierr,(*f)(PC,PC*); 4274b9ad928SBarry Smith 4284b9ad928SBarry Smith PetscFunctionBegin; 4294482741eSBarry Smith PetscValidHeaderSpecific(pc,PC_COOKIE,1); 4304482741eSBarry Smith PetscValidPointer(innerpc,2); 4314b9ad928SBarry Smith ierr = PetscObjectQueryFunction((PetscObject)pc,"PCRedundantGetPC_C",(void (**)(void))&f);CHKERRQ(ierr); 4324b9ad928SBarry Smith if (f) { 4334b9ad928SBarry Smith ierr = (*f)(pc,innerpc);CHKERRQ(ierr); 4344b9ad928SBarry Smith } 4354b9ad928SBarry Smith PetscFunctionReturn(0); 4364b9ad928SBarry Smith } 4374b9ad928SBarry Smith 4384b9ad928SBarry Smith EXTERN_C_BEGIN 4394b9ad928SBarry Smith #undef __FUNCT__ 4404b9ad928SBarry Smith #define __FUNCT__ "PCRedundantGetOperators_Redundant" 441dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCRedundantGetOperators_Redundant(PC pc,Mat *mat,Mat *pmat) 4424b9ad928SBarry Smith { 4434b9ad928SBarry Smith PC_Redundant *red = (PC_Redundant*)pc->data; 4444b9ad928SBarry Smith 4454b9ad928SBarry Smith PetscFunctionBegin; 446b3804887SHong Zhang if (mat) *mat = red->pmats; 447b3804887SHong Zhang if (pmat) *pmat = red->pmats; 4484b9ad928SBarry Smith PetscFunctionReturn(0); 4494b9ad928SBarry Smith } 4504b9ad928SBarry Smith EXTERN_C_END 4514b9ad928SBarry Smith 4524b9ad928SBarry Smith #undef __FUNCT__ 4534b9ad928SBarry Smith #define __FUNCT__ "PCRedundantGetOperators" 4544b9ad928SBarry Smith /*@ 4554b9ad928SBarry Smith PCRedundantGetOperators - gets the sequential matrix and preconditioner matrix 4564b9ad928SBarry Smith 4574b9ad928SBarry Smith Not Collective 4584b9ad928SBarry Smith 4594b9ad928SBarry Smith Input Parameter: 4604b9ad928SBarry Smith . pc - the preconditioner context 4614b9ad928SBarry Smith 4624b9ad928SBarry Smith Output Parameters: 4634b9ad928SBarry Smith + mat - the matrix 4644b9ad928SBarry Smith - pmat - the (possibly different) preconditioner matrix 4654b9ad928SBarry Smith 4664b9ad928SBarry Smith Level: advanced 4674b9ad928SBarry Smith 4684b9ad928SBarry Smith .keywords: PC, redundant solve 4694b9ad928SBarry Smith @*/ 470dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCRedundantGetOperators(PC pc,Mat *mat,Mat *pmat) 4714b9ad928SBarry Smith { 472dfbe8321SBarry Smith PetscErrorCode ierr,(*f)(PC,Mat*,Mat*); 4734b9ad928SBarry Smith 4744b9ad928SBarry Smith PetscFunctionBegin; 4754482741eSBarry Smith PetscValidHeaderSpecific(pc,PC_COOKIE,1); 4764482741eSBarry Smith if (mat) PetscValidPointer(mat,2); 4774482741eSBarry Smith if (pmat) PetscValidPointer(pmat,3); 4784b9ad928SBarry Smith ierr = PetscObjectQueryFunction((PetscObject)pc,"PCRedundantGetOperators_C",(void (**)(void))&f);CHKERRQ(ierr); 4794b9ad928SBarry Smith if (f) { 4804b9ad928SBarry Smith ierr = (*f)(pc,mat,pmat);CHKERRQ(ierr); 4814b9ad928SBarry Smith } 4824b9ad928SBarry Smith PetscFunctionReturn(0); 4834b9ad928SBarry Smith } 4844b9ad928SBarry Smith 4854b9ad928SBarry Smith /* -------------------------------------------------------------------------------------*/ 48637a17b4dSBarry Smith /*MC 487a98ce0f4SHong Zhang PCREDUNDANT - Runs a preconditioner for the entire problem on subgroups of processors 48837a17b4dSBarry Smith 48937a17b4dSBarry Smith Options for the redundant preconditioners can be set with -redundant_pc_xxx 49037a17b4dSBarry Smith 49109391456SBarry Smith Options Database: 49209391456SBarry Smith . -pc_redundant_number_comm - number of sub communicators to use 49309391456SBarry Smith 49437a17b4dSBarry Smith Level: intermediate 49537a17b4dSBarry Smith 49637a17b4dSBarry Smith .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PCRedundantSetScatter(), 49737a17b4dSBarry Smith PCRedundantGetPC(), PCRedundantGetOperators() 49837a17b4dSBarry Smith M*/ 49937a17b4dSBarry Smith 5004b9ad928SBarry Smith EXTERN_C_BEGIN 5014b9ad928SBarry Smith #undef __FUNCT__ 5024b9ad928SBarry Smith #define __FUNCT__ "PCCreate_Redundant" 503dba47a55SKris Buschelman PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_Redundant(PC pc) 5044b9ad928SBarry Smith { 505dfbe8321SBarry Smith PetscErrorCode ierr; 5064b9ad928SBarry Smith PC_Redundant *red; 507*69db28dcSHong Zhang PetscMPIInt size; 5083f457be1SHong Zhang 5094b9ad928SBarry Smith PetscFunctionBegin; 5104b9ad928SBarry Smith ierr = PetscNew(PC_Redundant,&red);CHKERRQ(ierr); 51152e6d16bSBarry Smith ierr = PetscLogObjectMemory(pc,sizeof(PC_Redundant));CHKERRQ(ierr); 512*69db28dcSHong Zhang ierr = MPI_Comm_size(pc->comm,&size);CHKERRQ(ierr); 513*69db28dcSHong Zhang red->nsubcomm = size; 5144b9ad928SBarry Smith red->useparallelmat = PETSC_TRUE; 5151fbd8f88SHong Zhang pc->data = (void*)red; 5164b9ad928SBarry Smith 5174b9ad928SBarry Smith pc->ops->apply = PCApply_Redundant; 5184b9ad928SBarry Smith pc->ops->applytranspose = 0; 5194b9ad928SBarry Smith pc->ops->setup = PCSetUp_Redundant; 5204b9ad928SBarry Smith pc->ops->destroy = PCDestroy_Redundant; 5214b9ad928SBarry Smith pc->ops->setfromoptions = PCSetFromOptions_Redundant; 5224b9ad928SBarry Smith pc->ops->view = PCView_Redundant; 5234b9ad928SBarry Smith ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCRedundantSetScatter_C","PCRedundantSetScatter_Redundant", 5244b9ad928SBarry Smith PCRedundantSetScatter_Redundant);CHKERRQ(ierr); 5254b9ad928SBarry Smith ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCRedundantGetPC_C","PCRedundantGetPC_Redundant", 5264b9ad928SBarry Smith PCRedundantGetPC_Redundant);CHKERRQ(ierr); 5274b9ad928SBarry Smith ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCRedundantGetOperators_C","PCRedundantGetOperators_Redundant", 5284b9ad928SBarry Smith PCRedundantGetOperators_Redundant);CHKERRQ(ierr); 5294b9ad928SBarry Smith PetscFunctionReturn(0); 5304b9ad928SBarry Smith } 5314b9ad928SBarry Smith EXTERN_C_END 532