14b9ad928SBarry Smith 24b9ad928SBarry Smith /* -------------------------------------------------------------------- 34b9ad928SBarry Smith 45a46d3faSBarry Smith This file implements a Jacobi preconditioner in PETSc as part of PC. 55a46d3faSBarry Smith You can use this as a starting point for implementing your own 65a46d3faSBarry Smith preconditioner that is not provided with PETSc. (You might also consider 75a46d3faSBarry Smith just using PCSHELL) 84b9ad928SBarry Smith 94b9ad928SBarry Smith The following basic routines are required for each preconditioner. 104b9ad928SBarry Smith PCCreate_XXX() - Creates a preconditioner context 114b9ad928SBarry Smith PCSetFromOptions_XXX() - Sets runtime options 124b9ad928SBarry Smith PCApply_XXX() - Applies the preconditioner 134b9ad928SBarry Smith PCDestroy_XXX() - Destroys the preconditioner context 144b9ad928SBarry Smith where the suffix "_XXX" denotes a particular implementation, in 154b9ad928SBarry Smith this case we use _Jacobi (e.g., PCCreate_Jacobi, PCApply_Jacobi). 164b9ad928SBarry Smith These routines are actually called via the common user interface 174b9ad928SBarry Smith routines PCCreate(), PCSetFromOptions(), PCApply(), and PCDestroy(), 184b9ad928SBarry Smith so the application code interface remains identical for all 194b9ad928SBarry Smith preconditioners. 204b9ad928SBarry Smith 214b9ad928SBarry Smith Another key routine is: 224b9ad928SBarry Smith PCSetUp_XXX() - Prepares for the use of a preconditioner 234b9ad928SBarry Smith by setting data structures and options. The interface routine PCSetUp() 244b9ad928SBarry Smith is not usually called directly by the user, but instead is called by 254b9ad928SBarry Smith PCApply() if necessary. 264b9ad928SBarry Smith 274b9ad928SBarry Smith Additional basic routines are: 284b9ad928SBarry Smith PCView_XXX() - Prints details of runtime options that 294b9ad928SBarry Smith have actually been used. 304b9ad928SBarry Smith These are called by application codes via the interface routines 314b9ad928SBarry Smith PCView(). 324b9ad928SBarry Smith 334b9ad928SBarry Smith The various types of solvers (preconditioners, Krylov subspace methods, 344b9ad928SBarry Smith nonlinear solvers, timesteppers) are all organized similarly, so the 354b9ad928SBarry Smith above description applies to these categories also. One exception is 364b9ad928SBarry Smith that the analogues of PCApply() for these components are KSPSolve(), 374b9ad928SBarry Smith SNESSolve(), and TSSolve(). 384b9ad928SBarry Smith 394b9ad928SBarry Smith Additional optional functionality unique to preconditioners is left and 404b9ad928SBarry Smith right symmetric preconditioner application via PCApplySymmetricLeft() 414b9ad928SBarry Smith and PCApplySymmetricRight(). The Jacobi implementation is 424b9ad928SBarry Smith PCApplySymmetricLeftOrRight_Jacobi(). 434b9ad928SBarry Smith 444b9ad928SBarry Smith -------------------------------------------------------------------- */ 454b9ad928SBarry Smith 464b9ad928SBarry Smith /* 474b9ad928SBarry Smith Include files needed for the Jacobi preconditioner: 484b9ad928SBarry Smith pcimpl.h - private include file intended for use by all preconditioners 494b9ad928SBarry Smith */ 504b9ad928SBarry Smith 51af0996ceSBarry Smith #include <petsc/private/pcimpl.h> /*I "petscpc.h" I*/ 524b9ad928SBarry Smith 530a545947SLisandro Dalcin const char *const PCJacobiTypes[] = {"DIAGONAL", "ROWMAX", "ROWSUM", "PCJacobiType", "PC_JACOBI_", NULL}; 54baa89ecbSBarry Smith 554b9ad928SBarry Smith /* 564b9ad928SBarry Smith Private context (data structure) for the Jacobi preconditioner. 574b9ad928SBarry Smith */ 584b9ad928SBarry Smith typedef struct { 592fa5cd67SKarl Rupp Vec diag; /* vector containing the reciprocals of the diagonal elements of the preconditioner matrix */ 604b9ad928SBarry Smith Vec diagsqrt; /* vector containing the reciprocals of the square roots of 614b9ad928SBarry Smith the diagonal elements of the preconditioner matrix (used 624b9ad928SBarry Smith only for symmetric preconditioner application) */ 63baa89ecbSBarry Smith PetscBool userowmax; /* set with PCJacobiSetType() */ 64ace3abfcSBarry Smith PetscBool userowsum; 65ace3abfcSBarry Smith PetscBool useabs; /* use the absolute values of the diagonal entries */ 6667ed0f3bSStefano Zampini PetscBool fixdiag; /* fix zero diagonal terms */ 674b9ad928SBarry Smith } PC_Jacobi; 684b9ad928SBarry Smith 697a660b1eSJed Brown static PetscErrorCode PCReset_Jacobi(PC); 707a660b1eSJed Brown 71d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCJacobiSetType_Jacobi(PC pc, PCJacobiType type) 72d71ae5a4SJacob Faibussowitsch { 73baa89ecbSBarry Smith PC_Jacobi *j = (PC_Jacobi *)pc->data; 747a660b1eSJed Brown PCJacobiType old_type; 754b9ad928SBarry Smith 764b9ad928SBarry Smith PetscFunctionBegin; 777a660b1eSJed Brown PetscCall(PCJacobiGetType(pc, &old_type)); 783ba16761SJacob Faibussowitsch if (old_type == type) PetscFunctionReturn(PETSC_SUCCESS); 797a660b1eSJed Brown PetscCall(PCReset_Jacobi(pc)); 80baa89ecbSBarry Smith j->userowmax = PETSC_FALSE; 81baa89ecbSBarry Smith j->userowsum = PETSC_FALSE; 82baa89ecbSBarry Smith if (type == PC_JACOBI_ROWMAX) { 834b9ad928SBarry Smith j->userowmax = PETSC_TRUE; 84baa89ecbSBarry Smith } else if (type == PC_JACOBI_ROWSUM) { 85baa89ecbSBarry Smith j->userowsum = PETSC_TRUE; 86baa89ecbSBarry Smith } 873ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 884b9ad928SBarry Smith } 894b9ad928SBarry Smith 90d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCJacobiGetType_Jacobi(PC pc, PCJacobiType *type) 91d71ae5a4SJacob Faibussowitsch { 92baa89ecbSBarry Smith PC_Jacobi *j = (PC_Jacobi *)pc->data; 9386697f06SMatthew Knepley 9486697f06SMatthew Knepley PetscFunctionBegin; 95baa89ecbSBarry Smith if (j->userowmax) { 96baa89ecbSBarry Smith *type = PC_JACOBI_ROWMAX; 97baa89ecbSBarry Smith } else if (j->userowsum) { 98baa89ecbSBarry Smith *type = PC_JACOBI_ROWSUM; 99baa89ecbSBarry Smith } else { 100baa89ecbSBarry Smith *type = PC_JACOBI_DIAGONAL; 101baa89ecbSBarry Smith } 1023ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 10386697f06SMatthew Knepley } 10486697f06SMatthew Knepley 105d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCJacobiSetUseAbs_Jacobi(PC pc, PetscBool flg) 106d71ae5a4SJacob Faibussowitsch { 107baa89ecbSBarry Smith PC_Jacobi *j = (PC_Jacobi *)pc->data; 108cd47f5d9SBarry Smith 109cd47f5d9SBarry Smith PetscFunctionBegin; 110baa89ecbSBarry Smith j->useabs = flg; 1113ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 112baa89ecbSBarry Smith } 113baa89ecbSBarry Smith 114d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCJacobiGetUseAbs_Jacobi(PC pc, PetscBool *flg) 115d71ae5a4SJacob Faibussowitsch { 116baa89ecbSBarry Smith PC_Jacobi *j = (PC_Jacobi *)pc->data; 117baa89ecbSBarry Smith 118baa89ecbSBarry Smith PetscFunctionBegin; 119baa89ecbSBarry Smith *flg = j->useabs; 1203ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 121cd47f5d9SBarry Smith } 122cd47f5d9SBarry Smith 123d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCJacobiSetFixDiagonal_Jacobi(PC pc, PetscBool flg) 124d71ae5a4SJacob Faibussowitsch { 12567ed0f3bSStefano Zampini PC_Jacobi *j = (PC_Jacobi *)pc->data; 12667ed0f3bSStefano Zampini 12767ed0f3bSStefano Zampini PetscFunctionBegin; 12867ed0f3bSStefano Zampini j->fixdiag = flg; 1293ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 13067ed0f3bSStefano Zampini } 13167ed0f3bSStefano Zampini 132d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCJacobiGetFixDiagonal_Jacobi(PC pc, PetscBool *flg) 133d71ae5a4SJacob Faibussowitsch { 13467ed0f3bSStefano Zampini PC_Jacobi *j = (PC_Jacobi *)pc->data; 13567ed0f3bSStefano Zampini 13667ed0f3bSStefano Zampini PetscFunctionBegin; 13767ed0f3bSStefano Zampini *flg = j->fixdiag; 1383ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 13967ed0f3bSStefano Zampini } 14067ed0f3bSStefano Zampini 1414b9ad928SBarry Smith /* 1424b9ad928SBarry Smith PCSetUp_Jacobi - Prepares for the use of the Jacobi preconditioner 1434b9ad928SBarry Smith by setting data structures and options. 1444b9ad928SBarry Smith 1454b9ad928SBarry Smith Input Parameter: 1464b9ad928SBarry Smith . pc - the preconditioner context 1474b9ad928SBarry Smith 1484b9ad928SBarry Smith Application Interface Routine: PCSetUp() 1494b9ad928SBarry Smith 150f1580f4eSBarry Smith Note: 1514b9ad928SBarry Smith The interface routine PCSetUp() is not usually called directly by 1524b9ad928SBarry Smith the user, but instead is called by PCApply() if necessary. 1534b9ad928SBarry Smith */ 154d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCSetUp_Jacobi(PC pc) 155d71ae5a4SJacob Faibussowitsch { 1564b9ad928SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 1574b9ad928SBarry Smith Vec diag, diagsqrt; 158cd47f5d9SBarry Smith PetscInt n, i; 1594b9ad928SBarry Smith PetscScalar *x; 160ace3abfcSBarry Smith PetscBool zeroflag = PETSC_FALSE; 1614b9ad928SBarry Smith 1624b9ad928SBarry Smith PetscFunctionBegin; 1634b9ad928SBarry Smith /* 1644b9ad928SBarry Smith For most preconditioners the code would begin here something like 1654b9ad928SBarry Smith 1664b9ad928SBarry Smith if (pc->setupcalled == 0) { allocate space the first time this is ever called 1679566063dSJacob Faibussowitsch PetscCall(MatCreateVecs(pc->mat,&jac->diag)); 1684b9ad928SBarry Smith } 1694b9ad928SBarry Smith 1704b9ad928SBarry Smith But for this preconditioner we want to support use of both the matrix' diagonal 1714b9ad928SBarry Smith elements (for left or right preconditioning) and square root of diagonal elements 1724b9ad928SBarry Smith (for symmetric preconditioning). Hence we do not allocate space here, since we 1734b9ad928SBarry Smith don't know at this point which will be needed (diag and/or diagsqrt) until the user 1744b9ad928SBarry Smith applies the preconditioner, and we don't want to allocate BOTH unless we need 1754b9ad928SBarry Smith them both. Thus, the diag and diagsqrt are allocated in PCSetUp_Jacobi_NonSymmetric() 1764b9ad928SBarry Smith and PCSetUp_Jacobi_Symmetric(), respectively. 1774b9ad928SBarry Smith */ 1784b9ad928SBarry Smith 1794b9ad928SBarry Smith /* 1804b9ad928SBarry Smith Here we set up the preconditioner; that is, we copy the diagonal values from 1814b9ad928SBarry Smith the matrix and put them into a format to make them quick to apply as a preconditioner. 1824b9ad928SBarry Smith */ 1834b9ad928SBarry Smith diag = jac->diag; 1844b9ad928SBarry Smith diagsqrt = jac->diagsqrt; 1854b9ad928SBarry Smith 1864b9ad928SBarry Smith if (diag) { 187b94d7dedSBarry Smith PetscBool isset, isspd; 18867ed0f3bSStefano Zampini 1894b9ad928SBarry Smith if (jac->userowmax) { 1909566063dSJacob Faibussowitsch PetscCall(MatGetRowMaxAbs(pc->pmat, diag, NULL)); 19186697f06SMatthew Knepley } else if (jac->userowsum) { 1929566063dSJacob Faibussowitsch PetscCall(MatGetRowSum(pc->pmat, diag)); 1934b9ad928SBarry Smith } else { 1949566063dSJacob Faibussowitsch PetscCall(MatGetDiagonal(pc->pmat, diag)); 1954b9ad928SBarry Smith } 1969566063dSJacob Faibussowitsch PetscCall(VecReciprocal(diag)); 1971baa6e33SBarry Smith if (jac->useabs) PetscCall(VecAbs(diag)); 198b94d7dedSBarry Smith PetscCall(MatIsSPDKnown(pc->pmat, &isset, &isspd)); 199b94d7dedSBarry Smith if (jac->fixdiag && (!isset || !isspd)) { 2009566063dSJacob Faibussowitsch PetscCall(VecGetLocalSize(diag, &n)); 2019566063dSJacob Faibussowitsch PetscCall(VecGetArray(diag, &x)); 2024b9ad928SBarry Smith for (i = 0; i < n; i++) { 2034b9ad928SBarry Smith if (x[i] == 0.0) { 2044b9ad928SBarry Smith x[i] = 1.0; 205cd47f5d9SBarry Smith zeroflag = PETSC_TRUE; 2064b9ad928SBarry Smith } 2074b9ad928SBarry Smith } 2089566063dSJacob Faibussowitsch PetscCall(VecRestoreArray(diag, &x)); 2094b9ad928SBarry Smith } 21067ed0f3bSStefano Zampini } 2114b9ad928SBarry Smith if (diagsqrt) { 2124b9ad928SBarry Smith if (jac->userowmax) { 2139566063dSJacob Faibussowitsch PetscCall(MatGetRowMaxAbs(pc->pmat, diagsqrt, NULL)); 21486697f06SMatthew Knepley } else if (jac->userowsum) { 2159566063dSJacob Faibussowitsch PetscCall(MatGetRowSum(pc->pmat, diagsqrt)); 2164b9ad928SBarry Smith } else { 2179566063dSJacob Faibussowitsch PetscCall(MatGetDiagonal(pc->pmat, diagsqrt)); 2184b9ad928SBarry Smith } 2199566063dSJacob Faibussowitsch PetscCall(VecGetLocalSize(diagsqrt, &n)); 2209566063dSJacob Faibussowitsch PetscCall(VecGetArray(diagsqrt, &x)); 2214b9ad928SBarry Smith for (i = 0; i < n; i++) { 2228f1a2a5eSBarry Smith if (x[i] != 0.0) x[i] = 1.0 / PetscSqrtReal(PetscAbsScalar(x[i])); 2234b9ad928SBarry Smith else { 2244b9ad928SBarry Smith x[i] = 1.0; 225cd47f5d9SBarry Smith zeroflag = PETSC_TRUE; 2264b9ad928SBarry Smith } 2274b9ad928SBarry Smith } 2289566063dSJacob Faibussowitsch PetscCall(VecRestoreArray(diagsqrt, &x)); 2294b9ad928SBarry Smith } 23048a46eb9SPierre Jolivet if (zeroflag) PetscCall(PetscInfo(pc, "Zero detected in diagonal of matrix, using 1 at those locations\n")); 2313ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 2324b9ad928SBarry Smith } 233f1580f4eSBarry Smith 2344b9ad928SBarry Smith /* 2354b9ad928SBarry Smith PCSetUp_Jacobi_Symmetric - Allocates the vector needed to store the 2364b9ad928SBarry Smith inverse of the square root of the diagonal entries of the matrix. This 2374b9ad928SBarry Smith is used for symmetric application of the Jacobi preconditioner. 2384b9ad928SBarry Smith 2394b9ad928SBarry Smith Input Parameter: 2404b9ad928SBarry Smith . pc - the preconditioner context 2414b9ad928SBarry Smith */ 242d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCSetUp_Jacobi_Symmetric(PC pc) 243d71ae5a4SJacob Faibussowitsch { 2444b9ad928SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 2454b9ad928SBarry Smith 2464b9ad928SBarry Smith PetscFunctionBegin; 2479566063dSJacob Faibussowitsch PetscCall(MatCreateVecs(pc->pmat, &jac->diagsqrt, NULL)); 2489566063dSJacob Faibussowitsch PetscCall(PCSetUp_Jacobi(pc)); 2493ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 2504b9ad928SBarry Smith } 251f1580f4eSBarry Smith 2524b9ad928SBarry Smith /* 2534b9ad928SBarry Smith PCSetUp_Jacobi_NonSymmetric - Allocates the vector needed to store the 2544b9ad928SBarry Smith inverse of the diagonal entries of the matrix. This is used for left of 2554b9ad928SBarry Smith right application of the Jacobi preconditioner. 2564b9ad928SBarry Smith 2574b9ad928SBarry Smith Input Parameter: 2584b9ad928SBarry Smith . pc - the preconditioner context 2594b9ad928SBarry Smith */ 260d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCSetUp_Jacobi_NonSymmetric(PC pc) 261d71ae5a4SJacob Faibussowitsch { 2624b9ad928SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 2634b9ad928SBarry Smith 2644b9ad928SBarry Smith PetscFunctionBegin; 2659566063dSJacob Faibussowitsch PetscCall(MatCreateVecs(pc->pmat, &jac->diag, NULL)); 2669566063dSJacob Faibussowitsch PetscCall(PCSetUp_Jacobi(pc)); 2673ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 2684b9ad928SBarry Smith } 269f1580f4eSBarry Smith 2704b9ad928SBarry Smith /* 2714b9ad928SBarry Smith PCApply_Jacobi - Applies the Jacobi preconditioner to a vector. 2724b9ad928SBarry Smith 2734b9ad928SBarry Smith Input Parameters: 2744b9ad928SBarry Smith . pc - the preconditioner context 2754b9ad928SBarry Smith . x - input vector 2764b9ad928SBarry Smith 2774b9ad928SBarry Smith Output Parameter: 2784b9ad928SBarry Smith . y - output vector 2794b9ad928SBarry Smith 2804b9ad928SBarry Smith Application Interface Routine: PCApply() 2814b9ad928SBarry Smith */ 282d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCApply_Jacobi(PC pc, Vec x, Vec y) 283d71ae5a4SJacob Faibussowitsch { 2844b9ad928SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 2854b9ad928SBarry Smith 2864b9ad928SBarry Smith PetscFunctionBegin; 28748a46eb9SPierre Jolivet if (!jac->diag) PetscCall(PCSetUp_Jacobi_NonSymmetric(pc)); 2889566063dSJacob Faibussowitsch PetscCall(VecPointwiseMult(y, x, jac->diag)); 2893ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 2904b9ad928SBarry Smith } 291f1580f4eSBarry Smith 2924b9ad928SBarry Smith /* 2934b9ad928SBarry Smith PCApplySymmetricLeftOrRight_Jacobi - Applies the left or right part of a 2944b9ad928SBarry Smith symmetric preconditioner to a vector. 2954b9ad928SBarry Smith 2964b9ad928SBarry Smith Input Parameters: 2974b9ad928SBarry Smith . pc - the preconditioner context 2984b9ad928SBarry Smith . x - input vector 2994b9ad928SBarry Smith 3004b9ad928SBarry Smith Output Parameter: 3014b9ad928SBarry Smith . y - output vector 3024b9ad928SBarry Smith 3034b9ad928SBarry Smith Application Interface Routines: PCApplySymmetricLeft(), PCApplySymmetricRight() 3044b9ad928SBarry Smith */ 305d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCApplySymmetricLeftOrRight_Jacobi(PC pc, Vec x, Vec y) 306d71ae5a4SJacob Faibussowitsch { 3074b9ad928SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 3084b9ad928SBarry Smith 3094b9ad928SBarry Smith PetscFunctionBegin; 31048a46eb9SPierre Jolivet if (!jac->diagsqrt) PetscCall(PCSetUp_Jacobi_Symmetric(pc)); 3119566063dSJacob Faibussowitsch PetscCall(VecPointwiseMult(y, x, jac->diagsqrt)); 3123ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 3134b9ad928SBarry Smith } 31467ed0f3bSStefano Zampini 315d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCReset_Jacobi(PC pc) 316d71ae5a4SJacob Faibussowitsch { 317a06653b4SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 318a06653b4SBarry Smith 319a06653b4SBarry Smith PetscFunctionBegin; 3209566063dSJacob Faibussowitsch PetscCall(VecDestroy(&jac->diag)); 3219566063dSJacob Faibussowitsch PetscCall(VecDestroy(&jac->diagsqrt)); 3223ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 323a06653b4SBarry Smith } 324a06653b4SBarry Smith 3254b9ad928SBarry Smith /* 3264b9ad928SBarry Smith PCDestroy_Jacobi - Destroys the private context for the Jacobi preconditioner 3274b9ad928SBarry Smith that was created with PCCreate_Jacobi(). 3284b9ad928SBarry Smith 3294b9ad928SBarry Smith Input Parameter: 3304b9ad928SBarry Smith . pc - the preconditioner context 3314b9ad928SBarry Smith 3324b9ad928SBarry Smith Application Interface Routine: PCDestroy() 3334b9ad928SBarry Smith */ 334d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCDestroy_Jacobi(PC pc) 335d71ae5a4SJacob Faibussowitsch { 3364b9ad928SBarry Smith PetscFunctionBegin; 3379566063dSJacob Faibussowitsch PetscCall(PCReset_Jacobi(pc)); 3389566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiSetType_C", NULL)); 3399566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiGetType_C", NULL)); 3409566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiSetUseAbs_C", NULL)); 3419566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiGetUseAbs_C", NULL)); 3429566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiSetFixDiagonal_C", NULL)); 3439566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiGetFixDiagonal_C", NULL)); 3444b9ad928SBarry Smith 3454b9ad928SBarry Smith /* 3464b9ad928SBarry Smith Free the private data structure that was hanging off the PC 3474b9ad928SBarry Smith */ 3489566063dSJacob Faibussowitsch PetscCall(PetscFree(pc->data)); 3493ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 3504b9ad928SBarry Smith } 3514b9ad928SBarry Smith 352d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCSetFromOptions_Jacobi(PC pc, PetscOptionItems *PetscOptionsObject) 353d71ae5a4SJacob Faibussowitsch { 3544b9ad928SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 355baa89ecbSBarry Smith PetscBool flg; 356baa89ecbSBarry Smith PCJacobiType deflt, type; 3574b9ad928SBarry Smith 3584b9ad928SBarry Smith PetscFunctionBegin; 3599566063dSJacob Faibussowitsch PetscCall(PCJacobiGetType(pc, &deflt)); 360d0609cedSBarry Smith PetscOptionsHeadBegin(PetscOptionsObject, "Jacobi options"); 3619566063dSJacob Faibussowitsch PetscCall(PetscOptionsEnum("-pc_jacobi_type", "How to construct diagonal matrix", "PCJacobiSetType", PCJacobiTypes, (PetscEnum)deflt, (PetscEnum *)&type, &flg)); 3621baa6e33SBarry Smith if (flg) PetscCall(PCJacobiSetType(pc, type)); 3639566063dSJacob Faibussowitsch PetscCall(PetscOptionsBool("-pc_jacobi_abs", "Use absolute values of diagonal entries", "PCJacobiSetUseAbs", jac->useabs, &jac->useabs, NULL)); 3649566063dSJacob Faibussowitsch PetscCall(PetscOptionsBool("-pc_jacobi_fixdiagonal", "Fix null terms on diagonal", "PCJacobiSetFixDiagonal", jac->fixdiag, &jac->fixdiag, NULL)); 365d0609cedSBarry Smith PetscOptionsHeadEnd(); 3663ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 3674b9ad928SBarry Smith } 3684b9ad928SBarry Smith 369d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCView_Jacobi(PC pc, PetscViewer viewer) 370d71ae5a4SJacob Faibussowitsch { 371569e28a7SMatthew G. Knepley PC_Jacobi *jac = (PC_Jacobi *)pc->data; 372569e28a7SMatthew G. Knepley PetscBool iascii; 373569e28a7SMatthew G. Knepley 374569e28a7SMatthew G. Knepley PetscFunctionBegin; 3759566063dSJacob Faibussowitsch PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii)); 376569e28a7SMatthew G. Knepley if (iascii) { 377569e28a7SMatthew G. Knepley PCJacobiType type; 37867ed0f3bSStefano Zampini PetscBool useAbs, fixdiag; 379569e28a7SMatthew G. Knepley PetscViewerFormat format; 380569e28a7SMatthew G. Knepley 3819566063dSJacob Faibussowitsch PetscCall(PCJacobiGetType(pc, &type)); 3829566063dSJacob Faibussowitsch PetscCall(PCJacobiGetUseAbs(pc, &useAbs)); 3839566063dSJacob Faibussowitsch PetscCall(PCJacobiGetFixDiagonal(pc, &fixdiag)); 3849566063dSJacob Faibussowitsch PetscCall(PetscViewerASCIIPrintf(viewer, " type %s%s%s\n", PCJacobiTypes[type], useAbs ? ", using absolute value of entries" : "", !fixdiag ? ", not checking null diagonal entries" : "")); 3859566063dSJacob Faibussowitsch PetscCall(PetscViewerGetFormat(viewer, &format)); 386b24842d1SPierre Jolivet if (format == PETSC_VIEWER_ASCII_INFO_DETAIL && jac->diag) PetscCall(VecView(jac->diag, viewer)); 387569e28a7SMatthew G. Knepley } 3883ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 389569e28a7SMatthew G. Knepley } 390569e28a7SMatthew G. Knepley 3914b9ad928SBarry Smith /* 3924b9ad928SBarry Smith PCCreate_Jacobi - Creates a Jacobi preconditioner context, PC_Jacobi, 3934b9ad928SBarry Smith and sets this as the private data within the generic preconditioning 3944b9ad928SBarry Smith context, PC, that was created within PCCreate(). 3954b9ad928SBarry Smith 3964b9ad928SBarry Smith Input Parameter: 3974b9ad928SBarry Smith . pc - the preconditioner context 3984b9ad928SBarry Smith 3994b9ad928SBarry Smith Application Interface Routine: PCCreate() 4004b9ad928SBarry Smith */ 4014b9ad928SBarry Smith 4024b9ad928SBarry Smith /*MC 4035a46d3faSBarry Smith PCJACOBI - Jacobi (i.e. diagonal scaling preconditioning) 4044b9ad928SBarry Smith 405f1580f4eSBarry Smith Options Database Keys: 406967c93d3SBarry Smith + -pc_jacobi_type <diagonal,rowmax,rowsum> - approach for forming the preconditioner 40767ed0f3bSStefano Zampini . -pc_jacobi_abs - use the absolute value of the diagonal entry 408147403d9SBarry Smith - -pc_jacobi_fixdiag - fix for zero diagonal terms by placing 1.0 in those locations 4094b9ad928SBarry Smith 4104b9ad928SBarry Smith Level: beginner 4114b9ad928SBarry Smith 41295452b02SPatrick Sanan Notes: 413f1580f4eSBarry Smith By using `KSPSetPCSide`(ksp,`PC_SYMMETRIC`) or -ksp_pc_side symmetric 4144b9ad928SBarry Smith can scale each side of the matrix by the square root of the diagonal entries. 4154b9ad928SBarry Smith 4164b9ad928SBarry Smith Zero entries along the diagonal are replaced with the value 1.0 4174b9ad928SBarry Smith 418f1580f4eSBarry Smith See `PCPBJACOBI` for fixed-size point block, `PCVPBJACOBI` for variable-sized point block, and `PCBJACOBI` for large size blocks 419422a814eSBarry Smith 420db781477SPatrick Sanan .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PC`, 421f1580f4eSBarry Smith `PCJacobiSetType()`, `PCJacobiSetUseAbs()`, `PCJacobiGetUseAbs()`, `PCASM`, 422db781477SPatrick Sanan `PCJacobiSetFixDiagonal()`, `PCJacobiGetFixDiagonal()` 423db781477SPatrick Sanan `PCJacobiSetType()`, `PCJacobiSetUseAbs()`, `PCJacobiGetUseAbs()`, `PCPBJACOBI`, `PCBJACOBI`, `PCVPBJACOBI` 4244b9ad928SBarry Smith M*/ 4254b9ad928SBarry Smith 426d71ae5a4SJacob Faibussowitsch PETSC_EXTERN PetscErrorCode PCCreate_Jacobi(PC pc) 427d71ae5a4SJacob Faibussowitsch { 4284b9ad928SBarry Smith PC_Jacobi *jac; 4294b9ad928SBarry Smith 4304b9ad928SBarry Smith PetscFunctionBegin; 4314b9ad928SBarry Smith /* 4324b9ad928SBarry Smith Creates the private data structure for this preconditioner and 4334b9ad928SBarry Smith attach it to the PC object. 4344b9ad928SBarry Smith */ 4354dfa11a4SJacob Faibussowitsch PetscCall(PetscNew(&jac)); 4364b9ad928SBarry Smith pc->data = (void *)jac; 4374b9ad928SBarry Smith 4384b9ad928SBarry Smith /* 4394b9ad928SBarry Smith Initialize the pointers to vectors to ZERO; these will be used to store 4404b9ad928SBarry Smith diagonal entries of the matrix for fast preconditioner application. 4414b9ad928SBarry Smith */ 4420a545947SLisandro Dalcin jac->diag = NULL; 4430a545947SLisandro Dalcin jac->diagsqrt = NULL; 4444b9ad928SBarry Smith jac->userowmax = PETSC_FALSE; 44586697f06SMatthew Knepley jac->userowsum = PETSC_FALSE; 446cd47f5d9SBarry Smith jac->useabs = PETSC_FALSE; 44767ed0f3bSStefano Zampini jac->fixdiag = PETSC_TRUE; 4484b9ad928SBarry Smith 4494b9ad928SBarry Smith /* 4504b9ad928SBarry Smith Set the pointers for the functions that are provided above. 4514b9ad928SBarry Smith Now when the user-level routines (such as PCApply(), PCDestroy(), etc.) 4524b9ad928SBarry Smith are called, they will automatically call these functions. Note we 4534b9ad928SBarry Smith choose not to provide a couple of these functions since they are 4544b9ad928SBarry Smith not needed. 4554b9ad928SBarry Smith */ 4564b9ad928SBarry Smith pc->ops->apply = PCApply_Jacobi; 4574b9ad928SBarry Smith pc->ops->applytranspose = PCApply_Jacobi; 4584b9ad928SBarry Smith pc->ops->setup = PCSetUp_Jacobi; 459a06653b4SBarry Smith pc->ops->reset = PCReset_Jacobi; 4604b9ad928SBarry Smith pc->ops->destroy = PCDestroy_Jacobi; 4614b9ad928SBarry Smith pc->ops->setfromoptions = PCSetFromOptions_Jacobi; 462569e28a7SMatthew G. Knepley pc->ops->view = PCView_Jacobi; 4630a545947SLisandro Dalcin pc->ops->applyrichardson = NULL; 4644b9ad928SBarry Smith pc->ops->applysymmetricleft = PCApplySymmetricLeftOrRight_Jacobi; 4654b9ad928SBarry Smith pc->ops->applysymmetricright = PCApplySymmetricLeftOrRight_Jacobi; 4662fa5cd67SKarl Rupp 4679566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiSetType_C", PCJacobiSetType_Jacobi)); 4689566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiGetType_C", PCJacobiGetType_Jacobi)); 4699566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiSetUseAbs_C", PCJacobiSetUseAbs_Jacobi)); 4709566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiGetUseAbs_C", PCJacobiGetUseAbs_Jacobi)); 4719566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiSetFixDiagonal_C", PCJacobiSetFixDiagonal_Jacobi)); 4729566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiGetFixDiagonal_C", PCJacobiGetFixDiagonal_Jacobi)); 4733ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 4744b9ad928SBarry Smith } 475cd47f5d9SBarry Smith 476cd47f5d9SBarry Smith /*@ 477f1580f4eSBarry Smith PCJacobiSetUseAbs - Causes the Jacobi preconditioner `PCJACOBI` to use the 47867ed0f3bSStefano Zampini absolute values of the diagonal divisors in the preconditioner 479cd47f5d9SBarry Smith 480c3339decSBarry Smith Logically Collective 481cd47f5d9SBarry Smith 482cd47f5d9SBarry Smith Input Parameters: 483baa89ecbSBarry Smith + pc - the preconditioner context 484baa89ecbSBarry Smith - flg - whether to use absolute values or not 485baa89ecbSBarry Smith 486baa89ecbSBarry Smith Options Database Key: 487147403d9SBarry Smith . -pc_jacobi_abs <bool> - use absolute values 488baa89ecbSBarry Smith 489f1580f4eSBarry Smith Note: 49095452b02SPatrick Sanan This takes affect at the next construction of the preconditioner 491baa89ecbSBarry Smith 492baa89ecbSBarry Smith Level: intermediate 493baa89ecbSBarry Smith 494f1580f4eSBarry Smith .seealso: `PCJACOBI`, `PCJacobiaSetType()`, `PCJacobiGetUseAbs()` 495baa89ecbSBarry Smith @*/ 496d71ae5a4SJacob Faibussowitsch PetscErrorCode PCJacobiSetUseAbs(PC pc, PetscBool flg) 497d71ae5a4SJacob Faibussowitsch { 498baa89ecbSBarry Smith PetscFunctionBegin; 499baa89ecbSBarry Smith PetscValidHeaderSpecific(pc, PC_CLASSID, 1); 500cac4c232SBarry Smith PetscTryMethod(pc, "PCJacobiSetUseAbs_C", (PC, PetscBool), (pc, flg)); 5013ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 502baa89ecbSBarry Smith } 503baa89ecbSBarry Smith 504baa89ecbSBarry Smith /*@ 505f1580f4eSBarry Smith PCJacobiGetUseAbs - Determines if the Jacobi preconditioner `PCJACOBI` uses the 50667ed0f3bSStefano Zampini absolute values of the diagonal divisors in the preconditioner 507baa89ecbSBarry Smith 508c3339decSBarry Smith Logically Collective 509baa89ecbSBarry Smith 510baa89ecbSBarry Smith Input Parameter: 511cd47f5d9SBarry Smith . pc - the preconditioner context 512cd47f5d9SBarry Smith 513baa89ecbSBarry Smith Output Parameter: 514baa89ecbSBarry Smith . flg - whether to use absolute values or not 515cd47f5d9SBarry Smith 516cd47f5d9SBarry Smith Level: intermediate 517cd47f5d9SBarry Smith 518f1580f4eSBarry Smith .seealso: `PCJACOBI`, `PCJacobiaSetType()`, `PCJacobiSetUseAbs()`, `PCJacobiGetType()` 519cd47f5d9SBarry Smith @*/ 520d71ae5a4SJacob Faibussowitsch PetscErrorCode PCJacobiGetUseAbs(PC pc, PetscBool *flg) 521d71ae5a4SJacob Faibussowitsch { 522cd47f5d9SBarry Smith PetscFunctionBegin; 5230700a824SBarry Smith PetscValidHeaderSpecific(pc, PC_CLASSID, 1); 524cac4c232SBarry Smith PetscUseMethod(pc, "PCJacobiGetUseAbs_C", (PC, PetscBool *), (pc, flg)); 5253ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 526cd47f5d9SBarry Smith } 527cd47f5d9SBarry Smith 5284b9ad928SBarry Smith /*@ 529147403d9SBarry Smith PCJacobiSetFixDiagonal - Check for zero values on the diagonal and replace them with 1.0 53067ed0f3bSStefano Zampini 531c3339decSBarry Smith Logically Collective 53267ed0f3bSStefano Zampini 53367ed0f3bSStefano Zampini Input Parameters: 53467ed0f3bSStefano Zampini + pc - the preconditioner context 53567ed0f3bSStefano Zampini - flg - the boolean flag 53667ed0f3bSStefano Zampini 53767ed0f3bSStefano Zampini Options Database Key: 538147403d9SBarry Smith . -pc_jacobi_fixdiagonal <bool> - check for zero values on the diagonal 53967ed0f3bSStefano Zampini 540f1580f4eSBarry Smith Note: 54167ed0f3bSStefano Zampini This takes affect at the next construction of the preconditioner 54267ed0f3bSStefano Zampini 54367ed0f3bSStefano Zampini Level: intermediate 54467ed0f3bSStefano Zampini 545f1580f4eSBarry Smith .seealso: `PCJACOBI`, `PCJacobiSetType()`, `PCJacobiGetFixDiagonal()`, `PCJacobiSetUseAbs()` 54667ed0f3bSStefano Zampini @*/ 547d71ae5a4SJacob Faibussowitsch PetscErrorCode PCJacobiSetFixDiagonal(PC pc, PetscBool flg) 548d71ae5a4SJacob Faibussowitsch { 54967ed0f3bSStefano Zampini PetscFunctionBegin; 55067ed0f3bSStefano Zampini PetscValidHeaderSpecific(pc, PC_CLASSID, 1); 551cac4c232SBarry Smith PetscTryMethod(pc, "PCJacobiSetFixDiagonal_C", (PC, PetscBool), (pc, flg)); 5523ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 55367ed0f3bSStefano Zampini } 55467ed0f3bSStefano Zampini 55567ed0f3bSStefano Zampini /*@ 556f1580f4eSBarry Smith PCJacobiGetFixDiagonal - Determines if the Jacobi preconditioner `PCJACOBI` checks for zero diagonal terms 55767ed0f3bSStefano Zampini 558c3339decSBarry Smith Logically Collective 55967ed0f3bSStefano Zampini 56067ed0f3bSStefano Zampini Input Parameter: 56167ed0f3bSStefano Zampini . pc - the preconditioner context 56267ed0f3bSStefano Zampini 56367ed0f3bSStefano Zampini Output Parameter: 56467ed0f3bSStefano Zampini . flg - the boolean flag 56567ed0f3bSStefano Zampini 56667ed0f3bSStefano Zampini Options Database Key: 56767b8a455SSatish Balay . -pc_jacobi_fixdiagonal <bool> - Fix 0 terms on diagonal by using 1 56867ed0f3bSStefano Zampini 56967ed0f3bSStefano Zampini Level: intermediate 57067ed0f3bSStefano Zampini 571f1580f4eSBarry Smith .seealso: `PCJACOBI`, `PCJacobiSetType()`, `PCJacobiSetFixDiagonal()` 57267ed0f3bSStefano Zampini @*/ 573d71ae5a4SJacob Faibussowitsch PetscErrorCode PCJacobiGetFixDiagonal(PC pc, PetscBool *flg) 574d71ae5a4SJacob Faibussowitsch { 57567ed0f3bSStefano Zampini PetscFunctionBegin; 57667ed0f3bSStefano Zampini PetscValidHeaderSpecific(pc, PC_CLASSID, 1); 577cac4c232SBarry Smith PetscUseMethod(pc, "PCJacobiGetFixDiagonal_C", (PC, PetscBool *), (pc, flg)); 5783ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 57967ed0f3bSStefano Zampini } 58067ed0f3bSStefano Zampini 58167ed0f3bSStefano Zampini /*@ 582baa89ecbSBarry Smith PCJacobiSetType - Causes the Jacobi preconditioner to use either the diagonal, the maximum entry in each row, 583baa89ecbSBarry Smith of the sum of rows entries for the diagonal preconditioner 5844b9ad928SBarry Smith 585c3339decSBarry Smith Logically Collective 5864b9ad928SBarry Smith 5874b9ad928SBarry Smith Input Parameters: 588baa89ecbSBarry Smith + pc - the preconditioner context 589f1580f4eSBarry Smith - type - `PC_JACOBI_DIAGONAL`, `PC_JACOBI_ROWMAX`, `PC_JACOBI_ROWSUM` 5904b9ad928SBarry Smith 5914b9ad928SBarry Smith Options Database Key: 592147403d9SBarry Smith . -pc_jacobi_type <diagonal,rowmax,rowsum> - the type of diagonal matrix to use for Jacobi 5934b9ad928SBarry Smith 5944b9ad928SBarry Smith Level: intermediate 5954b9ad928SBarry Smith 596*feefa0e1SJacob Faibussowitsch Developer Notes: 597f1580f4eSBarry Smith Why is there a separate function for using the absolute value? 598f1580f4eSBarry Smith 599f1580f4eSBarry Smith .seealso: `PCJACOBI`, `PCJacobiSetUseAbs()`, `PCJacobiGetType()` 6004b9ad928SBarry Smith @*/ 601d71ae5a4SJacob Faibussowitsch PetscErrorCode PCJacobiSetType(PC pc, PCJacobiType type) 602d71ae5a4SJacob Faibussowitsch { 6034b9ad928SBarry Smith PetscFunctionBegin; 6040700a824SBarry Smith PetscValidHeaderSpecific(pc, PC_CLASSID, 1); 605cac4c232SBarry Smith PetscTryMethod(pc, "PCJacobiSetType_C", (PC, PCJacobiType), (pc, type)); 6063ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 6074b9ad928SBarry Smith } 6084b9ad928SBarry Smith 60986697f06SMatthew Knepley /*@ 610baa89ecbSBarry Smith PCJacobiGetType - Gets how the diagonal matrix is produced for the preconditioner 61186697f06SMatthew Knepley 612c3339decSBarry Smith Not Collective 61386697f06SMatthew Knepley 614baa89ecbSBarry Smith Input Parameter: 61586697f06SMatthew Knepley . pc - the preconditioner context 61686697f06SMatthew Knepley 617baa89ecbSBarry Smith Output Parameter: 618f1580f4eSBarry Smith . type - `PC_JACOBI_DIAGONAL`, `PC_JACOBI_ROWMAX`, `PC_JACOBI_ROWSUM` 61986697f06SMatthew Knepley 62086697f06SMatthew Knepley Level: intermediate 62186697f06SMatthew Knepley 622f1580f4eSBarry Smith .seealso: `PCJACOBI`, `PCJacobiaUseAbs()`, `PCJacobiSetType()` 62386697f06SMatthew Knepley @*/ 624d71ae5a4SJacob Faibussowitsch PetscErrorCode PCJacobiGetType(PC pc, PCJacobiType *type) 625d71ae5a4SJacob Faibussowitsch { 62686697f06SMatthew Knepley PetscFunctionBegin; 6270700a824SBarry Smith PetscValidHeaderSpecific(pc, PC_CLASSID, 1); 628cac4c232SBarry Smith PetscUseMethod(pc, "PCJacobiGetType_C", (PC, PCJacobiType *), (pc, type)); 6293ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 63086697f06SMatthew Knepley } 631