[16222] | 1 | /* ----------------------------------------------------------------- |
---|
| 2 | * Programmer(s): Daniel R. Reynolds @ SMU |
---|
| 3 | * Radu Serban @ LLNL |
---|
| 4 | * ----------------------------------------------------------------- |
---|
| 5 | * LLNS/SMU Copyright Start |
---|
| 6 | * Copyright (c) 2017, Southern Methodist University and |
---|
| 7 | * Lawrence Livermore National Security |
---|
| 8 | * |
---|
| 9 | * This work was performed under the auspices of the U.S. Department |
---|
| 10 | * of Energy by Southern Methodist University and Lawrence Livermore |
---|
| 11 | * National Laboratory under Contract DE-AC52-07NA27344. |
---|
| 12 | * Produced at Southern Methodist University and the Lawrence |
---|
| 13 | * Livermore National Laboratory. |
---|
| 14 | * |
---|
| 15 | * All rights reserved. |
---|
| 16 | * For details, see the LICENSE file. |
---|
| 17 | * LLNS/SMU Copyright End |
---|
| 18 | * ----------------------------------------------------------------- |
---|
| 19 | * Common header file for the direct linear solver interface in |
---|
| 20 | * CVODES. |
---|
| 21 | * |
---|
| 22 | * Part I contains type definitions and function prototypes for |
---|
| 23 | * using a CVDLS linear solver on forward problems (IVP |
---|
| 24 | * integration and/or FSA) |
---|
| 25 | * |
---|
| 26 | * Part II contains type definitions and function prototypes for |
---|
| 27 | * using a CVDLS linear solver on adjoint (backward) problems |
---|
| 28 | * -----------------------------------------------------------------*/ |
---|
| 29 | |
---|
| 30 | #ifndef _CVSDLS_H |
---|
| 31 | #define _CVSDLS_H |
---|
| 32 | |
---|
| 33 | #include <sundials/sundials_direct.h> |
---|
| 34 | #include <sundials/sundials_matrix.h> |
---|
| 35 | #include <sundials/sundials_linearsolver.h> |
---|
| 36 | #include <sundials/sundials_nvector.h> |
---|
| 37 | |
---|
| 38 | #ifdef __cplusplus /* wrapper to enable C++ usage */ |
---|
| 39 | extern "C" { |
---|
| 40 | #endif |
---|
| 41 | |
---|
| 42 | /*================================================================= |
---|
| 43 | CVSDLS Constants |
---|
| 44 | =================================================================*/ |
---|
| 45 | |
---|
| 46 | |
---|
| 47 | /*----------------------------------------------------------------- |
---|
| 48 | return values |
---|
| 49 | -----------------------------------------------------------------*/ |
---|
| 50 | |
---|
| 51 | #define CVDLS_SUCCESS 0 |
---|
| 52 | #define CVDLS_MEM_NULL -1 |
---|
| 53 | #define CVDLS_LMEM_NULL -2 |
---|
| 54 | #define CVDLS_ILL_INPUT -3 |
---|
| 55 | #define CVDLS_MEM_FAIL -4 |
---|
| 56 | |
---|
| 57 | /* Additional last_flag values */ |
---|
| 58 | |
---|
| 59 | #define CVDLS_JACFUNC_UNRECVR -5 |
---|
| 60 | #define CVDLS_JACFUNC_RECVR -6 |
---|
| 61 | #define CVDLS_SUNMAT_FAIL -7 |
---|
| 62 | |
---|
| 63 | /* Return values for the adjoint module */ |
---|
| 64 | |
---|
| 65 | #define CVDLS_NO_ADJ -101 |
---|
| 66 | #define CVDLS_LMEMB_NULL -102 |
---|
| 67 | |
---|
| 68 | |
---|
| 69 | /*================================================================= |
---|
| 70 | PART I: Forward Problems |
---|
| 71 | =================================================================*/ |
---|
| 72 | |
---|
| 73 | /*----------------------------------------------------------------- |
---|
| 74 | FUNCTION TYPES |
---|
| 75 | -----------------------------------------------------------------*/ |
---|
| 76 | |
---|
| 77 | /*----------------------------------------------------------------- |
---|
| 78 | Type: CVDlsJacFn |
---|
| 79 | ----------------------------------------------------------------- |
---|
| 80 | |
---|
| 81 | A Jacobian approximation function Jac must be of type CVDlsJacFn. |
---|
| 82 | Its parameters are: |
---|
| 83 | |
---|
| 84 | Jac is the SUNMatrix that will be loaded by a CVDlsJacFn with an |
---|
| 85 | approximation to the Jacobian matrix J = (df_i/dy_j) at |
---|
| 86 | the point (t,y). |
---|
| 87 | |
---|
| 88 | t is the current value of the independent variable. |
---|
| 89 | |
---|
| 90 | y is the current value of the dependent variable vector, |
---|
| 91 | namely the predicted value of y(t). |
---|
| 92 | |
---|
| 93 | fy is the vector f(t,y). |
---|
| 94 | |
---|
| 95 | user_data is a pointer to user data - the same as the user_data |
---|
| 96 | parameter passed to CVodeSetUserdata. |
---|
| 97 | |
---|
| 98 | tmp1, tmp2, and tmp3 are pointers to memory allocated for |
---|
| 99 | vectors of length N which can be used by a CVDlsJacFn |
---|
| 100 | as temporary storage or work space. |
---|
| 101 | |
---|
| 102 | A CVDlsJacFn should return 0 if successful, a positive |
---|
| 103 | value if a recoverable error occurred, and a negative value if |
---|
| 104 | an unrecoverable error occurred. |
---|
| 105 | |
---|
| 106 | ----------------------------------------------------------------- |
---|
| 107 | |
---|
| 108 | NOTE: See the relevant SUNMatrix implementation header files and |
---|
| 109 | documentation for mechanisms to inquire about matrix |
---|
| 110 | dimensions, and for efficient ways to set matrix entries. |
---|
| 111 | |
---|
| 112 | NOTE: If the user's Jacobian routine needs other quantities, |
---|
| 113 | they are accessible as follows: hcur (the current stepsize) |
---|
| 114 | and ewt (the error weight vector) are accessible through |
---|
| 115 | CVodeGetCurrentStep and CVodeGetErrWeights, respectively |
---|
| 116 | (see cvode.h). The unit roundoff is available as |
---|
| 117 | UNIT_ROUNDOFF defined in sundials_types.h. |
---|
| 118 | |
---|
| 119 | -----------------------------------------------------------------*/ |
---|
| 120 | typedef int (*CVDlsJacFn)(realtype t, N_Vector y, N_Vector fy, |
---|
| 121 | SUNMatrix Jac, void *user_data, |
---|
| 122 | N_Vector tmp1, N_Vector tmp2, N_Vector tmp3); |
---|
| 123 | |
---|
| 124 | |
---|
| 125 | /*----------------------------------------------------------------- |
---|
| 126 | EXPORTED FUNCTIONS |
---|
| 127 | -----------------------------------------------------------------*/ |
---|
| 128 | |
---|
| 129 | /*----------------------------------------------------------------- |
---|
| 130 | Required inputs to the CVSDLS linear solver interface |
---|
| 131 | ----------------------------------------------------------------- |
---|
| 132 | |
---|
| 133 | CVDlsSetLinearSolver specifies the direct SUNLinearSolver object |
---|
| 134 | that should be used. This is required if CVodes is solving a |
---|
| 135 | problem using the Newton nonlinear solver (not the function |
---|
| 136 | iteration). |
---|
| 137 | |
---|
| 138 | The return value is one of: |
---|
| 139 | CVDLS_SUCCESS if successful |
---|
| 140 | CVDLS_MEM_NULL if the CVODE memory was NULL |
---|
| 141 | CVDLS_ILL_INPUT if the arguments are incompatible |
---|
| 142 | -----------------------------------------------------------------*/ |
---|
| 143 | SUNDIALS_EXPORT int CVDlsSetLinearSolver(void *cvode_mem, |
---|
| 144 | SUNLinearSolver LS, |
---|
| 145 | SUNMatrix A); |
---|
| 146 | |
---|
| 147 | |
---|
| 148 | /*----------------------------------------------------------------- |
---|
| 149 | Optional inputs to the CVSDLS linear solver |
---|
| 150 | ----------------------------------------------------------------- |
---|
| 151 | |
---|
| 152 | CVDlsSetJacFn specifies the dense/band/sparse Jacobian |
---|
| 153 | approximation routine to be used for a direct linear solver. |
---|
| 154 | By default, a difference quotient approximation is used for |
---|
| 155 | dense/band; no default exists for sparse (so this must be |
---|
| 156 | user-supplied). |
---|
| 157 | |
---|
| 158 | The return value is one of: |
---|
| 159 | CVDLS_SUCCESS if successful |
---|
| 160 | CVDLS_MEM_NULL if the CVODE memory was NULL |
---|
| 161 | CVDLS_LMEM_NULL if the linear solver memory was NULL |
---|
| 162 | -----------------------------------------------------------------*/ |
---|
| 163 | SUNDIALS_EXPORT int CVDlsSetJacFn(void *cvode_mem, CVDlsJacFn jac); |
---|
| 164 | |
---|
| 165 | |
---|
| 166 | /*----------------------------------------------------------------- |
---|
| 167 | Optional outputs from the CVSSDLS linear solver interface |
---|
| 168 | ----------------------------------------------------------------- |
---|
| 169 | |
---|
| 170 | CVDlsGetWorkSpace returns the real and integer workspace used |
---|
| 171 | by the direct linear solver. |
---|
| 172 | CVDlsGetNumJacEvals returns the number of calls made to the |
---|
| 173 | Jacobian evaluation routine jac. |
---|
| 174 | CVDlsGetNumRhsEvals returns the number of calls to the user |
---|
| 175 | f routine due to finite difference Jacobian |
---|
| 176 | evaluation. |
---|
| 177 | CVDlsGetLastFlag returns the last error flag set by any of |
---|
| 178 | the CVSDIRECT interface functions. |
---|
| 179 | |
---|
| 180 | The return value of CVDlsGet* is one of: |
---|
| 181 | CVDLS_SUCCESS if successful |
---|
| 182 | CVDLS_MEM_NULL if the CVODES memory was NULL |
---|
| 183 | CVDLS_LMEM_NULL if the linear solver memory was NULL |
---|
| 184 | -----------------------------------------------------------------*/ |
---|
| 185 | |
---|
| 186 | SUNDIALS_EXPORT int CVDlsGetWorkSpace(void *cvode_mem, |
---|
| 187 | long int *lenrwLS, |
---|
| 188 | long int *leniwLS); |
---|
| 189 | SUNDIALS_EXPORT int CVDlsGetNumJacEvals(void *cvode_mem, |
---|
| 190 | long int *njevals); |
---|
| 191 | SUNDIALS_EXPORT int CVDlsGetNumRhsEvals(void *cvode_mem, |
---|
| 192 | long int *nfevalsLS); |
---|
| 193 | SUNDIALS_EXPORT int CVDlsGetLastFlag(void *cvode_mem, |
---|
| 194 | long int *flag); |
---|
| 195 | |
---|
| 196 | /*----------------------------------------------------------------- |
---|
| 197 | The following function returns the name of the constant |
---|
| 198 | associated with a CVSDLS return flag |
---|
| 199 | -----------------------------------------------------------------*/ |
---|
| 200 | SUNDIALS_EXPORT char *CVDlsGetReturnFlagName(long int flag); |
---|
| 201 | |
---|
| 202 | /*================================================================= |
---|
| 203 | PART II: Backward Problems |
---|
| 204 | =================================================================*/ |
---|
| 205 | |
---|
| 206 | /*----------------------------------------------------------------- |
---|
| 207 | FUNCTION TYPES |
---|
| 208 | -----------------------------------------------------------------*/ |
---|
| 209 | |
---|
| 210 | /*----------------------------------------------------------------- |
---|
| 211 | Type: CVDlsJacFnB |
---|
| 212 | ----------------------------------------------------------------- |
---|
| 213 | A Jacobian approximation function jacB for the adjoint |
---|
| 214 | (backward) problem must have the prototype given below. |
---|
| 215 | -----------------------------------------------------------------*/ |
---|
| 216 | typedef int (*CVDlsJacFnB)(realtype t, N_Vector y, N_Vector yB, |
---|
| 217 | N_Vector fyB, SUNMatrix JB, |
---|
| 218 | void *user_dataB, N_Vector tmp1B, |
---|
| 219 | N_Vector tmp2B, N_Vector tmp3B); |
---|
| 220 | |
---|
| 221 | /*----------------------------------------------------------------- |
---|
| 222 | Type: CVDlsJacFnBS |
---|
| 223 | ----------------------------------------------------------------- |
---|
| 224 | A Jacobian approximation function jacBS for the adjoint |
---|
| 225 | (backward) problem, sensitivity-dependent case, must have the |
---|
| 226 | prototype given below. |
---|
| 227 | -----------------------------------------------------------------*/ |
---|
| 228 | typedef int (*CVDlsJacFnBS)(realtype t, N_Vector y, N_Vector *yS, |
---|
| 229 | N_Vector yB, N_Vector fyB, SUNMatrix JB, |
---|
| 230 | void *user_dataB, N_Vector tmp1B, |
---|
| 231 | N_Vector tmp2B, N_Vector tmp3B); |
---|
| 232 | |
---|
| 233 | |
---|
| 234 | /*----------------------------------------------------------------- |
---|
| 235 | EXPORTED FUNCTIONS |
---|
| 236 | -----------------------------------------------------------------*/ |
---|
| 237 | |
---|
| 238 | /*----------------------------------------------------------------- |
---|
| 239 | Required inputs for the CVSDLS linear solver interface: |
---|
| 240 | |
---|
| 241 | CVDlsSetLinearSolverB specifies the direct SUNLinearSolver |
---|
| 242 | object that should be used for backward integration. The |
---|
| 243 | 'which' argument is the int returned by CVodeCreateB. |
---|
| 244 | |
---|
| 245 | The return value is one of: |
---|
| 246 | CVDLS_SUCCESS if successful |
---|
| 247 | CVDLS_MEM_NULL if the cvode memory was NULL |
---|
| 248 | CVDLS_ILL_INPUT if the arguments are incompatible |
---|
| 249 | ---------------------------------------------------------------*/ |
---|
| 250 | SUNDIALS_EXPORT int CVDlsSetLinearSolverB(void *cvode_mem, |
---|
| 251 | int which, |
---|
| 252 | SUNLinearSolver LS, |
---|
| 253 | SUNMatrix A); |
---|
| 254 | |
---|
| 255 | /*----------------------------------------------------------------- |
---|
| 256 | Functions: CVDlsSetJacFnB and CVDlsSetJacFnBS |
---|
| 257 | ----------------------------------------------------------------- |
---|
| 258 | CVDlsSetJacFnB specifies the Jacobian function to be used by a |
---|
| 259 | CVSDLS linear solver for the backward integration phase, when |
---|
| 260 | the backward problem does not depend on forward sensitivities. |
---|
| 261 | CVDlsSetJacFnBS specifies the Jacobian function when the backward |
---|
| 262 | problem does depend on sensitivities. |
---|
| 263 | |
---|
| 264 | The 'which' argument is the int returned by CVodeCreateB. |
---|
| 265 | |
---|
| 266 | The return value is one of: |
---|
| 267 | CVDLS_SUCCESS if successful |
---|
| 268 | CVDLS_MEM_NULL if the CVODE memory was NULL |
---|
| 269 | CVDLS_LMEM_NULL if the linear solver memory was NULL |
---|
| 270 | -----------------------------------------------------------------*/ |
---|
| 271 | SUNDIALS_EXPORT int CVDlsSetJacFnB(void *cvode_mem, int which, |
---|
| 272 | CVDlsJacFnB jacB); |
---|
| 273 | SUNDIALS_EXPORT int CVDlsSetJacFnBS(void *cvode_mem, int which, |
---|
| 274 | CVDlsJacFnBS jacBS); |
---|
| 275 | |
---|
| 276 | #ifdef __cplusplus |
---|
| 277 | } |
---|
| 278 | #endif |
---|
| 279 | |
---|
| 280 | #endif |
---|