Chaste  Release::3.4
BackwardEulerIvpOdeSolver.hpp
1 /*
2 
3 Copyright (c) 2005-2016, University of Oxford.
4 All rights reserved.
5 
6 University of Oxford means the Chancellor, Masters and Scholars of the
7 University of Oxford, having an administrative office at Wellington
8 Square, Oxford OX1 2JD, UK.
9 
10 This file is part of Chaste.
11 
12 Redistribution and use in source and binary forms, with or without
13 modification, 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 
23 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
24 AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
25 IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
26 ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
27 LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
28 CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE
29 GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
30 HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
31 LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
32 OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
33 
34 */
35 
36 #ifndef BACKWARDEULERIVPODESOLVER_HPP_
37 #define BACKWARDEULERIVPODESOLVER_HPP_
38 
39 #include "AbstractOneStepIvpOdeSolver.hpp"
40 #include "AbstractOdeSystemWithAnalyticJacobian.hpp"
41 
42 #include "ChasteSerialization.hpp"
43 #include <boost/serialization/base_object.hpp>
44 
45 #include "AbstractOneStepIvpOdeSolver.hpp"
46 
52 {
53 private:
62  template<class Archive>
63  void serialize(Archive & archive, const unsigned int version)
64  {
65  // This calls serialize on the base class.
66  archive & boost::serialization::base_object<AbstractOneStepIvpOdeSolver>(*this);
67  //archive & mSizeOfOdeSystem; - this done in save and load construct now.
68  archive & mNumericalJacobianEpsilon;
70  }
71 
73  unsigned mSizeOfOdeSystem;
74 
77 
83 
84  /*
85  * NOTE: we use (unsafe) double pointers here rather than
86  * std::vectors because using std::vectors would lead to a
87  * slow down by a factor of about 4.
88  */
89 
91  double* mResidual;
92 
94  double** mJacobian;
95 
97  double* mUpdate;
98 
108  void ComputeResidual(AbstractOdeSystem* pAbstractOdeSystem,
109  double timeStep,
110  double time,
111  std::vector<double>& rCurrentYValues,
112  std::vector<double>& rCurrentGuess);
113 
123  void ComputeJacobian(AbstractOdeSystem* pAbstractOdeSystem,
124  double timeStep,
125  double time,
126  std::vector<double>& rCurrentYValues,
127  std::vector<double>& rCurrentGuess);
128 
135  void SolveLinearSystem();
136 
144  double ComputeNorm(double* pVector);
145 
155  void ComputeNumericalJacobian(AbstractOdeSystem* pAbstractOdeSystem,
156  double timeStep,
157  double time,
158  std::vector<double>& rCurrentYValues,
159  std::vector<double>& rCurrentGuess);
160 
161 protected:
162 
176  void CalculateNextYValue(AbstractOdeSystem* pAbstractOdeSystem,
177  double timeStep,
178  double time,
179  std::vector<double>& rCurrentYValues,
180  std::vector<double>& rNextYValues);
181 
182 public:
183 
189  BackwardEulerIvpOdeSolver(unsigned sizeOfOdeSystem);
190 
195 
202  void SetEpsilonForNumericalJacobian(double epsilon);
203 
209 
215  unsigned GetSystemSize() const {return mSizeOfOdeSystem;};
216 };
217 
220 
221 namespace boost
222 {
223 namespace serialization
224 {
229 template<class Archive>
230 inline void save_construct_data(
231  Archive & ar, const BackwardEulerIvpOdeSolver * t, const unsigned int file_version)
232 {
233  const unsigned system_size = t->GetSystemSize();
234  ar & system_size;
235 }
236 
243 template<class Archive>
244 inline void load_construct_data(
245  Archive & ar, BackwardEulerIvpOdeSolver * t, const unsigned int file_version)
246 {
247  unsigned ode_system_size;
248  ar >> ode_system_size;
249  ::new(t)BackwardEulerIvpOdeSolver(ode_system_size);
250 }
251 }
252 } // namespace ...
253 
254 #endif /*BACKWARDEULERIVPODESOLVER_HPP_*/
void SetEpsilonForNumericalJacobian(double epsilon)
BackwardEulerIvpOdeSolver(unsigned sizeOfOdeSystem)
void serialize(Archive &archive, const unsigned int version)
void ComputeJacobian(AbstractOdeSystem *pAbstractOdeSystem, double timeStep, double time, std::vector< double > &rCurrentYValues, std::vector< double > &rCurrentGuess)
void ComputeNumericalJacobian(AbstractOdeSystem *pAbstractOdeSystem, double timeStep, double time, std::vector< double > &rCurrentYValues, std::vector< double > &rCurrentGuess)
#define CHASTE_CLASS_EXPORT(T)
void ComputeResidual(AbstractOdeSystem *pAbstractOdeSystem, double timeStep, double time, std::vector< double > &rCurrentYValues, std::vector< double > &rCurrentGuess)
void CalculateNextYValue(AbstractOdeSystem *pAbstractOdeSystem, double timeStep, double time, std::vector< double > &rCurrentYValues, std::vector< double > &rNextYValues)
friend class boost::serialization::access