18 #ifndef __DOLFIN_C_VODE_H
19 #define __DOLFIN_C_VODE_H
23 #include <dolfin/la/SUNDIALSNVector.h>
25 #include <cvode/cvode.h>
26 #include <cvode/cvode_impl.h>
27 #include <sunlinsol/sunlinsol_spgmr.h>
28 #include <sundials/sundials_dense.h>
29 #include <sundials/sundials_types.h>
30 #include <sundials/sundials_iterative.h>
40 enum LMM { cv_bdf = CV_BDF, cv_adams = CV_ADAMS };
42 enum ITER { cv_functional = CV_FUNCTIONAL, cv_newton = CV_NEWTON };
49 CVode(LMM cv_lmm, ITER cv_iter);
63 void init(std::shared_ptr<GenericVector> u0,
double atol,
double rtol,
long int mxsteps = 0);
70 double step(
double dt);
90 virtual void derivs(
double t, std::shared_ptr<GenericVector> u,
91 std::shared_ptr<GenericVector> udot);
106 virtual int jacobian(std::shared_ptr<const GenericVector> v,
107 std::shared_ptr<GenericVector> Jv,
108 double t, std::shared_ptr<const GenericVector> y,
109 std::shared_ptr<const GenericVector> fy);
123 std::shared_ptr<GenericVector> Jv,
124 std::shared_ptr<GenericVector> y);
145 virtual int psolve(
double tn, std::shared_ptr<const GenericVector>y,
146 std::shared_ptr<const GenericVector> fy,
147 std::shared_ptr<const GenericVector> r,
148 std::shared_ptr<GenericVector> z,
149 double gamma,
double delta,
int lr);
160 static int f(realtype t, N_Vector u, N_Vector udot,
void *user_data);
164 static int f_jac_setup(
double t, N_Vector y, N_Vector fy,
void *user_data);
168 static int f_jac(N_Vector u, N_Vector fu,
double t, N_Vector y, N_Vector fy,
void* , N_Vector tmp);
172 static int prec_solve(
double, N_Vector, N_Vector, N_Vector, N_Vector,
double,
double,
int,
void*);
175 std::shared_ptr<SUNDIALSNVector> _u;