1 /*----------------------------------------------------------------- 2 * Programmer(s): Daniel R. Reynolds @ SMU 3 * David J. Gardner, Radu Serban and Aaron Collier @ LLNL 4 *----------------------------------------------------------------- 5 * SUNDIALS Copyright Start 6 * Copyright (c) 2002-2020, Lawrence Livermore National Security 7 * and Southern Methodist University. 8 * All rights reserved. 9 * 10 * See the top-level LICENSE and NOTICE files for details. 11 * 12 * SPDX-License-Identifier: BSD-3-Clause 13 * SUNDIALS Copyright End 14 *----------------------------------------------------------------- 15 * Implementation header file for KINSOL's linear solver interface. 16 *-----------------------------------------------------------------*/ 17 18 #ifndef _KINLS_IMPL_H 19 #define _KINLS_IMPL_H 20 21 #include <kinsol/kinsol_ls.h> 22 #include "kinsol_impl.h" 23 24 #ifdef __cplusplus /* wrapper to enable C++ usage */ 25 extern "C" { 26 #endif 27 28 29 /*------------------------------------------------------------------ 30 keys for KINPrintInfo (do not use 1 -> conflict with PRNT_RETVAL) 31 ------------------------------------------------------------------*/ 32 #define PRNT_NLI 101 33 #define PRNT_EPS 102 34 35 36 /*------------------------------------------------------------------ 37 Types : struct KINLsMemRec, struct *KINLsMem 38 39 The type KINLsMem is a pointer to a KINLsMemRec, which is a 40 structure containing fields that must be accessible by LS module 41 routines. 42 ------------------------------------------------------------------*/ 43 typedef struct KINLsMemRec { 44 45 /* Linear solver type information */ 46 booleantype iterative; /* is the solver iterative? */ 47 booleantype matrixbased; /* is a matrix structure used? */ 48 49 /* Jacobian construction & storage */ 50 booleantype jacDQ; /* SUNTRUE if using internal DQ Jacobian approx. */ 51 KINLsJacFn jac; /* Jacobian routine to be called */ 52 void *J_data; /* J_data is passed to jac */ 53 54 /* Linear solver, matrix and vector objects/pointers */ 55 SUNLinearSolver LS; /* generic iterative linear solver object */ 56 SUNMatrix J; /* problem Jacobian */ 57 58 /* Solver tolerance adjustment factor (if needed, see kinLsSolve) */ 59 realtype tol_fac; 60 61 /* Statistics and associated parameters */ 62 long int nje; /* no. of calls to jac */ 63 long int nfeDQ; /* no. of calls to F due to DQ Jacobian or J*v 64 approximations */ 65 long int npe; /* npe = total number of precond calls */ 66 long int nli; /* nli = total number of linear iterations */ 67 long int nps; /* nps = total number of psolve calls */ 68 long int ncfl; /* ncfl = total number of convergence failures */ 69 long int njtimes; /* njtimes = total number of calls to jtimes */ 70 71 booleantype new_uu; /* flag indicating if the iterate has been 72 updated - the Jacobian must be updated or 73 reevaluated (meant to be used by a 74 user-supplied jtimes function */ 75 76 int last_flag; /* last error return flag */ 77 78 /* Preconditioner computation 79 (a) user-provided: 80 - pdata == user_data 81 - pfree == NULL (the user dealocates memory) 82 (b) internal preconditioner module 83 - pdata == kin_mem 84 - pfree == set by the prec. module and called in kinLsFree */ 85 KINLsPrecSetupFn pset; 86 KINLsPrecSolveFn psolve; 87 int (*pfree)(KINMem kin_mem); 88 void *pdata; 89 90 /* Jacobian times vector compuation 91 (a) jtimes function provided by the user: 92 - jt_data == user_data 93 - jtimesDQ == SUNFALSE 94 (b) internal jtimes 95 - jt_data == kin_mem 96 - jtimesDQ == SUNTRUE */ 97 booleantype jtimesDQ; 98 KINLsJacTimesVecFn jtimes; 99 KINSysFn jt_func; 100 void *jt_data; 101 102 } *KINLsMem; 103 104 105 /*------------------------------------------------------------------ 106 Prototypes of internal functions 107 ------------------------------------------------------------------*/ 108 109 /* Interface routines called by system SUNLinearSolvers */ 110 int kinLsATimes(void *kinmem, N_Vector v, N_Vector z); 111 int kinLsPSetup(void *kinmem); 112 int kinLsPSolve(void *kinmem, N_Vector r, N_Vector z, 113 realtype tol, int lr); 114 115 /* Difference quotient approximation for Jacobian times vector */ 116 int kinLsDQJtimes(N_Vector v, N_Vector Jv, N_Vector u, 117 booleantype *new_u, void *data); 118 119 /* Difference-quotient Jacobian approximation routines */ 120 int kinLsDQJac(N_Vector u, N_Vector fu, SUNMatrix Jac, 121 void *data, N_Vector tmp1, N_Vector tmp2); 122 123 int kinLsDenseDQJac(N_Vector u, N_Vector fu, SUNMatrix Jac, 124 KINMem kin_mem, N_Vector tmp1, N_Vector tmp2); 125 126 int kinLsBandDQJac(N_Vector u, N_Vector fu, SUNMatrix Jac, 127 KINMem kin_mem, N_Vector tmp1, N_Vector tmp2); 128 129 /* Generic linit/lsetup/lsolve/lfree interface routines for KINSOL to call */ 130 int kinLsInitialize(KINMem kin_mem); 131 int kinLsSetup(KINMem kin_mem); 132 int kinLsSolve(KINMem kin_mem, N_Vector x, N_Vector b, 133 realtype *sJpnorm, realtype *sFdotJp); 134 int kinLsFree(KINMem kin_mem); 135 136 /* Auxilliary functions */ 137 int kinLsInitializeCounters(KINLsMem kinls_mem); 138 int kinLs_AccessLMem(void* kinmem, const char *fname, 139 KINMem* kin_mem, KINLsMem *kinls_mem); 140 141 142 /*------------------------------------------------------------------ 143 Error messages 144 ------------------------------------------------------------------*/ 145 146 #define MSG_LS_KINMEM_NULL "KINSOL memory is NULL." 147 #define MSG_LS_MEM_FAIL "A memory request failed." 148 #define MSG_LS_BAD_NVECTOR "A required vector operation is not implemented." 149 #define MSG_LS_LMEM_NULL "Linear solver memory is NULL." 150 #define MSG_LS_NEG_MAXRS "maxrs < 0 illegal." 151 #define MSG_LS_BAD_SIZES "Illegal bandwidth parameter(s). Must have 0 <= ml, mu <= N-1." 152 153 #define MSG_LS_JACFUNC_FAILED "The Jacobian routine failed in an unrecoverable manner." 154 #define MSG_LS_PSET_FAILED "The preconditioner setup routine failed in an unrecoverable manner." 155 #define MSG_LS_PSOLVE_FAILED "The preconditioner solve routine failed in an unrecoverable manner." 156 #define MSG_LS_JTIMES_FAILED "The Jacobian x vector routine failed in an unrecoverable manner." 157 #define MSG_LS_MATZERO_FAILED "The SUNMatZero routine failed in an unrecoverable manner." 158 159 160 /*------------------------------------------------------------------ 161 Info messages 162 ------------------------------------------------------------------*/ 163 164 #define INFO_NLI "nli_inc = %d" 165 166 #if defined(SUNDIALS_EXTENDED_PRECISION) 167 168 #define INFO_EPS "residual norm = %12.3Lg eps = %12.3Lg" 169 170 #elif defined(SUNDIALS_DOUBLE_PRECISION) 171 172 #define INFO_EPS "residual norm = %12.3lg eps = %12.3lg" 173 174 #else 175 176 #define INFO_EPS "residual norm = %12.3g eps = %12.3g" 177 178 #endif 179 180 181 #ifdef __cplusplus 182 } 183 #endif 184 185 #endif 186