/******************************************************************* * File : cvbbdpre.h * * Programmers : Michael Wittman, Alan C. Hindmarsh, and * * Radu Serban @ LLNL * * Version of : 26 June 2002 * *-----------------------------------------------------------------* * Copyright (c) 2002, The Regents of the University of California * * Produced at the Lawrence Livermore National Laboratory * * All rights reserved * * For details, see sundials/cvode/LICENSE * *-----------------------------------------------------------------* * This is the header file for the CVBBDPRE module, for a * * band-block-diagonal preconditioner, i.e. a block-diagonal * * matrix with banded blocks, for use with CVODE, CVSpgmr, and * * the parallel implementation of the NVECTOR module. * * * * Summary: * * * * These routines provide a preconditioner matrix for CVODE that * * is block-diagonal with banded blocks. The blocking corresponds * * to the distribution of the dependent variable vector y among * * the processors. Each preconditioner block is generated from * * the Jacobian of the local part (on the current processor) of a * * given function g(t,y) approximating f(t,y). The blocks are * * generated by a difference quotient scheme on each processor * * independently. This scheme utilizes an assumed banded * * structure with given half-bandwidths, mudq and mldq. * * However, the banded Jacobian block kept by the scheme has * * half-bandwiths mukeep and mlkeep, which may be smaller. * * * * The user's calling program should have the following form: * * * * #include "nvector_parallel.h" * * #include "cvbbdpre.h" * * ... * * CVBBDData p_data; * * ... * * machEnv = M_EnvInit_Parallel(...); * * ... * * cvode_mem = CVodeMalloc(...); * * ... * * p_data = CVBBDAlloc(Nlocal, mudq ,mldq, mukeep, mlkeep, * * dqrely, gloc,cfn, f_data, cvode_mem); * * ... * * flag = CVSpgmr(cvode_mem, pretype, gstype, maxl, delt, * * CVBBDPrecon, CVBBDPSol, p_data); * * ... * * ier = CVode(...); * * ... * * CVBBDFree(p_data); * * ... * * CVodeFree(...); * * * * M_EnvFree_Parallel(machEnv); * * * * * * The user-supplied routines required are: * * * * f = function defining the ODE right-hand side f(t,y). * * * * gloc = function defining the approximation g(t,y). * * * * cfn = function to perform communication need for gloc. * * * * * * Notes: * * * * 1) This header file is included by the user for the definition * * of the CVBBDData type and for needed function prototypes. * * * * 2) The CVBBDAlloc call includes half-bandwiths mudq and mldq * * to be used in the difference-quotient calculation of the * * approximate Jacobian. They need not be the true * * half-bandwidths of the Jacobian of the local block of g, * * when smaller values may provide a greater efficiency. * * Also, the half-bandwidths mukeep and mlkeep of the retained * * banded approximate Jacobian block may be even smaller, * * to reduce storage and computation costs further. * * For all four half-bandwidths, the values need not be the * * same on every processor. * * * * 3) The actual name of the user's f function is passed to * * CVodeMalloc, and the names of the user's gloc and cfn * * functions are passed to CVBBDAlloc. * * * * 4) The pointer to the user-defined data block f_data, which is * * passed to CVodeMalloc, is also passed to CVBBDAlloc, and * * is available to the user in gloc and cfn. * * * * 5) For the CVSpgmr call, the preconditioning and Gram-Schmidt * * types, pretype and gstype, are left to the user to specify. * * * * 6) Functions CVBBDPrecon and CVBBDPSol are never called by the * * user explicitly; their names are simply passed to CVSpgmr. * * * * 7) Optional outputs specific to this module are available by * * way of macros listed below. These include work space sizes * * and the cumulative number of gloc calls. The costs * * associated with this module also include nsetups banded LU * * factorizations, nsetups cfn calls, and nps banded * * backsolve calls, where nsetups and nps are CVODE optional * * outputs. * *******************************************************************/ #ifdef __cplusplus /* wrapper to enable C++ usage */ extern "C" { #endif #ifndef _cvbbdpre_h #define _cvbbdpre_h #include "cvode.h" #include "sundialstypes.h" #include "nvector.h" #include "band.h" /****************************************************************** * Type : CVLocalFn * *----------------------------------------------------------------* * The user must supply a function g(t,y) which approximates the * * right-hand side function f for the system y'=f(t,y), and which * * is computed locally (without inter-processor communication). * * (The case where g is mathematically identical to f is allowed.)* * The implementation of this function must have type CVLocalFn. * * * * This function takes as input the local vector size Nlocal, the * * independent variable value t, the local realtype dependent * * variable array ylocal, and a pointer to the user-defined data * * block f_data. It is to compute the local part of g(t,y) and * * store this in the realtype array glocal. * * (Allocation of memory for ylocal and glocal is handled within * * the preconditioner module.) * * The f_data parameter is the same as that passed by the user * * to the CVodeMalloc routine. * * A CVLocalFn gloc does not have a return value. * ******************************************************************/ typedef void (*CVLocalFn)(integertype Nlocal, realtype t, realtype *ylocal, realtype *glocal, void *f_data); /****************************************************************** * Type : CVCommFn * *----------------------------------------------------------------* * The user must supply a function of type CVCommFn which performs* * all inter-processor communication necessary to evaluate the * * approximate right-hand side function described above. * * * * This function takes as input the local vector size Nlocal, * * the independent variable value t, the dependent variable * * vector y, and a pointer to the user-defined data block f_data. * * The f_data parameter is the same as that passed by the user to * * the CVodeMalloc routine. The CVCommFn cfn is expected to save * * communicated data in space defined with the structure *f_data. * * A CVCommFn cfn does not have a return value. * * * * Each call to the CVCommFn cfn is preceded by a call to the * * RhsFn f with the same (t,y) arguments. Thus cfn can omit any * * communications done by f if relevant to the evaluation of g. * ******************************************************************/ typedef void (*CVCommFn)(integertype Nlocal, realtype t, N_Vector y, void *f_data); /*********************** Definition of CVBBDData *****************/ typedef struct { /* passed by user to CVBBDAlloc, used by Precond/Psolve functions: */ void *f_data; integertype mudq, mldq, mukeep, mlkeep; realtype dqrely; CVLocalFn gloc; CVCommFn cfn; /* set by CVBBDPrecon and used by CVBBDPSol: */ BandMat savedJ; BandMat savedP; integertype *pivots; /* set by CVBBDAlloc and used by CVBBDPrecon */ integertype n_local; /* available for optional output: */ integertype rpwsize; integertype ipwsize; integertype nge; } *CVBBDData; /*************** Macros for optional outputs ********************** * * * CVBBD_RPWSIZE(pdata) returns the size of the real work space, * * in realtype words, used by this preconditioner module. * * This size is local to the current processor. * * * * CVBBD_IPWSIZE(pdata) returns the size of the integer work * * space, in integertype words, used by this preconditioner * * module. This size is local to the current processor. * * * * CVBBD_NGE(pdata) returns the number of g(t,y) evaluations, * * i.e. the number of calls to the gloc function, so far. * ******************************************************************/ #define CVBBD_RPWSIZE(pdata) (pdata->rpwsize) #define CVBBD_IPWSIZE(pdata) (pdata->ipwsize) #define CVBBD_NGE(pdata) (pdata->nge) /****************************************************************** * Function : CVBBDAlloc * *----------------------------------------------------------------* * CVBBDAlloc allocates and initializes a CVBBDData structure * * to be passed to CVSpgmr (and subsequently used by CVBBDPrecon * * and CVBBDPSol. * * * * The parameters of CVBBDAlloc are as follows: * * * * Nlocal is the length of the local block of the vectors y etc. * * on the current processor. * * * * mudq, mldq are the upper and lower half-bandwidths to be used * * in the difference-quotient computation of the local * * Jacobian block. * * * * mukeep, mlkeep are the upper and lower half-bandwidths of the * * retained banded approximation to the local Jacobian * * block. * * * * dqrely is an optional input. It is the relative increment * * in components of y used in the difference quotient * * approximations. To specify the default, pass 0. * * The default is dqrely = sqrt(unit roundoff). * * * * gloc is the name of the user-supplied function g(t,y) that * * approximates f and whose local Jacobian blocks are * * to form the preconditioner. * * * * cfn is the name of the user-defined function that performs * * necessary inter-processor communication for the * * execution of gloc. * * * * f_data is a pointer to the optional user data block. * * * * CVBBDAlloc returns the storage allocated (type CVBBDData), * * or NULL if the request for storage cannot be satisfied. * ******************************************************************/ CVBBDData CVBBDAlloc(integertype Nlocal, integertype mudq, integertype mldq, integertype mukeep, integertype mlkeep, realtype dqrely, CVLocalFn gloc, CVCommFn cfn, void *f_data, void *cvode_mem ); /****************************************************************** * Function : CVReInitBBD * *----------------------------------------------------------------* * CVReInitBBD re-initializes the BBDPRE module when solving a * * sequence of problems of the same size with CVSPGMR/CVBBDPRE, * * provided there is no change in Nlocal, mukeep, or mlkeep. * * After solving one problem, and after calling CVReInit to * * re-initialize CVODE for a subsequent problem, call CVReInitBBD.* * Then call CVReInitSpgmr or CVSpgmr if necessary, depending on * * changes made in the CVSpgmr parameters, before calling CVode. * * * * The first argument to CVReInitBBD must be the pointer pdata * * that was returned by CVBBDAlloc. All other arguments have * * the same names and meanings as those of CVBBDAlloc. * * * * The return value of CVReInitBBD is 0, indicating success. * ******************************************************************/ int CVReInitBBD(CVBBDData pdata, integertype Nlocal, integertype mudq, integertype mldq, integertype mukeep, integertype mlkeep, realtype dqrely, CVLocalFn gloc, CVCommFn cfn, void *f_data); /****************************************************************** * Function : CVBBDFree * *----------------------------------------------------------------* * CVBBDFree frees the memory block p_data allocated by the call * * to CVBBDAlloc. * ******************************************************************/ void CVBBDFree(CVBBDData pdata); /****** Prototypes of functions CVBBDPrecon and CVBBDPSol *********/ int CVBBDPrecon(integertype N, realtype t, N_Vector y, N_Vector fy, booleantype jok, booleantype *jcurPtr, realtype gamma, N_Vector ewt, realtype h, realtype uround, long int *nfePtr, void *P_data, N_Vector vtemp1, N_Vector vtemp2, N_Vector vtemp3); int CVBBDPSol(integertype N, realtype t, N_Vector y,N_Vector fy, N_Vector vtemp, realtype gamma, N_Vector ewt, realtype delta, long int *nfePtr, N_Vector r, int lr, void *P_data, N_Vector z); #endif #ifdef __cplusplus } #endif