1 | /* |
---|
2 | * ----------------------------------------------------------------- |
---|
3 | * Programmer(s): Daniel R. Reynolds @ SMU |
---|
4 | * Radu Serban @ LLNL |
---|
5 | * ----------------------------------------------------------------- |
---|
6 | * LLNS/SMU Copyright Start |
---|
7 | * Copyright (c) 2017, Southern Methodist University and |
---|
8 | * Lawrence Livermore National Security |
---|
9 | * |
---|
10 | * This work was performed under the auspices of the U.S. Department |
---|
11 | * of Energy by Southern Methodist University and Lawrence Livermore |
---|
12 | * National Laboratory under Contract DE-AC52-07NA27344. |
---|
13 | * Produced at Southern Methodist University and the Lawrence |
---|
14 | * Livermore National Laboratory. |
---|
15 | * |
---|
16 | * All rights reserved. |
---|
17 | * For details, see the LICENSE file. |
---|
18 | * LLNS/SMU Copyright End |
---|
19 | * ----------------------------------------------------------------- |
---|
20 | * This is the header file for the CVBBDPRE module, for a |
---|
21 | * band-block-diagonal preconditioner, i.e. a block-diagonal |
---|
22 | * matrix with banded blocks, for use with the CVSPILS interface, |
---|
23 | * and the MPI-parallel implementation of the NVECTOR module. |
---|
24 | * |
---|
25 | * Part I contains type definitions and function prototypes for using |
---|
26 | * CVBBDPRE on forward problems (IVP integration and/or FSA) |
---|
27 | * |
---|
28 | * Part II contains type definitions and function prototypes for using |
---|
29 | * CVBBDPRE on adjopint (backward) problems |
---|
30 | * ----------------------------------------------------------------- |
---|
31 | */ |
---|
32 | |
---|
33 | #ifndef _CVSBBDPRE_H |
---|
34 | #define _CVSBBDPRE_H |
---|
35 | |
---|
36 | #include <sundials/sundials_nvector.h> |
---|
37 | |
---|
38 | #ifdef __cplusplus /* wrapper to enable C++ usage */ |
---|
39 | extern "C" { |
---|
40 | #endif |
---|
41 | |
---|
42 | /*================================================================= |
---|
43 | PART I - forward problems |
---|
44 | =================================================================*/ |
---|
45 | |
---|
46 | /*----------------------------------------------------------------- |
---|
47 | SUMMARY |
---|
48 | |
---|
49 | These routines provide a preconditioner matrix that is |
---|
50 | block-diagonal with banded blocks. The blocking corresponds |
---|
51 | to the distribution of the dependent variable vector y among |
---|
52 | the processors. Each preconditioner block is generated from |
---|
53 | the Jacobian of the local part (on the current processor) of a |
---|
54 | given function g(t,y) approximating f(t,y). The blocks are |
---|
55 | generated by a difference quotient scheme on each processor |
---|
56 | independently. This scheme utilizes an assumed banded |
---|
57 | structure with given half-bandwidths, mudq and mldq. |
---|
58 | However, the banded Jacobian block kept by the scheme has |
---|
59 | half-bandwiths mukeep and mlkeep, which may be smaller. |
---|
60 | |
---|
61 | The user's calling program should have the following form: |
---|
62 | |
---|
63 | #include <cvodes/cvodes_spils.h> |
---|
64 | #include <cvodes/cvodes_bbdpre.h> |
---|
65 | #include <nvector_parallel.h> |
---|
66 | ... |
---|
67 | void *cvode_mem; |
---|
68 | ... |
---|
69 | Set y0 |
---|
70 | ... |
---|
71 | SUNLinearSolver LS = SUNSPBCGS(y0, pretype, maxl); |
---|
72 | -or- |
---|
73 | SUNLinearSolver LS = SUNSPFGMR(y0, pretype, maxl); |
---|
74 | -or- |
---|
75 | SUNLinearSolver LS = SUNSPGMR(y0, pretype, maxl); |
---|
76 | -or- |
---|
77 | SUNLinearSolver LS = SUNSPTFQMR(y0, pretype, maxl); |
---|
78 | -or- |
---|
79 | SUNLinearSolver LS = SUNPCG(y0, pretype, maxl); |
---|
80 | ... |
---|
81 | cvode_mem = CVodeCreate(...); |
---|
82 | flag = CVodeInit(...); |
---|
83 | ... |
---|
84 | flag = CVSpilsSetLinearSolver(cvode_mem, LS); |
---|
85 | ... |
---|
86 | flag = CVBBDPrecInit(cvode_mem, Nlocal, mudq ,mldq, |
---|
87 | mukeep, mlkeep, dqrely, gloc, cfn); |
---|
88 | ... |
---|
89 | flag = CVode(...); |
---|
90 | ... |
---|
91 | CVodeFree(&cvode_mem); |
---|
92 | ... |
---|
93 | Free y0 |
---|
94 | ... |
---|
95 | CVodeFree(&cvode_mem); |
---|
96 | ... |
---|
97 | SUNLinSolFree(LS); |
---|
98 | ... |
---|
99 | |
---|
100 | The user-supplied routines required are: |
---|
101 | |
---|
102 | f = function defining the ODE right-hand side f(t,y). |
---|
103 | |
---|
104 | gloc = function defining the approximation g(t,y). |
---|
105 | |
---|
106 | cfn = function to perform communication need for gloc. |
---|
107 | |
---|
108 | Notes: |
---|
109 | |
---|
110 | 1) This header file is included by the user for the definition |
---|
111 | of the CVBBDData type and for needed function prototypes. |
---|
112 | |
---|
113 | 2) The CVBBDPrecInit call includes half-bandwiths mudq and mldq |
---|
114 | to be used in the difference quotient calculation of the |
---|
115 | approximate Jacobian. They need not be the true |
---|
116 | half-bandwidths of the Jacobian of the local block of g, |
---|
117 | when smaller values may provide a greater efficiency. |
---|
118 | Also, the half-bandwidths mukeep and mlkeep of the retained |
---|
119 | banded approximate Jacobian block may be even smaller, |
---|
120 | to reduce storage and computation costs further. |
---|
121 | For all four half-bandwidths, the values need not be the |
---|
122 | same on every processor. |
---|
123 | |
---|
124 | 3) The actual name of the user's f function is passed to |
---|
125 | CVodeInit, and the names of the user's gloc and cfn |
---|
126 | functions are passed to CVBBDPrecInit. |
---|
127 | |
---|
128 | 4) The pointer to the user-defined data block user_data, which is |
---|
129 | set through CVodeSetUserData is also available to the user in |
---|
130 | gloc and cfn. |
---|
131 | |
---|
132 | 5) Optional outputs specific to this module are available by |
---|
133 | way of routines listed below. These include work space sizes |
---|
134 | and the cumulative number of gloc calls. The costs |
---|
135 | associated with this module also include nsetups banded LU |
---|
136 | factorizations, nlinsetups cfn calls, and npsolves banded |
---|
137 | backsolve calls, where nlinsetups and npsolves are |
---|
138 | integrator/CVSPILS optional outputs. |
---|
139 | -----------------------------------------------------------------*/ |
---|
140 | |
---|
141 | /*----------------------------------------------------------------- |
---|
142 | Type : CVLocalFn |
---|
143 | ----------------------------------------------------------------- |
---|
144 | The user must supply a function g(t,y) which approximates the |
---|
145 | right-hand side function f for the system y'=f(t,y), and which |
---|
146 | is computed locally (without interprocess communication). |
---|
147 | (The case where g is mathematically identical to f is allowed.) |
---|
148 | The implementation of this function must have type CVLocalFn. |
---|
149 | |
---|
150 | This function takes as input the local vector size Nlocal, the |
---|
151 | independent variable value t, the local real dependent |
---|
152 | variable vector y, and a pointer to the user-defined data |
---|
153 | block user_data. It is to compute the local part of g(t,y) and |
---|
154 | store this in the vector g. |
---|
155 | (Allocation of memory for y and g is handled within the |
---|
156 | preconditioner module.) |
---|
157 | The user_data parameter is the same as that specified by the user |
---|
158 | through the CVodeSetUserdata routine. |
---|
159 | |
---|
160 | A CVLocalFn should return 0 if successful, a positive value if |
---|
161 | a recoverable error occurred, and a negative value if an |
---|
162 | unrecoverable error occurred. |
---|
163 | -----------------------------------------------------------------*/ |
---|
164 | typedef int (*CVLocalFn)(sunindextype Nlocal, realtype t, |
---|
165 | N_Vector y, N_Vector g, void *user_data); |
---|
166 | |
---|
167 | |
---|
168 | /*----------------------------------------------------------------- |
---|
169 | Type : CVCommFn |
---|
170 | ----------------------------------------------------------------- |
---|
171 | The user may supply a function of type CVCommFn which performs |
---|
172 | all interprocess communication necessary to evaluate the |
---|
173 | approximate right-hand side function described above. |
---|
174 | |
---|
175 | This function takes as input the local vector size Nlocal, |
---|
176 | the independent variable value t, the dependent variable |
---|
177 | vector y, and a pointer to the user-defined data block |
---|
178 | user_data. The user_data parameter is the same as that |
---|
179 | specified by the user through the CVodeSetUserData routine. |
---|
180 | The CVCommFn cfn is expected to save communicated data in |
---|
181 | space defined within the structure user_data. |
---|
182 | Note: A CVCommFn cfn does not have a return value. |
---|
183 | |
---|
184 | Each call to the CVCommFn cfn is preceded by a call to the |
---|
185 | CVRhsFn f with the same (t,y) arguments. Thus cfn can omit any |
---|
186 | communications done by f if relevant to the evaluation of g. |
---|
187 | If all necessary communication was done by f, the user can |
---|
188 | pass NULL for cfn in CVBBDPrecInit (see below). |
---|
189 | |
---|
190 | A CVCommFn should return 0 if successful, a positive value if |
---|
191 | a recoverable error occurred, and a negative value if an |
---|
192 | unrecoverable error occurred. |
---|
193 | -----------------------------------------------------------------*/ |
---|
194 | typedef int (*CVCommFn)(sunindextype Nlocal, realtype t, |
---|
195 | N_Vector y, void *user_data); |
---|
196 | |
---|
197 | |
---|
198 | /*----------------------------------------------------------------- |
---|
199 | Function : CVBBDPrecInit |
---|
200 | ----------------------------------------------------------------- |
---|
201 | CVBBDPrecInit allocates and initializes the BBD preconditioner. |
---|
202 | |
---|
203 | The parameters of CVBBDPrecInit are as follows: |
---|
204 | |
---|
205 | cvode_mem is the pointer to the integrator memory. |
---|
206 | |
---|
207 | Nlocal is the length of the local block of the vectors y etc. |
---|
208 | on the current processor. |
---|
209 | |
---|
210 | mudq, mldq are the upper and lower half-bandwidths to be used |
---|
211 | in the difference quotient computation of the local |
---|
212 | Jacobian block. |
---|
213 | |
---|
214 | mukeep, mlkeep are the upper and lower half-bandwidths of the |
---|
215 | retained banded approximation to the local Jacobian |
---|
216 | block. |
---|
217 | |
---|
218 | dqrely is an optional input. It is the relative increment |
---|
219 | in components of y used in the difference quotient |
---|
220 | approximations. To specify the default, pass 0. |
---|
221 | The default is dqrely = sqrt(unit roundoff). |
---|
222 | |
---|
223 | gloc is the name of the user-supplied function g(t,y) that |
---|
224 | approximates f and whose local Jacobian blocks are |
---|
225 | to form the preconditioner. |
---|
226 | |
---|
227 | cfn is the name of the user-defined function that performs |
---|
228 | necessary interprocess communication for the |
---|
229 | execution of gloc. |
---|
230 | |
---|
231 | The return value of CVBBDPrecInit is one of: |
---|
232 | CVSPILS_SUCCESS if no errors occurred |
---|
233 | CVSPILS_MEM_NULL if the integrator memory is NULL |
---|
234 | CVSPILS_LMEM_NULL if the linear solver memory is NULL |
---|
235 | CVSPILS_ILL_INPUT if an input has an illegal value |
---|
236 | CVSPILS_MEM_FAIL if a memory allocation request failed |
---|
237 | -----------------------------------------------------------------*/ |
---|
238 | SUNDIALS_EXPORT int CVBBDPrecInit(void *cvode_mem, |
---|
239 | sunindextype Nlocal, |
---|
240 | sunindextype mudq, |
---|
241 | sunindextype mldq, |
---|
242 | sunindextype mukeep, |
---|
243 | sunindextype mlkeep, |
---|
244 | realtype dqrely, |
---|
245 | CVLocalFn gloc, |
---|
246 | CVCommFn cfn); |
---|
247 | |
---|
248 | |
---|
249 | /*----------------------------------------------------------------- |
---|
250 | Function : CVBBDPrecReInit |
---|
251 | ----------------------------------------------------------------- |
---|
252 | CVBBDPrecReInit re-initializes the BBDPRE module when solving a |
---|
253 | sequence of problems of the same size with CVSPILS/CVBBDPRE |
---|
254 | provided there is no change in Nlocal, mukeep, or mlkeep. After |
---|
255 | solving one problem, and after calling CVodeReInit to |
---|
256 | re-initialize the integrator for a subsequent problem, call |
---|
257 | CVBBDPrecReInit. |
---|
258 | |
---|
259 | All arguments have the same names and meanings as those |
---|
260 | of CVBBDPrecInit. |
---|
261 | |
---|
262 | The return value of CVBBDPrecReInit is one of: |
---|
263 | CVSPILS_SUCCESS if no errors occurred |
---|
264 | CVSPILS_MEM_NULL if the integrator memory is NULL |
---|
265 | CVSPILS_LMEM_NULL if the linear solver memory is NULL |
---|
266 | CVSPILS_PMEM_NULL if the preconditioner memory is NULL |
---|
267 | -----------------------------------------------------------------*/ |
---|
268 | SUNDIALS_EXPORT int CVBBDPrecReInit(void *cvode_mem, |
---|
269 | sunindextype mudq, |
---|
270 | sunindextype mldq, |
---|
271 | realtype dqrely); |
---|
272 | |
---|
273 | |
---|
274 | /*----------------------------------------------------------------- |
---|
275 | BBDPRE optional output extraction routines |
---|
276 | ----------------------------------------------------------------- |
---|
277 | CVBBDPrecGetWorkSpace returns the BBDPRE real and integer work |
---|
278 | space sizes. |
---|
279 | CVBBDPrecGetNumGfnEvals returns the number of calls to gfn. |
---|
280 | |
---|
281 | The return value of CVBBDPrecGet* is one of: |
---|
282 | CVSPILS_SUCCESS if no errors occurred |
---|
283 | CVSPILS_MEM_NULL if the integrator memory is NULL |
---|
284 | CVSPILS_LMEM_NULL if the linear solver memory is NULL |
---|
285 | CVSPILS_PMEM_NULL if the preconditioner memory is NULL |
---|
286 | -----------------------------------------------------------------*/ |
---|
287 | |
---|
288 | SUNDIALS_EXPORT int CVBBDPrecGetWorkSpace(void *cvode_mem, |
---|
289 | long int *lenrwBBDP, |
---|
290 | long int *leniwBBDP); |
---|
291 | SUNDIALS_EXPORT int CVBBDPrecGetNumGfnEvals(void *cvode_mem, |
---|
292 | long int *ngevalsBBDP); |
---|
293 | |
---|
294 | |
---|
295 | /*================================================================= |
---|
296 | PART II - backward problems |
---|
297 | =================================================================*/ |
---|
298 | |
---|
299 | /*----------------------------------------------------------------- |
---|
300 | Types: CVLocalFnB and CVCommFnB |
---|
301 | ----------------------------------------------------------------- |
---|
302 | Local approximation function and inter-process communication |
---|
303 | function for the BBD preconditioner on the backward phase. |
---|
304 | -----------------------------------------------------------------*/ |
---|
305 | typedef int (*CVLocalFnB)(sunindextype NlocalB, realtype t, |
---|
306 | N_Vector y, N_Vector yB, N_Vector gB, |
---|
307 | void *user_dataB); |
---|
308 | typedef int (*CVCommFnB)(sunindextype NlocalB, realtype t, |
---|
309 | N_Vector y, N_Vector yB, void *user_dataB); |
---|
310 | |
---|
311 | |
---|
312 | /*----------------------------------------------------------------- |
---|
313 | Functions: CVBBDPrecInitB, CVBBDPrecReInit |
---|
314 | ----------------------------------------------------------------- |
---|
315 | Interface functions for the CVBBDPRE preconditioner to be used on |
---|
316 | the backward phase. |
---|
317 | The 'which' argument is the int returned by CVodeCreateB. |
---|
318 | -----------------------------------------------------------------*/ |
---|
319 | SUNDIALS_EXPORT int CVBBDPrecInitB(void *cvode_mem, int which, |
---|
320 | sunindextype NlocalB, |
---|
321 | sunindextype mudqB, |
---|
322 | sunindextype mldqB, |
---|
323 | sunindextype mukeepB, |
---|
324 | sunindextype mlkeepB, |
---|
325 | realtype dqrelyB, |
---|
326 | CVLocalFnB glocB, |
---|
327 | CVCommFnB cfnB); |
---|
328 | SUNDIALS_EXPORT int CVBBDPrecReInitB(void *cvode_mem, int which, |
---|
329 | sunindextype mudqB, |
---|
330 | sunindextype mldqB, |
---|
331 | realtype dqrelyB); |
---|
332 | |
---|
333 | #ifdef __cplusplus |
---|
334 | } |
---|
335 | #endif |
---|
336 | |
---|
337 | #endif |
---|