Sacado Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Sacado_LFad_LogicalSparseOps.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Sacado Package
4 //
5 // Copyright 2006 NTESS and the Sacado contributors.
6 // SPDX-License-Identifier: LGPL-2.1-or-later
7 //
8 // ***********************************************************************
9 //
10 // The forward-mode AD classes in Sacado are a derivative work of the
11 // expression template classes in the Fad package by Nicolas Di Cesare.
12 // The following banner is included in the original Fad source code:
13 //
14 // ************ DO NOT REMOVE THIS BANNER ****************
15 //
16 // Nicolas Di Cesare <Nicolas.Dicesare@ann.jussieu.fr>
17 // http://www.ann.jussieu.fr/~dicesare
18 //
19 // CEMRACS 98 : C++ courses,
20 // templates : new C++ techniques
21 // for scientific computing
22 //
23 //********************************************************
24 //
25 // A short implementation ( not all operators and
26 // functions are overloaded ) of 1st order Automatic
27 // Differentiation in forward mode (FAD) using
28 // EXPRESSION TEMPLATES.
29 //
30 //********************************************************
31 // @HEADER
32 
33 #ifndef SACADO_LFAD_LOGICALSPARSEOPS_HPP
34 #define SACADO_LFAD_LOGICALSPARSEOPS_HPP
35 
36 #include <cmath>
37 #include <algorithm> // for std::min and std::max
38 #include <ostream> // for std::ostream
39 
40 #define FAD_UNARYOP_MACRO(OPNAME,OP,VALUE) \
41 namespace Sacado { \
42  namespace LFad { \
43  \
44  template <typename ExprT> \
45  class OP {}; \
46  \
47  template <typename ExprT> \
48  class Expr< OP<ExprT> > { \
49  public: \
50  \
51  typedef typename ExprT::value_type value_type; \
52  typedef typename ExprT::scalar_type scalar_type; \
53  typedef typename ExprT::base_expr_type base_expr_type; \
54  typedef typename ExprT::logical_type logical_type; \
55  \
56  Expr(const ExprT& expr_) : expr(expr_) {} \
57  \
58  int size() const { return expr.size(); } \
59  \
60  bool hasFastAccess() const { return expr.hasFastAccess(); } \
61  \
62  bool isPassive() const { return expr.isPassive();} \
63  \
64  value_type val() const { \
65  return VALUE; \
66  } \
67  \
68  logical_type dx(int i) const { \
69  return expr.dx(i); \
70  } \
71  \
72  logical_type fastAccessDx(int i) const { \
73  return expr.fastAccessDx(i); \
74  } \
75  \
76  protected: \
77  \
78  const ExprT& expr; \
79  }; \
80  \
81  template <typename T> \
82  inline Expr< OP< Expr<T> > > \
83  OPNAME (const Expr<T>& expr) \
84  { \
85  typedef OP< Expr<T> > expr_t; \
86  \
87  return Expr<expr_t>(expr); \
88  } \
89  } \
90 }
91 
92 FAD_UNARYOP_MACRO(operator+,
93  UnaryPlusOp,
94  expr.val())
95 FAD_UNARYOP_MACRO(operator-,
97  -expr.val())
100  std::exp(expr.val()))
103  std::log(expr.val()))
106  std::log10(expr.val()))
109  std::sqrt(expr.val()))
112  std::cos(expr.val()))
115  std::sin(expr.val()))
118  std::tan(expr.val()))
121  std::acos(expr.val()))
124  std::asin(expr.val()))
127  std::atan(expr.val()))
130  std::cosh(expr.val()))
133  std::sinh(expr.val()))
136  std::tanh(expr.val()))
139  acosh(expr.val()))
142  asinh(expr.val()))
145  atanh(expr.val()))
148  std::abs(expr.val()))
151  std::fabs(expr.val()))
154  std::cbrt(expr.val()))
155 
156 #undef FAD_UNARYOP_MACRO
157 
158 #define FAD_BINARYOP_MACRO(OPNAME,OP,VALUE,DX,FASTACCESSDX,VAL_CONST_DX_1,VAL_CONST_DX_2,CONST_DX_1,CONST_DX_2,CONST_FASTACCESSDX_1,CONST_FASTACCESSDX_2) \
159 namespace Sacado { \
160  namespace LFad { \
161  \
162  template <typename ExprT1, typename ExprT2> \
163  class OP {}; \
164  \
165  template <typename T1, typename T2> \
166  class Expr< OP< Expr<T1>, Expr<T2> > > { \
167  \
168  public: \
169  \
170  typedef Expr<T1> ExprT1; \
171  typedef Expr<T2> ExprT2; \
172  typedef typename ExprT1::value_type value_type_1; \
173  typedef typename ExprT2::value_type value_type_2; \
174  typedef typename Sacado::Promote<value_type_1, \
175  value_type_2>::type value_type; \
176  \
177  typedef typename ExprT1::scalar_type scalar_type_1; \
178  typedef typename ExprT2::scalar_type scalar_type_2; \
179  typedef typename Sacado::Promote<scalar_type_1, \
180  scalar_type_2>::type scalar_type; \
181  \
182  typedef typename ExprT1::base_expr_type base_expr_type_1; \
183  typedef typename ExprT2::base_expr_type base_expr_type_2; \
184  typedef typename Sacado::Promote<base_expr_type_1, \
185  base_expr_type_2>::type base_expr_type; \
186  \
187  typedef typename ExprT1::logical_type logical_type; \
188  \
189  Expr(const ExprT1& expr1_, const ExprT2& expr2_) : \
190  expr1(expr1_), expr2(expr2_) {} \
191  \
192  int size() const { \
193  int sz1 = expr1.size(), sz2 = expr2.size(); \
194  return sz1 > sz2 ? sz1 : sz2; \
195  } \
196  \
197  bool hasFastAccess() const { \
198  return expr1.hasFastAccess() && expr2.hasFastAccess(); \
199  } \
200  \
201  bool isPassive() const { \
202  return expr1.isPassive() && expr2.isPassive(); \
203  } \
204  \
205  const value_type val() const { \
206  return VALUE; \
207  } \
208  \
209  const logical_type dx(int i) const { \
210  return DX; \
211  } \
212  \
213  const logical_type fastAccessDx(int i) const { \
214  return FASTACCESSDX; \
215  } \
216  \
217  protected: \
218  \
219  const ExprT1& expr1; \
220  const ExprT2& expr2; \
221  \
222  }; \
223  \
224  template <typename T1> \
225  class Expr< OP< Expr<T1>, typename Expr<T1>::value_type> > { \
226  \
227  public: \
228  \
229  typedef Expr<T1> ExprT1; \
230  typedef typename ExprT1::value_type value_type; \
231  typedef typename ExprT1::scalar_type scalar_type; \
232  typedef typename ExprT1::base_expr_type base_expr_type; \
233  typedef typename ExprT1::value_type ConstT; \
234  typedef typename ExprT1::logical_type logical_type; \
235  \
236  Expr(const ExprT1& expr1_, const ConstT& c_) : \
237  expr1(expr1_), c(c_) {} \
238  \
239  int size() const { \
240  return expr1.size(); \
241  } \
242  \
243  bool hasFastAccess() const { \
244  return expr1.hasFastAccess(); \
245  } \
246  \
247  bool isPassive() const { \
248  return expr1.isPassive(); \
249  } \
250  \
251  const value_type val() const { \
252  return VAL_CONST_DX_2; \
253  } \
254  \
255  const logical_type dx(int i) const { \
256  return CONST_DX_2; \
257  } \
258  \
259  const logical_type fastAccessDx(int i) const { \
260  return CONST_FASTACCESSDX_2; \
261  } \
262  \
263  protected: \
264  \
265  const ExprT1& expr1; \
266  const ConstT& c; \
267  }; \
268  \
269  template <typename T2> \
270  class Expr< OP< typename Expr<T2>::value_type, Expr<T2> > > { \
271  \
272  public: \
273  \
274  typedef Expr<T2> ExprT2; \
275  typedef typename ExprT2::value_type value_type; \
276  typedef typename ExprT2::scalar_type scalar_type; \
277  typedef typename ExprT2::base_expr_type base_expr_type; \
278  typedef typename ExprT2::value_type ConstT; \
279  typedef typename ExprT2::logical_type logical_type; \
280  \
281  Expr(const ConstT& c_, const ExprT2& expr2_) : \
282  c(c_), expr2(expr2_) {} \
283  \
284  int size() const { \
285  return expr2.size(); \
286  } \
287  \
288  bool hasFastAccess() const { \
289  return expr2.hasFastAccess(); \
290  } \
291  \
292  bool isPassive() const { \
293  return expr2.isPassive(); \
294  } \
295  \
296  const value_type val() const { \
297  return VAL_CONST_DX_1; \
298  } \
299  \
300  const logical_type dx(int i) const { \
301  return CONST_DX_1; \
302  } \
303  \
304  const logical_type fastAccessDx(int i) const { \
305  return CONST_FASTACCESSDX_1; \
306  } \
307  \
308  protected: \
309  \
310  const ConstT& c; \
311  const ExprT2& expr2; \
312  }; \
313  \
314  template <typename T1, typename T2> \
315  inline Expr< OP< Expr<T1>, Expr<T2> > > \
316  OPNAME (const Expr<T1>& expr1, const Expr<T2>& expr2) \
317  { \
318  typedef OP< Expr<T1>, Expr<T2> > expr_t; \
319  \
320  return Expr<expr_t>(expr1, expr2); \
321  } \
322  \
323  template <typename T> \
324  inline Expr< OP< Expr<T>, Expr<T> > > \
325  OPNAME (const Expr<T>& expr1, const Expr<T>& expr2) \
326  { \
327  typedef OP< Expr<T>, Expr<T> > expr_t; \
328  \
329  return Expr<expr_t>(expr1, expr2); \
330  } \
331  \
332  template <typename T> \
333  inline Expr< OP< typename Expr<T>::value_type, Expr<T> > > \
334  OPNAME (const typename Expr<T>::value_type& c, \
335  const Expr<T>& expr) \
336  { \
337  typedef typename Expr<T>::value_type ConstT; \
338  typedef OP< ConstT, Expr<T> > expr_t; \
339  \
340  return Expr<expr_t>(c, expr); \
341  } \
342  \
343  template <typename T> \
344  inline Expr< OP< Expr<T>, typename Expr<T>::value_type > > \
345  OPNAME (const Expr<T>& expr, \
346  const typename Expr<T>::value_type& c) \
347  { \
348  typedef typename Expr<T>::value_type ConstT; \
349  typedef OP< Expr<T>, ConstT > expr_t; \
350  \
351  return Expr<expr_t>(expr, c); \
352  } \
353  } \
354 }
355 
356 
357 FAD_BINARYOP_MACRO(operator+,
359  expr1.val() + expr2.val(),
360  expr1.dx(i) || expr2.dx(i),
361  expr1.fastAccessDx(i) || expr2.fastAccessDx(i),
362  c + expr2.val(),
363  expr1.val() + c,
364  expr2.dx(i),
365  expr1.dx(i),
366  expr2.fastAccessDx(i),
367  expr1.fastAccessDx(i))
368 FAD_BINARYOP_MACRO(operator-,
370  expr1.val() - expr2.val(),
371  expr1.dx(i) || expr2.dx(i),
372  expr1.fastAccessDx(i) || expr2.fastAccessDx(i),
373  c - expr2.val(),
374  expr1.val() - c,
375  expr2.dx(i),
376  expr1.dx(i),
377  expr2.fastAccessDx(i),
378  expr1.fastAccessDx(i))
379 FAD_BINARYOP_MACRO(operator*,
381  expr1.val() * expr2.val(),
382  expr1.dx(i) || expr2.dx(i),
383  expr1.fastAccessDx(i) || expr2.fastAccessDx(i),
384  c * expr2.val(),
385  expr1.val() * c,
386  expr2.dx(i),
387  expr1.dx(i),
388  expr2.fastAccessDx(i),
389  expr1.fastAccessDx(i))
390 FAD_BINARYOP_MACRO(operator/,
392  expr1.val() / expr2.val(),
393  expr1.dx(i) || expr2.dx(i),
394  expr1.fastAccessDx(i) || expr2.fastAccessDx(i),
395  c / expr2.val(),
396  expr1.val() / c,
397  expr2.dx(i),
398  expr1.dx(i),
399  expr2.fastAccessDx(i),
400  expr1.fastAccessDx(i))
403  std::atan2(expr1.val(), expr2.val()),
404  expr1.dx(i) || expr2.dx(i),
405  expr1.fastAccessDx(i) || expr2.fastAccessDx(i),
406  std::atan2(c, expr2.val()),
407  std::atan2(expr1.val(), c),
408  expr2.dx(i),
409  expr1.dx(i),
410  expr2.fastAccessDx(i),
411  expr1.fastAccessDx(i))
414  std::pow(expr1.val(), expr2.val()),
415  expr1.dx(i) || expr2.dx(i),
416  expr1.fastAccessDx(i) || expr2.fastAccessDx(i),
417  std::pow(c, expr2.val()),
418  std::pow(expr1.val(), c),
419  expr2.dx(i),
420  expr1.dx(i),
421  expr2.fastAccessDx(i),
422  expr1.fastAccessDx(i))
425  std::max(expr1.val(), expr2.val()),
426  expr1.val() >= expr2.val() ? expr1.dx(i) : expr2.dx(i),
427  expr1.val() >= expr2.val() ? expr1.fastAccessDx(i) :
428  expr2.fastAccessDx(i),
429  std::max(c, expr2.val()),
430  std::max(expr1.val(), c),
431  c >= expr2.val() ? logical_type(0) : expr2.dx(i),
432  expr1.val() >= c ? expr1.dx(i) : logical_type(0),
433  c >= expr2.val() ? logical_type(0) : expr2.fastAccessDx(i),
434  expr1.val() >= c ? expr1.fastAccessDx(i) : logical_type(0))
436  MinOp,
437  std::min(expr1.val(), expr2.val()),
438  expr1.val() <= expr2.val() ? expr1.dx(i) : expr2.dx(i),
439  expr1.val() <= expr2.val() ? expr1.fastAccessDx(i) :
440  expr2.fastAccessDx(i),
441  std::min(c, expr2.val()),
442  std::min(expr1.val(), c),
443  c <= expr2.val() ? logical_type(0) : expr2.dx(i),
444  expr1.val() <= c ? expr1.dx(i) : logical_type(0),
445  c <= expr2.val() ? logical_type(0) : expr2.fastAccessDx(i),
446  expr1.val() <= c ? expr1.fastAccessDx(i) : logical_type(0))
447 
448 #undef FAD_BINARYOP_MACRO
449 
450 //-------------------------- Relational Operators -----------------------
451 
452 #define FAD_RELOP_MACRO(OP) \
453 namespace Sacado { \
454  namespace LFad { \
455  template <typename ExprT1, typename ExprT2> \
456  inline bool \
457  operator OP (const Expr<ExprT1>& expr1, \
458  const Expr<ExprT2>& expr2) \
459  { \
460  return expr1.val() OP expr2.val(); \
461  } \
462  \
463  template <typename ExprT2> \
464  inline bool \
465  operator OP (const typename Expr<ExprT2>::value_type& a, \
466  const Expr<ExprT2>& expr2) \
467  { \
468  return a OP expr2.val(); \
469  } \
470  \
471  template <typename ExprT1> \
472  inline bool \
473  operator OP (const Expr<ExprT1>& expr1, \
474  const typename Expr<ExprT1>::value_type& b) \
475  { \
476  return expr1.val() OP b; \
477  } \
478  } \
479 }
480 
481 FAD_RELOP_MACRO(==)
482 FAD_RELOP_MACRO(!=)
485 FAD_RELOP_MACRO(<=)
486 FAD_RELOP_MACRO(>=)
487 FAD_RELOP_MACRO(<<=)
488 FAD_RELOP_MACRO(>>=)
491 
492 #undef FAD_RELOP_MACRO
493 
494 namespace Sacado {
495 
496  namespace LFad {
497 
498  template <typename ExprT>
499  inline bool operator ! (const Expr<ExprT>& expr)
500  {
501  return ! expr.val();
502  }
503 
504  } // namespace LFad
505 
506 } // namespace Sacado
507 
508 //-------------------------- Boolean Operators -----------------------
509 namespace Sacado {
510 
511  namespace LFad {
512 
513  template <typename ExprT>
514  bool toBool(const Expr<ExprT>& x) {
515  bool is_zero = (x.val() == 0.0);
516  for (int i=0; i<x.size(); i++)
517  is_zero = is_zero && (x.dx(i) == 0.0);
518  return !is_zero;
519  }
520 
521  } // namespace Fad
522 
523 } // namespace Sacado
524 
525 #define FAD_BOOL_MACRO(OP) \
526 namespace Sacado { \
527  namespace LFad { \
528  template <typename ExprT1, typename ExprT2> \
529  inline bool \
530  operator OP (const Expr<ExprT1>& expr1, \
531  const Expr<ExprT2>& expr2) \
532  { \
533  return toBool(expr1) OP toBool(expr2); \
534  } \
535  \
536  template <typename ExprT2> \
537  inline bool \
538  operator OP (const typename Expr<ExprT2>::value_type& a, \
539  const Expr<ExprT2>& expr2) \
540  { \
541  return a OP toBool(expr2); \
542  } \
543  \
544  template <typename ExprT1> \
545  inline bool \
546  operator OP (const Expr<ExprT1>& expr1, \
547  const typename Expr<ExprT1>::value_type& b) \
548  { \
549  return toBool(expr1) OP b; \
550  } \
551  } \
552 }
553 
554 FAD_BOOL_MACRO(&&)
555 FAD_BOOL_MACRO(||)
556 
557 #undef FAD_BOOL_MACRO
558 
559 //-------------------------- I/O Operators -----------------------
560 
561 namespace Sacado {
562 
563  namespace LFad {
564 
565  template <typename ExprT>
566  std::ostream& operator << (std::ostream& os, const Expr<ExprT>& x) {
567  os << x.val() << " [";
568 
569  for (int i=0; i< x.size(); i++) {
570  os << " " << x.dx(i);
571  }
572 
573  os << " ]";
574  return os;
575  }
576 
577  } // namespace LFad
578 
579 } // namespace Sacado
580 
581 
582 #endif // SACADO_LFAD_LOGICALSPARSEOPS_HPP
cbrt(expr.val())
expr expr SinOp
expr2 expr1 expr2 expr2 c *expr2 c *expr1 c *expr2 c *expr1 MaxOp
asinh(expr.val())
SACADO_INLINE_FUNCTION bool operator!(const Expr< ExprT > &expr)
#define FAD_UNARYOP_MACRO(OPNAME, OP, USING, VALUE, DX, FASTACCESSDX)
asin(expr.val())
cosh(expr.val())
expr expr dx(i)
abs(expr.val())
atanh(expr.val())
expr expr CoshOp
expr expr ATanhOp
expr expr TanhOp
expr expr SqrtOp
expr expr ASinhOp
atan(expr.val())
expr1 expr1 expr2 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 MultiplicationOp
expr1 expr1 expr2 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 expr1 expr1 c *expr2 expr1 c *expr2 expr1 c *expr2 expr1 expr1 expr1 expr2 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 Atan2Op
expr val()
tanh(expr.val())
expr expr CosOp
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
#define FAD_BINARYOP_MACRO(OPNAME, OP, USING, VALUE, DX, CDX1, CDX2, FASTACCESSDX, VAL_CONST_DX_1, VAL_CONST_DX_2, CONST_DX_1, CONST_DX_2, CONST_FASTACCESSDX_1, CONST_FASTACCESSDX_2)
expr expr ATanOp
SimpleFad< ValueT > min(const SimpleFad< ValueT > &a, const SimpleFad< ValueT > &b)
expr expr ACosOp
sqrt(expr.val())
sinh(expr.val())
tan(expr.val())
expr1 expr1 expr2 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 expr1 expr1 c *expr2 expr1 c *expr2 expr1 c *expr2 expr1 DivisionOp
#define FAD_RELOP_MACRO(OP)
atan2(expr1.val(), expr2.val())
sin(expr.val())
expr expr expr fastAccessDx(i)) FAD_UNARYOP_MACRO(exp
log(expr.val())
expr expr ACoshOp
expr expr Log10Op
expr expr SinhOp
SACADO_INLINE_FUNCTION bool toBool(const Expr< T > &xx)
acosh(expr.val())
acos(expr.val())
SimpleFad< ValueT > max(const SimpleFad< ValueT > &a, const SimpleFad< ValueT > &b)
expr expr ASinOp
SACADO_INLINE_FUNCTION mpl::enable_if_c< ExprLevel< Expr< T1 > >::value==ExprLevel< Expr< T2 > >::value, Expr< PowerOp< Expr< T1 >, Expr< T2 > > > >::type pow(const Expr< T1 > &expr1, const Expr< T2 > &expr2)
#define FAD_BOOL_MACRO(OP)
exp(expr.val())
expr expr expr ExpOp
fabs(expr.val())
expr expr AbsOp
expr expr TanOp
log10(expr.val())
cos(expr.val())
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 PowerOp