104c3f3b8SBarry Smith /* 25a46d3faSBarry Smith This file implements a Jacobi preconditioner in PETSc as part of PC. 35a46d3faSBarry Smith You can use this as a starting point for implementing your own 45a46d3faSBarry Smith preconditioner that is not provided with PETSc. (You might also consider 55a46d3faSBarry Smith just using PCSHELL) 64b9ad928SBarry Smith 74b9ad928SBarry Smith The following basic routines are required for each preconditioner. 84b9ad928SBarry Smith PCCreate_XXX() - Creates a preconditioner context 94b9ad928SBarry Smith PCSetFromOptions_XXX() - Sets runtime options 104b9ad928SBarry Smith PCApply_XXX() - Applies the preconditioner 114b9ad928SBarry Smith PCDestroy_XXX() - Destroys the preconditioner context 124b9ad928SBarry Smith where the suffix "_XXX" denotes a particular implementation, in 134b9ad928SBarry Smith this case we use _Jacobi (e.g., PCCreate_Jacobi, PCApply_Jacobi). 144b9ad928SBarry Smith These routines are actually called via the common user interface 154b9ad928SBarry Smith routines PCCreate(), PCSetFromOptions(), PCApply(), and PCDestroy(), 164b9ad928SBarry Smith so the application code interface remains identical for all 174b9ad928SBarry Smith preconditioners. 184b9ad928SBarry Smith 194b9ad928SBarry Smith Another key routine is: 204b9ad928SBarry Smith PCSetUp_XXX() - Prepares for the use of a preconditioner 214b9ad928SBarry Smith by setting data structures and options. The interface routine PCSetUp() 224b9ad928SBarry Smith is not usually called directly by the user, but instead is called by 234b9ad928SBarry Smith PCApply() if necessary. 244b9ad928SBarry Smith 254b9ad928SBarry Smith Additional basic routines are: 264b9ad928SBarry Smith PCView_XXX() - Prints details of runtime options that 274b9ad928SBarry Smith have actually been used. 284b9ad928SBarry Smith These are called by application codes via the interface routines 294b9ad928SBarry Smith PCView(). 304b9ad928SBarry Smith 314b9ad928SBarry Smith The various types of solvers (preconditioners, Krylov subspace methods, 324b9ad928SBarry Smith nonlinear solvers, timesteppers) are all organized similarly, so the 334b9ad928SBarry Smith above description applies to these categories also. One exception is 344b9ad928SBarry Smith that the analogues of PCApply() for these components are KSPSolve(), 354b9ad928SBarry Smith SNESSolve(), and TSSolve(). 364b9ad928SBarry Smith 374b9ad928SBarry Smith Additional optional functionality unique to preconditioners is left and 384b9ad928SBarry Smith right symmetric preconditioner application via PCApplySymmetricLeft() 394b9ad928SBarry Smith and PCApplySymmetricRight(). The Jacobi implementation is 404b9ad928SBarry Smith PCApplySymmetricLeftOrRight_Jacobi(). 4104c3f3b8SBarry Smith */ 424b9ad928SBarry Smith 434b9ad928SBarry Smith /* 444b9ad928SBarry Smith Include files needed for the Jacobi preconditioner: 454b9ad928SBarry Smith pcimpl.h - private include file intended for use by all preconditioners 464b9ad928SBarry Smith */ 474b9ad928SBarry Smith 48af0996ceSBarry Smith #include <petsc/private/pcimpl.h> /*I "petscpc.h" I*/ 494b9ad928SBarry Smith 500a545947SLisandro Dalcin const char *const PCJacobiTypes[] = {"DIAGONAL", "ROWMAX", "ROWSUM", "PCJacobiType", "PC_JACOBI_", NULL}; 51baa89ecbSBarry Smith 524b9ad928SBarry Smith /* 534b9ad928SBarry Smith Private context (data structure) for the Jacobi preconditioner. 544b9ad928SBarry Smith */ 554b9ad928SBarry Smith typedef struct { 562fa5cd67SKarl Rupp Vec diag; /* vector containing the reciprocals of the diagonal elements of the preconditioner matrix */ 574b9ad928SBarry Smith Vec diagsqrt; /* vector containing the reciprocals of the square roots of 584b9ad928SBarry Smith the diagonal elements of the preconditioner matrix (used 594b9ad928SBarry Smith only for symmetric preconditioner application) */ 60baa89ecbSBarry Smith PetscBool userowmax; /* set with PCJacobiSetType() */ 61ace3abfcSBarry Smith PetscBool userowsum; 62ace3abfcSBarry Smith PetscBool useabs; /* use the absolute values of the diagonal entries */ 6367ed0f3bSStefano Zampini PetscBool fixdiag; /* fix zero diagonal terms */ 644b9ad928SBarry Smith } PC_Jacobi; 654b9ad928SBarry Smith 667a660b1eSJed Brown static PetscErrorCode PCReset_Jacobi(PC); 677a660b1eSJed Brown 68d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCJacobiSetType_Jacobi(PC pc, PCJacobiType type) 69d71ae5a4SJacob Faibussowitsch { 70baa89ecbSBarry Smith PC_Jacobi *j = (PC_Jacobi *)pc->data; 717a660b1eSJed Brown PCJacobiType old_type; 724b9ad928SBarry Smith 734b9ad928SBarry Smith PetscFunctionBegin; 747a660b1eSJed Brown PetscCall(PCJacobiGetType(pc, &old_type)); 753ba16761SJacob Faibussowitsch if (old_type == type) PetscFunctionReturn(PETSC_SUCCESS); 767a660b1eSJed Brown PetscCall(PCReset_Jacobi(pc)); 77baa89ecbSBarry Smith j->userowmax = PETSC_FALSE; 78baa89ecbSBarry Smith j->userowsum = PETSC_FALSE; 79baa89ecbSBarry Smith if (type == PC_JACOBI_ROWMAX) { 804b9ad928SBarry Smith j->userowmax = PETSC_TRUE; 81baa89ecbSBarry Smith } else if (type == PC_JACOBI_ROWSUM) { 82baa89ecbSBarry Smith j->userowsum = PETSC_TRUE; 83baa89ecbSBarry Smith } 843ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 854b9ad928SBarry Smith } 864b9ad928SBarry Smith 87d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCJacobiGetType_Jacobi(PC pc, PCJacobiType *type) 88d71ae5a4SJacob Faibussowitsch { 89baa89ecbSBarry Smith PC_Jacobi *j = (PC_Jacobi *)pc->data; 9086697f06SMatthew Knepley 9186697f06SMatthew Knepley PetscFunctionBegin; 92baa89ecbSBarry Smith if (j->userowmax) { 93baa89ecbSBarry Smith *type = PC_JACOBI_ROWMAX; 94baa89ecbSBarry Smith } else if (j->userowsum) { 95baa89ecbSBarry Smith *type = PC_JACOBI_ROWSUM; 96baa89ecbSBarry Smith } else { 97baa89ecbSBarry Smith *type = PC_JACOBI_DIAGONAL; 98baa89ecbSBarry Smith } 993ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 10086697f06SMatthew Knepley } 10186697f06SMatthew Knepley 102d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCJacobiSetUseAbs_Jacobi(PC pc, PetscBool flg) 103d71ae5a4SJacob Faibussowitsch { 104baa89ecbSBarry Smith PC_Jacobi *j = (PC_Jacobi *)pc->data; 105cd47f5d9SBarry Smith 106cd47f5d9SBarry Smith PetscFunctionBegin; 107baa89ecbSBarry Smith j->useabs = flg; 1083ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 109baa89ecbSBarry Smith } 110baa89ecbSBarry Smith 111d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCJacobiGetUseAbs_Jacobi(PC pc, PetscBool *flg) 112d71ae5a4SJacob Faibussowitsch { 113baa89ecbSBarry Smith PC_Jacobi *j = (PC_Jacobi *)pc->data; 114baa89ecbSBarry Smith 115baa89ecbSBarry Smith PetscFunctionBegin; 116baa89ecbSBarry Smith *flg = j->useabs; 1173ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 118cd47f5d9SBarry Smith } 119cd47f5d9SBarry Smith 120d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCJacobiSetFixDiagonal_Jacobi(PC pc, PetscBool flg) 121d71ae5a4SJacob Faibussowitsch { 12267ed0f3bSStefano Zampini PC_Jacobi *j = (PC_Jacobi *)pc->data; 12367ed0f3bSStefano Zampini 12467ed0f3bSStefano Zampini PetscFunctionBegin; 12567ed0f3bSStefano Zampini j->fixdiag = flg; 1263ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 12767ed0f3bSStefano Zampini } 12867ed0f3bSStefano Zampini 129d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCJacobiGetFixDiagonal_Jacobi(PC pc, PetscBool *flg) 130d71ae5a4SJacob Faibussowitsch { 13167ed0f3bSStefano Zampini PC_Jacobi *j = (PC_Jacobi *)pc->data; 13267ed0f3bSStefano Zampini 13367ed0f3bSStefano Zampini PetscFunctionBegin; 13467ed0f3bSStefano Zampini *flg = j->fixdiag; 1353ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 13667ed0f3bSStefano Zampini } 13767ed0f3bSStefano Zampini 1384b9ad928SBarry Smith /* 1394b9ad928SBarry Smith PCSetUp_Jacobi - Prepares for the use of the Jacobi preconditioner 1404b9ad928SBarry Smith by setting data structures and options. 1414b9ad928SBarry Smith 1424b9ad928SBarry Smith Input Parameter: 1434b9ad928SBarry Smith . pc - the preconditioner context 1444b9ad928SBarry Smith 1454b9ad928SBarry Smith Application Interface Routine: PCSetUp() 1464b9ad928SBarry Smith 147f1580f4eSBarry Smith Note: 1484b9ad928SBarry Smith The interface routine PCSetUp() is not usually called directly by 1494b9ad928SBarry Smith the user, but instead is called by PCApply() if necessary. 1504b9ad928SBarry Smith */ 151d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCSetUp_Jacobi(PC pc) 152d71ae5a4SJacob Faibussowitsch { 1534b9ad928SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 1544b9ad928SBarry Smith Vec diag, diagsqrt; 155cd47f5d9SBarry Smith PetscInt n, i; 1564b9ad928SBarry Smith PetscScalar *x; 157ace3abfcSBarry Smith PetscBool zeroflag = PETSC_FALSE; 1584b9ad928SBarry Smith 1594b9ad928SBarry Smith PetscFunctionBegin; 1604b9ad928SBarry Smith /* 1614b9ad928SBarry Smith For most preconditioners the code would begin here something like 1624b9ad928SBarry Smith 1634b9ad928SBarry Smith if (pc->setupcalled == 0) { allocate space the first time this is ever called 1649566063dSJacob Faibussowitsch PetscCall(MatCreateVecs(pc->mat,&jac->diag)); 1654b9ad928SBarry Smith } 1664b9ad928SBarry Smith 1674b9ad928SBarry Smith But for this preconditioner we want to support use of both the matrix' diagonal 1684b9ad928SBarry Smith elements (for left or right preconditioning) and square root of diagonal elements 1694b9ad928SBarry Smith (for symmetric preconditioning). Hence we do not allocate space here, since we 1704b9ad928SBarry Smith don't know at this point which will be needed (diag and/or diagsqrt) until the user 1714b9ad928SBarry Smith applies the preconditioner, and we don't want to allocate BOTH unless we need 1724b9ad928SBarry Smith them both. Thus, the diag and diagsqrt are allocated in PCSetUp_Jacobi_NonSymmetric() 1734b9ad928SBarry Smith and PCSetUp_Jacobi_Symmetric(), respectively. 1744b9ad928SBarry Smith */ 1754b9ad928SBarry Smith 1764b9ad928SBarry Smith /* 1774b9ad928SBarry Smith Here we set up the preconditioner; that is, we copy the diagonal values from 1784b9ad928SBarry Smith the matrix and put them into a format to make them quick to apply as a preconditioner. 1794b9ad928SBarry Smith */ 1804b9ad928SBarry Smith diag = jac->diag; 1814b9ad928SBarry Smith diagsqrt = jac->diagsqrt; 1824b9ad928SBarry Smith 1834b9ad928SBarry Smith if (diag) { 184b94d7dedSBarry Smith PetscBool isset, isspd; 18567ed0f3bSStefano Zampini 1864b9ad928SBarry Smith if (jac->userowmax) { 1879566063dSJacob Faibussowitsch PetscCall(MatGetRowMaxAbs(pc->pmat, diag, NULL)); 18886697f06SMatthew Knepley } else if (jac->userowsum) { 1899566063dSJacob Faibussowitsch PetscCall(MatGetRowSum(pc->pmat, diag)); 1904b9ad928SBarry Smith } else { 1919566063dSJacob Faibussowitsch PetscCall(MatGetDiagonal(pc->pmat, diag)); 1924b9ad928SBarry Smith } 1939566063dSJacob Faibussowitsch PetscCall(VecReciprocal(diag)); 1941baa6e33SBarry Smith if (jac->useabs) PetscCall(VecAbs(diag)); 195b94d7dedSBarry Smith PetscCall(MatIsSPDKnown(pc->pmat, &isset, &isspd)); 196b94d7dedSBarry Smith if (jac->fixdiag && (!isset || !isspd)) { 1979566063dSJacob Faibussowitsch PetscCall(VecGetLocalSize(diag, &n)); 1989566063dSJacob Faibussowitsch PetscCall(VecGetArray(diag, &x)); 1994b9ad928SBarry Smith for (i = 0; i < n; i++) { 2004b9ad928SBarry Smith if (x[i] == 0.0) { 2014b9ad928SBarry Smith x[i] = 1.0; 202cd47f5d9SBarry Smith zeroflag = PETSC_TRUE; 2034b9ad928SBarry Smith } 2044b9ad928SBarry Smith } 2059566063dSJacob Faibussowitsch PetscCall(VecRestoreArray(diag, &x)); 2064b9ad928SBarry Smith } 20767ed0f3bSStefano Zampini } 2084b9ad928SBarry Smith if (diagsqrt) { 2094b9ad928SBarry Smith if (jac->userowmax) { 2109566063dSJacob Faibussowitsch PetscCall(MatGetRowMaxAbs(pc->pmat, diagsqrt, NULL)); 21186697f06SMatthew Knepley } else if (jac->userowsum) { 2129566063dSJacob Faibussowitsch PetscCall(MatGetRowSum(pc->pmat, diagsqrt)); 2134b9ad928SBarry Smith } else { 2149566063dSJacob Faibussowitsch PetscCall(MatGetDiagonal(pc->pmat, diagsqrt)); 2154b9ad928SBarry Smith } 2169566063dSJacob Faibussowitsch PetscCall(VecGetLocalSize(diagsqrt, &n)); 2179566063dSJacob Faibussowitsch PetscCall(VecGetArray(diagsqrt, &x)); 2184b9ad928SBarry Smith for (i = 0; i < n; i++) { 2198f1a2a5eSBarry Smith if (x[i] != 0.0) x[i] = 1.0 / PetscSqrtReal(PetscAbsScalar(x[i])); 2204b9ad928SBarry Smith else { 2214b9ad928SBarry Smith x[i] = 1.0; 222cd47f5d9SBarry Smith zeroflag = PETSC_TRUE; 2234b9ad928SBarry Smith } 2244b9ad928SBarry Smith } 2259566063dSJacob Faibussowitsch PetscCall(VecRestoreArray(diagsqrt, &x)); 2264b9ad928SBarry Smith } 22748a46eb9SPierre Jolivet if (zeroflag) PetscCall(PetscInfo(pc, "Zero detected in diagonal of matrix, using 1 at those locations\n")); 2283ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 2294b9ad928SBarry Smith } 230f1580f4eSBarry Smith 2314b9ad928SBarry Smith /* 2324b9ad928SBarry Smith PCSetUp_Jacobi_Symmetric - Allocates the vector needed to store the 2334b9ad928SBarry Smith inverse of the square root of the diagonal entries of the matrix. This 2344b9ad928SBarry Smith is used for symmetric application of the Jacobi preconditioner. 2354b9ad928SBarry Smith 2364b9ad928SBarry Smith Input Parameter: 2374b9ad928SBarry Smith . pc - the preconditioner context 2384b9ad928SBarry Smith */ 239d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCSetUp_Jacobi_Symmetric(PC pc) 240d71ae5a4SJacob Faibussowitsch { 2414b9ad928SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 2424b9ad928SBarry Smith 2434b9ad928SBarry Smith PetscFunctionBegin; 2449566063dSJacob Faibussowitsch PetscCall(MatCreateVecs(pc->pmat, &jac->diagsqrt, NULL)); 2459566063dSJacob Faibussowitsch PetscCall(PCSetUp_Jacobi(pc)); 2463ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 2474b9ad928SBarry Smith } 248f1580f4eSBarry Smith 2494b9ad928SBarry Smith /* 2504b9ad928SBarry Smith PCSetUp_Jacobi_NonSymmetric - Allocates the vector needed to store the 2514b9ad928SBarry Smith inverse of the diagonal entries of the matrix. This is used for left of 2524b9ad928SBarry Smith right application of the Jacobi preconditioner. 2534b9ad928SBarry Smith 2544b9ad928SBarry Smith Input Parameter: 2554b9ad928SBarry Smith . pc - the preconditioner context 2564b9ad928SBarry Smith */ 257d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCSetUp_Jacobi_NonSymmetric(PC pc) 258d71ae5a4SJacob Faibussowitsch { 2594b9ad928SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 2604b9ad928SBarry Smith 2614b9ad928SBarry Smith PetscFunctionBegin; 2629566063dSJacob Faibussowitsch PetscCall(MatCreateVecs(pc->pmat, &jac->diag, NULL)); 2639566063dSJacob Faibussowitsch PetscCall(PCSetUp_Jacobi(pc)); 2643ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 2654b9ad928SBarry Smith } 266f1580f4eSBarry Smith 2674b9ad928SBarry Smith /* 2684b9ad928SBarry Smith PCApply_Jacobi - Applies the Jacobi preconditioner to a vector. 2694b9ad928SBarry Smith 2704b9ad928SBarry Smith Input Parameters: 2714b9ad928SBarry Smith . pc - the preconditioner context 2724b9ad928SBarry Smith . x - input vector 2734b9ad928SBarry Smith 2744b9ad928SBarry Smith Output Parameter: 2754b9ad928SBarry Smith . y - output vector 2764b9ad928SBarry Smith 2774b9ad928SBarry Smith Application Interface Routine: PCApply() 2784b9ad928SBarry Smith */ 279d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCApply_Jacobi(PC pc, Vec x, Vec y) 280d71ae5a4SJacob Faibussowitsch { 2814b9ad928SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 2824b9ad928SBarry Smith 2834b9ad928SBarry Smith PetscFunctionBegin; 28448a46eb9SPierre Jolivet if (!jac->diag) PetscCall(PCSetUp_Jacobi_NonSymmetric(pc)); 2859566063dSJacob Faibussowitsch PetscCall(VecPointwiseMult(y, x, jac->diag)); 2863ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 2874b9ad928SBarry Smith } 288f1580f4eSBarry Smith 2894b9ad928SBarry Smith /* 2904b9ad928SBarry Smith PCApplySymmetricLeftOrRight_Jacobi - Applies the left or right part of a 2914b9ad928SBarry Smith symmetric preconditioner to a vector. 2924b9ad928SBarry Smith 2934b9ad928SBarry Smith Input Parameters: 2944b9ad928SBarry Smith . pc - the preconditioner context 2954b9ad928SBarry Smith . x - input vector 2964b9ad928SBarry Smith 2974b9ad928SBarry Smith Output Parameter: 2984b9ad928SBarry Smith . y - output vector 2994b9ad928SBarry Smith 3004b9ad928SBarry Smith Application Interface Routines: PCApplySymmetricLeft(), PCApplySymmetricRight() 3014b9ad928SBarry Smith */ 302d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCApplySymmetricLeftOrRight_Jacobi(PC pc, Vec x, Vec y) 303d71ae5a4SJacob Faibussowitsch { 3044b9ad928SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 3054b9ad928SBarry Smith 3064b9ad928SBarry Smith PetscFunctionBegin; 30748a46eb9SPierre Jolivet if (!jac->diagsqrt) PetscCall(PCSetUp_Jacobi_Symmetric(pc)); 3089566063dSJacob Faibussowitsch PetscCall(VecPointwiseMult(y, x, jac->diagsqrt)); 3093ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 3104b9ad928SBarry Smith } 31167ed0f3bSStefano Zampini 312d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCReset_Jacobi(PC pc) 313d71ae5a4SJacob Faibussowitsch { 314a06653b4SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 315a06653b4SBarry Smith 316a06653b4SBarry Smith PetscFunctionBegin; 3179566063dSJacob Faibussowitsch PetscCall(VecDestroy(&jac->diag)); 3189566063dSJacob Faibussowitsch PetscCall(VecDestroy(&jac->diagsqrt)); 3193ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 320a06653b4SBarry Smith } 321a06653b4SBarry Smith 3224b9ad928SBarry Smith /* 3234b9ad928SBarry Smith PCDestroy_Jacobi - Destroys the private context for the Jacobi preconditioner 3244b9ad928SBarry Smith that was created with PCCreate_Jacobi(). 3254b9ad928SBarry Smith 3264b9ad928SBarry Smith Input Parameter: 3274b9ad928SBarry Smith . pc - the preconditioner context 3284b9ad928SBarry Smith 3294b9ad928SBarry Smith Application Interface Routine: PCDestroy() 3304b9ad928SBarry Smith */ 331d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCDestroy_Jacobi(PC pc) 332d71ae5a4SJacob Faibussowitsch { 3334b9ad928SBarry Smith PetscFunctionBegin; 3349566063dSJacob Faibussowitsch PetscCall(PCReset_Jacobi(pc)); 3359566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiSetType_C", NULL)); 3369566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiGetType_C", NULL)); 3379566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiSetUseAbs_C", NULL)); 3389566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiGetUseAbs_C", NULL)); 3399566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiSetFixDiagonal_C", NULL)); 3409566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiGetFixDiagonal_C", NULL)); 3414b9ad928SBarry Smith 3424b9ad928SBarry Smith /* 3434b9ad928SBarry Smith Free the private data structure that was hanging off the PC 3444b9ad928SBarry Smith */ 3459566063dSJacob Faibussowitsch PetscCall(PetscFree(pc->data)); 3463ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 3474b9ad928SBarry Smith } 3484b9ad928SBarry Smith 349d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCSetFromOptions_Jacobi(PC pc, PetscOptionItems *PetscOptionsObject) 350d71ae5a4SJacob Faibussowitsch { 3514b9ad928SBarry Smith PC_Jacobi *jac = (PC_Jacobi *)pc->data; 352baa89ecbSBarry Smith PetscBool flg; 353baa89ecbSBarry Smith PCJacobiType deflt, type; 3544b9ad928SBarry Smith 3554b9ad928SBarry Smith PetscFunctionBegin; 3569566063dSJacob Faibussowitsch PetscCall(PCJacobiGetType(pc, &deflt)); 357d0609cedSBarry Smith PetscOptionsHeadBegin(PetscOptionsObject, "Jacobi options"); 3589566063dSJacob Faibussowitsch PetscCall(PetscOptionsEnum("-pc_jacobi_type", "How to construct diagonal matrix", "PCJacobiSetType", PCJacobiTypes, (PetscEnum)deflt, (PetscEnum *)&type, &flg)); 3591baa6e33SBarry Smith if (flg) PetscCall(PCJacobiSetType(pc, type)); 3609566063dSJacob Faibussowitsch PetscCall(PetscOptionsBool("-pc_jacobi_abs", "Use absolute values of diagonal entries", "PCJacobiSetUseAbs", jac->useabs, &jac->useabs, NULL)); 3619566063dSJacob Faibussowitsch PetscCall(PetscOptionsBool("-pc_jacobi_fixdiagonal", "Fix null terms on diagonal", "PCJacobiSetFixDiagonal", jac->fixdiag, &jac->fixdiag, NULL)); 362d0609cedSBarry Smith PetscOptionsHeadEnd(); 3633ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 3644b9ad928SBarry Smith } 3654b9ad928SBarry Smith 366d71ae5a4SJacob Faibussowitsch static PetscErrorCode PCView_Jacobi(PC pc, PetscViewer viewer) 367d71ae5a4SJacob Faibussowitsch { 368569e28a7SMatthew G. Knepley PC_Jacobi *jac = (PC_Jacobi *)pc->data; 369569e28a7SMatthew G. Knepley PetscBool iascii; 370569e28a7SMatthew G. Knepley 371569e28a7SMatthew G. Knepley PetscFunctionBegin; 3729566063dSJacob Faibussowitsch PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii)); 373569e28a7SMatthew G. Knepley if (iascii) { 374569e28a7SMatthew G. Knepley PCJacobiType type; 37567ed0f3bSStefano Zampini PetscBool useAbs, fixdiag; 376569e28a7SMatthew G. Knepley PetscViewerFormat format; 377569e28a7SMatthew G. Knepley 3789566063dSJacob Faibussowitsch PetscCall(PCJacobiGetType(pc, &type)); 3799566063dSJacob Faibussowitsch PetscCall(PCJacobiGetUseAbs(pc, &useAbs)); 3809566063dSJacob Faibussowitsch PetscCall(PCJacobiGetFixDiagonal(pc, &fixdiag)); 3819566063dSJacob Faibussowitsch PetscCall(PetscViewerASCIIPrintf(viewer, " type %s%s%s\n", PCJacobiTypes[type], useAbs ? ", using absolute value of entries" : "", !fixdiag ? ", not checking null diagonal entries" : "")); 3829566063dSJacob Faibussowitsch PetscCall(PetscViewerGetFormat(viewer, &format)); 383b24842d1SPierre Jolivet if (format == PETSC_VIEWER_ASCII_INFO_DETAIL && jac->diag) PetscCall(VecView(jac->diag, viewer)); 384569e28a7SMatthew G. Knepley } 3853ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 386569e28a7SMatthew G. Knepley } 387569e28a7SMatthew G. Knepley 3884b9ad928SBarry Smith /* 3894b9ad928SBarry Smith PCCreate_Jacobi - Creates a Jacobi preconditioner context, PC_Jacobi, 3904b9ad928SBarry Smith and sets this as the private data within the generic preconditioning 3914b9ad928SBarry Smith context, PC, that was created within PCCreate(). 3924b9ad928SBarry Smith 3934b9ad928SBarry Smith Input Parameter: 3944b9ad928SBarry Smith . pc - the preconditioner context 3954b9ad928SBarry Smith 3964b9ad928SBarry Smith Application Interface Routine: PCCreate() 3974b9ad928SBarry Smith */ 3984b9ad928SBarry Smith 3994b9ad928SBarry Smith /*MC 4005a46d3faSBarry Smith PCJACOBI - Jacobi (i.e. diagonal scaling preconditioning) 4014b9ad928SBarry Smith 402f1580f4eSBarry Smith Options Database Keys: 403967c93d3SBarry Smith + -pc_jacobi_type <diagonal,rowmax,rowsum> - approach for forming the preconditioner 40467ed0f3bSStefano Zampini . -pc_jacobi_abs - use the absolute value of the diagonal entry 405147403d9SBarry Smith - -pc_jacobi_fixdiag - fix for zero diagonal terms by placing 1.0 in those locations 4064b9ad928SBarry Smith 4074b9ad928SBarry Smith Level: beginner 4084b9ad928SBarry Smith 40995452b02SPatrick Sanan Notes: 410f1580f4eSBarry Smith By using `KSPSetPCSide`(ksp,`PC_SYMMETRIC`) or -ksp_pc_side symmetric 4114b9ad928SBarry Smith can scale each side of the matrix by the square root of the diagonal entries. 4124b9ad928SBarry Smith 4134b9ad928SBarry Smith Zero entries along the diagonal are replaced with the value 1.0 4144b9ad928SBarry Smith 415f1580f4eSBarry Smith See `PCPBJACOBI` for fixed-size point block, `PCVPBJACOBI` for variable-sized point block, and `PCBJACOBI` for large size blocks 416422a814eSBarry Smith 417db781477SPatrick Sanan .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PC`, 418f1580f4eSBarry Smith `PCJacobiSetType()`, `PCJacobiSetUseAbs()`, `PCJacobiGetUseAbs()`, `PCASM`, 419db781477SPatrick Sanan `PCJacobiSetFixDiagonal()`, `PCJacobiGetFixDiagonal()` 420db781477SPatrick Sanan `PCJacobiSetType()`, `PCJacobiSetUseAbs()`, `PCJacobiGetUseAbs()`, `PCPBJACOBI`, `PCBJACOBI`, `PCVPBJACOBI` 4214b9ad928SBarry Smith M*/ 4224b9ad928SBarry Smith 423d71ae5a4SJacob Faibussowitsch PETSC_EXTERN PetscErrorCode PCCreate_Jacobi(PC pc) 424d71ae5a4SJacob Faibussowitsch { 4254b9ad928SBarry Smith PC_Jacobi *jac; 4264b9ad928SBarry Smith 4274b9ad928SBarry Smith PetscFunctionBegin; 4284b9ad928SBarry Smith /* 4294b9ad928SBarry Smith Creates the private data structure for this preconditioner and 4304b9ad928SBarry Smith attach it to the PC object. 4314b9ad928SBarry Smith */ 4324dfa11a4SJacob Faibussowitsch PetscCall(PetscNew(&jac)); 4334b9ad928SBarry Smith pc->data = (void *)jac; 4344b9ad928SBarry Smith 4354b9ad928SBarry Smith /* 4364b9ad928SBarry Smith Initialize the pointers to vectors to ZERO; these will be used to store 4374b9ad928SBarry Smith diagonal entries of the matrix for fast preconditioner application. 4384b9ad928SBarry Smith */ 4390a545947SLisandro Dalcin jac->diag = NULL; 4400a545947SLisandro Dalcin jac->diagsqrt = NULL; 4414b9ad928SBarry Smith jac->userowmax = PETSC_FALSE; 44286697f06SMatthew Knepley jac->userowsum = PETSC_FALSE; 443cd47f5d9SBarry Smith jac->useabs = PETSC_FALSE; 44467ed0f3bSStefano Zampini jac->fixdiag = PETSC_TRUE; 4454b9ad928SBarry Smith 4464b9ad928SBarry Smith /* 4474b9ad928SBarry Smith Set the pointers for the functions that are provided above. 4484b9ad928SBarry Smith Now when the user-level routines (such as PCApply(), PCDestroy(), etc.) 4494b9ad928SBarry Smith are called, they will automatically call these functions. Note we 4504b9ad928SBarry Smith choose not to provide a couple of these functions since they are 4514b9ad928SBarry Smith not needed. 4524b9ad928SBarry Smith */ 4534b9ad928SBarry Smith pc->ops->apply = PCApply_Jacobi; 4544b9ad928SBarry Smith pc->ops->applytranspose = PCApply_Jacobi; 4554b9ad928SBarry Smith pc->ops->setup = PCSetUp_Jacobi; 456a06653b4SBarry Smith pc->ops->reset = PCReset_Jacobi; 4574b9ad928SBarry Smith pc->ops->destroy = PCDestroy_Jacobi; 4584b9ad928SBarry Smith pc->ops->setfromoptions = PCSetFromOptions_Jacobi; 459569e28a7SMatthew G. Knepley pc->ops->view = PCView_Jacobi; 4600a545947SLisandro Dalcin pc->ops->applyrichardson = NULL; 4614b9ad928SBarry Smith pc->ops->applysymmetricleft = PCApplySymmetricLeftOrRight_Jacobi; 4624b9ad928SBarry Smith pc->ops->applysymmetricright = PCApplySymmetricLeftOrRight_Jacobi; 4632fa5cd67SKarl Rupp 4649566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiSetType_C", PCJacobiSetType_Jacobi)); 4659566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiGetType_C", PCJacobiGetType_Jacobi)); 4669566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiSetUseAbs_C", PCJacobiSetUseAbs_Jacobi)); 4679566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiGetUseAbs_C", PCJacobiGetUseAbs_Jacobi)); 4689566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiSetFixDiagonal_C", PCJacobiSetFixDiagonal_Jacobi)); 4699566063dSJacob Faibussowitsch PetscCall(PetscObjectComposeFunction((PetscObject)pc, "PCJacobiGetFixDiagonal_C", PCJacobiGetFixDiagonal_Jacobi)); 4703ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 4714b9ad928SBarry Smith } 472cd47f5d9SBarry Smith 473cd47f5d9SBarry Smith /*@ 474f1580f4eSBarry Smith PCJacobiSetUseAbs - Causes the Jacobi preconditioner `PCJACOBI` to use the 47567ed0f3bSStefano Zampini absolute values of the diagonal divisors in the preconditioner 476cd47f5d9SBarry Smith 477c3339decSBarry Smith Logically Collective 478cd47f5d9SBarry Smith 479cd47f5d9SBarry Smith Input Parameters: 480baa89ecbSBarry Smith + pc - the preconditioner context 481baa89ecbSBarry Smith - flg - whether to use absolute values or not 482baa89ecbSBarry Smith 483baa89ecbSBarry Smith Options Database Key: 484147403d9SBarry Smith . -pc_jacobi_abs <bool> - use absolute values 485baa89ecbSBarry Smith 486f1580f4eSBarry Smith Note: 48795452b02SPatrick Sanan This takes affect at the next construction of the preconditioner 488baa89ecbSBarry Smith 489baa89ecbSBarry Smith Level: intermediate 490baa89ecbSBarry Smith 491*562efe2eSBarry Smith .seealso: [](ch_ksp), `PCJACOBI`, `PCJacobiaSetType()`, `PCJacobiGetUseAbs()` 492baa89ecbSBarry Smith @*/ 493d71ae5a4SJacob Faibussowitsch PetscErrorCode PCJacobiSetUseAbs(PC pc, PetscBool flg) 494d71ae5a4SJacob Faibussowitsch { 495baa89ecbSBarry Smith PetscFunctionBegin; 496baa89ecbSBarry Smith PetscValidHeaderSpecific(pc, PC_CLASSID, 1); 497cac4c232SBarry Smith PetscTryMethod(pc, "PCJacobiSetUseAbs_C", (PC, PetscBool), (pc, flg)); 4983ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 499baa89ecbSBarry Smith } 500baa89ecbSBarry Smith 501baa89ecbSBarry Smith /*@ 502f1580f4eSBarry Smith PCJacobiGetUseAbs - Determines if the Jacobi preconditioner `PCJACOBI` uses the 50367ed0f3bSStefano Zampini absolute values of the diagonal divisors in the preconditioner 504baa89ecbSBarry Smith 505c3339decSBarry Smith Logically Collective 506baa89ecbSBarry Smith 507baa89ecbSBarry Smith Input Parameter: 508cd47f5d9SBarry Smith . pc - the preconditioner context 509cd47f5d9SBarry Smith 510baa89ecbSBarry Smith Output Parameter: 511baa89ecbSBarry Smith . flg - whether to use absolute values or not 512cd47f5d9SBarry Smith 513cd47f5d9SBarry Smith Level: intermediate 514cd47f5d9SBarry Smith 515*562efe2eSBarry Smith .seealso: [](ch_ksp), `PCJACOBI`, `PCJacobiaSetType()`, `PCJacobiSetUseAbs()`, `PCJacobiGetType()` 516cd47f5d9SBarry Smith @*/ 517d71ae5a4SJacob Faibussowitsch PetscErrorCode PCJacobiGetUseAbs(PC pc, PetscBool *flg) 518d71ae5a4SJacob Faibussowitsch { 519cd47f5d9SBarry Smith PetscFunctionBegin; 5200700a824SBarry Smith PetscValidHeaderSpecific(pc, PC_CLASSID, 1); 521cac4c232SBarry Smith PetscUseMethod(pc, "PCJacobiGetUseAbs_C", (PC, PetscBool *), (pc, flg)); 5223ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 523cd47f5d9SBarry Smith } 524cd47f5d9SBarry Smith 5254b9ad928SBarry Smith /*@ 526147403d9SBarry Smith PCJacobiSetFixDiagonal - Check for zero values on the diagonal and replace them with 1.0 52767ed0f3bSStefano Zampini 528c3339decSBarry Smith Logically Collective 52967ed0f3bSStefano Zampini 53067ed0f3bSStefano Zampini Input Parameters: 53167ed0f3bSStefano Zampini + pc - the preconditioner context 53267ed0f3bSStefano Zampini - flg - the boolean flag 53367ed0f3bSStefano Zampini 53467ed0f3bSStefano Zampini Options Database Key: 535147403d9SBarry Smith . -pc_jacobi_fixdiagonal <bool> - check for zero values on the diagonal 53667ed0f3bSStefano Zampini 537f1580f4eSBarry Smith Note: 53867ed0f3bSStefano Zampini This takes affect at the next construction of the preconditioner 53967ed0f3bSStefano Zampini 54067ed0f3bSStefano Zampini Level: intermediate 54167ed0f3bSStefano Zampini 542*562efe2eSBarry Smith .seealso: [](ch_ksp), `PCJACOBI`, `PCJacobiSetType()`, `PCJacobiGetFixDiagonal()`, `PCJacobiSetUseAbs()` 54367ed0f3bSStefano Zampini @*/ 544d71ae5a4SJacob Faibussowitsch PetscErrorCode PCJacobiSetFixDiagonal(PC pc, PetscBool flg) 545d71ae5a4SJacob Faibussowitsch { 54667ed0f3bSStefano Zampini PetscFunctionBegin; 54767ed0f3bSStefano Zampini PetscValidHeaderSpecific(pc, PC_CLASSID, 1); 548cac4c232SBarry Smith PetscTryMethod(pc, "PCJacobiSetFixDiagonal_C", (PC, PetscBool), (pc, flg)); 5493ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 55067ed0f3bSStefano Zampini } 55167ed0f3bSStefano Zampini 55267ed0f3bSStefano Zampini /*@ 553f1580f4eSBarry Smith PCJacobiGetFixDiagonal - Determines if the Jacobi preconditioner `PCJACOBI` checks for zero diagonal terms 55467ed0f3bSStefano Zampini 555c3339decSBarry Smith Logically Collective 55667ed0f3bSStefano Zampini 55767ed0f3bSStefano Zampini Input Parameter: 55867ed0f3bSStefano Zampini . pc - the preconditioner context 55967ed0f3bSStefano Zampini 56067ed0f3bSStefano Zampini Output Parameter: 56167ed0f3bSStefano Zampini . flg - the boolean flag 56267ed0f3bSStefano Zampini 56367ed0f3bSStefano Zampini Options Database Key: 56467b8a455SSatish Balay . -pc_jacobi_fixdiagonal <bool> - Fix 0 terms on diagonal by using 1 56567ed0f3bSStefano Zampini 56667ed0f3bSStefano Zampini Level: intermediate 56767ed0f3bSStefano Zampini 568*562efe2eSBarry Smith .seealso: [](ch_ksp), `PCJACOBI`, `PCJacobiSetType()`, `PCJacobiSetFixDiagonal()` 56967ed0f3bSStefano Zampini @*/ 570d71ae5a4SJacob Faibussowitsch PetscErrorCode PCJacobiGetFixDiagonal(PC pc, PetscBool *flg) 571d71ae5a4SJacob Faibussowitsch { 57267ed0f3bSStefano Zampini PetscFunctionBegin; 57367ed0f3bSStefano Zampini PetscValidHeaderSpecific(pc, PC_CLASSID, 1); 574cac4c232SBarry Smith PetscUseMethod(pc, "PCJacobiGetFixDiagonal_C", (PC, PetscBool *), (pc, flg)); 5753ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 57667ed0f3bSStefano Zampini } 57767ed0f3bSStefano Zampini 57867ed0f3bSStefano Zampini /*@ 579baa89ecbSBarry Smith PCJacobiSetType - Causes the Jacobi preconditioner to use either the diagonal, the maximum entry in each row, 580baa89ecbSBarry Smith of the sum of rows entries for the diagonal preconditioner 5814b9ad928SBarry Smith 582c3339decSBarry Smith Logically Collective 5834b9ad928SBarry Smith 5844b9ad928SBarry Smith Input Parameters: 585baa89ecbSBarry Smith + pc - the preconditioner context 586f1580f4eSBarry Smith - type - `PC_JACOBI_DIAGONAL`, `PC_JACOBI_ROWMAX`, `PC_JACOBI_ROWSUM` 5874b9ad928SBarry Smith 5884b9ad928SBarry Smith Options Database Key: 589147403d9SBarry Smith . -pc_jacobi_type <diagonal,rowmax,rowsum> - the type of diagonal matrix to use for Jacobi 5904b9ad928SBarry Smith 5914b9ad928SBarry Smith Level: intermediate 5924b9ad928SBarry Smith 593feefa0e1SJacob Faibussowitsch Developer Notes: 594f1580f4eSBarry Smith Why is there a separate function for using the absolute value? 595f1580f4eSBarry Smith 596*562efe2eSBarry Smith .seealso: [](ch_ksp), `PCJACOBI`, `PCJacobiSetUseAbs()`, `PCJacobiGetType()` 5974b9ad928SBarry Smith @*/ 598d71ae5a4SJacob Faibussowitsch PetscErrorCode PCJacobiSetType(PC pc, PCJacobiType type) 599d71ae5a4SJacob Faibussowitsch { 6004b9ad928SBarry Smith PetscFunctionBegin; 6010700a824SBarry Smith PetscValidHeaderSpecific(pc, PC_CLASSID, 1); 602cac4c232SBarry Smith PetscTryMethod(pc, "PCJacobiSetType_C", (PC, PCJacobiType), (pc, type)); 6033ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 6044b9ad928SBarry Smith } 6054b9ad928SBarry Smith 60686697f06SMatthew Knepley /*@ 607baa89ecbSBarry Smith PCJacobiGetType - Gets how the diagonal matrix is produced for the preconditioner 60886697f06SMatthew Knepley 609c3339decSBarry Smith Not Collective 61086697f06SMatthew Knepley 611baa89ecbSBarry Smith Input Parameter: 61286697f06SMatthew Knepley . pc - the preconditioner context 61386697f06SMatthew Knepley 614baa89ecbSBarry Smith Output Parameter: 615f1580f4eSBarry Smith . type - `PC_JACOBI_DIAGONAL`, `PC_JACOBI_ROWMAX`, `PC_JACOBI_ROWSUM` 61686697f06SMatthew Knepley 61786697f06SMatthew Knepley Level: intermediate 61886697f06SMatthew Knepley 619*562efe2eSBarry Smith .seealso: [](ch_ksp), `PCJACOBI`, `PCJacobiaUseAbs()`, `PCJacobiSetType()` 62086697f06SMatthew Knepley @*/ 621d71ae5a4SJacob Faibussowitsch PetscErrorCode PCJacobiGetType(PC pc, PCJacobiType *type) 622d71ae5a4SJacob Faibussowitsch { 62386697f06SMatthew Knepley PetscFunctionBegin; 6240700a824SBarry Smith PetscValidHeaderSpecific(pc, PC_CLASSID, 1); 625cac4c232SBarry Smith PetscUseMethod(pc, "PCJacobiGetType_C", (PC, PCJacobiType *), (pc, type)); 6263ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS); 62786697f06SMatthew Knepley } 628