Chaste Commit::baa90ac2819b962188b7562f2326be23c47859a7
CvodeAdaptor.cpp
1/*
2
3Copyright (c) 2005-2024, University of Oxford.
4All rights reserved.
5
6University of Oxford means the Chancellor, Masters and Scholars of the
7University of Oxford, having an administrative office at Wellington
8Square, Oxford OX1 2JD, UK.
9
10This file is part of Chaste.
11
12Redistribution and use in source and binary forms, with or without
13modification, are permitted provided that the following conditions are met:
14 * Redistributions of source code must retain the above copyright notice,
15 this list of conditions and the following disclaimer.
16 * Redistributions in binary form must reproduce the above copyright notice,
17 this list of conditions and the following disclaimer in the documentation
18 and/or other materials provided with the distribution.
19 * Neither the name of the University of Oxford nor the names of its
20 contributors may be used to endorse or promote products derived from this
21 software without specific prior written permission.
22
23THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
24AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
25IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
26ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
27LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
28CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE
29GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
30HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
31LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
32OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
33
34*/
35
36#ifdef CHASTE_CVODE
37
38#include "CvodeAdaptor.hpp"
39
40#include "Exception.hpp"
41#include "MathsCustomFunctions.hpp" // For tolerance check.
42#include "TimeStepper.hpp"
44
45#include <iostream>
46#include <sstream>
47
48// CVODE headers
49#include <sundials/sundials_nvector.h>
50
51#if CHASTE_SUNDIALS_VERSION >= 30000
52#include <cvode/cvode_direct.h> /* access to CVDls interface */
53#include <sundials/sundials_types.h> /* defs. of realtype, sunindextype */
54#include <sunlinsol/sunlinsol_dense.h> /* access to dense SUNLinearSolver */
55#include <sunmatrix/sunmatrix_dense.h> /* access to dense SUNMatrix */
56#else
57#include <cvode/cvode_dense.h>
58#endif
59
60#if CHASTE_SUNDIALS_VERSION >= 60000
61#include "CvodeContextManager.hpp" // access to shared SUNContext object required by Sundials 6.0+
62#endif
63
80int CvodeRhsAdaptor(realtype t, N_Vector y, N_Vector ydot, void* pData)
81{
82 assert(pData != nullptr);
83 CvodeData* p_data = (CvodeData*)pData;
84 // Get y, ydot into std::vector<>s
85 static std::vector<realtype> ydot_vec;
86 CopyToStdVector(y, *p_data->pY);
87 CopyToStdVector(ydot, ydot_vec);
88 // Call our function
89 try
90 {
91 p_data->pSystem->EvaluateYDerivatives(t, *(p_data->pY), ydot_vec);
92 }
93 catch (const Exception& e)
94 {
95 std::cerr << "CVODE RHS Exception: " << e.GetMessage() << std::endl
96 << std::flush;
97 return -1;
98 }
99 // Copy derivative back
100 CopyFromStdVector(ydot_vec, ydot);
101 return 0;
102}
103
124int CvodeRootAdaptor(realtype t, N_Vector y, realtype* pGOut, void* pData)
125{
126 assert(pData != nullptr);
127 CvodeData* p_data = (CvodeData*)pData;
128 // Get y into a std::vector
129 CopyToStdVector(y, *p_data->pY);
130 // Call our function
131 try
132 {
133 *pGOut = p_data->pSystem->CalculateRootFunction(t, *p_data->pY);
134 }
135 catch (const Exception& e)
136 {
137 std::cerr << "CVODE Root Exception: " << e.GetMessage() << std::endl
138 << std::flush;
139 return -1;
140 }
141 return 0;
142}
143
144// /**
145// * Jacobian computation adaptor function.
146// *
147// * If solving an AbstractOdeSystemWithAnalyticJacobian, this function
148// * can be used to allow CVODE to compute the Jacobian analytically.
149// *
150// * Note to self: can test using pSystem->GetUseAnalyticJacobian().
151// */
152// int CvodeDenseJacobianAdaptor(long int numberOfStateVariables, DenseMat J,
153// realtype t, N_Vector y, N_Vector fy,
154// void* pData,
155// N_Vector tmp1, N_Vector tmp2, N_Vector tmp3)
156// {
157// AbstractOdeSystemWithAnalyticJacobian* pSystem
158// = (AbstractOdeSystemWithAnalyticJacobian*) pData;
159// // Get the current time step
160// double dt;
161// CVodeGetCurrentStep(CvodeMem, &dt);
162// // Get std::vector<> for y and double** for J
163// std::vector<realtype>& y_vec = *NV_DATA_STL(y);
164// double** ppJ = J->data; // organised column-wise: J_{i,j} = ppJ[j][i]
165// // Call our function
166// try
167// {
168// pSystem->AnalyticJacobian(y_vec, ppJ, t, dt);
169// }
170// catch (const Exception &e)
171// {
172// std::cerr << "CVODE J Exception: " << e.GetMessage() << std::endl << std::flush;
173// return -1;
174// }
175// // Update J (if needed)
176// return 0;
177// }
178
179void CvodeErrorHandler(int errorCode, const char* module, const char* function,
180 char* message, void* pData)
181{
182 std::stringstream err;
183 err << "CVODE Error " << errorCode << " in module " << module
184 << " function " << function << ": " << message;
185 std::cerr << "*" << err.str() << std::endl
186 << std::flush;
187 // Throwing the exception here causes termination on Maverick (g++ 4.4)
188 //EXCEPTION(err.str());
189}
190
192 std::vector<double>& rInitialY,
193 double startTime,
194 double maxStep)
195{
196 assert(rInitialY.size() == pOdeSystem->GetNumberOfStateVariables());
197 assert(maxStep > 0.0);
198
200#if CHASTE_SUNDIALS_VERSION >= 60000
201 N_Vector initial_values = N_VMake_Serial(rInitialY.size(), &(rInitialY[0]), CvodeContextManager::Instance()->GetSundialsContext());
202#else
203 N_Vector initial_values = N_VMake_Serial(rInitialY.size(), &(rInitialY[0]));
204#endif
205 assert(NV_DATA_S(initial_values) == &(rInitialY[0]));
206 assert(!NV_OWN_DATA_S(initial_values));
207 // std::cout << " Initial values: "; N_VPrint_Serial(initial_values);
208 // std::cout << " Rtol: " << mRelTol << ", Atol: " << mAbsTol << std::endl;
209 // std::cout << " Start: " << startTime << " max dt=" << maxStep << std::endl << std::flush;
210
211 // Find out if we need to (re-)initialise
213 if (!reinit && !mForceMinimalReset)
214 {
215 const unsigned size = GetVectorSize(rInitialY);
216 for (unsigned i = 0; i < size; i++)
217 {
219 {
220 reinit = true;
221 break;
222 }
223 }
224 }
225
226 if (!mpCvodeMem) // First run of this solver, set up CVODE memory
227 {
228 // Set up CVODE's memory.
229#if CHASTE_SUNDIALS_VERSION >= 60000
230 mpCvodeMem = CVodeCreate(CV_BDF, CvodeContextManager::Instance()->GetSundialsContext());
231#elif CHASTE_SUNDIALS_VERSION >= 40000
232 // v4.0.0 release notes: instead of specifying the nonlinear iteration type when creating the CVODE(S) memory structure,
233 // CVODE(S) uses the SUNNONLINSOL_NEWTON module implementation of a Newton iteration by default.
234 mpCvodeMem = CVodeCreate(CV_BDF);
235#else
236 mpCvodeMem = CVodeCreate(CV_BDF, CV_NEWTON);
237#endif
238
239 if (mpCvodeMem == nullptr)
240 EXCEPTION("Failed to SetupCvode CVODE"); // LCOV_EXCL_LINE
241 // Set error handler
242 CVodeSetErrHandlerFn(mpCvodeMem, CvodeErrorHandler, nullptr);
243 // Set the user data
244 mData.pSystem = pOdeSystem;
245 mData.pY = &rInitialY;
246#if CHASTE_SUNDIALS_VERSION >= 20400
247 CVodeSetUserData(mpCvodeMem, (void*)(&mData));
248#else
249 CVodeSetFdata(mpCvodeMem, (void*)(&mData));
250#endif
251
252// Setup CVODE
253#if CHASTE_SUNDIALS_VERSION >= 20400
254 CVodeInit(mpCvodeMem, CvodeRhsAdaptor, startTime, initial_values);
255 CVodeSStolerances(mpCvodeMem, mRelTol, mAbsTol);
256#else
257 CVodeMalloc(mpCvodeMem, CvodeRhsAdaptor, startTime, initial_values,
258 CV_SS, mRelTol, &mAbsTol);
259#endif
260
261 // Set the rootfinder function if wanted
262 if (mCheckForRoots)
263 {
264#if CHASTE_SUNDIALS_VERSION >= 20400
265 CVodeRootInit(mpCvodeMem, 1, CvodeRootAdaptor);
266#else
267 CVodeRootInit(mpCvodeMem, 1, CvodeRootAdaptor, (void*)(&mData));
268#endif
269 }
270
271 /* Create dense SUNMatrix for use in linear solves */
272#if CHASTE_SUNDIALS_VERSION >= 60000
273 mpSundialsDenseMatrix = SUNDenseMatrix(rInitialY.size(), rInitialY.size(), CvodeContextManager::Instance()->GetSundialsContext());
274#elif CHASTE_SUNDIALS_VERSION >= 30000
275 mpSundialsDenseMatrix = SUNDenseMatrix(rInitialY.size(), rInitialY.size());
276#endif
277
278#if CHASTE_SUNDIALS_VERSION >= 60000
279 /* Create dense SUNLinearSolver object for use by CVode */
280 mpSundialsLinearSolver = SUNLinSol_Dense(initial_values, mpSundialsDenseMatrix, CvodeContextManager::Instance()->GetSundialsContext());
281
282 /* Call CVodeSetLinearSolver to attach the matrix and linear solver to CVode */
283 CVodeSetLinearSolver(mpCvodeMem, mpSundialsLinearSolver, mpSundialsDenseMatrix);
284#elif CHASTE_SUNDIALS_VERSION >= 40000
285 /* Create dense SUNLinearSolver object for use by CVode */
286 mpSundialsLinearSolver = SUNLinSol_Dense(initial_values, mpSundialsDenseMatrix);
287
288 /* Call CVodeSetLinearSolver to attach the matrix and linear solver to CVode */
289 CVodeSetLinearSolver(mpCvodeMem, mpSundialsLinearSolver, mpSundialsDenseMatrix);
290#elif CHASTE_SUNDIALS_VERSION >= 30000
291 /* Create dense SUNLinearSolver object for use by CVode */
292 mpSundialsLinearSolver = SUNDenseLinearSolver(initial_values, mpSundialsDenseMatrix);
293
294 /* Call CVDlsSetLinearSolver to attach the matrix and linear solver to CVode */
295 CVDlsSetLinearSolver(mpCvodeMem, mpSundialsLinearSolver, mpSundialsDenseMatrix);
296#else
297 // CVODE < v3.0.0
298 // Attach a linear solver for Newton iteration
299 CVDense(mpCvodeMem, rInitialY.size());
300#endif
301 }
302 else if (reinit) // Could be new ODE system, or new Y values
303 {
304 // Set the user data
305 mData.pSystem = pOdeSystem; // stays the same on a re-initialize
306 mData.pY = &rInitialY; // changes on a re-initialize
307#if CHASTE_SUNDIALS_VERSION >= 20400
308 CVodeSetUserData(mpCvodeMem, (void*)(&mData));
309#else
310 CVodeSetFdata(mpCvodeMem, (void*)(&mData));
311#endif
312
313#if CHASTE_SUNDIALS_VERSION >= 20400
314 CVodeReInit(mpCvodeMem, startTime, initial_values);
315 CVodeSStolerances(mpCvodeMem, mRelTol, mAbsTol);
316#else
317 CVodeReInit(mpCvodeMem, CvodeRhsAdaptor, startTime, initial_values,
318 CV_SS, mRelTol, &mAbsTol);
319#endif
320
321#if CHASTE_SUNDIALS_VERSION >= 30000
322 if (mpSundialsLinearSolver)
323 {
324 /* Free the linear solver memory */
325 SUNLinSolFree(mpSundialsLinearSolver);
326 }
327 if (mpSundialsDenseMatrix)
328 {
329 /* Free the matrix memory */
330 SUNMatDestroy(mpSundialsDenseMatrix);
331 }
332
333 /* Create dense matrix of type SUNDenseMatrix for use in linear solves */
334#if CHASTE_SUNDIALS_VERSION >= 60000
335 mpSundialsDenseMatrix = SUNDenseMatrix(rInitialY.size(), rInitialY.size(), CvodeContextManager::Instance()->GetSundialsContext());
336#else
337 mpSundialsDenseMatrix = SUNDenseMatrix(rInitialY.size(), rInitialY.size());
338#endif
339
340#if CHASTE_SUNDIALS_VERSION >= 60000
341 /* Create dense SUNLinSol_Dense object for use by CVode */
342 mpSundialsLinearSolver = SUNLinSol_Dense(initial_values, mpSundialsDenseMatrix, CvodeContextManager::Instance()->GetSundialsContext());
343
344 /* Call CVodeSetLinearSolver to attach the matrix and linear solver to CVode */
345 CVodeSetLinearSolver(mpCvodeMem, mpSundialsLinearSolver, mpSundialsDenseMatrix);
346#elif CHASTE_SUNDIALS_VERSION >= 40000
347 /* Create dense SUNLinSol_Dense object for use by CVode */
348 mpSundialsLinearSolver = SUNLinSol_Dense(initial_values, mpSundialsDenseMatrix);
349
350 /* Call CVodeSetLinearSolver to attach the matrix and linear solver to CVode */
351 CVodeSetLinearSolver(mpCvodeMem, mpSundialsLinearSolver, mpSundialsDenseMatrix);
352#else
353 /* Create dense SUNLinearSolver object for use by CVode */
354 mpSundialsLinearSolver = SUNDenseLinearSolver(initial_values, mpSundialsDenseMatrix);
355
356 /* Call CVDlsSetLinearSolver to attach the matrix and linear solver to CVode */
357 CVDlsSetLinearSolver(mpCvodeMem, mpSundialsLinearSolver, mpSundialsDenseMatrix);
358#endif
359#else
360 // Attach a linear solver for Newton iteration
361 CVDense(mpCvodeMem, rInitialY.size());
362#endif
363 }
364
365 CVodeSetMaxStep(mpCvodeMem, maxStep);
366 // Change max steps if wanted
367 if (mMaxSteps > 0)
368 {
369 CVodeSetMaxNumSteps(mpCvodeMem, mMaxSteps);
370 CVodeSetMaxErrTestFails(mpCvodeMem, 15);
371 }
372 DeleteVector(initial_values);
373}
374
376{
377 if (mpCvodeMem)
378 {
379 CVodeFree(&mpCvodeMem);
380 }
381 mpCvodeMem = nullptr;
382
383#if CHASTE_SUNDIALS_VERSION >= 30000
384 if (mpSundialsLinearSolver)
385 {
386 /* Free the linear solver memory */
387 SUNLinSolFree(mpSundialsLinearSolver);
388 }
389 mpSundialsLinearSolver = nullptr;
390
391 if (mpSundialsDenseMatrix)
392 {
393 /* Free the matrix memory */
394 SUNMatDestroy(mpSundialsDenseMatrix);
395 }
396 mpSundialsDenseMatrix = nullptr;
397#endif
398}
399
400void CvodeAdaptor::SetForceReset(bool autoReset)
401{
402 mForceReset = autoReset;
403 if (mForceReset)
404 {
405 SetMinimalReset(false);
406 ResetSolver();
407 }
408}
409
410void CvodeAdaptor::SetMinimalReset(bool minimalReset)
411{
412 mForceMinimalReset = minimalReset;
414 {
415 SetForceReset(false);
416 }
417}
418
420{
421 // Doing this makes the SetupCvode think it needs to reset things.
423}
424
425void CvodeAdaptor::CvodeError(int flag, const char* msg)
426{
427 std::stringstream err;
428 char* p_flag_name = CVodeGetReturnFlagName(flag);
429 err << msg << ": " << p_flag_name;
430 free(p_flag_name);
431 std::cerr << err.str() << std::endl
432 << std::flush;
433 EXCEPTION(err.str());
434}
435
437 std::vector<double>& rYValues,
438 double startTime,
439 double endTime,
440 double maxStep,
441 double timeSampling)
442{
443 assert(endTime > startTime);
444 assert(timeSampling > 0.0);
445
447 if (mCheckForRoots && pOdeSystem->CalculateStoppingEvent(startTime, rYValues) == true)
448 {
449 EXCEPTION("(Solve with sampling) Stopping event is true for initial condition");
450 }
451
452 SetupCvode(pOdeSystem, rYValues, startTime, maxStep);
453
454 TimeStepper stepper(startTime, endTime, timeSampling);
455#if CHASTE_SUNDIALS_VERSION >= 60000
456 N_Vector yout = N_VMake_Serial(rYValues.size(), &(rYValues[0]), CvodeContextManager::Instance()->GetSundialsContext());
457#else
458 N_Vector yout = N_VMake_Serial(rYValues.size(), &(rYValues[0]));
459#endif
460
461 // Set up ODE solution
462 OdeSolution solutions;
463 solutions.SetNumberOfTimeSteps(stepper.EstimateTimeSteps());
464 solutions.rGetSolutions().push_back(rYValues);
465 solutions.rGetTimes().push_back(startTime);
466 solutions.SetOdeSystemInformation(pOdeSystem->GetSystemInformation());
467
468 // Main time sampling loop
469 while (!stepper.IsTimeAtEnd() && !mStoppingEventOccurred)
470 {
471 // This should stop CVODE going past the end of where we wanted and interpolating back.
472 int ierr = CVodeSetStopTime(mpCvodeMem, stepper.GetNextTime());
473 assert(ierr == CV_SUCCESS);
474 UNUSED_OPT(ierr); // avoid unused var warning
475
476 double tend;
477 ierr = CVode(mpCvodeMem, stepper.GetNextTime(), yout, &tend, CV_NORMAL);
478 if (ierr < 0)
479 {
481 DeleteVector(yout);
482 CvodeError(ierr, "CVODE failed to solve system");
483 }
484 // Store solution
485 solutions.rGetSolutions().push_back(rYValues);
486 solutions.rGetTimes().push_back(tend);
487 if (ierr == CV_ROOT_RETURN)
488 {
489 // Stopping event occurred
491 mStoppingTime = tend;
492 }
493 mLastSolutionTime = tend;
494 stepper.AdvanceOneTimeStep();
495 }
496
497 // stepper.EstimateTimeSteps may have been an overestimate...
498 solutions.SetNumberOfTimeSteps(stepper.GetTotalTimeStepsTaken());
499
500 int ierr = CVodeGetLastStep(mpCvodeMem, &mLastInternalStepSize);
501 assert(ierr == CV_SUCCESS);
502 UNUSED_OPT(ierr); // avoid unused var warning
504 DeleteVector(yout);
505
506 return solutions;
507}
508
510 std::vector<double>& rYValues,
511 double startTime,
512 double endTime,
513 double maxStep)
514{
515 assert(endTime > startTime);
516
518 if (mCheckForRoots && pOdeSystem->CalculateStoppingEvent(startTime, rYValues) == true)
519 {
520 EXCEPTION("(Solve) Stopping event is true for initial condition");
521 }
522
523 SetupCvode(pOdeSystem, rYValues, startTime, maxStep);
524
525#if CHASTE_SUNDIALS_VERSION >= 60000
526 N_Vector yout = N_VMake_Serial(rYValues.size(), &(rYValues[0]), CvodeContextManager::Instance()->GetSundialsContext());
527#else
528 N_Vector yout = N_VMake_Serial(rYValues.size(), &(rYValues[0]));
529#endif
530
531 // This should stop CVODE going past the end of where we wanted and interpolating back.
532 int ierr = CVodeSetStopTime(mpCvodeMem, endTime);
533 assert(ierr == CV_SUCCESS);
534 UNUSED_OPT(ierr); // avoid unused var warning
535
536 double tend;
537 ierr = CVode(mpCvodeMem, endTime, yout, &tend, CV_NORMAL);
538 if (ierr < 0)
539 {
541 DeleteVector(yout);
542 CvodeError(ierr, "CVODE failed to solve system");
543 }
544 if (ierr == CV_ROOT_RETURN)
545 {
546 // Stopping event occurred
548 mStoppingTime = tend;
549 }
550 assert(NV_DATA_S(yout) == &(rYValues[0]));
551 assert(!NV_OWN_DATA_S(yout));
552
553 // long int steps;
554 // CVodeGetNumSteps(mpCvodeMem, &steps);
555 // std::cout << " Solved to " << endTime << " in " << steps << " steps.\n";
556
557 ierr = CVodeGetLastStep(mpCvodeMem, &mLastInternalStepSize);
558 assert(ierr == CV_SUCCESS);
559 UNUSED_OPT(ierr); // avoid unused var warning
560 RecordStoppingPoint(tend, yout);
561 DeleteVector(yout);
562}
563
564CvodeAdaptor::CvodeAdaptor(double relTol, double absTol)
566 mpCvodeMem(nullptr),
567 mRelTol(relTol),
568 mAbsTol(absTol),
569 mLastInternalStepSize(-0.0),
570 mMaxSteps(0),
571 mCheckForRoots(false),
572 mLastSolutionState(nullptr),
573 mLastSolutionTime(0.0),
574#if CHASTE_SUNDIALS_VERSION >= 20400
575 mForceReset(false),
576#else
577 mForceReset(true),
578#endif
579 mForceMinimalReset(false)
580#if CHASTE_SUNDIALS_VERSION >= 30000
581 ,
582 mpSundialsDenseMatrix(nullptr),
583 mpSundialsLinearSolver(nullptr)
584#endif
585{
586}
587
593
595{
596 if (!mForceReset)
597 {
598 const unsigned size = GetVectorSize(yEnd);
600 for (unsigned i = 0; i < size; i++)
601 {
603 }
604 mLastSolutionTime = stopTime;
605 }
606}
607
608void CvodeAdaptor::SetTolerances(double relTol, double absTol)
609{
610 mRelTol = relTol;
611 mAbsTol = absTol;
612}
613
615{
616 return mRelTol;
617}
618
620{
621 return mAbsTol;
622}
623
628
633
634void CvodeAdaptor::SetMaxSteps(long int numSteps)
635{
636 mMaxSteps = numSteps;
637}
638
640{
641 return mMaxSteps;
642}
643
644// Serialization for Boost >= 1.36
647
648#endif // CHASTE_CVODE
#define EXCEPTION(message)
#define UNUSED_OPT(var)
#define CHASTE_CLASS_EXPORT(T)
void DeleteVector(VECTOR &rVec)
void CreateVectorIfEmpty(VECTOR &rVec, unsigned size)
void SetVectorComponent(VECTOR &rVec, unsigned index, double value)
double GetVectorComponent(const VECTOR &rVec, unsigned index)
unsigned GetVectorSize(const VECTOR &rVec)
void CopyFromStdVector(const std::vector< double > &rSrc, VECTOR &rDest)
void CopyToStdVector(const VECTOR &rSrc, std::vector< double > &rDest)
virtual double CalculateRootFunction(double time, const std::vector< double > &rY)
virtual void EvaluateYDerivatives(double time, const std::vector< double > &rY, std::vector< double > &rDY)=0
virtual bool CalculateStoppingEvent(double time, const std::vector< double > &rY)
boost::shared_ptr< const AbstractOdeSystemInformation > GetSystemInformation() const
static bool WithinAnyTolerance(double number1, double number2, double relTol=DBL_EPSILON, double absTol=DBL_EPSILON, bool printError=false)
double GetAbsoluteTolerance()
long int GetMaxSteps()
void CvodeError(int flag, const char *msg)
long int mMaxSteps
void SetMinimalReset(bool minimalReset)
N_Vector mLastSolutionState
void SetMaxSteps(long int numSteps)
void FreeCvodeMemory()
void SetForceReset(bool autoReset)
void CheckForStoppingEvents()
CvodeData mData
OdeSolution Solve(AbstractOdeSystem *pOdeSystem, std::vector< double > &rYValues, double startTime, double endTime, double maxStep, double timeSampling)
double GetLastStepSize()
double mLastInternalStepSize
void RecordStoppingPoint(double stopTime, N_Vector yEnd)
void SetupCvode(AbstractOdeSystem *pOdeSystem, std::vector< double > &rInitialY, double startTime, double maxStep)
double mLastSolutionTime
double GetRelativeTolerance()
void SetTolerances(double relTol=1e-4, double absTol=1e-6)
CvodeAdaptor(double relTol=1e-4, double absTol=1e-6)
void SetNumberOfTimeSteps(unsigned numTimeSteps)
std::vector< std::vector< double > > & rGetSolutions()
std::vector< double > & rGetTimes()
void SetOdeSystemInformation(boost::shared_ptr< const AbstractOdeSystemInformation > pOdeSystemInfo)
bool IsTimeAtEnd() const
unsigned GetTotalTimeStepsTaken() const
void AdvanceOneTimeStep()
double GetNextTime() const
unsigned EstimateTimeSteps() const
AbstractOdeSystem * pSystem
std::vector< realtype > * pY