Rythmos - Transient Integration for Differential Equations  Version of the Day
 All Classes Functions Variables Typedefs Pages
Rythmos_CubicSplineInterpolator_decl.hpp
1 //@HEADER
2 // ***********************************************************************
3 //
4 // Rythmos Package
5 // Copyright (2006) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // This library is free software; you can redistribute it and/or modify
11 // it under the terms of the GNU Lesser General Public License as
12 // published by the Free Software Foundation; either version 2.1 of the
13 // License, or (at your option) any later version.
14 //
15 // This library is distributed in the hope that it will be useful, but
16 // WITHOUT ANY WARRANTY; without even the implied warranty of
17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18 // Lesser General Public License for more details.
19 //
20 // You should have received a copy of the GNU Lesser General Public
21 // License along with this library; if not, write to the Free Software
22 // Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23 // USA
24 // Questions? Contact Todd S. Coffey (tscoffe@sandia.gov)
25 //
26 // ***********************************************************************
27 //@HEADER
28 
29 #ifndef Rythmos_CUBIC_SPLINE_INTERPOLATOR_DECL_H
30 #define Rythmos_CUBIC_SPLINE_INTERPOLATOR_DECL_H
31 
32 #include "Rythmos_InterpolatorBase.hpp"
33 #include "Rythmos_Types.hpp"
34 
35 
36 namespace Rythmos {
37 
38 using Teuchos::Ptr;
39 
40 template<class Scalar>
41 class CubicSplineCoeff {
42  public:
43  Array<Scalar> t;
44  Array<RCP<Thyra::VectorBase<Scalar> > > a;
45  Array<RCP<Thyra::VectorBase<Scalar> > > b;
46  Array<RCP<Thyra::VectorBase<Scalar> > > c;
47  Array<RCP<Thyra::VectorBase<Scalar> > > d;
48 };
49 
50 
54 template<class Scalar>
55 class CubicSplineInterpolator : virtual public InterpolatorBase<Scalar>
56 {
57 public:
58 
61 
64 
66  bool supportsCloning() const;
67 
69  RCP<InterpolatorBase<Scalar> > cloneInterpolator() const;
70 
72  void setNodes(
73  const RCP<const typename DataStore<Scalar>::DataStoreVector_t> & nodes
74  );
75 
77  void interpolate(
78  const Array<Scalar> &t_values,
79  typename DataStore<Scalar>::DataStoreVector_t *data_out
80  ) const;
81 
83  int order() const;
84 
86  std::string description() const;
87 
89  void describe(
90  FancyOStream &out,
91  const Teuchos::EVerbosityLevel verbLevel
92  ) const;
93 
95  void setParameterList(RCP<ParameterList> const& paramList);
96 
98  RCP<ParameterList> getNonconstParameterList();
99 
101  RCP<ParameterList> unsetParameterList();
102 
104  RCP<const Teuchos::ParameterList> getValidParameters() const;
105 
106 private:
107 
108 
109  RCP<const typename DataStore<Scalar>::DataStoreVector_t> nodes_;
110 #ifdef HAVE_RYTHMOS_DEBUG
111  RCP<typename DataStore<Scalar>::DataStoreVector_t> nodes_copy_;
112 #endif // HAVE_RYTHMOS_DEBUG
113 
114  mutable CubicSplineCoeff<Scalar> splineCoeff_;
115  mutable bool splineCoeffComputed_;
116  bool nodesSet_;
117 
118  RCP<ParameterList> parameterList_;
119 
120 };
121 
122 // non-member constructor
123 template<class Scalar>
124 RCP<CubicSplineInterpolator<Scalar> > cubicSplineInterpolator();
125 
126 // Non-member helper function
127 // Algorithm from "CRC Standard Mathematical Tables and Formulae" by Daniel
128 // Zwillinger, 31st Edition, pg 738, natural cubic splines
129 // input: n, {x_0, x_1, ..., x_n}
130 // a_0 = f(x_0), a_1 = f(x_1), ... a_n = f(x_n)
131 // output: {a_j,b_j,c_j,d_j}, j=0..n-1
132 template<class Scalar>
133 void computeCubicSplineCoeff(
134  const typename DataStore<Scalar>::DataStoreVector_t & data,
135  const Ptr<CubicSplineCoeff<Scalar> > & coeffPtr
136  );
137 
138 template<class Scalar>
139 void validateCubicSplineCoeff(const CubicSplineCoeff<Scalar>& coeff);
140 
141 // Non-member helper function
142 // Evaluate cubic spline
143 template<class Scalar>
144 void evaluateCubicSpline(
145  const CubicSplineCoeff<Scalar>& coeff,
146  Teuchos::Ordinal j,
147  const Scalar& t,
148  const Ptr<Thyra::VectorBase<Scalar> >& S,
149  const Ptr<Thyra::VectorBase<Scalar> >& Sp = Teuchos::null,
150  const Ptr<Thyra::VectorBase<Scalar> >& Spp = Teuchos::null
151  );
152 
153 
154 } // namespace Rythmos
155 
156 
157 #endif // Rythmos_CUBIC_SPLINE_INTERPOLATOR_DECL_H
Base strategy class for interpolation functionality.
void setNodes(const RCP< const typename DataStore< Scalar >::DataStoreVector_t > &nodes)
RCP< const Teuchos::ParameterList > getValidParameters() const
RCP< InterpolatorBase< Scalar > > cloneInterpolator() const
Concrete implemenation of InterpolatorBase that implements cubic spline interpolation.
void describe(FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
void setParameterList(RCP< ParameterList > const &paramList)
void interpolate(const Array< Scalar > &t_values, typename DataStore< Scalar >::DataStoreVector_t *data_out) const