Sacado Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Sacado_Tay_Taylor.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Sacado 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 David M. Gay (dmgay@sandia.gov) or Eric T. Phipps
25 // (etphipp@sandia.gov).
26 //
27 // ***********************************************************************
28 // @HEADER
29 
30 #ifndef SACADO_TAY_TAYLOR_HPP
31 #define SACADO_TAY_TAYLOR_HPP
32 
33 #include "Sacado_ConfigDefs.h"
34 #include "Sacado_Base.hpp"
35 #include "Sacado_Handle.hpp"
36 #include <cmath>
37 #include <algorithm> // for std::min and std::max
38 #include <ostream> // for std::ostream
39 #include "Sacado_dummy_arg.hpp"
40 
41 namespace Sacado {
42 
44  namespace Tay {
45 
47 
51  template <typename T>
52  class Taylor : public Base< Taylor<T> > {
53  public:
54 
56  template <typename U>
57  struct apply {
58  typedef Taylor<U> type;
59  };
60 
62  typedef T value_type;
63 
66 
68  Taylor();
69 
71 
74  Taylor(const T& x);
75 
77 
82 
84 
87  Taylor(int d, const T & x);
88 
90 
93  explicit Taylor(int d);
94 
96  Taylor(const Taylor& x);
97 
99  ~Taylor();
100 
102 
106  void resize(int d, bool keep_coeffs);
107 
109 
112  void reserve(int d);
113 
115 
124  void copyForWrite() { th.makeOwnCopy(); }
125 
127  bool isEqualTo(const Taylor& x) const {
128  typedef IsEqual<value_type> IE;
129  if (x.degree() != this->degree()) return false;
130  bool eq = true;
131  for (int i=0; i<=this->degree(); i++)
132  eq = eq && IE::eval(x.coeff(i), this->coeff(i));
133  return eq;
134  }
135 
140 
142  Taylor<T>& operator=(const T& val);
143 
145 
150 
152  Taylor<T>& operator=(const Taylor<T>& x);
153 
155 
160 
162  const T& val() const { return th->coeff_[0];}
163 
165  T& val() { return th->coeff_[0];}
166 
168 
173 
175  int degree() const { return th->deg_;}
176 
178  bool hasFastAccess(int d) const { return th->deg_>=d;}
179 
181  const T* coeff() const { return th->coeff_;}
182 
184  T* coeff() { return th->coeff_;}
185 
187  T coeff(int i) const {
188  T tmp= i<=th->deg_ ? th->coeff_[i]:T(0.); return tmp;}
189 
191  T& fastAccessCoeff(int i) { return th->coeff_[i];}
192 
194  const T& fastAccessCoeff(int i) const { return th->coeff_[i];}
195 
197 
202 
204  Taylor<T> operator + () const;
205 
207  Taylor<T> operator - () const;
208 
210  Taylor<T>& operator += (const T& x);
211 
213  Taylor<T>& operator -= (const T& x);
214 
216  Taylor<T>& operator *= (const T& x);
217 
219  Taylor<T>& operator /= (const T& x);
220 
222  Taylor<T>& operator += (const Taylor<T>& x);
223 
225  Taylor<T>& operator -= (const Taylor<T>& x);
226 
228  Taylor<T>& operator *= (const Taylor<T>& x);
229 
231  Taylor<T>& operator /= (const Taylor<T>& x);
232 
234 
235  protected:
236 
238  int length() const { return th->len_; }
239 
241  void resizeCoeffs(int len);
242 
243  protected:
244 
245  struct TaylorData {
246 
249 
251  int deg_;
252 
254  int len_;
255 
257  TaylorData();
258 
260  TaylorData(const T& x);
261 
263  TaylorData(int d, const T & x);
264 
266  TaylorData(int d);
267 
269  TaylorData(int d, int l);
270 
272  TaylorData(const TaylorData& x);
273 
275  ~TaylorData();
276 
278  TaylorData& operator=(const TaylorData& x);
279 
280  };
281 
283 
284  }; // class Taylor
285 
287  template <typename T>
288  Taylor<T> diff(const Taylor<T>& x, int n = 1) {
289  const int d = x.degree();
290  if (n <= 0)
291  return x;
292  else if (n > d) {
293  Taylor<T> y(0);
294  return y;
295  }
296  Taylor<T> y(d-n);
297  int c = 1;
298  for (int i=1; i<=n; ++i)
299  c *= i;
300  for (int i=n; i<=d; ++i) {
301  y.fastAccessCoeff(i-n) = x.fastAccessCoeff(i) * T(c);
302  c = (c / (i-n+1)) * (i+1);
303  }
304  return y;
305  }
306 
307  // Operations
308  template <typename T> Taylor<T> operator+(const Base< Taylor<T> >& a,
309  const Base< Taylor<T> >& b);
310  template <typename T> Taylor<T> operator+(const typename Taylor<T>::value_type& a,
311  const Base< Taylor<T> >& b);
312  template <typename T> Taylor<T> operator+(const Base< Taylor<T> >& a,
313  const typename Taylor<T>::value_type& b);
314  template <typename T> Taylor<T> operator-(const Base< Taylor<T> >& a,
315  const Base< Taylor<T> >& b);
316  template <typename T> Taylor<T> operator-(const typename Taylor<T>::value_type& a,
317  const Base< Taylor<T> >& b);
318  template <typename T> Taylor<T> operator-(const Base< Taylor<T> >& a,
319  const typename Taylor<T>::value_type& b);
320  template <typename T> Taylor<T> operator*(const Base< Taylor<T> >& a,
321  const Base< Taylor<T> >& b);
322  template <typename T> Taylor<T> operator*(const typename Taylor<T>::value_type& a,
323  const Base< Taylor<T> >& b);
324  template <typename T> Taylor<T> operator*(const Base< Taylor<T> >& a,
325  const typename Taylor<T>::value_type& b);
326  template <typename T> Taylor<T> operator/(const Base< Taylor<T> >& a,
327  const Base< Taylor<T> >& b);
328  template <typename T> Taylor<T> operator/(const typename Taylor<T>::value_type& a,
329  const Base< Taylor<T> >& b);
330  template <typename T> Taylor<T> operator/(const Base< Taylor<T> >& a,
331  const typename Taylor<T>::value_type& b);
332  template <typename T> Taylor<T> exp(const Base< Taylor<T> >& a);
333  template <typename T> Taylor<T> log(const Base< Taylor<T> >& a);
334  template <typename T> Taylor<T> log10(const Base< Taylor<T> >& a);
335  template <typename T> Taylor<T> sqrt(const Base< Taylor<T> >& a);
336 #ifdef HAVE_SACADO_CXX11
337  template <typename T> Taylor<T> cbrt(const Base< Taylor<T> >& a);
338 #endif
339  template <typename T> Taylor<T> pow(const Base< Taylor<T> >& a,
340  const Base< Taylor<T> >& b);
341  template <typename T> Taylor<T> pow(const typename Taylor<T>::value_type& a,
342  const Base< Taylor<T> >& b);
343  template <typename T> Taylor<T> pow(const Base< Taylor<T> >& a,
344  const typename Taylor<T>::value_type& b);
345  template <typename T> void sincos(const Base< Taylor<T> >& a,
346  Taylor<T>& s, Taylor<T>& c);
347  template <typename T> Taylor<T> cos(const Base< Taylor<T> >& a);
348  template <typename T> Taylor<T> sin(const Base< Taylor<T> >& a);
349  template <typename T> Taylor<T> tan(const Base< Taylor<T> >& a);
350  template <typename T> void sinhcosh(const Base< Taylor<T> >& a,
351  Taylor<T>& s, Taylor<T>& c);
352  template <typename T> Taylor<T> cosh(const Base< Taylor<T> >& a);
353  template <typename T> Taylor<T> sinh(const Base< Taylor<T> >& a);
354  template <typename T> Taylor<T> tanh(const Base< Taylor<T> >& a);
355  template <typename T> Taylor<T> quad(const typename Taylor<T>::value_type& c0,
356  const Base< Taylor<T> >& a,
357  const Base< Taylor<T> >& b);
358  template <typename T> Taylor<T> acos(const Base< Taylor<T> >& a);
359  template <typename T> Taylor<T> asin(const Base< Taylor<T> >& a);
360  template <typename T> Taylor<T> atan(const Base< Taylor<T> >& a);
361  template <typename T> Taylor<T> atan2(const Base< Taylor<T> >& a,
362  const Base< Taylor<T> >& b);
363  template <typename T> Taylor<T> atan2(const typename Taylor<T>::value_type& a,
364  const Base< Taylor<T> >& b);
365  template <typename T> Taylor<T> atan2(const Base< Taylor<T> >& a,
366  const typename Taylor<T>::value_type& b);
367  template <typename T> Taylor<T> acosh(const Base< Taylor<T> >& a);
368  template <typename T> Taylor<T> asinh(const Base< Taylor<T> >& a);
369  template <typename T> Taylor<T> atanh(const Base< Taylor<T> >& a);
370  template <typename T> Taylor<T> abs(const Base< Taylor<T> >& a);
371  template <typename T> Taylor<T> fabs(const Base< Taylor<T> >& a);
372  template <typename T> Taylor<T> max(const Base< Taylor<T> >& a,
373  const Base< Taylor<T> >& b);
374  template <typename T> Taylor<T> max(const typename Taylor<T>::value_type& a,
375  const Base< Taylor<T> >& b);
376  template <typename T> Taylor<T> max(const Base< Taylor<T> >& a,
377  const typename Taylor<T>::value_type& b);
378  template <typename T> Taylor<T> min(const Base< Taylor<T> >& a,
379  const Base< Taylor<T> >& b);
380  template <typename T> Taylor<T> min(const typename Taylor<T>::value_type& a,
381  const Base< Taylor<T> >& b);
382  template <typename T> Taylor<T> min(const Base< Taylor<T> >& a,
383  const typename Taylor<T>::value_type& b);
384  template <typename T> bool operator==(const Base< Taylor<T> >& a,
385  const Base< Taylor<T> >& b);
386  template <typename T> bool operator==(const typename Taylor<T>::value_type& a,
387  const Base< Taylor<T> >& b);
388  template <typename T> bool operator==(const Base< Taylor<T> >& a,
389  const typename Taylor<T>::value_type& b);
390  template <typename T> bool operator!=(const Base< Taylor<T> >& a,
391  const Base< Taylor<T> >& b);
392  template <typename T> bool operator!=(const typename Taylor<T>::value_type& a,
393  const Base< Taylor<T> >& b);
394  template <typename T> bool operator!=(const Base< Taylor<T> >& a,
395  const typename Taylor<T>::value_type& b);
396  template <typename T> bool operator<=(const Base< Taylor<T> >& a,
397  const Base< Taylor<T> >& b);
398  template <typename T> bool operator<=(const typename Taylor<T>::value_type& a,
399  const Base< Taylor<T> >& b);
400  template <typename T> bool operator<=(const Base< Taylor<T> >& a,
401  const typename Taylor<T>::value_type& b);
402  template <typename T> bool operator>=(const Base< Taylor<T> >& a,
403  const Base< Taylor<T> >& b);
404  template <typename T> bool operator>=(const typename Taylor<T>::value_type& a,
405  const Base< Taylor<T> >& b);
406  template <typename T> bool operator>=(const Base< Taylor<T> >& a,
407  const typename Taylor<T>::value_type& b);
408  template <typename T> bool operator<(const Base< Taylor<T> >& a,
409  const Base< Taylor<T> >& b);
410  template <typename T> bool operator<(const typename Taylor<T>::value_type& a,
411  const Base< Taylor<T> >& b);
412  template <typename T> bool operator<(const Base< Taylor<T> >& a,
413  const typename Taylor<T>::value_type& b);
414  template <typename T> bool operator>(const Base< Taylor<T> >& a,
415  const Base< Taylor<T> >& b);
416  template <typename T> bool operator>(const typename Taylor<T>::value_type& a,
417  const Base< Taylor<T> >& b);
418  template <typename T> bool operator>(const Base< Taylor<T> >& a,
419  const typename Taylor<T>::value_type& b);
420  template <typename T> std::ostream& operator << (std::ostream& os,
421  const Base< Taylor<T> >& a);
422 
423  } // namespace Tay
424 
425 } // namespace Sacado
426 
428 #include "Sacado_Tay_TaylorImp.hpp"
429 
430 #endif // SACADO_TAY_TAYLOR_HPP
cbrt(expr.val())
Taylor< T > operator+(const Base< Taylor< T > > &a, const Base< Taylor< T > > &b)
A generic handle class.
ASinExprType< T >::expr_type asin(const Expr< T > &expr)
void resizeCoeffs(int len)
Resize coefficient array to new size.
void resize(int d, bool keep_coeffs)
Resize polynomial to degree d.
void sinhcosh(const Base< Taylor< T > > &a, Taylor< T > &s, Taylor< T > &c)
T value_type
Typename of values.
int deg_
Degree of polynomial.
Taylor< T > log(const Base< Taylor< T > > &a)
bool hasFastAccess(int d) const
Returns true if polynomial has degree &gt;= d.
Taylor< T > & operator-=(const T &x)
Subtraction-assignment operator with constant right-hand-side.
Taylor< T > asinh(const Base< Taylor< T > > &a)
TanhExprType< T >::expr_type tanh(const Expr< T > &expr)
T & val()
Returns value.
PowExprType< Expr< T1 >, Expr< T2 > >::expr_type pow(const Expr< T1 > &expr1, const Expr< T2 > &expr2)
TanExprType< T >::expr_type tan(const Expr< T > &expr)
T coeff(int i) const
Returns degree i term with bounds checking.
bool operator>(const Base< Taylor< T > > &a, const Base< Taylor< T > > &b)
Taylor< T > operator+() const
Unary-plus operator.
int len_
Length of allocated polynomial array.
TaylorData & operator=(const TaylorData &x)
Assignment operator.
Taylor< T > & operator*=(const T &x)
Multiplication-assignment operator with constant right-hand-side.
Taylor< T > operator/(const Base< Taylor< T > > &a, const Base< Taylor< T > > &b)
Taylor< T > & operator=(const T &val)
Assignment operator with constant right-hand-side.
ACosExprType< T >::expr_type acos(const Expr< T > &expr)
bool operator>=(const Base< Taylor< T > > &a, const Base< Taylor< T > > &b)
Turn Taylor into a meta-function class usable with mpl::apply.
Base class for Sacado types to control overload resolution.
Definition: Sacado_Base.hpp:46
Taylor< T > sin(const Base< Taylor< T > > &a)
#define T
Definition: Sacado_rad.hpp:573
void copyForWrite()
Prepare polynomial for writing.
const T & fastAccessCoeff(int i) const
Returns degree i term without bounds checking.
CacheTaylor< T > diff(const CacheTaylor< T > &x, int n=1)
Compute Taylor series of n-th derivative of x.
expr expr1 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 c *expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr2 expr1 expr2 expr1 expr1 expr1 c
Taylor< T > cos(const Base< Taylor< T > > &a)
T & fastAccessCoeff(int i)
Returns degree i term without bounds checking.
Taylor< T > operator-() const
Unary-minus operator.
Taylor< T > quad(const typename Taylor< T >::value_type &c0, const Base< Taylor< T > > &a, const Base< Taylor< T > > &b)
Base template specification for testing equivalence.
bool isEqualTo(const Taylor &x) const
Returns whether two Taylor objects have the same values.
Taylor< T > operator*(const Base< Taylor< T > > &a, const Base< Taylor< T > > &b)
Taylor< T > sinh(const Base< Taylor< T > > &a)
Taylor< T > sqrt(const Base< Taylor< T > > &a)
int length() const
Return length of array.
Log10ExprType< T >::expr_type log10(const Expr< T > &expr)
Taylor< T > max(const Base< Taylor< T > > &a, const Base< Taylor< T > > &b)
Taylor< T > fabs(const Base< Taylor< T > > &a)
const T & val() const
Returns value.
const T * coeff() const
Returns Taylor coefficient array.
Sacado::Handle< TaylorData > th
ATanExprType< T >::expr_type atan(const Expr< T > &expr)
Taylor< T > abs(const Base< Taylor< T > > &a)
void reserve(int d)
Reserve space for a degree d polynomial.
std::ostream & operator<<(std::ostream &os, const Expr< ExprT > &x)
Taylor< T > atan2(const Base< Taylor< T > > &a, const Base< Taylor< T > > &b)
Taylor< T > exp(const Base< Taylor< T > > &a)
T * coeff_
Taylor polynomial coefficients.
ScalarType< value_type >::type scalar_type
Typename of scalar&#39;s (which may be different from value_type)
void sincos(const Base< Taylor< T > > &a, Taylor< T > &s, Taylor< T > &c)
Taylor< T > atanh(const Base< Taylor< T > > &a)
Taylor polynomial class.
Taylor< T > cosh(const Base< Taylor< T > > &a)
bool operator==(const Base< Taylor< T > > &a, const Base< Taylor< T > > &b)
Taylor< T > min(const Base< Taylor< T > > &a, const Base< Taylor< T > > &b)
Taylor< T > operator-(const Base< Taylor< T > > &a, const Base< Taylor< T > > &b)
Taylor< T > & operator/=(const T &x)
Division-assignment operator with constant right-hand-side.
Taylor< T > acosh(const Base< Taylor< T > > &a)
Taylor< T > & operator+=(const T &x)
Addition-assignment operator with constant right-hand-side.
int degree() const
Returns degree of polynomial.
int n
bool operator!=(const Base< Taylor< T > > &a, const Base< Taylor< T > > &b)
Taylor()
Default constructor.
T * coeff()
Returns Taylor coefficient array.