Sacado Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Sacado_Fad_SFad.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Sacado Package
5 // Copyright (2006) Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
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 //
29 // The forward-mode AD classes in Sacado are a derivative work of the
30 // expression template classes in the Fad package by Nicolas Di Cesare.
31 // The following banner is included in the original Fad source code:
32 //
33 // ************ DO NOT REMOVE THIS BANNER ****************
34 //
35 // Nicolas Di Cesare <Nicolas.Dicesare@ann.jussieu.fr>
36 // http://www.ann.jussieu.fr/~dicesare
37 //
38 // CEMRACS 98 : C++ courses,
39 // templates : new C++ techniques
40 // for scientific computing
41 //
42 //********************************************************
43 //
44 // A short implementation ( not all operators and
45 // functions are overloaded ) of 1st order Automatic
46 // Differentiation in forward mode (FAD) using
47 // EXPRESSION TEMPLATES.
48 //
49 //********************************************************
50 // @HEADER
51 
52 #ifndef SACADO_FAD_SFAD_HPP
53 #define SACADO_FAD_SFAD_HPP
54 
55 #include "Sacado_ConfigDefs.h"
56 
57 #ifdef SACADO_NEW_FAD_DESIGN_IS_DEFAULT
58 
59 #include "Sacado_Fad_Exp_SFad.hpp"
60 
61 namespace Sacado {
62  namespace Fad {
63  template <typename T, int Num>
64  using SFad = Exp::GeneralFad< Exp::StaticFixedStorage<T,Num> >;
65  }
66 }
67 
68 #else
69 
73 
74 namespace Sacado {
75 
77  namespace Fad {
78 
79 #ifndef SACADO_FAD_DERIV_LOOP
80 #if defined(SACADO_VIEW_CUDA_HIERARCHICAL_DFAD) && !defined(SACADO_DISABLE_CUDA_IN_KOKKOS) && defined(__CUDA_ARCH__)
81 #define SACADO_FAD_DERIV_LOOP(I,SZ) for (int I=threadIdx.x; I<SZ; I+=blockDim.x)
82 #else
83 #define SACADO_FAD_DERIV_LOOP(I,SZ) for (int I=0; I<SZ; ++I)
84 #endif
85 #endif
86 
87 #ifndef SACADO_FAD_THREAD_SINGLE
88 #if (defined(SACADO_VIEW_CUDA_HIERARCHICAL) || defined(SACADO_VIEW_CUDA_HIERARCHICAL_DFAD)) && !defined(SACADO_DISABLE_CUDA_IN_KOKKOS) && defined(__CUDA_ARCH__)
89 #define SACADO_FAD_THREAD_SINGLE if (threadIdx.x == 0)
90 #else
91 #define SACADO_FAD_THREAD_SINGLE /* */
92 #endif
93 #endif
94 
96  template <typename T, int Num>
97  struct SFadExprTag {};
98 
99  // Forward declaration
100  template <typename T, int Num> class SFad;
101 
109  template <typename T, int Num>
111 
112  public:
113 
116 
119 
122 
127 
130  Expr() : val_( T(0.)) {
131  ss_array<T>::zero(dx_, Num); }
132 
134 
137  template <typename S>
140  val_(x) {
141  ss_array<T>::zero(dx_, Num);
142  }
143 
145 
149  Expr(const int sz, const T & x, const DerivInit zero_out = InitDerivArray) : val_(x) {
150 #if defined(SACADO_DEBUG) && !defined(__CUDA_ARCH__ )
151  if (sz != Num)
152  throw "SFad::SFad() Error: Supplied derivative dimension does not match compile time length.";
153 #endif
154  if (zero_out == InitDerivArray)
155  ss_array<T>::zero(dx_, Num);
156  }
157 
159 
165  Expr(const int sz, const int i, const T & x) :
166  val_(x) {
167 #if defined(SACADO_DEBUG) && !defined(__CUDA_ARCH__ )
168  if (sz != Num)
169  throw "SFad::SFad() Error: Supplied derivative dimension does not match compile time length.";
170  if (i >= Num)
171  throw "SFad::SFad() Error: Invalid derivative index.";
172 #endif
173 
174  ss_array<T>::zero(dx_, Num);
175  dx_[i]=1.;
176  }
177 
180  Expr(const Expr& x) :
181  val_(x.val()) {
182  for (int i=0; i<Num; i++)
183  dx_[i] = x.dx_[i];
184  }
185 
187  template <typename S>
190 #if defined(SACADO_DEBUG) && !defined(__CUDA_ARCH__ )
191  if (x.size() != Num)
192  throw "SFad::SFad() Error: Attempt to assign with incompatible sizes";
193 #endif
194 
195  for(int i=0; i<Num; ++i)
196  dx_[i] = x.fastAccessDx(i);
197 
198  this->val() = x.val();
199  }
200 
203  ~Expr() {}
204 
206 
213  void diff(const int ith, const int n) {
214 #if defined(SACADO_DEBUG) && !defined(__CUDA_ARCH__ )
215  if (n != Num)
216  throw "SFad::diff() Error: Supplied derivative dimension does not match compile time length.";
217 #endif
218 
219  ss_array<T>::zero(dx_, Num);
220  dx_[ith] = T(1.);
221  }
222 
224 
229  void resize(int sz) {
230 #if defined(SACADO_DEBUG) && !defined(__CUDA_ARCH__ )
231  if (sz != Num)
232  throw "SFad::resize() Error: Cannot resize fixed derivative array dimension";
233 #endif
234  }
235 
237 
242  void expand(int sz) { resize(sz); }
243 
246  void zero() { ss_array<T>::zero(dx_, Num); }
247 
250  void setUpdateValue(bool update_val) { }
251 
254  bool updateValue() const { return true; }
255 
258  void cache() const {}
259 
261  template <typename S>
263  SACADO_ENABLE_EXPR_FUNC(bool) isEqualTo(const Expr<S>& x) const {
264  typedef IsEqual<value_type> IE;
265  if (x.size() != this->size()) return false;
266  bool eq = IE::eval(x.val(), this->val());
267  for (int i=0; i<this->size(); i++)
268  eq = eq && IE::eval(x.dx(i), this->dx(i));
269  return eq;
270  }
271 
273 
278 
281  const T& val() const { return val_;}
282 
285  T& val() { return val_;}
286 
288 
293 
296  int size() const { return Num;}
297 
303  int availableSize() const { return Num; }
304 
307  bool hasFastAccess() const { return true; }
308 
311  bool isPassive() const { return false; }
312 
315  void setIsConstant(bool is_const) {}
316 
319  const T* dx() const { return &(dx_[0]);}
320 
323  const T& dx(int i) const { return dx_[i]; }
324 
327  T& fastAccessDx(int i) { return dx_[i];}
328 
331  const T& fastAccessDx(int i) const { return dx_[i];}
332 
334 
339 
341  template <typename S>
343  SACADO_ENABLE_VALUE_FUNC(Expr&) operator=(const S& v) {
344  val_ = v;
345  ss_array<T>::zero(dx_, Num);
346  return *this;
347  }
348 
351  Expr& operator=(const Expr& x) {
352  if (this != &x) {
353  // Copy value
354  val_ = x.val_;
355 
356  // Copy dx_
357  for (int i=0; i<Num; i++)
358  dx_[i] = x.dx_[i];
359  }
360  return *this;
361  }
362 
364  template <typename S>
366  SACADO_ENABLE_EXPR_FUNC(Expr&) operator=(const Expr<S>& x) {
367 #if defined(SACADO_DEBUG) && !defined(__CUDA_ARCH__ )
368  if (x.size() != Num)
369  throw "SFad::operator=() Error: Attempt to assign with incompatible sizes";
370 #endif
371 
372  for(int i=0; i<Num; ++i)
373  dx_[i] = x.fastAccessDx(i);
374 
375  val_ = x.val();
376 
377  return *this;
378  }
379 
381 
386 
388  template <typename S>
390  SACADO_ENABLE_VALUE_FUNC(Expr&) operator += (const S& v) {
391  this->val() += v;
392  return *this;
393  }
394 
396  template <typename S>
398  SACADO_ENABLE_VALUE_FUNC(Expr&) operator -= (const S& v) {
399  this->val() -= v;
400  return *this;
401  }
402 
404  template <typename S>
406  SACADO_ENABLE_VALUE_FUNC(Expr&) operator *= (const S& v) {
407  this->val() *= v;
408  for (int i=0; i<Num; ++i)
409  dx_[i] *= v;
410  return *this;
411  }
412 
414  template <typename S>
416  SACADO_ENABLE_VALUE_FUNC(Expr&) operator /= (const S& v) {
417  this->val() /= v;
418  for (int i=0; i<Num; ++i)
419  dx_[i] /= v;
420  return *this;
421  }
422 
424  template <typename S>
426  SACADO_ENABLE_EXPR_FUNC(Expr&) operator += (const Expr<S>& x) {
427 #if defined(SACADO_DEBUG) && !defined(__CUDA_ARCH__ )
428  if (x.size() != Num)
429  throw "SFad::operator+=() Error: Attempt to assign with incompatible sizes";
430 #endif
431 
432  for (int i=0; i<Num; ++i)
433  dx_[i] += x.fastAccessDx(i);
434 
435  val_ += x.val();
436 
437  return *this;
438  }
439 
441  template <typename S>
443  SACADO_ENABLE_EXPR_FUNC(Expr&) operator -= (const Expr<S>& x) {
444 #if defined(SACADO_DEBUG) && !defined(__CUDA_ARCH__ )
445  if (x.size() != Num)
446  throw "SFad::operator-=() Error: Attempt to assign with incompatible sizes";
447 #endif
448 
449  for(int i=0; i<Num; ++i)
450  dx_[i] -= x.fastAccessDx(i);
451 
452  val_ -= x.val();
453 
454  return *this;
455  }
456 
458  template <typename S>
460  SACADO_ENABLE_EXPR_FUNC(Expr&) operator *= (const Expr<S>& x) {
461  T xval = x.val();
462 
463 #if defined(SACADO_DEBUG) && !defined(__CUDA_ARCH__ )
464  if (x.size() != Num)
465  throw "SFad::operator*=() Error: Attempt to assign with incompatible sizes";
466 #endif
467 
468  for(int i=0; i<Num; ++i)
469  dx_[i] = val_ * x.fastAccessDx(i) + dx_[i] * xval;
470 
471  val_ *= xval;
472 
473  return *this;
474  }
475 
477  template <typename S>
479  SACADO_ENABLE_EXPR_FUNC(Expr&) operator /= (const Expr<S>& x) {
480  T xval = x.val();
481 
482 #if defined(SACADO_DEBUG) && !defined(__CUDA_ARCH__ )
483  if (x.size() != Num)
484  throw "SFad::operator/=() Error: Attempt to assign with incompatible sizes";
485 #endif
486 
487  for(int i=0; i<Num; ++i)
488  dx_[i] = ( dx_[i]*xval - val_*x.fastAccessDx(i) )/ (xval*xval);
489 
490  val_ /= xval;
491 
492  return *this;
493  }
494 
496 
497  protected:
498 
500  T dx_[Num];
501 
504 
505  }; // class Expr<SFadExprTag>
506 
507  } // namespace Fad
508 
509 } // namespace Sacado
510 
511 #define FAD_NS Fad
512 #include "Sacado_Fad_SFad_tmpl.hpp"
513 #undef FAD_NS
514 
515 #include "Sacado_Fad_Ops.hpp"
516 
517 #endif // SACADO_NEW_FAD_DESIGN_IS_DEFAULT
518 
519 #include "Sacado_Fad_ViewFad.hpp"
520 
521 #endif // SACADO_FAD_SFAD_HPP
Wrapper for a generic expression template.
expr expr dx(i)
KOKKOS_INLINE_FUNCTION SACADO_ENABLE_EXPR_FUNC(bool) isEqualTo(const Expr< S > &x) const
Returns whether two Fad objects have the same values.
KOKKOS_INLINE_FUNCTION const T & fastAccessDx(int i) const
Returns derivative component i without bounds checking.
#define SACADO_ENABLE_VALUE_CTOR_DECL
KOKKOS_INLINE_FUNCTION void zero()
Zero out the derivative array.
KOKKOS_INLINE_FUNCTION bool hasFastAccess() const
Returns true if derivative array is not empty.
KOKKOS_INLINE_FUNCTION Expr()
Default constructor.
KOKKOS_INLINE_FUNCTION int size() const
Returns number of derivative components.
RemoveConst< T >::type value_type
Typename of values.
KOKKOS_INLINE_FUNCTION const T * dx() const
Returns derivative array.
#define SACADO_ENABLE_EXPR_CTOR_DECL
KOKKOS_INLINE_FUNCTION Expr(const Expr &x)
Copy constructor.
static KOKKOS_INLINE_FUNCTION void zero(T *dest, int sz)
Zero out array dest of length sz.
KOKKOS_INLINE_FUNCTION void diff(const int ith, const int n)
Set Fad object as the ith independent variable.
KOKKOS_INLINE_FUNCTION T & fastAccessDx(int i)
Returns derivative component i without bounds checking.
expr val()
#define KOKKOS_INLINE_FUNCTION
#define T
Definition: Sacado_rad.hpp:573
KOKKOS_INLINE_FUNCTION bool isPassive() const
Returns true if derivative array is empty.
KOKKOS_INLINE_FUNCTION int availableSize() const
Returns number of derivative components that can be stored without reallocation.
KOKKOS_INLINE_FUNCTION const T & val() const
Returns value.
#define SACADO_ENABLE_VALUE_FUNC(RETURN_TYPE)
KOKKOS_INLINE_FUNCTION Expr(const S &x, SACADO_ENABLE_VALUE_CTOR_DECL)
Constructor with supplied value x.
Base template specification for testing equivalence.
KOKKOS_INLINE_FUNCTION bool updateValue() const
Return whether this Fad object has an updated value.
KOKKOS_INLINE_FUNCTION void setIsConstant(bool is_const)
Set whether variable is constant.
KOKKOS_INLINE_FUNCTION void setUpdateValue(bool update_val)
Set whether this Fad object should update values.
KOKKOS_INLINE_FUNCTION void resize(int sz)
Resize derivative array to length sz.
#define SACADO_ENABLE_EXPR_FUNC(RETURN_TYPE)
KOKKOS_INLINE_FUNCTION void expand(int sz)
Expand derivative array to size sz.
DerivInit
Enum use to signal whether the derivative array should be initialized in AD object constructors...
KOKKOS_INLINE_FUNCTION T & val()
Returns value.
ScalarType< value_type >::type scalar_type
Typename of scalar&#39;s (which may be different from T)
KOKKOS_INLINE_FUNCTION Expr(const int sz, const T &x, const DerivInit zero_out=InitDerivArray)
Constructor with size sz and value x.
KOKKOS_INLINE_FUNCTION Expr(const Expr< S > &x, SACADO_ENABLE_EXPR_CTOR_DECL)
Copy constructor from any Expression object.
Initialize the derivative array.
KOKKOS_INLINE_FUNCTION Expr(const int sz, const int i, const T &x)
Constructor with size sz, index i, and value x.
KOKKOS_INLINE_FUNCTION const T & dx(int i) const
Returns derivative component i with bounds checking.
KOKKOS_INLINE_FUNCTION void cache() const
Cache values.
A tag for specializing Expr for SFad expressions.
SFad< value_type, Num > base_expr_type
Typename of base-expressions.
GeneralFad< StaticFixedStorage< T, Num > > SFad