MBDyn-1.7.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups
matvec.h
Go to the documentation of this file.
1 /* $Header: /var/cvs/mbdyn/mbdyn/mbdyn-1.0/libraries/libmbmath/matvec.h,v 1.11 2017/01/12 14:43:54 masarati Exp $ */
2 /*
3  * MBDyn (C) is a multibody analysis code.
4  * http://www.mbdyn.org
5  *
6  * Copyright (C) 1996-2017
7  *
8  * Pierangelo Masarati <masarati@aero.polimi.it>
9  * Paolo Mantegazza <mantegazza@aero.polimi.it>
10  *
11  * Dipartimento di Ingegneria Aerospaziale - Politecnico di Milano
12  * via La Masa, 34 - 20156 Milano, Italy
13  * http://www.aero.polimi.it
14  *
15  * Changing this copyright notice is forbidden.
16  *
17  * This program is free software; you can redistribute it and/or modify
18  * it under the terms of the GNU General Public License as published by
19  * the Free Software Foundation (version 2 of the License).
20  *
21  *
22  * This program is distributed in the hope that it will be useful,
23  * but WITHOUT ANY WARRANTY; without even the implied warranty of
24  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25  * GNU General Public License for more details.
26  *
27  * You should have received a copy of the GNU General Public License
28  * along with this program; if not, write to the Free Software
29  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
30  */
31 
32 /*
33  AUTHOR: Reinhard Resch <r.resch@secop.com>
34  Copyright (C) 2013(-2017) all rights reserved.
35 
36  The copyright of this code is transferred
37  to Pierangelo Masarati and Paolo Mantegazza
38  for use in the software MBDyn as described
39  in the GNU Public License version 2.1
40 */
41 
42 #ifndef ___MAT_VEC_H__INCLUDED___
43 #define ___MAT_VEC_H__INCLUDED___
44 
45 #include <cassert>
46 #include <iostream>
47 #include "myassert.h"
48 #include "gradient.h"
49 #include "matvec3.h" // FIXME: We need that for compatibility reasons
50 #include "matvec6.h"
51 #include "matvec3n.h"
52 #include "RotCoeff.hh"
53 #include "Rot.hh"
54 
55 #ifndef MATVEC_DEBUG
56  #ifdef DEBUG
57  #define MATVEC_DEBUG 1
58  #else
59  #define MATVEC_DEBUG 0
60  #endif
61 #endif
62 
63 #if MATVEC_DEBUG == 0 || defined(DEBUG)
64  #define MATVEC_ASSERT(expr) ASSERT(expr)
65 #elif MATVEC_DEBUG > 0
66  #define MATVEC_ASSERT(expr) assert(expr)
67 #endif
68 
69 namespace grad {
70 
71 template <typename T>
72  struct VectorSize {
73  static const int N = -1; // Note must not be zero because that would mean variable size!
74 };
75 
76 template <>
78  static const int N = 1;
79 };
80 
81 template <>
82  struct VectorSize<Vec3> {
83  static const int N = 3;
84 };
85 
86 template <>
87  struct VectorSize<Vec6> {
88  static const int N = 6;
89 };
90 
91 template <typename T, index_type N_rows, index_type N_cols>
92 class Matrix;
93 
94 template <typename T, index_type N_rows>
95 class Vector;
96 
101 template <typename ScalarBinaryFunction, typename T, typename ScalarLhsExpr, typename ScalarRhsExpr>
103 public:
104  static const bool bAlias = false;
105  static const index_type iMaxDerivatives = 0;
106  static const bool bVectorize = true;
107  static const index_type iDimension = -1;
108  typedef char vector_deriv_type;
109  typedef T ScalarType;
110  typedef T ExpressionType;
111 
112  GenericBinaryExpression(const ScalarLhsExpr& lhs, const ScalarRhsExpr& rhs)
113  :a(ScalarBinaryFunction::f(lhs, rhs)) {
114  };
115 
116  operator ExpressionType() const {
117  return a;
118  }
119 
121  return a;
122  }
123 
125  return 0;
126  }
127 
129  return vector_deriv_type();
130  }
131 
133  return std::numeric_limits<index_type>::max();
134  }
135 
137  return 0;
138  }
139 
141  return std::numeric_limits<index_type>::max();
142  }
143 
145  return 0;
146  }
147 
149  return 0;
150  }
151 
152  bool bHaveReferenceTo(const void* p) const {
153  return false;
154  }
155 
157  return iMaxDerivatives;
158  }
159 
160  void Compute() const {}
161 
162 private:
164 };
165 
166 template <typename ScalarUnaryFunction, typename T, typename ScalarExpr>
168 public:
169  static const bool bAlias = false;
170  static const index_type iMaxDerivatives = 0;
171  static const bool bVectorize = true;
172  static const index_type iDimension = -1;
173  typedef char vector_deriv_type;
174  typedef T ScalarType;
175  typedef T ExpressionType;
176 
177  GenericUnaryExpression(const ScalarExpr& expr)
178  :a(ScalarUnaryFunction::f(expr)) {
179  };
180 
181  operator ExpressionType() const {
182  return a;
183  }
184 
186  return a;
187  }
188 
190  return 0;
191  }
192 
194  return vector_deriv_type();
195  }
196 
198  return std::numeric_limits<index_type>::max();
199  }
200 
202  return 0;
203  }
204 
206  return std::numeric_limits<index_type>::max();
207  }
208 
210  return 0;
211  }
212 
214  return 0;
215  }
216 
217  bool bHaveReferenceTo(const void* p) const {
218  return false;
219  }
220 
222  return iMaxDerivatives;
223  }
224 
225  void Compute() const {}
226 
227 private:
229 };
230 
231 template <typename ScalarTypeLhs, typename ScalarTypeRhs>
233 private:
234  typedef void ScalarType;
235 };
236 
237 template <index_type N_SIZE>
238 struct CommonScalarType<Gradient<N_SIZE>, Gradient<N_SIZE> > {
240 };
241 
242 template <index_type N_SIZE>
245 };
246 
247 template <index_type N_SIZE>
250 };
251 
252 template <>
255 };
256 
257 template <typename T>
259  typedef T ScalarType;
260 };
261 
262 template <typename BinFunc, typename LhsExpr, typename RhsExpr>
263 struct BasicScalarType<BinaryExpr<BinFunc, LhsExpr, RhsExpr> > {
265 };
266 
267 template <typename UnFunc, typename Expr>
268 struct BasicScalarType<UnaryExpr<UnFunc, Expr> > {
270 };
271 
272 template <typename Expression>
275 };
276 
277 template <index_type N_SIZE, bool ALIAS>
278 struct BasicScalarType<DirectExpr<Gradient<N_SIZE>, ALIAS> > {
280 };
281 
282 template <index_type N_SIZE>
283 struct BasicScalarType<ConstExpr<Gradient<N_SIZE> > > {
285 };
286 
287 template <typename ScalarBinaryFunction, typename T, typename ScalarLhsExpr, typename ScalarRhsExpr>
288 struct BasicScalarType<GenericBinaryExpression<ScalarBinaryFunction, T, ScalarLhsExpr, ScalarRhsExpr> > {
290 };
291 
292 template <typename ScalarUnaryFunction, typename T, typename ScalarExpr>
293 struct BasicScalarType<GenericUnaryExpression<ScalarUnaryFunction, T, ScalarExpr> > {
295 };
296 
297 template <typename T>
299  typedef T ScalarType;
301 
302  template <typename Expression, typename PointerType>
303  static bool bHaveReferenceTo(const Expression&, const PointerType*, const PointerType*) {
304  return false;
305  }
306 };
307 
308 template <index_type N_SIZE>
309 struct ScalarTypeTraits<Gradient<N_SIZE> > {
312 
313  template <typename Expression, typename PointerType>
314  static bool bHaveReferenceTo(const GradientExpression<Expression>& g, const PointerType* pFirst, const PointerType* pLast) {
315  for (const PointerType* p = pFirst; p <= pLast; ++p) {
316  if (g.bHaveReferenceTo(p)) {
317  return true;
318  }
319  }
320 
321  return false;
322  }
323 };
324 
325 template <typename ScalarBinaryFunction, typename T, typename ScalarLhsExpr, typename ScalarRhsExpr>
328 };
329 
330 template <typename ScalarBinaryFunction, index_type N_SIZE, typename ScalarLhsExpr, typename ScalarRhsExpr>
331 struct ScalarBinaryExpressionTraits<ScalarBinaryFunction, Gradient<N_SIZE>, ScalarLhsExpr, ScalarRhsExpr>: ScalarTypeTraits<Gradient<N_SIZE> > {
333 };
334 
335 template <typename ScalarBinaryFunction, index_type N_SIZE, typename ScalarLhsExpr>
336 struct ScalarBinaryExpressionTraits<ScalarBinaryFunction, Gradient<N_SIZE>, ScalarLhsExpr, Gradient<N_SIZE> >: ScalarTypeTraits<Gradient<N_SIZE> > {
338 };
339 
340 template <typename ScalarBinaryFunction, index_type N_SIZE, typename ScalarLhsExpr>
341 struct ScalarBinaryExpressionTraits<ScalarBinaryFunction, Gradient<N_SIZE>, ScalarLhsExpr, scalar_func_type>: ScalarTypeTraits<Gradient<N_SIZE> > {
343 };
344 
345 template <typename ScalarBinaryFunction, index_type N_SIZE, typename ScalarRhsExpr>
346 struct ScalarBinaryExpressionTraits<ScalarBinaryFunction, Gradient<N_SIZE>, Gradient<N_SIZE>, ScalarRhsExpr>: ScalarTypeTraits<Gradient<N_SIZE> > {
348 };
349 
350 template <typename ScalarBinaryFunction, index_type N_SIZE, typename ScalarRhsExpr>
351 struct ScalarBinaryExpressionTraits<ScalarBinaryFunction, Gradient<N_SIZE>, scalar_func_type, ScalarRhsExpr>: ScalarTypeTraits<Gradient<N_SIZE> > {
353 };
354 
355 template <typename ScalarBinaryFunction, index_type N_SIZE>
356 struct ScalarBinaryExpressionTraits<ScalarBinaryFunction, Gradient<N_SIZE>, Gradient<N_SIZE>, Gradient<N_SIZE> >: ScalarTypeTraits<Gradient<N_SIZE> > {
358 };
359 
360 template <typename ScalarBinaryFunction, index_type N_SIZE>
361 struct ScalarBinaryExpressionTraits<ScalarBinaryFunction, Gradient<N_SIZE>, Gradient<N_SIZE>, scalar_func_type>: ScalarTypeTraits<Gradient<N_SIZE> > {
363 };
364 
365 template <typename ScalarBinaryFunction, index_type N_SIZE>
366 struct ScalarBinaryExpressionTraits<ScalarBinaryFunction, Gradient<N_SIZE>, scalar_func_type, Gradient<N_SIZE> >: ScalarTypeTraits<Gradient<N_SIZE> > {
368 };
369 
370 template <typename ScalarUnaryFunction, typename T, typename ScalarExpr>
373 };
374 
375 template <typename ScalarUnaryFunction, index_type N_SIZE, typename ScalarExpr>
376 struct ScalarUnaryExpressionTraits<ScalarUnaryFunction, Gradient<N_SIZE>, ScalarExpr>: ScalarTypeTraits<Gradient<N_SIZE> > {
378 };
379 
380 template <typename ScalarUnaryFunction, index_type N_SIZE>
381 struct ScalarUnaryExpressionTraits<ScalarUnaryFunction, Gradient<N_SIZE>, Gradient<N_SIZE> >: ScalarTypeTraits<Gradient<N_SIZE> > {
383 };
384 
385 template <typename ScalarBinaryFunction, typename ScalarLhsExpr, typename ScalarRhsExpr>
387 private:
390 
391 public:
393 
394  typedef typename ScalarBinaryExpressionTraits<ScalarBinaryFunction,
395  ScalarType,
396  ScalarLhsExpr,
398 
399  static ExpressionType f(const ScalarLhsExpr& u, const ScalarRhsExpr& v) {
400  return ExpressionType(u, v);
401  }
402 };
403 
404 template <typename ScalarUnaryFunction, typename ScalarExpr>
406 public:
408 
409  typedef typename ScalarUnaryExpressionTraits<ScalarUnaryFunction,
410  ScalarType,
412 
413  static ExpressionType f(const ScalarExpr& u) {
414  return ExpressionType(u);
415  }
416 };
417 
429 template <typename ScalarType>
430 inline bool bArrayOverlap(const ScalarType* pFirstArray1,
431  const ScalarType* pLastArray1,
432  const ScalarType* pFirstArray2,
433  const ScalarType* pLastArray2) {
434  MATVEC_ASSERT(pLastArray1 >= pFirstArray1);
435  MATVEC_ASSERT(pLastArray2 >= pFirstArray2);
436 
437  return (pFirstArray1 >= pFirstArray2 && pFirstArray1 <= pLastArray2)
438  || (pLastArray1 >= pFirstArray2 && pLastArray1 <= pLastArray2)
439  || (pFirstArray2 >= pFirstArray1 && pFirstArray2 <= pLastArray1)
440  || (pLastArray2 >= pFirstArray1 && pLastArray2 <= pLastArray1);
441 }
442 
443 template <typename ScalarType1, typename ScalarType2>
444 inline bool bArrayOverlap(const ScalarType1*, const ScalarType1*, const ScalarType2*, const ScalarType2*) {
445  // Self reference with different data types is not possible since unions are not supported
446  return false;
447 }
448 
449 template <typename T>
450 inline void ZeroInit(T* first, T* last) {
451 
452 }
453 
454 template <>
455 inline void ZeroInit<float>(float* first, float* last) {
456  array_fill(first, last, 0.0F);
457 }
458 
459 template <>
460 inline void ZeroInit<double>(double* first, double* last) {
461  array_fill(first, last, 0.0);
462 }
463 
464 template <>
465 inline void ZeroInit<long double>(long double* first, long double* last) {
466  array_fill(first, last, 0.0L);
467 }
468 
478 template <long DIFF>
479 struct IndexCheck {
480 private:
481  /* Whenever we get an compilation error here,
482  * the value of iNumRows or iNumCols in VectorExpression or MatrixExpression
483  * is not consistent with their template argument Expression
484  */
486 };
487 
488 template <>
489 struct IndexCheck<0L> {
491 };
492 
493 namespace MatVecHelp
494 {
495  template <typename T>
497 
498  };
499 
500  template <>
502  static const bool bAlias = false;
503  };
504 
505  template <typename Expression>
507  static const bool bAlias = GradientExpression<Expression>::bAlias;
508  };
509 
510  template <index_type N_SIZE, bool ALIAS>
511  struct AliasTypeHelper<DirectExpr<Gradient<N_SIZE>, ALIAS> > {
512  static const bool bAlias = DirectExpr<Gradient<N_SIZE>, ALIAS>::bAlias;
513  };
514 
515  template <index_type N_SIZE>
516  struct AliasTypeHelper<Gradient<N_SIZE> > {
517  static const bool bAlias = false;
518  };
519 
520  template <typename BinFunc, typename LhsExpr, typename RhsExpr>
521  struct AliasTypeHelper<BinaryExpr<BinFunc, LhsExpr, RhsExpr> > {
523  };
524 
525  template <typename UnFunc, typename Expr>
526  struct AliasTypeHelper<UnaryExpr<UnFunc, Expr> > {
527  static const bool bAlias = UnaryExpr<UnFunc, Expr>::bAlias;
528  };
529 
530  template <bool bAlias>
532 
533  template <>
535  template <typename MatrixType, typename Func, typename Expression>
536  static inline void ApplyMatrixFunc(MatrixType& A, const Expression& B, const Func& f) {
537  A.ApplyMatrixFuncNoAlias(B, f);
538  }
539  };
540 
541  template <>
543  template <typename MatrixType, typename Func, typename Expression>
544  static inline void ApplyMatrixFunc(MatrixType& A, const Expression& B, const Func& f) {
545  A.ApplyMatrixFuncAlias(B, f);
546  }
547  };
548 
549  struct Assign {
550  template <typename T, typename U>
551  static inline void Eval(T& b, const U& a) {
552  b = a;
553  }
554  };
555 
556  struct Add {
557  template <typename T, typename U>
558  static inline void Eval(T& b, const U& a) {
559  b += a;
560  }
561  };
562 
563  struct Sub {
564  template <typename T, typename U>
565  static inline void Eval(T& b, const U& a) {
566  b -= a;
567  }
568  };
569 
570  struct Mul {
571  template <typename T, typename U>
572  static inline void Eval(T& b, const U& a) {
573  b *= a;
574  }
575  };
576 
577  struct Div {
578  template <typename T, typename U>
579  static inline void Eval(T& b, const U& a) {
580  b /= a;
581  }
582  };
583 }
584 
585 template <typename Expression, index_type N_rows>
587 public:
588  static const index_type iNumRows = N_rows;
589  typedef typename Expression::ScalarType ScalarType;
590  typedef typename Expression::ExpressionType ExpressionType;
591 
592  explicit VectorExpression(const Expression& e)
593  :Expression(e) {
594 #if MATVEC_DEBUG > 0
595  AssertValid();
596 #endif
597  }
598 
599  template <typename Expr>
600  explicit VectorExpression(const Expr& u)
601  :Expression(u) {
602 #if MATVEC_DEBUG > 0
603  AssertValid();
604 #endif
605  }
606 
607  template <typename LhsExpr, typename RhsExpr>
608  explicit VectorExpression(const LhsExpr& u, const RhsExpr& v)
609  :Expression(u, v) {
610 #if MATVEC_DEBUG > 0
611  AssertValid();
612 #endif
613  }
614 
615  template <typename ScalarType>
616  explicit VectorExpression(const ScalarType* pArray, index_type iRows, index_type iOffset)
617  :Expression(pArray, iRows, iOffset) {
618 #if MATVEC_DEBUG > 0
619  AssertValid();
620 #endif
621  }
622 
623  index_type iGetNumRows() const { return Expression::iGetNumRows(); }
624 
625 private:
626 #if MATVEC_DEBUG > 0
627  void AssertValid() const {
628  MATVEC_ASSERT((Expression::iGetNumRows() == iNumRows) || (iNumRows == DYNAMIC_SIZE && Expression::iGetNumRows() >= 0));
629  }
630 #endif
631  typedef typename IndexCheck<iNumRows - Expression::iNumRows>::CheckType check_iNumRows;
632 };
633 
634 template <typename Expression, index_type N_rows, index_type N_cols, bool CLEAR_ALIAS=false>
636 public:
637  static const bool bAlias = CLEAR_ALIAS ? false : Expression::bAlias;
638  static const index_type iNumRows = N_rows;
639  static const index_type iNumCols = N_cols;
640  typedef typename Expression::ScalarType ScalarType;
641  typedef typename Expression::ExpressionType ExpressionType;
642 
643  explicit MatrixExpression(const Expression& e)
644  :Expression(e) {
645 #if MATVEC_DEBUG > 0
646  AssertValid();
647 #endif
648  }
649 
650  template <typename Expr>
651  explicit MatrixExpression(const Expr& u)
652  :Expression(u) {
653 #if MATVEC_DEBUG > 0
654  AssertValid();
655 #endif
656  }
657 
658  template <typename LhsExpr, typename RhsExpr>
659  explicit MatrixExpression(const LhsExpr& u, const RhsExpr& v)
660  :Expression(u, v) {
661 #if MATVEC_DEBUG > 0
662  AssertValid();
663 #endif
664  }
665 
667  return Expression::iGetNumRows();
668  }
670  return Expression::iGetNumCols();
671  }
672 
673 private:
674 
675 #if MATVEC_DEBUG > 0
676  void AssertValid() const {
677  MATVEC_ASSERT((Expression::iGetNumRows() == iNumRows) || (iNumRows == DYNAMIC_SIZE && Expression::iGetNumRows() >= 0));
678  MATVEC_ASSERT((Expression::iGetNumCols() == iNumCols) || (iNumCols == DYNAMIC_SIZE && Expression::iGetNumCols() >= 0));
679  }
680 #endif
681 
682  typedef typename IndexCheck<iNumRows - Expression::iNumRows>::CheckType check_iNumRows;
683  typedef typename IndexCheck<iNumCols - Expression::iNumCols>::CheckType check_iNumCols;
684 };
685 
690 template <typename ScalarBinFunc, typename VectorLhsExpr, typename VectorRhsExpr>
692 public:
693  static const bool bAlias = VectorLhsExpr::bAlias || VectorRhsExpr::bAlias;
694  static const index_type iNumRows = VectorLhsExpr::iNumRows;
695  typedef typename ScalarBinFunc::ScalarType ScalarType;
696  typedef typename ScalarBinFunc::ExpressionType ExpressionType;
697 
698  VectorVectorVectorBinaryExpr(const VectorLhsExpr& u, const VectorRhsExpr& v)
699  :oU(u), oV(v) {
700 
701  }
702 
704  MATVEC_ASSERT(i >= 1);
705  MATVEC_ASSERT(i <= iGetNumRows());
706  return ScalarBinFunc::f(oU(i), oV(i));
707  }
708 
710  MATVEC_ASSERT(oU.iGetNumRows() == oV.iGetNumRows());
711  MATVEC_ASSERT((oU.iGetNumRows() == iNumRows) || (iNumRows == DYNAMIC_SIZE && oU.iGetNumRows() >= 0));
712  return oU.iGetNumRows();
713  }
714 
715  template <typename ScalarType2>
716  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
717  return oU.bHaveReferenceTo(pFirst, pLast) || oV.bHaveReferenceTo(pFirst, pLast);
718  }
719 
720 private:
721  const VectorLhsExpr oU;
722  const VectorRhsExpr oV;
723 };
724 
729 template <typename ScalarBinFunc, typename VectorLhsExpr, typename ScalarRhsExpr>
731 public:
732  static const bool bAlias = VectorLhsExpr::bAlias || MatVecHelp::AliasTypeHelper<ScalarRhsExpr>::bAlias;
733  static const index_type iNumRows = VectorLhsExpr::iNumRows;
734  typedef typename ScalarBinFunc::ScalarType ScalarType;
735  typedef typename ScalarBinFunc::ExpressionType ExpressionType;
736 
737  VectorScalarVectorBinaryExpr(const VectorLhsExpr& u, const ScalarRhsExpr& v)
738  :oU(u), oV(v) {
739 
740  }
741 
743  MATVEC_ASSERT(i >= 1);
744  MATVEC_ASSERT(i <= iGetNumRows());
745  return ScalarBinFunc::f(oU(i), oV);
746  }
747 
749  MATVEC_ASSERT((oU.iGetNumRows() == iNumRows) || (iNumRows == DYNAMIC_SIZE && oU.iGetNumRows() >= 0));
750  return oU.iGetNumRows();
751  }
752 
753  template <typename ScalarType2>
754  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
756  return oU.bHaveReferenceTo(pFirst, pLast)
757  || ScalarTraits::bHaveReferenceTo(oV, pFirst, pLast);
758  }
759 
760 private:
761  const VectorLhsExpr oU;
762  const ScalarRhsExpr oV;
763 };
764 
769 template <typename ScalarUnaryFunc, typename VectorExpr>
771 public:
772  static const bool bAlias = VectorExpr::bAlias;
773  static const index_type iNumRows = VectorExpr::iNumRows;
774  typedef typename ScalarUnaryFunc::ScalarType ScalarType;
775  typedef typename ScalarUnaryFunc::ExpressionType ExpressionType;
776 
777  VectorVectorUnaryExpr(const VectorExpr& u)
778  :oU(u) {
779 
780  }
781 
783  MATVEC_ASSERT(i >= 1);
784  MATVEC_ASSERT(i <= iGetNumRows());
785  return ScalarUnaryFunc::f(oU(i));
786  }
787 
789  MATVEC_ASSERT((oU.iGetNumRows() == iNumRows) || (iNumRows == DYNAMIC_SIZE && oU.iGetNumRows() >= 0));
790  return oU.iGetNumRows();
791  }
792 
793  template <typename ScalarType2>
794  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
795  return oU.bHaveReferenceTo(pFirst, pLast);
796  }
797 
798 private:
799  const VectorExpr oU;
800 };
801 
802 template <index_type iStartIndex, index_type iEndIndex, typename VectorExpr>
804 public:
805  static const bool bAlias = VectorExpr::bAlias;
806  static const index_type iNumRows = iEndIndex - iStartIndex + 1;
807  typedef typename VectorExpr::ScalarType ScalarType;
808  typedef typename VectorExpr::ExpressionType ExpressionType;
809 
810  SubVectorExpr(const VectorExpr& u)
811  :oU(u) {
812 
813  }
814 
816  MATVEC_ASSERT(i >= 1);
817  MATVEC_ASSERT(i <= iNumRows);
818  return oU(i + iStartIndex - 1);
819  }
820 
822  MATVEC_ASSERT(oU.iGetNumRows() == VectorExpr::iNumRows);
823  return iNumRows;
824  }
825 
826  template <typename ScalarType2>
827  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
828  return oU.bHaveReferenceTo(pFirst, pLast);
829  }
830 
831 private:
832  typedef typename MaxSizeCheck<iEndIndex <= VectorExpr::iNumRows>::CheckType check_iEndIndex;
833  typedef typename MaxSizeCheck<iStartIndex >= 1>::CheckType check_iStartIndex;
835  const VectorExpr oU;
836 };
837 
838 template <typename VectorType, bool ALIAS=false>
840 public:
841  static const bool bAlias = ALIAS;
842  static const index_type iNumRows = VectorType::iNumRows;
843  typedef typename VectorType::ScalarType ScalarType;
844  typedef typename VectorType::ExpressionType ExpressionType;
845 
846  VectorDirectExpr(const VectorType& u)
847  :oU(u) {
848  }
849 
850  const ScalarType& operator()(index_type i) const {
851  MATVEC_ASSERT(i >= 1);
852  MATVEC_ASSERT(i <= iGetNumRows());
853  return oU(i);
854  }
855 
857  MATVEC_ASSERT((oU.iGetNumRows() == iNumRows) || (iNumRows == DYNAMIC_SIZE && oU.iGetNumRows() >= 0));
858  return oU.iGetNumRows();
859  }
860 
861  template <typename ScalarType2>
862  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
863  return oU.bHaveReferenceTo(pFirst, pLast);
864  }
865 
866 private:
867  const VectorType& oU;
868 };
869 
871 public:
872  static const bool bAlias = false;
873  static const index_type iNumRows = 3;
876 
878  :oU(u) {
879 
880  }
881 
882  const ScalarType& operator()(index_type i) const {
883  MATVEC_ASSERT(i >= 1);
884  MATVEC_ASSERT(i <= iNumRows);
885  return oU(i);
886  }
887 
889  return iNumRows;
890  }
891 
892  template <typename ScalarType2>
893  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
894  return bArrayOverlap(oU.pGetVec(), oU.pGetVec() + 2, pFirst, pLast);
895  }
896 
897 private:
898  const Vec3& oU;
899 };
900 
901 template <typename T, index_type N_rows, index_type N_offset>
902 class SliceVector {
903 public:
905  static const index_type iNumRows = N_rows;
908 
909  SliceVector(const ScalarType* p, index_type iRows, index_type iOffset)
910  :pVec(p) {
911  MATVEC_ASSERT(iRows == iNumRows);
912  MATVEC_ASSERT(iOffset = N_offset);
913  }
914 
915  const ScalarType& operator()(index_type iRow) const {
916  --iRow; // row index is 1-based
917  MATVEC_ASSERT(iRow >= 0);
918  MATVEC_ASSERT(iRow < iGetNumRows());
919  return *(pVec + iRow * N_offset);
920  }
921 
922  index_type iGetNumRows() const { return iNumRows; }
923 
924  template <typename ScalarType2>
925  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
926  return bArrayOverlap(pVec, pVec + (iNumRows - 1) * N_offset, pFirst, pLast);
927  }
928 
929  private:
930  const ScalarType* const pVec;
931 
934  };
935 
936  template <typename T, index_type N_rows>
937  class SliceVector<T, N_rows, DYNAMIC_SIZE> {
938  public:
940  static const index_type iNumRows = N_rows;
943 
944  SliceVector(const ScalarType* p, index_type iRows, index_type iOffset)
945  :pVec(p), iOffset(iOffset) {
946 
947  MATVEC_ASSERT(iRows == iNumRows);
948  }
949 
950  const ScalarType& operator()(index_type iRow) const {
951  --iRow; // row index is 1-based
952  MATVEC_ASSERT(iRow >= 0);
953  MATVEC_ASSERT(iRow < iGetNumRows());
954  return *(pVec + iRow * iOffset);
955  }
956 
957  index_type iGetNumRows() const { return iNumRows; }
958 
959  template <typename ScalarType2>
960  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
961  return bArrayOverlap(pVec, pVec + (iNumRows - 1) * iOffset, pFirst, pLast);
962  }
963 
964  private:
965  const ScalarType* const pVec;
967 
969  };
970 
971  template <typename T, index_type N_offset>
972  class SliceVector<T, DYNAMIC_SIZE, N_offset> {
973  public:
978 
979  SliceVector(const ScalarType* p, index_type iRows, index_type iOffset)
980  :pVec(p), iCurrRows(iRows) {
981 
982  MATVEC_ASSERT(iOffset == N_offset);
983  }
984 
985  const ScalarType& operator()(index_type iRow) const {
986  --iRow; // row index is 1-based
987  MATVEC_ASSERT(iRow >= 0);
988  MATVEC_ASSERT(iRow < iGetNumRows());
989  return *(pVec + iRow * N_offset);
990  }
991 
992  index_type iGetNumRows() const { return iCurrRows; }
993 
994  template <typename ScalarType2>
995  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
996  return bArrayOverlap(pVec, pVec + (iCurrRows - 1) * N_offset, pFirst, pLast);
997  }
998 
999  private:
1000  const ScalarType* const pVec;
1002 
1004  };
1005 
1006 
1007  template <typename T>
1009  public:
1014 
1015  SliceVector(const ScalarType* p, index_type iRows, index_type iOffset)
1016  :pVec(p), iCurrRows(iRows), iOffset(iOffset) {
1017 
1018  }
1019 
1020  const ScalarType& operator()(index_type iRow) const {
1021  --iRow; // row index is 1-based
1022  MATVEC_ASSERT(iRow >= 0);
1023  MATVEC_ASSERT(iRow < iGetNumRows());
1024  return *(pVec + iRow * iOffset);
1025  }
1026 
1027  index_type iGetNumRows() const { return iCurrRows; }
1028 
1029  template <typename ScalarType2>
1030  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
1031  return bArrayOverlap(pVec, pVec + (iCurrRows - 1) * iOffset, pFirst, pLast);
1032  }
1033 
1034 private:
1035  const ScalarType* const pVec;
1038 
1040 };
1041 
1042 template <typename MatrixExpr>
1044 public:
1045  static const bool bAlias = MatrixExpr::bAlias;
1046  static const index_type iNumRows = MatrixExpr::iNumRows;
1047  typedef typename MatrixExpr::ScalarType ScalarType;
1048  typedef typename MatrixExpr::ExpressionType ExpressionType;
1049 
1050  ColumnVectorExpr(const MatrixExpr& A, index_type iCol)
1051  :A(A), iCol(iCol) {
1052  MATVEC_ASSERT((iNumRows == A.iGetNumRows()) || (iNumRows == DYNAMIC_SIZE && A.iGetNumRows() >= 0));
1053  }
1054 
1056  MATVEC_ASSERT(iRow >= 1);
1057  MATVEC_ASSERT(iRow <= iGetNumRows());
1058  return A(iRow, iCol);
1059  }
1060 
1061  index_type iGetNumRows() const { return A.iGetNumRows(); }
1062 
1063  template <typename ScalarType2>
1064  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
1065  return A.bHaveReferenceTo(pFirst, pLast);
1066  }
1067 
1068 private:
1069  const MatrixExpr A;
1071 };
1072 
1073 template <typename MatrixExpr>
1075 public:
1076  static const bool bAlias = MatrixExpr::bAlias;
1077  static const index_type iNumRows = MatrixExpr::iNumCols;
1078  typedef typename MatrixExpr::ScalarType ScalarType;
1079  typedef typename MatrixExpr::ExpressionType ExpressionType;
1080 
1081  RowVectorExpr(const MatrixExpr& A, index_type iRow)
1082  :A(A), iRow(iRow) {
1083  MATVEC_ASSERT((iNumRows == A.iGetNumCols()) || (iNumRows == DYNAMIC_SIZE && A.iGetNumCols() >= 0));
1084  }
1085 
1087  MATVEC_ASSERT(iCol >= 1);
1088  MATVEC_ASSERT(iCol <= iGetNumRows());
1089  return A(iRow, iCol);
1090  }
1091 
1092  index_type iGetNumRows() const { return A.iGetNumCols(); }
1093 
1094  template <typename ScalarType2>
1095  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
1096  return A.bHaveReferenceTo(pFirst, pLast);
1097  }
1098 
1099 private:
1100  const MatrixExpr A;
1102 };
1103 
1104 template <typename MatrixExpr>
1106 public:
1107  static const bool bAlias = MatrixExpr::bAlias;
1108  static const index_type iNumRows = MatrixExpr::iNumCols;
1109  static const index_type iNumCols = MatrixExpr::iNumRows;
1110  typedef typename MatrixExpr::RowVectorType ColumnVectorType;
1111  typedef typename MatrixExpr::ColumnVectorType RowVectorType;
1112  typedef typename MatrixExpr::ScalarType ScalarType;
1113  typedef typename MatrixExpr::ExpressionType ExpressionType;
1114 
1115  TransposedMatrix(const MatrixExpr& A)
1116  :A(A) {
1117  MATVEC_ASSERT((iNumRows == A.iGetNumCols()) || (iNumRows == DYNAMIC_SIZE && A.iGetNumCols() >= 0));
1118  MATVEC_ASSERT((iNumCols == A.iGetNumRows()) || (iNumCols == DYNAMIC_SIZE && A.iGetNumRows() >= 0));
1119  }
1120 
1122  MATVEC_ASSERT(i >= 1);
1123  MATVEC_ASSERT(i <= iGetNumRows());
1124  MATVEC_ASSERT(j >= 1);
1125  MATVEC_ASSERT(j <= iGetNumCols());
1126  return A(j, i);
1127  }
1128 
1130  MATVEC_ASSERT((iNumRows == A.iGetNumCols()) || (iNumRows == DYNAMIC_SIZE && A.iGetNumCols() >= 0));
1131  return A.iGetNumCols();
1132  }
1133 
1135  MATVEC_ASSERT((iNumCols == A.iGetNumRows()) || (iNumCols == DYNAMIC_SIZE && A.iGetNumRows() >= 0));
1136  return A.iGetNumRows();
1137  }
1138 
1140  MATVEC_ASSERT(iRow >= 1);
1141  MATVEC_ASSERT(iRow <= iGetNumRows());
1142  return A.GetCol(iRow);
1143  }
1144 
1146  MATVEC_ASSERT(iCol >= 1);
1147  MATVEC_ASSERT(iCol <= iGetNumCols());
1148  return A.GetRow(iCol);
1149  }
1150 
1151  template <typename ScalarType2>
1152  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
1153  return A.bHaveReferenceTo(pFirst, pLast);
1154  }
1155 
1156 private:
1157  const MatrixExpr A;
1158 };
1159 
1160 template <index_type iRowStart,
1161  index_type iRowEnd,
1162  index_type iColStart,
1163  index_type iColEnd,
1164  typename MatrixExpr>
1166 public:
1167  static const bool bAlias = MatrixExpr::bAlias;
1168  static const index_type iNumRows = iRowEnd - iRowStart + 1;
1169  static const index_type iNumCols = iColEnd - iColStart + 1;
1170 
1173  typedef typename MatrixExpr::ScalarType ScalarType;
1174  typedef typename MatrixExpr::ExpressionType ExpressionType;
1175 
1176  SubMatrixExpr(const MatrixExpr& A)
1177  :A(A) {
1178  MATVEC_ASSERT(MatrixExpr::iNumCols == A.iGetNumCols());
1179  MATVEC_ASSERT(MatrixExpr::iNumRows == A.iGetNumRows());
1180  }
1181 
1183  MATVEC_ASSERT(i >= 1);
1184  MATVEC_ASSERT(i <= iGetNumRows());
1185  MATVEC_ASSERT(j >= 1);
1186  MATVEC_ASSERT(j <= iGetNumCols());
1187  return A(i + iRowStart - 1, j + iColStart - 1);
1188  }
1189 
1191  return iNumRows;
1192  }
1193 
1195  return iNumCols;
1196  }
1197 
1199  MATVEC_ASSERT(iRow >= 1);
1200  MATVEC_ASSERT(iRow <= iGetNumRows());
1201  return RowVectorType(A.GetRow(iRow + iRowStart - 1));
1202  }
1203 
1205  MATVEC_ASSERT(iCol >= 1);
1206  MATVEC_ASSERT(iCol <= iGetNumCols());
1207  return ColumnVectorType(A.GetCol(iCol + iColStart - 1));
1208  }
1209 
1210  template <typename ScalarType2>
1211  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
1212  return A.bHaveReferenceTo(pFirst, pLast);
1213  }
1214 
1215 private:
1216  typedef typename MaxSizeCheck<iRowStart >= 1>::CheckType check_iRowStart;
1217  typedef typename MaxSizeCheck<iRowEnd <= MatrixExpr::iNumRows>::CheckType check_iRowEnd;
1218  typedef typename MaxSizeCheck<iColStart >= 1>::CheckType check_iColStart;
1219  typedef typename MaxSizeCheck<iColEnd <= MatrixExpr::iNumCols>::CheckType check_iColEnd;
1220  typedef typename MaxSizeCheck<iNumRows != DYNAMIC_SIZE>::CheckType check_iNumRows;
1221  typedef typename MaxSizeCheck<iNumCols != DYNAMIC_SIZE>::CheckType check_iNumCols;
1222  const MatrixExpr A;
1223 };
1224 
1225 template <typename ScalarBinFunc, typename MatrixLhsExpr, typename MatrixRhsExpr>
1227 public:
1228  static const bool bAlias = MatrixLhsExpr::bAlias || MatrixRhsExpr::bAlias;
1229  static const index_type iNumRows = MatrixLhsExpr::iNumRows;
1230  static const index_type iNumCols = MatrixLhsExpr::iNumCols;
1231  typedef typename ScalarBinFunc::ScalarType ScalarType;
1232  typedef typename ScalarBinFunc::ExpressionType ExpressionType;
1235 
1236  MatrixMatrixMatrixBinaryExpr(const MatrixLhsExpr& u, const MatrixRhsExpr& v)
1237  :oU(u), oV(v) {
1238 
1239  MATVEC_ASSERT(iNumRows == oU.iGetNumRows());
1240  MATVEC_ASSERT(iNumCols == oU.iGetNumCols());
1241  MATVEC_ASSERT(iNumRows == oV.iGetNumRows());
1242  MATVEC_ASSERT(iNumCols == oV.iGetNumCols());
1243  }
1244 
1246  MATVEC_ASSERT(i >= 1);
1247  MATVEC_ASSERT(i <= iGetNumRows());
1248  MATVEC_ASSERT(j >= 1);
1249  MATVEC_ASSERT(j <= iGetNumCols());
1250  return ScalarBinFunc::f(oU(i, j), oV(i, j));
1251  }
1252 
1254  MATVEC_ASSERT(oU.iGetNumRows() == oV.iGetNumRows());
1255  MATVEC_ASSERT(iNumRows == oU.iGetNumRows());
1256  return oU.iGetNumRows();
1257  }
1258 
1260  MATVEC_ASSERT(oU.iGetNumCols() == oV.iGetNumCols());
1261  MATVEC_ASSERT(iNumCols == oU.iGetNumCols());
1262  return oU.iGetNumCols();
1263  }
1264 
1266  return RowVectorType(*this, iRow);
1267  }
1268 
1270  return ColumnVectorType(*this, iCol);
1271  }
1272 
1273  template <typename ScalarType2>
1274  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
1275  return oU.bHaveReferenceTo(pFirst, pLast) || oV.bHaveReferenceTo(pFirst, pLast);
1276  }
1277 
1278 private:
1279  const MatrixLhsExpr oU;
1280  const MatrixRhsExpr oV;
1281 
1282  // check if the dimensions of both matrices are the same
1283  typedef typename IndexCheck<MatrixLhsExpr::iNumRows - MatrixRhsExpr::iNumRows>::CheckType check_iNumRows;
1284  typedef typename IndexCheck<MatrixLhsExpr::iNumCols - MatrixRhsExpr::iNumCols>::CheckType check_iNumCols;
1285 };
1286 
1287 
1288 
1293 template <typename ScalarBinFunc, typename MatrixLhsExpr, typename ScalarRhsExpr>
1295 public:
1296  static const bool bAlias = MatrixLhsExpr::bAlias || MatVecHelp::AliasTypeHelper<ScalarRhsExpr>::bAlias;
1297  static const index_type iNumRows = MatrixLhsExpr::iNumRows;
1298  static const index_type iNumCols = MatrixLhsExpr::iNumCols;
1299  typedef typename ScalarBinFunc::ScalarType ScalarType;
1300  typedef typename ScalarBinFunc::ExpressionType ExpressionType;
1303 
1304  MatrixScalarMatrixBinaryExpr(const MatrixLhsExpr& u, const ScalarRhsExpr& v)
1305  :oU(u), oV(v) {
1306 
1307  MATVEC_ASSERT((iNumRows == oU.iGetNumRows()) || (iNumRows == DYNAMIC_SIZE && oU.iGetNumRows() >= 0));
1308  MATVEC_ASSERT((iNumCols == oU.iGetNumCols()) || (iNumCols == DYNAMIC_SIZE && oU.iGetNumCols() >= 0)) ;
1309  }
1310 
1312  MATVEC_ASSERT(i >= 1);
1313  MATVEC_ASSERT(i <= iGetNumRows());
1314  MATVEC_ASSERT(j >= 1);
1315  MATVEC_ASSERT(j <= iGetNumCols());
1316  return ScalarBinFunc::f(oU(i, j), oV);
1317  }
1318 
1320  MATVEC_ASSERT((iNumRows == oU.iGetNumRows()) || (iNumRows == DYNAMIC_SIZE && oU.iGetNumRows() >= 0));
1321  return oU.iGetNumRows();
1322  }
1323 
1325  MATVEC_ASSERT((iNumCols == oU.iGetNumCols()) || (iNumCols == DYNAMIC_SIZE && oU.iGetNumCols() >= 0));
1326  return oU.iGetNumCols();
1327  }
1328 
1330  MATVEC_ASSERT(iRow >= 1);
1331  MATVEC_ASSERT(iRow <= iGetNumRows());
1332  return RowVectorType(*this, iRow);
1333  }
1334 
1336  MATVEC_ASSERT(iCol >= 1);
1337  MATVEC_ASSERT(iCol <= iGetNumCols());
1338  return ColumnVectorType(*this, iCol);
1339  }
1340 
1341  template <typename ScalarType2>
1342  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
1343  typedef ScalarTypeTraits<BasicScalarType<ScalarRhsExpr> > ScalarTraits;
1344  return oU.bHaveReferenceTo(pFirst, pLast)
1345  || ScalarTraits::bHaveReferenceTo(oV, pFirst, pLast);
1346  }
1347 
1348 private:
1349  const MatrixLhsExpr oU;
1350  const ScalarRhsExpr oV;
1351 };
1352 
1353 template <typename ScalarUnaryFunc, typename MatrixExpr>
1355 public:
1356  static const bool bAlias = MatrixExpr::bAlias;
1357  static const index_type iNumRows = MatrixExpr::iNumRows;
1358  static const index_type iNumCols = MatrixExpr::iNumCols;
1359  typedef typename ScalarUnaryFunc::ScalarType ScalarType;
1360  typedef typename ScalarUnaryFunc::ExpressionType ExpressionType;
1363 
1364  MatrixMatrixUnaryExpr(const MatrixExpr& u)
1365  :oU(u) {
1366  MATVEC_ASSERT((iNumRows == oU.iGetNumRows()) || (iNumRows == DYNAMIC_SIZE && oU.iGetNumRows() >= 0));
1367  MATVEC_ASSERT((iNumCols == oU.iGetNumCols()) || (iNumCols == DYNAMIC_SIZE && oU.iGetNumCols() >= 0));
1368  }
1369 
1371  MATVEC_ASSERT(i >= 1);
1372  MATVEC_ASSERT(i <= iGetNumRows());
1373  MATVEC_ASSERT(j >= 1);
1374  MATVEC_ASSERT(j <= iGetNumCols());
1375  return ScalarUnaryFunc::f(oU(i, j));
1376  }
1377 
1379  MATVEC_ASSERT((iNumRows == oU.iGetNumRows()) || (iNumRows == DYNAMIC_SIZE && oU.iGetNumRows() >= 0));
1380  return oU.iGetNumRows();
1381  }
1382 
1384  MATVEC_ASSERT(iNumCols == oU.iGetNumCols());
1385  return oU.iGetNumCols();
1386  }
1387 
1389  MATVEC_ASSERT(iRow >= 1);
1390  MATVEC_ASSERT(iRow <= iGetNumRows());
1391  return RowVectorType(*this, iRow);
1392  }
1393 
1395  MATVEC_ASSERT(iCol >= 1);
1396  MATVEC_ASSERT(iCol <= iGetNumCols());
1397  return ColumnVectorType(*this, iCol);
1398  }
1399 
1400  template <typename ScalarType2>
1401  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
1402  return oU.bHaveReferenceTo(pFirst, pLast);
1403  }
1404 
1405 private:
1406  const MatrixExpr oU;
1407 };
1408 
1409 template <typename MatrixType, bool ALIAS=false>
1411 public:
1412  static const bool bAlias = ALIAS;
1413  static const index_type iNumRows = MatrixType::iNumRows;
1414  static const index_type iNumCols = MatrixType::iNumCols;
1415  typedef typename MatrixType::ScalarType ScalarType;
1416  typedef typename MatrixType::ExpressionType ExpressionType;
1417  typedef typename MatrixType::RowVectorType RowVectorType;
1418  typedef typename MatrixType::ColumnVectorType ColumnVectorType;
1419 
1420  MatrixDirectExpr(const MatrixType& A)
1421  :A(A) {
1422  MATVEC_ASSERT((iNumRows == A.iGetNumRows()) || (iNumRows == DYNAMIC_SIZE && A.iGetNumRows() >= 0));
1423  MATVEC_ASSERT((iNumCols == A.iGetNumCols()) || (iNumCols == DYNAMIC_SIZE && A.iGetNumCols() >= 0));
1424  }
1425 
1427  MATVEC_ASSERT(i >= 1);
1428  MATVEC_ASSERT(i <= iGetNumRows());
1429  MATVEC_ASSERT(j >= 1);
1430  MATVEC_ASSERT(j <= iGetNumCols());
1431  return A(i, j);
1432  }
1433 
1435  MATVEC_ASSERT((iNumRows == A.iGetNumRows()) || (iNumRows == DYNAMIC_SIZE && A.iGetNumRows() >= 0));
1436  return A.iGetNumRows();
1437  }
1438 
1440  MATVEC_ASSERT((iNumCols == A.iGetNumCols()) || (iNumCols == DYNAMIC_SIZE && A.iGetNumCols() >= 0));
1441  return A.iGetNumCols();
1442  }
1443 
1445  MATVEC_ASSERT(iRow >= 1);
1446  MATVEC_ASSERT(iRow <= iGetNumRows());
1447  return A.GetRow(iRow);
1448  }
1449 
1451  MATVEC_ASSERT(iCol >= 1);
1452  MATVEC_ASSERT(iCol <= iGetNumCols());
1453  return A.GetCol(iCol);
1454  }
1455 
1456  template <typename ScalarType2>
1457  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
1458  return A.bHaveReferenceTo(pFirst, pLast);
1459  }
1460 
1461 private:
1462  const MatrixType& A;
1463 };
1464 
1466 public:
1467  static const bool bAlias = false;
1468  static const index_type iNumRows = 3;
1469  static const index_type iNumCols = 3;
1470 
1473 
1476 
1478  :A(A) {
1479 
1480  }
1481 
1483  MATVEC_ASSERT(i >= 1);
1484  MATVEC_ASSERT(i <= iGetNumRows());
1485  MATVEC_ASSERT(j >= 1);
1486  MATVEC_ASSERT(j <= iGetNumCols());
1487  return A(i, j);
1488  }
1489 
1491  return iNumRows;
1492  }
1493 
1495  return iNumCols;
1496  }
1497 
1499  MATVEC_ASSERT(iRow >= 1);
1500  MATVEC_ASSERT(iRow <= iGetNumRows());
1501  --iRow;
1502  return RowVectorType(A.pGetMat() + iRow, iGetNumCols(), iGetNumRows());
1503  }
1504 
1506  MATVEC_ASSERT(iCol >= 1);
1507  MATVEC_ASSERT(iCol <= iGetNumCols());
1508  --iCol;
1509  return ColumnVectorType(A.pGetMat() + iCol * iNumRows, iGetNumRows(), 1);
1510  }
1511 
1512  template <typename ScalarType2>
1513  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
1514  return bArrayOverlap(A.pGetMat(), A.pGetMat() + (iNumRows - 1) * (iNumCols - 1), pFirst, pLast);
1515  }
1516 
1517 private:
1518  const Mat3x3& A;
1519 };
1520 
1522 public:
1523  static const bool bAlias = false;
1524  static const index_type iNumRows = 3;
1526 
1527  typedef void RowVectorType;
1528  typedef void ColumnVectorType;
1531 
1533  :A(A) {
1534 
1535  }
1536 
1538  MATVEC_ASSERT(i >= 1);
1539  MATVEC_ASSERT(i <= iGetNumRows());
1540  MATVEC_ASSERT(j >= 1);
1541  MATVEC_ASSERT(j <= iGetNumCols());
1542  return A(i, j);
1543  }
1544 
1546  return iNumRows;
1547  }
1548 
1550  return A.iGetNumCols();
1551  }
1552 
1553  template <typename ScalarType2>
1554  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
1555  return bArrayOverlap(&A(1, 1), &A(A.iGetNumRows(), A.iGetNumCols()), pFirst, pLast);
1556  }
1557 
1558  private:
1559  const Mat3xN& A;
1560 };
1561 
1563 public:
1564  static const bool bAlias = false;
1567 
1568  typedef void RowVectorType;
1569  typedef void ColumnVectorType;
1572 
1574  :A(A) {
1575 
1576  }
1577 
1579  MATVEC_ASSERT(i >= 1);
1580  MATVEC_ASSERT(i <= iGetNumRows());
1581  MATVEC_ASSERT(j >= 1);
1582  MATVEC_ASSERT(j <= iGetNumCols());
1583  return A(i, j);
1584  }
1585 
1587  return A.iGetNumRows();
1588  }
1589 
1591  return A.iGetNumCols();
1592  }
1593 
1594  template <typename ScalarType2>
1595  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
1596  return bArrayOverlap(&A(1, 1), &A(A.iGetNumRows(), A.iGetNumCols()), pFirst, pLast);
1597  }
1598 
1599  private:
1600  const MatNxN& A;
1601  };
1602 
1603  template <typename InitClass, typename T, index_type N_rows, index_type N_cols>
1604  class MatrixInit: public InitClass {
1605  public:
1606  template <typename InitArg>
1607  explicit MatrixInit(const InitArg& v)
1608  :InitClass(v) {
1609 
1610  }
1611 
1612  template <typename InitArg1, typename InitArg2>
1613  explicit MatrixInit(const InitArg1& v1, const InitArg2& a2)
1614  :InitClass(v1, a2) {
1615 
1616  }
1617  };
1618 
1619  template <typename InitClass, typename T, index_type N_rows>
1620  class VectorInit: public InitClass {
1621  public:
1622  template <typename InitArg>
1623  explicit VectorInit(const InitArg& R)
1624  :InitClass(R) {
1625 
1626  }
1627  };
1628 
1629  template <typename T, index_type N_rows>
1630  class VectorData {
1631  public:
1632  explicit VectorData(index_type N, bool bZeroInit) {
1633  MATVEC_ASSERT(N == N_rows);
1634 
1635  if (bZeroInit) {
1636  ZeroInit(&rgData[0], &rgData[N_rows]);
1637  }
1638  }
1639 
1640  VectorData(const VectorData& oData) {
1641  std::copy(&oData.rgData[0], &oData.rgData[N_rows], &rgData[0]);
1642  }
1643 
1644  template <typename U>
1646  std::copy(&oData.rgData[0], &oData.rgData[N_rows], &rgData[0]);
1647  }
1648 
1650  MATVEC_ASSERT(i >= 0);
1651  MATVEC_ASSERT(i < N_rows);
1652  return rgData[i];
1653  }
1654 
1655  const T& operator[](integer i) const {
1656  MATVEC_ASSERT(i >= 0);
1657  MATVEC_ASSERT(i < N_rows);
1658  return rgData[i];
1659  }
1660 
1662  return N_rows;
1663  }
1664 
1665  void Resize(index_type iNumRows) {
1666  MATVEC_ASSERT(iNumRows == N_rows);
1667  }
1668 
1669  const T* pGetFirstElem() const {
1670  return &rgData[0];
1671  }
1672 
1673  const T* pGetLastElem() const {
1674  return &rgData[N_rows - 1];
1675  }
1676 
1677  private:
1678  T rgData[N_rows];
1679  };
1680 
1681  template <typename T>
1683  public:
1684  explicit VectorData(index_type N, bool)
1685  :rgData(N) {
1686  }
1687 
1688  VectorData(const VectorData& oVec)
1689  :rgData(oVec.rgData) {
1690  }
1691 
1693  MATVEC_ASSERT(i >= 0);
1694  MATVEC_ASSERT(i < index_type(rgData.size()));
1695  return rgData[i];
1696  }
1697 
1698  const T& operator[](integer i) const {
1699  MATVEC_ASSERT(i >= 0);
1700  MATVEC_ASSERT(i < index_type(rgData.size()));
1701  return rgData[i];
1702  }
1703 
1705  return rgData.size();
1706  }
1707 
1708  void Resize(index_type iNumRows) {
1709  if (iNumRows != iGetNumRows()) {
1710  rgData.resize(iNumRows);
1711  }
1712  }
1713 
1714  const T* pGetFirstElem() const {
1715  return &rgData.front();
1716  }
1717 
1718  const T* pGetLastElem() const {
1719  return &rgData.back();
1720  }
1721 
1722  private:
1723  std::vector<T, GradientAllocator<T> > rgData;
1724  };
1725 
1726 
1727  template <typename T, index_type N_rows, index_type N_cols>
1729  {
1730  public:
1731  MatrixData(index_type iNumRows, index_type iNumCols, bool bZeroInit) {
1732  MATVEC_ASSERT(iNumRows == N_rows);
1733  MATVEC_ASSERT(iNumCols == N_cols);
1734 
1735  if (bZeroInit) {
1736  ZeroInit(&rgData[0][0], &rgData[iNumRows - 1][iNumCols - 1] + 1);
1737  }
1738  }
1739 
1740  index_type iGetNumRows() const { return N_rows; }
1741  index_type iGetNumCols() const { return N_cols; }
1742 
1744  MATVEC_ASSERT(i >= 1 && i <= N_rows);
1745  MATVEC_ASSERT(j >= 1 && j <= N_cols);
1746  return rgData[i - 1][j - 1];
1747  }
1748 
1749  const T& operator() (index_type i, index_type j) const {
1750  MATVEC_ASSERT(i >= 1 && i <= N_rows);
1751  MATVEC_ASSERT(j >= 1 && j <= N_cols);
1752  return rgData[i - 1][j - 1];
1753  }
1754 
1755  void Resize(index_type iRows, index_type iCols) {
1756  MATVEC_ASSERT(iRows == N_rows);
1757  MATVEC_ASSERT(iCols == N_cols);
1758  }
1759 
1760  const T* pGetFirstElem() const {
1761  return &rgData[0][0];
1762  }
1763 
1764  const T* pGetLastElem() const {
1765  return &rgData[N_rows - 1][N_cols - 1];
1766  }
1767 
1768  private:
1769  T rgData[N_rows][N_cols];
1770  };
1771 
1772  template <typename T>
1774  {
1775  public:
1776  MatrixData(index_type iNumRows, index_type iNumCols, bool bZeroInit)
1777  :iNumRows(iNumRows),
1778  iNumCols(iNumCols),
1779  rgData(iNumRows * iNumCols, bZeroInit) {
1780  }
1781 
1782  index_type iGetNumRows() const { return iNumRows; }
1783  index_type iGetNumCols() const { return iNumCols; }
1784 
1786  MATVEC_ASSERT(i >= 1 && i <= iNumRows);
1787  MATVEC_ASSERT(j >= 1 && j <= iNumCols);
1788  return rgData[(i - 1) * iNumCols + (j - 1)];
1789  }
1790 
1791  const T& operator() (index_type i, index_type j) const {
1792  MATVEC_ASSERT(i >= 1 && i <= iNumRows);
1793  MATVEC_ASSERT(j >= 1 && j <= iNumCols);
1794  return rgData[(i - 1) * iNumCols + (j - 1)];
1795  }
1796 
1797  void Resize(index_type iRows, index_type iCols) {
1798  if (iRows != iNumRows || iCols != iCols) {
1799  rgData.Resize(iRows * iCols);
1800  iNumRows = iRows;
1801  iNumCols = iCols;
1802  }
1803  }
1804 
1805  const T* pGetFirstElem() const {
1806  return rgData.pGetFirstElem();
1807  }
1808 
1809  const T* pGetLastElem() const {
1810  return rgData.pGetLastElem();
1811  }
1812  private:
1813  integer iNumRows, iNumCols;
1815  };
1816 
1817  template <typename T, index_type N_rows>
1818  class MatrixData<T, N_rows, DYNAMIC_SIZE>
1819  {
1820  public:
1821  MatrixData(index_type iNumRows, index_type iNumCols, bool bZeroInit)
1822  :iNumCols(iNumCols),
1823  rgData(N_rows * iNumCols, bZeroInit) {
1824  MATVEC_ASSERT(iNumRows == N_rows);
1825  }
1826 
1827  index_type iGetNumRows() const { return N_rows; }
1828  index_type iGetNumCols() const { return iNumCols; }
1829 
1831  MATVEC_ASSERT(i >= 1 && i <= N_rows);
1832  MATVEC_ASSERT(j >= 1 && j <= iNumCols);
1833  return rgData[(i - 1) * iNumCols + (j - 1)];
1834  }
1835 
1836  const T& operator() (index_type i, index_type j) const {
1837  MATVEC_ASSERT(i >= 1 && i <= N_rows);
1838  MATVEC_ASSERT(j >= 1 && j <= iNumCols);
1839  return rgData[(i - 1) * iNumCols + (j - 1)];
1840  }
1841 
1842  void Resize(index_type iRows, index_type iCols) {
1843  MATVEC_ASSERT(iRows == N_rows);
1844 
1845  if (iCols != iNumCols) {
1846  rgData.Resize(iRows * iCols);
1847  iNumCols = iCols;
1848  }
1849  }
1850 
1851  const T* pGetFirstElem() const {
1852  return rgData.pGetFirstElem();
1853  }
1854 
1855  const T* pGetLastElem() const {
1856  return rgData.pGetLastElem();
1857  }
1858  private:
1861  };
1862 
1863  template <typename T, index_type N_cols>
1864  class MatrixData<T, DYNAMIC_SIZE, N_cols>
1865  {
1866  public:
1867  MatrixData(index_type iNumRows, index_type iNumCols, bool bZeroInit)
1868  :iNumRows(iNumRows),
1869  rgData(iNumRows * N_cols, bZeroInit) {
1870  MATVEC_ASSERT(iNumCols == N_cols);
1871  }
1872 
1873  index_type iGetNumRows() const { return iNumRows; }
1874  index_type iGetNumCols() const { return N_cols; }
1875 
1877  MATVEC_ASSERT(i >= 1 && i <= iNumRows);
1878  MATVEC_ASSERT(j >= 1 && j <= N_cols);
1879  return rgData[(i - 1) * N_cols + (j - 1)];
1880  }
1881 
1882  const T& operator() (index_type i, index_type j) const {
1883  MATVEC_ASSERT(i >= 1 && i <= iNumRows);
1884  MATVEC_ASSERT(j >= 1 && j <= N_cols);
1885  return rgData[(i - 1) * N_cols + (j - 1)];
1886  }
1887 
1888  void Resize(index_type iRows, index_type iCols) {
1889  MATVEC_ASSERT(iCols == N_cols);
1890 
1891  if (iRows != iNumRows) {
1892  rgData.Resize(iRows * iCols);
1893  iNumRows = iRows;
1894  }
1895  }
1896 
1897  const T* pGetFirstElem() const {
1898  return rgData.pGetFirstElem();
1899  }
1900 
1901  const T* pGetLastElem() const {
1902  return rgData.pGetLastElem();
1903  }
1904  private:
1907  };
1908 
1909  template <typename T, index_type N_rows, index_type N_cols>
1910  class Matrix {
1911  public:
1912  static const index_type iNumRows = N_rows;
1913  static const index_type iNumCols = N_cols;
1918 
1921 
1924  }
1925 
1927  :rgMat(iRows, iCols, true) {
1928  }
1929 
1930  Matrix(const T& A11, const T& A21, const T& A12, const T& A22)
1932  typedef typename IndexCheck<iNumRows - 2>::CheckType check_iNumRows;
1933  typedef typename IndexCheck<iNumCols - 2>::CheckType check_iNumCols;
1934 
1935  (*this)(1, 1) = A11;
1936  (*this)(2, 1) = A21;
1937  (*this)(1, 2) = A12;
1938  (*this)(2, 2) = A22;
1939  }
1940 
1941  Matrix(const Matrix& A)
1942  :rgMat(A.rgMat) {
1943  }
1944 
1945  explicit inline Matrix(const Mat3x3& A);
1946 
1947  template <typename InitClass>
1949  :rgMat(func.iGetNumRows(), func.iGetNumCols(), false) {
1950  func.Initialize(*this);
1951  }
1952 
1953  template <typename Expression>
1955  :rgMat(A.iGetNumRows(), A.iGetNumCols(), false) {
1956  // No aliases are possible because the object did not exist before
1957  using namespace MatVecHelp;
1958  ApplyMatrixFuncNoAlias(A, Assign());
1959  }
1960 
1961  template <typename T2>
1963  :rgMat(A.iGetNumRows(), A.iGetNumCols(), false) {
1964  Copy(A, pDofMap);
1965  }
1966 
1967  template <typename T2>
1968  void Copy(const Matrix<T2, N_rows, N_cols>& A, LocalDofMap* pDofMap) {
1969  rgMat.Resize(A.iGetNumRows(), A.iGetNumCols());
1970 
1973 
1974  for (index_type i = 1; i <= A.iGetNumRows(); ++i) {
1975  for (index_type j = 1; j <= A.iGetNumCols(); ++j) {
1976  grad::Copy((*this)(i, j), A(i, j), pDofMap);
1977  }
1978  }
1979  }
1980 
1981  void Reset() {
1982  for (index_type i = 1; i <= iGetNumRows(); ++i) {
1983  for (index_type j = 1; j <= iGetNumCols(); ++j) {
1984  grad::Reset((*this)(i, j));
1985  }
1986  }
1987  }
1988 
1989  void Resize(index_type iRows, index_type iCols) {
1990  rgMat.Resize(iRows, iCols);
1991  }
1992 
1993  template <typename Expression>
1995 
1996  rgMat.Resize(A.iGetNumRows(), A.iGetNumCols());
1997 
1998  using namespace MatVecHelp;
1999 
2000  ApplyMatrixFunc<Assign>(A);
2001 
2002  return *this;
2003  }
2004 
2005  template <typename T_Rhs>
2009  using namespace MatVecHelp;
2010 
2011  ApplyMatrixFunc<Add>(MatrixExpression<MatrixDirectExpr<Matrix<T_Rhs, N_rows, N_cols> >, N_rows, N_cols>(A));
2012 
2013  return *this;
2014  }
2015 
2016  template <typename T_Rhs>
2020 
2021  using namespace MatVecHelp;
2022 
2023  ApplyMatrixFunc<Sub>(MatrixExpression<MatrixDirectExpr<Matrix<T_Rhs, N_rows, N_cols> >, N_rows, N_cols>(A));
2024 
2025  return *this;
2026  }
2027 
2028  template <typename Expression>
2032  using namespace MatVecHelp;
2033 
2034  ApplyMatrixFunc<Add>(A);
2035 
2036  return *this;
2037  }
2038 
2039  template <typename Expression>
2043  using namespace MatVecHelp;
2044 
2045  ApplyMatrixFunc<Sub>(A);
2046 
2047  return *this;
2048  }
2049 
2050  template <typename T_Rhs>
2051  Matrix& operator*=(const T_Rhs& a) {
2052  using namespace MatVecHelp;
2053 
2054  ApplyScalarFunc<Mul>(a);
2055 
2056  return *this;
2057  }
2058 
2059  template <typename T_Rhs>
2060  Matrix& operator/=(const T_Rhs& a) {
2061  using namespace MatVecHelp;
2062 
2063  ApplyScalarFunc<Div>(a);
2064 
2065  return *this;
2066  }
2067 
2068  template <typename ScalarExpression>
2070  using namespace MatVecHelp;
2071 
2072  ApplyScalarFunc<Mul>(ScalarType(a));
2073 
2074  return *this;
2075  }
2076 
2077  template <typename ScalarExpression>
2079  using namespace MatVecHelp;
2080 
2081  ApplyScalarFunc<Div>(ScalarType(a));
2082 
2083  return *this;
2084  }
2085 
2086  inline Matrix& operator=(const Mat3x3& A);
2087 
2088  const ScalarType& operator()(index_type iRow, index_type iCol) const {
2089  MATVEC_ASSERT(iRow >= 1);
2090  MATVEC_ASSERT(iRow <= iGetNumRows());
2091  MATVEC_ASSERT(iCol >= 1);
2092  MATVEC_ASSERT(iCol <= iGetNumCols());
2093  return rgMat(iRow, iCol);
2094  }
2095 
2097  MATVEC_ASSERT(iRow >= 1);
2098  MATVEC_ASSERT(iRow <= iGetNumRows());
2099  MATVEC_ASSERT(iCol >= 1);
2100  MATVEC_ASSERT(iCol <= iGetNumCols());
2101  return rgMat(iRow, iCol);
2102  }
2103 
2105  GetRow(index_type iRow) const {
2106  MATVEC_ASSERT(iRow >= 1);
2107  MATVEC_ASSERT(iRow <= iGetNumRows());
2108  return RowVectorType(&rgMat(iRow, 1), iGetNumCols(), 1);
2109  }
2110 
2112  GetCol(index_type iCol) const {
2113  MATVEC_ASSERT(iCol >= 1);
2114  MATVEC_ASSERT(iCol <= iGetNumCols());
2115  return ColumnVectorType(&rgMat(1, iCol), iGetNumRows(), iGetNumCols());
2116  }
2117 
2118  index_type iGetNumRows() const { return rgMat.iGetNumRows(); }
2119  index_type iGetNumCols() const { return rgMat.iGetNumCols(); }
2120  const ScalarType* pGetMat() const { return pGetFirstElem(); }
2121 
2122  template <typename ScalarType2>
2123  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
2124  return bArrayOverlap(pGetFirstElem(),
2125  pGetLastElem(),
2126  pFirst,
2127  pLast);
2128  }
2129 
2130 private:
2131  const ScalarType* pGetFirstElem() const {
2132  return rgMat.pGetFirstElem();
2133  }
2134 
2135  const ScalarType* pGetLastElem() const {
2136  return rgMat.pGetLastElem();
2137  }
2138 
2141 
2142  template <typename Func, typename U>
2143  void ApplyScalarFunc(const U& a) {
2144  for (index_type i = 1; i <= iGetNumRows(); ++i) {
2145  for (index_type j = 1; j <= iGetNumCols(); ++j) {
2146  Func::Eval((*this)(i, j), a);
2147  }
2148  }
2149  }
2150 
2151  template <typename Func, typename Expression>
2153  using namespace MatVecHelp;
2154 
2155  ApplyAliasHelperMatrix<MatrixExpression<Expression, N_rows, N_cols>::bAlias>::ApplyMatrixFunc(*this, A, Func());
2156  }
2157 
2158  template <typename Func, typename Expression>
2160  MATVEC_ASSERT(!A.bHaveReferenceTo(pGetFirstElem(), pGetLastElem()));
2161 
2162  rgMat.Resize(A.iGetNumRows(), A.iGetNumCols());
2163 
2166 
2167  for (index_type i = 1; i <= iGetNumRows(); ++i) {
2168  for (index_type j = 1; j <= iGetNumCols(); ++j) {
2169  Func::Eval((*this)(i, j), A(i, j));
2170  }
2171  }
2172  }
2173 
2174  template <typename Func, typename Expression>
2177  }
2178 
2179 private:
2181 };
2182 
2183 template <typename T, index_type N_rows, index_type N_cols>
2187  }
2188 
2189  inline MatrixExpression<Mat3x3DirectExpr, 3, 3>
2190  Direct(const Mat3x3& A) {
2192  }
2193 
2194  inline MatrixExpression<Mat3xNDirectExpr, 3, DYNAMIC_SIZE>
2195  Direct(const Mat3xN& A) {
2197  }
2198 
2199  inline MatrixExpression<MatNxNDirectExpr, DYNAMIC_SIZE, DYNAMIC_SIZE>
2200  Direct(const MatNxN& A) {
2202  }
2203 
2204  template <typename T, index_type N_rows, index_type N_cols>
2205 inline MatrixExpression<TransposedMatrix<MatrixDirectExpr<Matrix<T, N_rows, N_cols> > >, N_cols, N_rows>
2208 }
2209 
2210  inline MatrixExpression<TransposedMatrix<Mat3xNDirectExpr>, DYNAMIC_SIZE, 3>
2211  Transpose(const Mat3xN& A) {
2213 }
2214 
2215  inline MatrixExpression<TransposedMatrix<MatNxNDirectExpr>, DYNAMIC_SIZE, DYNAMIC_SIZE>
2216  Transpose(const MatNxN& A) {
2218 }
2219 
2220  template <typename MatrixExpr, index_type N_rows, index_type N_cols>
2221  inline MatrixExpression<TransposedMatrix<MatrixExpr>, N_cols, N_rows>
2223  return MatrixExpression<TransposedMatrix<MatrixExpr>, N_cols, N_rows>(A);
2224 }
2225 
2226 template <typename T, index_type N_rows, index_type N_cols>
2227 inline MatrixExpression<MatrixDirectExpr<Matrix<T, N_rows, N_cols>, true>, N_rows, N_cols>
2229  return MatrixExpression<MatrixDirectExpr<Matrix<T, N_rows, N_cols>, true>, N_rows, N_cols>(A);
2230 }
2231 
2232 template <typename T, index_type N_rows, index_type N_cols>
2234  :rgMat(N_rows, N_cols, false) {
2235  using namespace MatVecHelp;
2236 
2237  ApplyMatrixFunc<Assign>(Direct(A));
2238 }
2239 
2240 template <typename T, index_type N_rows, index_type N_cols>
2245 
2246  rgMat.Resize(3, 3);
2247 
2248  using namespace MatVecHelp;
2249 
2250  ApplyMatrixFunc<Assign>(Direct(A));
2251 
2252  return *this;
2253 }
2254 
2255 template <typename T, index_type N_rows>
2256 class Vector {
2257 public:
2258  static const index_type iNumRows = N_rows;
2259  static const index_type iInitNumRows = iNumRows == DYNAMIC_SIZE ? 0 : iNumRows;
2262 
2264  :rgVec(iInitNumRows, true) {
2265  }
2266 
2267  explicit Vector(index_type N)
2268  :rgVec(N, true) {
2269  }
2270 
2271  Vector(const Vector& oVec)
2272  :rgVec(oVec.rgVec) {
2273  }
2274 
2275  Vector(const T& v1, const T& v2)
2276  :rgVec(iInitNumRows, false) {
2277  typedef typename IndexCheck<iNumRows - 2>::CheckType check_iNumRows;
2278 
2279  (*this)(1) = v1;
2280  (*this)(2) = v2;
2281  }
2282 
2283  template <typename Expr1, typename Expr2>
2285  :rgVec(iNumRows, false) {
2286  typedef typename IndexCheck<iNumRows - 2>::CheckType check_iNumRows;
2287 
2288  (*this)(1) = v1;
2289  (*this)(2) = v2;
2290  }
2291 
2292  Vector(const T& v1, const T& v2, const T& v3)
2293  :rgVec(iNumRows, false) {
2294  typedef typename IndexCheck<iNumRows - 3>::CheckType check_iNumRows;
2295 
2296  (*this)(1) = v1;
2297  (*this)(2) = v2;
2298  (*this)(3) = v3;
2299  }
2300 
2301  template <typename Expr1, typename Expr2, typename Expr3>
2303  :rgVec(iNumRows, false) {
2304  typedef typename IndexCheck<iNumRows - 3>::CheckType check_iNumRows;
2305 
2306  (*this)(1) = v1;
2307  (*this)(2) = v2;
2308  (*this)(3) = v3;
2309  }
2310 
2311  explicit inline Vector(const Vec3& v);
2312 
2313  template <typename Expression>
2315  :rgVec(f.iGetNumRows(), false) {
2316  using namespace MatVecHelp;
2317 
2318  ApplyMatrixFuncNoAlias(f, Assign());
2319  }
2320 
2321  template <typename InitClass>
2323  :rgVec(func.iGetNumRows(), false) {
2324  func.Initialize(*this);
2325  }
2326 
2327  template <typename T2>
2328  explicit Vector(const Vector<T2, N_rows>& v)
2329  :rgVec(v.iGetNumRows(), false) {
2330 
2332 
2333  for (index_type i = 1; i <= v.iGetNumRows(); ++i) {
2334  grad::Convert((*this)(i), v(i));
2335  }
2336  }
2337 
2338  template <typename T2>
2340  :rgVec(v.iGetNumRows(), false) {
2341 
2342  Copy(v, pDofMap);
2343  }
2344 
2345  template <typename T2>
2346  void Copy(const Vector<T2, N_rows>& v, LocalDofMap* pDofMap) {
2347  rgVec.Resize(v.iGetNumRows());
2348 
2350 
2351  for (index_type i = 1; i <= iGetNumRows(); ++i) {
2352  grad::Copy((*this)(i), v(i), pDofMap);
2353  }
2354  }
2355 
2356  void Reset() {
2357  for (index_type i = 1; i <= iGetNumRows(); ++i) {
2358  grad::Reset((*this)(i));
2359  }
2360  }
2361 
2363  rgVec.Resize(iNumRows);
2364  }
2365 
2366  template <typename Expression>
2368  rgVec.Resize(f.iGetNumRows());
2369 
2370  using namespace MatVecHelp;
2371 
2372  ApplyMatrixFunc<Assign>(f);
2373 
2374  return *this;
2375  }
2376 
2377  template <typename T_Rhs>
2380  using namespace MatVecHelp;
2381 
2382  ApplyMatrixFunc<Add>(VectorExpression<VectorDirectExpr<Vector<T_Rhs, N_rows> >, N_rows>(v));
2383 
2384  return *this;
2385  }
2386 
2387  Vector& operator+=(const Vec3& v) {
2388  MATVEC_ASSERT(iGetNumRows() == 3);
2390 
2391  using namespace MatVecHelp;
2392 
2393  ApplyMatrixFunc<Add>(VectorExpression<Vec3DirectExpr, N_rows>(v));
2394 
2395  return *this;
2396  }
2397 
2398  template <typename T_Rhs>
2401  using namespace MatVecHelp;
2402 
2403  ApplyMatrixFunc<Sub>(VectorExpression<VectorDirectExpr<Vector<T_Rhs, N_rows> >, N_rows>(v));
2404 
2405  return *this;
2406  }
2407 
2408  template <typename Expression>
2411  using namespace MatVecHelp;
2412 
2413  ApplyMatrixFunc<Add>(f);
2414 
2415  return *this;
2416  }
2417 
2418  template <typename Expression>
2421  using namespace MatVecHelp;
2422 
2423  ApplyMatrixFunc<Sub>(f);
2424 
2425  return *this;
2426  }
2427 
2428  template <typename T_Rhs>
2429  Vector& operator*=(const T_Rhs& g) {
2430  using namespace MatVecHelp;
2431 
2432  ApplyScalarFunc<Mul>(g);
2433 
2434  return *this;
2435  }
2436 
2437  template <typename T_Rhs>
2438  Vector& operator/=(const T_Rhs& g) {
2439  using namespace MatVecHelp;
2440 
2441  ApplyScalarFunc<Div>(g);
2442 
2443  return *this;
2444  }
2445 
2446  template <typename ScalarExpression>
2448  using namespace MatVecHelp;
2449 
2450  ApplyScalarFunc<Mul>(ScalarType(f));
2451 
2452  return *this;
2453  }
2454 
2455  template <typename ScalarExpression>
2457  using namespace MatVecHelp;
2458 
2459  ApplyScalarFunc<Div>(ScalarType(f));
2460 
2461  return *this;
2462  }
2463 
2464  inline Vector& operator=(const Vec3& v);
2465 
2466  const ScalarType& operator()(index_type iRow) const {
2467  --iRow; // Row index is 1-based for compatibility reasons
2468  MATVEC_ASSERT(iRow >= 0);
2469  MATVEC_ASSERT(iRow < iGetNumRows());
2470  return rgVec[iRow];
2471  }
2472 
2474  --iRow; // Row index is 1-based for compatibility reasons
2475  MATVEC_ASSERT(iRow >= 0);
2476  MATVEC_ASSERT(iRow < iGetNumRows());
2477  return rgVec[iRow];
2478  }
2479 
2481  const index_type iRows = rgVec.iGetNumRows();
2482 
2483  MATVEC_ASSERT((N_rows == DYNAMIC_SIZE) || (iRows == N_rows));
2484 
2485  return iRows;
2486  }
2487 
2488  ScalarType* pGetVec() { return &rgVec[0]; }
2489  const ScalarType* pGetVec() const { return pGetFirstElem(); }
2490 
2491  template <typename ScalarType2>
2492  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
2493  return bArrayOverlap(pGetFirstElem(),
2494  pGetLastElem(),
2495  pFirst,
2496  pLast);
2497  }
2498 
2499 private:
2500  const ScalarType* pGetFirstElem() const {
2501  return rgVec.pGetFirstElem();
2502  }
2503 
2504  const ScalarType* pGetLastElem() const {
2505  return rgVec.pGetLastElem();
2506  }
2507 
2510 
2511  template <typename Func, typename U>
2512  void ApplyScalarFunc(const U& a) {
2513  for (index_type i = 1; i <= iGetNumRows(); ++i) {
2514  Func::Eval((*this)(i), a);
2515  }
2516  }
2517 
2518  template <typename Func, typename Expression>
2520  using namespace MatVecHelp;
2521 
2522  ApplyAliasHelperMatrix<VectorExpression<Expression, N_rows>::bAlias>::ApplyMatrixFunc(*this, A, Func());
2523  }
2524 
2525  template <typename Func, typename Expression>
2527  MATVEC_ASSERT((N_rows == A.iGetNumRows()) || (N_rows == DYNAMIC_SIZE && A.iGetNumRows() >= 0));
2529  MATVEC_ASSERT(!A.bHaveReferenceTo(pGetFirstElem(), pGetLastElem()));
2530 
2531  for (index_type i = 1; i <= iGetNumRows(); ++i) {
2532  Func::Eval((*this)(i), A(i));
2533  }
2534  }
2535 
2536  template <typename Func, typename Expression>
2539  }
2540 
2541 private:
2543 };
2544 
2545 template <typename T, index_type N_rows>
2549 }
2550 
2551 inline VectorExpression<Vec3DirectExpr, 3>
2552 Direct(const Vec3& v) {
2554 }
2555 
2556 template <typename T, index_type N_rows>
2557 inline VectorExpression<VectorDirectExpr<Vector<T, N_rows>, true>, N_rows>
2559  return VectorExpression<VectorDirectExpr<Vector<T, N_rows>, true>, N_rows>(v);
2560 }
2561 
2562 template <typename T, index_type N_rows>
2563 inline
2565  :rgVec(iNumRows, false) {
2566  using namespace MatVecHelp;
2567 
2568  typedef typename IndexCheck<iNumRows - 3>::CheckType check_iNumRows;
2569 
2570  ApplyMatrixFunc<Assign>(Direct(v));
2571 }
2572 
2573 template <typename T, index_type N_rows>
2574 inline Vector<T, N_rows>&
2576  using namespace MatVecHelp;
2577 
2578  ApplyMatrixFunc<Assign>(Direct(v));
2579 
2580  return *this;
2581 }
2582 
2583 template <index_type iStartIndex, index_type iEndIndex, typename T, index_type N_rows>
2586  MATVEC_ASSERT(iStartIndex >= 1);
2587  MATVEC_ASSERT(iEndIndex <= v.iGetNumRows());
2589 }
2590 
2591 template <index_type iStartIndex, index_type iEndIndex, typename VectorExpr, index_type N_rows>
2592 VectorExpression<SubVectorExpr<iStartIndex, iEndIndex, VectorExpr> , iEndIndex - iStartIndex + 1>
2594  MATVEC_ASSERT(iStartIndex >= 1);
2595  MATVEC_ASSERT(iEndIndex <= v.iGetNumRows());
2596  return VectorExpression<SubVectorExpr<iStartIndex, iEndIndex, VectorExpr> , iEndIndex - iStartIndex + 1>(v);
2597 }
2598 
2599 template <index_type iRowStart, index_type iRowEnd, index_type iColStart, index_type iColEnd, typename T, index_type N_rows, index_type N_cols>
2600 MatrixExpression<SubMatrixExpr<iRowStart, iRowEnd, iColStart, iColEnd, MatrixExpression<MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, N_rows, N_cols> >, iRowEnd - iRowStart + 1, iColEnd - iColStart + 1>
2602  return MatrixExpression<SubMatrixExpr<iRowStart, iRowEnd, iColStart, iColEnd, MatrixExpression<MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, N_rows, N_cols> >, iRowEnd - iRowStart + 1, iColEnd - iColStart + 1>(Direct(A));
2603 }
2604 
2605 template <index_type iRowStart, index_type iRowEnd, index_type iColStart, index_type iColEnd, typename MatrixExpr, index_type N_rows, index_type N_cols>
2606 MatrixExpression<SubMatrixExpr<iRowStart, iRowEnd, iColStart, iColEnd, MatrixExpression<MatrixExpr, N_rows, N_cols> >, iRowEnd - iRowStart + 1, iColEnd - iColStart + 1>
2609 }
2610 
2611 template <typename T, typename VectorExpr>
2613 public:
2614  static const index_type iNumRows = 3;
2615  static const index_type iNumCols = 3;
2616 
2618  :v(v), d(d) {
2619 
2620  }
2621 
2622  void Initialize(Matrix<T, 3, 3>& A) const {
2623  const index_type x = 1, y = 2, z = 3;
2624  /*
2625  * A = [ 0, -z, y;
2626  * z, 0, -x;
2627  * -y, x, 0];
2628  */
2629  A(x, x) = A(y, y) = A(z, z) = d;
2630  A(x, y) = -v(z);
2631  A(x, z) = v(y);
2632  A(y, x) = v(z);
2633  A(y, z) = -v(x);
2634  A(z, x) = -v(y);
2635  A(z, y) = v(x);
2636  }
2637 
2638  index_type iGetNumRows() const { return iNumRows; }
2639  index_type iGetNumCols() const { return iNumCols; }
2640 private:
2642  const doublereal d;
2643  };
2644 
2645  template <typename T, typename VectorExpr>
2646  class MatRInit {
2647  public:
2648  static const index_type iNumRows = 3;
2649  static const index_type iNumCols = 3;
2650 
2652  :g(g) {
2653  }
2654 
2655  void Initialize(Matrix<T, 3, 3>& RDelta) const {
2656  const T d = 4. / (4. + Dot(g, g));
2657 
2658  const T tmp1 = -g(3) * g(3);
2659  const T tmp2 = -g(2) * g(2);
2660  const T tmp3 = -g(1) * g(1);
2661  const T tmp4 = g(1) * g(2) * 0.5;
2662  const T tmp5 = g(2) * g(3) * 0.5;
2663  const T tmp6 = g(1) * g(3) * 0.5;
2664 
2665  RDelta(1,1) = (tmp1 + tmp2) * d * 0.5 + 1;
2666  RDelta(1,2) = (tmp4 - g(3)) * d;
2667  RDelta(1,3) = (tmp6 + g(2)) * d;
2668  RDelta(2,1) = (g(3) + tmp4) * d;
2669  RDelta(2,2) = (tmp1 + tmp3) * d * 0.5 + 1.;
2670  RDelta(2,3) = (tmp5 - g(1)) * d;
2671  RDelta(3,1) = (tmp6 - g(2)) * d;
2672  RDelta(3,2) = (tmp5 + g(1)) * d;
2673  RDelta(3,3) = (tmp2 + tmp3) * d * 0.5 + 1.;
2674  }
2675 
2676  index_type iGetNumRows() const { return iNumRows; }
2677  index_type iGetNumCols() const { return iNumCols; }
2678 
2679  private:
2681 };
2682 
2683 template <typename T, typename VectorExpr>
2685 public:
2686  static const index_type iNumRows = 3;
2687  static const index_type iNumCols = 3;
2688 
2690  :v(v) {
2691 
2692  }
2693 
2694  void Initialize(Matrix<T, 3, 3>& A) const {
2695  const index_type x = 1, y = 2, z = 3;
2696 
2697  const T vxvx = v(x) * v(x);
2698  const T vyvy = v(y) * v(y);
2699  const T vzvz = v(z) * v(z);
2700  const T vxvy = v(x) * v(y);
2701  const T vxvz = v(x) * v(z);
2702  const T vyvz = v(y) * v(z);
2703 
2704  A(x, x) = -vzvz - vyvy;
2705  A(x, y) = vxvy;
2706  A(x, z) = vxvz;
2707  A(y, x) = vxvy;
2708  A(y, y) = -vzvz - vxvx;
2709  A(y, z) = vyvz;
2710  A(z, x) = vxvz;
2711  A(z, y) = vyvz;
2712  A(z, z) = -vyvy - vxvx;
2713  }
2714 
2715  index_type iGetNumRows() const { return iNumRows; }
2716  index_type iGetNumCols() const { return iNumCols; }
2717 private:
2719 };
2720 
2721 template <typename T, typename VectorExpr>
2722 class MatGInit {
2723 public:
2724  static const index_type iNumRows = 3;
2725  static const index_type iNumCols = 3;
2726 
2728  :g(g) {
2729 
2730  }
2731 
2732  void Initialize(Matrix<T, 3, 3>& G) const {
2733 /*
2734  * G = 4/(4 + g' * g) * (eye(3) + 1/2*skew(g))
2735  *
2736  * skew(g) = [ 0, -g(3), g(2);
2737  * g(3), 0, -g(1);
2738  * -g(2), g(1), 0]
2739  */
2740  const T d = (4./(4.+Dot(g, g)));
2741 
2742  G(1, 1) = d;
2743  G(2, 1) = g(3) * d / 2.;
2744  G(3, 1) = -g(2) * d / 2.;
2745  G(1, 2) = -g(3) * d / 2.;
2746  G(2, 2) = d;
2747  G(3, 2) = g(1) * d / 2.;
2748  G(1, 3) = g(2) * d / 2.;
2749  G(2, 3) = -g(1) * d / 2.;
2750  G(3, 3) = d;
2751  }
2752 
2753  index_type iGetNumRows() const { return iNumRows; }
2754  index_type iGetNumCols() const { return iNumCols; }
2755 private:
2757 };
2758 
2759 template <typename T>
2763  }
2764 
2765  template <typename T>
2766  inline MatrixInit<MatCrossInit<T, Vec3DirectExpr>, T, 3, 3>
2767  MatCrossVec(const Vec3& v, doublereal d=0.) {
2768  return MatrixInit<MatCrossInit<T, Vec3DirectExpr>, T, 3, 3>(Direct(v), d);
2769 }
2770 
2771 template <typename VectorExpr>
2772 inline MatrixInit<MatCrossInit<typename VectorExpr::ScalarType, VectorExpr>, typename VectorExpr::ScalarType, 3, 3>
2774  return MatrixInit<MatCrossInit<typename VectorExpr::ScalarType, VectorExpr>, typename VectorExpr::ScalarType, 3, 3>(v, d);
2775 }
2776 
2777 template <typename T>
2778 inline MatrixInit<MatCrossCrossInit<T, VectorDirectExpr<Vector<T, 3> > >, T, 3, 3>
2781 }
2782 
2783 template <typename VectorExpr>
2784 inline MatrixInit<MatCrossCrossInit<typename VectorExpr::ScalarType, VectorExpr>, typename VectorExpr::ScalarType, 3, 3>
2786  return MatrixInit<MatCrossCrossInit<typename VectorExpr::ScalarType, VectorExpr>, typename VectorExpr::ScalarType, 3, 3>(v);
2787 }
2788 
2789 template <typename T>
2790 inline MatrixInit<MatGInit<T, VectorDirectExpr<Vector<T, 3> > >, T, 3, 3>
2793 }
2794 
2795 template <typename VectorExpr>
2796 inline MatrixInit<MatGInit<typename VectorExpr::ScalarType, VectorExpr>, typename VectorExpr::ScalarType, 3, 3>
2798  return MatrixInit<MatGInit<typename VectorExpr::ScalarType, VectorExpr>, typename VectorExpr::ScalarType, 3, 3>(g);
2799 }
2800 
2801  template <typename T>
2802  inline MatrixInit<MatRInit<T, VectorDirectExpr<Vector<T, 3> > >, T, 3, 3>
2803  MatRVec(const Vector<T, 3>& g) {
2805  }
2806 
2807  template <typename VectorExpr>
2808  inline MatrixInit<MatRInit<typename VectorExpr::ScalarType, VectorExpr>, typename VectorExpr::ScalarType, 3, 3>
2810  return MatrixInit<MatRInit<typename VectorExpr::ScalarType, VectorExpr>, typename VectorExpr::ScalarType, 3, 3>(g);
2811  }
2812 
2813  template <typename T, typename MatrixExpr>
2814  class VecRotInit {
2815  public:
2816  static const index_type iNumRows = 3;
2817 
2819  :R(R) {
2820  }
2821 
2822  index_type iGetNumRows() const { return iNumRows; }
2823 #if 0
2824  void Initialize(Vector<T, 3>& Phi) const {
2825  using std::sqrt;
2826  using std::fabs;
2827  using std::acos;
2828  /*
2829  This algorithm was taken from
2830  Rebecca M. Brannon and coauthors
2831  Computational Physics and Mechanics
2832  Sandia National Laboratories
2833  Albuquerque, NM 87185-0820
2834  http://www.me.unm.edu/~rmbrann/gobag.html
2835  */
2836  T cth = (R(1, 1) + R(2, 2) + R(3, 3) - 1.) * 0.5;
2837  T sth(0.);
2838 
2839  if (cth < 1.) {
2840  sth = sqrt(1. - cth * cth);
2841  }
2842 
2843  const scalar_func_type puny = 1e-12;
2844 
2845  if (fabs(sth) > puny) {
2846  const T angle = acos(cth);
2847  const T a1 = 0.5 * angle / sth;
2848  Phi(1) = (R(3, 2) - R(2, 3)) * a1;
2849  Phi(2) = (R(1, 3) - R(3, 1)) * a1;
2850  Phi(3) = (R(2, 1) - R(1, 2)) * a1;
2851  } else if (cth > 0.) {
2852  Phi(1) = R(3, 2);
2853  Phi(2) = R(1, 3);
2854  Phi(3) = R(2, 1);
2855  } else {
2856  const T angle = acos(cth);
2857  Matrix<T, 3, 3> scr = R;
2858 
2859  cth = 0.;
2860 
2861  index_type j = -1;
2862 
2863  for (index_type k = 1; k <= 3; ++k) {
2864  scr(k, k) += 1.;
2865  sth = Dot(scr.GetCol(k), scr.GetCol(k));
2866 
2867  if (sth > cth) {
2868  cth = sqrt(sth);
2869  j = k;
2870 
2871  for (index_type i = 1; i <= 3; ++i) {
2872  scr(i, j) /= cth;
2873  }
2874  }
2875  }
2876 
2877  if (j < 1) {
2878  MATVEC_ASSERT(false);
2880  }
2881 
2882  for (index_type i = 1; i <= 3; ++i) {
2883  Phi(i) = scr(i, j) * angle;
2884  }
2885  }
2886  }
2887 #else
2888  void Initialize(Vector<T, 3>& unit) const {
2889  // Modified from Appendix 2.4 of
2890  //
2891  // author = {Marco Borri and Lorenzo Trainelli and Carlo L. Bottasso},
2892  // title = {On Representations and Parameterizations of Motion},
2893  // journal = {Multibody System Dynamics},
2894  // volume = {4},
2895  // pages = {129--193},
2896  // year = {2000}
2897  using std::atan2;
2898  using std::sqrt;
2899 
2900  const T cosphi = 0.5 * (R(1, 1) + R(2, 2) + R(3, 3) - 1.);
2901 
2902  if (cosphi > 0.) {
2903  unit(1) = 0.5*(R(3, 2) - R(2, 3));
2904  unit(2) = 0.5*(R(1, 3) - R(3, 1));
2905  unit(3) = 0.5*(R(2, 1) - R(1, 2));
2906 
2907  const T sinphi2 = Dot(unit, unit);
2908  T sinphi;
2909 
2910  if (sinphi2 != 0) {
2911  sinphi = sqrt(sinphi2);
2912  } else {
2913  sinphi = unit(1);
2914  }
2915 
2916  const T phi = atan2(sinphi, cosphi);
2917  T a;
2918  RotCo(phi, a);
2919  unit /= a;
2920  } else {
2921  // -1 <= cosphi <= 0
2922  Matrix<T, 3, 3> eet = (R + Transpose(R)) * 0.5;
2923  eet(1, 1) -= cosphi;
2924  eet(2, 2) -= cosphi;
2925  eet(3, 3) -= cosphi;
2926  // largest (abs) component of unit vector phi/|phi|
2927  index_type maxcol = 1;
2928  if (eet(2, 2) > eet(1, 1)) {
2929  maxcol = 2;
2930  }
2931  if (eet(3, 3) > eet(maxcol, maxcol)) {
2932  maxcol = 3;
2933  }
2934  unit = (eet.GetCol(maxcol)/sqrt(eet(maxcol, maxcol)*(1. - cosphi)));
2935  T sinphi(0.);
2936  for (index_type i = 1; i <= 3; ++i) {
2937  sinphi -= Cross(unit, R.GetCol(i))(i) * 0.5;
2938  }
2939 
2940  unit *= atan2(sinphi, cosphi);
2941  }
2942  }
2943 
2944  private:
2945  static void RotCo(const T& phi, T& cf) {
2946  // This algorithm is a simplified version of RotCo in RotCoeff.hc
2947  // from Marco Morandini <morandini@aero.polimi.it>
2948  // and Teodoro Merlini <merlini@aero.polimi.it>
2949  using std::sin;
2950  using std::cos;
2951  using std::sqrt;
2952  using std::fabs;
2953 
2954  T phip[10];
2955  T phi2(phi * phi);
2956 
2957  if (fabs(phi) < RotCoeff::SerThrsh[0]) {
2958  phip[0] = 1.;
2959 
2960  for (index_type j = 1; j <= 9; j++) {
2961  phip[j] = phip[j - 1] * phi2;
2962  }
2963 
2964  cf = 0.;
2965 
2966  for (index_type j = 0; j < RotCoeff::SerTrunc[0]; j++) {
2967  cf += phip[j] / RotCoeff::SerCoeff[0][j];
2968  }
2969 
2970  return;
2971  }
2972 
2973  const T pd(sqrt(phi2));
2974  cf = sin(pd) / pd; // a = sin(phi)/phi
2975  };
2976 #endif
2977 
2978  private:
2980  };
2981 
2982  template <typename T>
2983  inline VectorInit<VecRotInit<T, MatrixDirectExpr<Matrix<T, 3, 3> > >, T, 3>
2986  }
2987 
2988 template <typename T, index_type N_rows, index_type N_cols>
2990 public:
2992  :A(A), iColWidth(iColWidth) {
2993 
2994  }
2995 
2996  void Print(std::ostream& os) const {
2997  for (index_type i = 1; i <= A.iGetNumRows(); ++i) {
2998  for (index_type j = 1; j <= A.iGetNumCols(); ++j) {
2999  os << std::setw(iColWidth) << A(i, j) << ' ';
3000  }
3001  os << std::endl;
3002  }
3003  }
3004 
3005 private:
3007  const int iColWidth;
3008 };
3009 
3010 template <typename T, index_type N_rows, index_type N_cols>
3012 Tabular(const Matrix<T, N_rows, N_cols>& A, int iColWidth=10) {
3013  return TabularMatrixView<T, N_rows, N_cols>(A, iColWidth);
3014 }
3015 
3016 template <typename T, index_type N_rows, index_type N_cols>
3017 inline std::ostream& operator<<(std::ostream& os, const Matrix<T, N_rows, N_cols>& A) {
3018  for (index_type i = 1; i <= A.iGetNumRows(); ++i) {
3019  for (index_type j = 1; j <= A.iGetNumCols(); ++j) {
3020  os << A(i, j);
3021 
3022  if (i < A.iGetNumRows() || j < A.iGetNumCols())
3023  os << ' ';
3024  }
3025  }
3026 
3027  return os;
3028 }
3029 
3030 template <typename T, index_type N_rows, index_type N_cols>
3031 inline std::ostream& operator<<(std::ostream& os, const TabularMatrixView<T, N_rows, N_cols>& tabA) {
3032  tabA.Print(os);
3033 
3034  return os;
3035 }
3036 
3037 template <typename T, index_type N_rows>
3038 inline std::ostream& operator<<(std::ostream& os, const Vector<T, N_rows>& x) {
3039  for (index_type i = 1; i <= x.iGetNumRows(); ++i) {
3040  os << x(i);
3041 
3042  if (i < x.iGetNumRows())
3043  os << ' ';
3044  }
3045 
3046  return os;
3047 }
3048 
3049 template <typename VectorExpressionType, index_type N_rows, index_type N_index>
3050 struct SumTraits {
3053  typedef SumTraits<VectorExpressionType, N_rows, N_index - 1> SumTraitsN_minus_1;
3056 
3057  static ExpressionType
3059  return ExpressionType(SumTraitsN_minus_1::Sum(v), v(N_index));
3060  }
3061 };
3062 
3063 template <typename VectorExpressionType, index_type N_rows>
3064 struct SumTraits<VectorExpressionType, N_rows, 1> {
3067 
3068  static ExpressionType
3070  return v(1);
3071  }
3072 };
3073 
3074 template <typename VectorExpressionType, index_type N_rows>
3075 inline typename SumTraits<VectorExpressionType, N_rows, N_rows>::ExpressionType
3078 }
3079 
3080 template <typename T, index_type N_rows>
3081 inline typename SumTraits<VectorDirectExpr<Vector<T, N_rows> >, N_rows, N_rows>::ExpressionType
3083  return SumTraits<VectorDirectExpr<Vector<T, N_rows> >, N_rows, N_rows>::Sum(Direct(v));
3084 }
3085 
3086 template <typename VectorExprLhs, typename VectorExprRhs, index_type N_rows, index_type N_index>
3087 struct DotTraits {
3093 
3094  typedef typename ScalarBinaryExpressionTraits<FuncMult,
3095  ScalarType,
3096  ScalarExprLhs,
3099 
3100  typedef typename ScalarBinaryExpressionTraits<FuncPlus,
3101  ScalarType,
3102  typename DotTraits<VectorExprLhs, VectorExprRhs, N_rows, N_index - 1>::ExpressionType,
3105 
3106  static ExpressionType
3109  }
3110 };
3111 
3112 template <typename VectorExprLhs, typename VectorExprRhs, index_type N_rows>
3113 struct DotTraits<VectorExprLhs, VectorExprRhs, N_rows, 1> {
3119  typedef typename ScalarBinaryExpressionTraits<FuncMult,
3120  ScalarType,
3121  ScalarExprLhs,
3124 
3125  static ExpressionType
3127  return ExpressionType(u(1), v(1));
3128  }
3129 };
3130 
3131 template <typename VectorExprLhs, typename VectorExprRhs, index_type N_rows>
3132 inline typename DotTraits<VectorExprLhs, VectorExprRhs, N_rows, N_rows>::ExpressionType
3135 }
3136 
3137 template <typename VectorExprLhs, typename T, index_type N_rows>
3138 inline typename DotTraits<VectorExprLhs, VectorDirectExpr<Vector<T, N_rows> >, N_rows, N_rows>::ExpressionType
3141 }
3142 
3143 template <typename T, typename VectorExprRhs, index_type N_rows>
3144 inline typename DotTraits<VectorDirectExpr<Vector<T, N_rows> >, VectorExprRhs, N_rows, N_rows>::ExpressionType
3146  return DotTraits<VectorDirectExpr<Vector<T, N_rows> >, VectorExprRhs, N_rows, N_rows>::Dot(Direct(u), v);
3147 }
3148 
3149 template <typename T_Lhs, typename T_Rhs, index_type N_rows>
3150 inline typename DotTraits<VectorDirectExpr<Vector<T_Lhs, N_rows> >, VectorDirectExpr<Vector<T_Rhs, N_rows> >, N_rows, N_rows>::ExpressionType
3152  return DotTraits<VectorDirectExpr<Vector<T_Lhs, N_rows> >, VectorDirectExpr<Vector<T_Rhs, N_rows> >, N_rows, N_rows>::Dot(Direct(u), Direct(v));
3153 }
3154 
3155  template <typename T_Lhs, index_type N_rows>
3156  inline typename DotTraits<VectorDirectExpr<Vector<T_Lhs, N_rows> >, Vec3DirectExpr, N_rows, N_rows>::ExpressionType
3157  Dot(const Vector<T_Lhs, N_rows>& u, const Vec3& v) {
3158  return DotTraits<VectorDirectExpr<Vector<T_Lhs, N_rows> >, Vec3DirectExpr, N_rows, N_rows>::Dot(Direct(u), Direct(v));
3159  }
3160 
3161 template <typename VectorExpr, index_type N_rows>
3162 inline typename VectorExpression<VectorExpr, N_rows>::ScalarType
3164  // avoid double evaluation
3166  using std::sqrt; // needed for g++-4.8
3167  return sqrt(Dot(u1, u1));
3168 }
3169 
3170 template <typename T, index_type N_rows>
3171 inline T
3173  return Norm(Direct(u));
3174 }
3175 
3176 template <typename VectorLhsExpr, typename VectorRhsExpr>
3177 struct CrossTraits {
3178  typedef typename VectorLhsExpr::ScalarType ScalarTypeLhs;
3179  typedef typename VectorLhsExpr::ExpressionType ExpressionTypeLhs;
3180  typedef typename VectorRhsExpr::ScalarType ScalarTypeRhs;
3181  typedef typename VectorRhsExpr::ExpressionType ExpressionTypeRhs;
3183  typedef typename ScalarBinaryExpressionTraits<FuncMult,
3184  ScalarType,
3188  typedef typename ScalarBinaryExpressionTraits<FuncMinus,
3189  ScalarType,
3190  ExprMult,
3191  ExprMult
3193 private:
3194  typedef typename IndexCheck<VectorLhsExpr::iNumRows - 3>::CheckType check_iNumRowsLhs;
3195  typedef typename IndexCheck<VectorRhsExpr::iNumRows - 3>::CheckType check_iNumRowsRhs;
3196 };
3197 
3198 template <typename VectorLhsExpr, typename VectorRhsExpr>
3200 public:
3201  static const bool bAlias = VectorLhsExpr::bAlias || VectorRhsExpr::bAlias;
3202  static const index_type iNumRows = 3;
3206 
3207  VectorCrossExpr(const VectorLhsExpr& u, const VectorRhsExpr& v)
3208  :oU(u), oV(v) {
3209 
3210  MATVEC_ASSERT(iNumRows == oU.iGetNumRows());
3211  MATVEC_ASSERT(iNumRows == oV.iGetNumRows());
3212  }
3213 
3215  --i;
3216  MATVEC_ASSERT(i >= 0);
3217  MATVEC_ASSERT(i < iNumRows);
3218 
3219  static const index_type x = 1, y = 2, z = 3;
3220  static const index_type e1[3] = {y, z, x},
3221  e2[3] = {z, x, y},
3222  e3[3] = {z, x, y},
3223  e4[3] = {y, z, x};
3224 
3225  return ExpressionType(ExprMult(oU(e1[i]), oV(e2[i])),
3226  ExprMult(oU(e3[i]), oV(e4[i])));
3227  }
3228 
3230  return iNumRows;
3231  }
3232 
3233  template <typename ScalarType2>
3234  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
3235  return oU.bHaveReferenceTo(pFirst, pLast) || oV.bHaveReferenceTo(pFirst, pLast);
3236  }
3237 
3238 private:
3239  const VectorLhsExpr oU;
3240  const VectorRhsExpr oV;
3241 
3242  typedef typename IndexCheck<iNumRows - VectorLhsExpr::iNumRows>::CheckType check_VectorLhsExpr;
3243  typedef typename IndexCheck<iNumRows - VectorRhsExpr::iNumRows>::CheckType check_VectorRhsExpr;
3244 };
3245 
3246 template <typename VectorLhsExpr, typename VectorRhsExpr>
3250 }
3251 
3252 template <typename VectorLhsExpr, typename T_Rhs>
3253 VectorExpression<VectorCrossExpr<VectorLhsExpr, VectorDirectExpr<Vector<T_Rhs, 3> > >, 3>
3256 }
3257 
3258 template <typename T_Lhs, typename VectorRhsExpr>
3259 VectorExpression<VectorCrossExpr<VectorDirectExpr<Vector<T_Lhs, 3> >, VectorRhsExpr>, 3>
3262 }
3263 
3264 template <typename T_Lhs, typename T_Rhs>
3265 VectorExpression<VectorCrossExpr<VectorDirectExpr<Vector<T_Lhs, 3> >, VectorDirectExpr<Vector<T_Rhs, 3> > >, 3>
3266 inline Cross(const Vector<T_Lhs, 3>& u, const Vector<T_Rhs, 3>& v) {
3267  return VectorExpression<VectorCrossExpr<VectorDirectExpr<Vector<T_Lhs, 3> >, VectorDirectExpr<Vector<T_Rhs, 3> > >, 3>(Direct(u), Direct(v));
3268 }
3269 
3270  template <typename T_Lhs>
3271  VectorExpression<VectorCrossExpr<VectorDirectExpr<Vector<T_Lhs, 3> >, Vec3DirectExpr>, 3>
3272  inline Cross(const Vector<T_Lhs, 3>& u, const Vec3& v) {
3274  }
3275 
3276 template <typename T>
3277 inline T Det(const Matrix<T, 2, 2>& A) {
3278  return A(1, 1) * A(2, 2) - A(1, 2) * A(2, 1);
3279 }
3280 
3281 template <typename T>
3283  const T detA = Det(A);
3284 
3285  return Matrix<T, 2, 2>(A(2, 2) / detA,
3286  -A(2, 1) / detA,
3287  -A(1, 2) / detA,
3288  A(1, 1) / detA);
3289 }
3290 
3291  template <index_type N_rows, index_type N_cols, typename MatrixLhsExpr, typename VectorRhsExpr>
3293 public:
3294  static const bool bAlias = MatrixLhsExpr::bAlias || VectorRhsExpr::bAlias;
3295  static const index_type iNumRows = MatrixLhsExpr::iNumRows;
3296  typedef typename MatrixLhsExpr::ScalarType MatrixLhsScalarExpr;
3297  typedef typename VectorRhsExpr::ScalarType VectorRhsScalarExpr;
3298  typedef typename MatrixLhsExpr::RowVectorType MatrixLhsRowVector;
3302 
3303  MatrixVectorProduct(const MatrixLhsExpr& A, const VectorRhsExpr& x)
3304  :A(A), x(x) {
3305  MATVEC_ASSERT(iNumRows == A.iGetNumRows());
3306  }
3307 
3309  MATVEC_ASSERT(i >= 1);
3310  MATVEC_ASSERT(i <= iGetNumRows());
3311  return Dot(A.GetRow(i), x);
3312  }
3313 
3315  MATVEC_ASSERT(iNumRows == A.iGetNumRows());
3316  return A.iGetNumRows();
3317  }
3318 
3319  template <typename ScalarType2>
3320  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
3321  return A.bHaveReferenceTo(pFirst, pLast) || x.bHaveReferenceTo(pFirst, pLast);
3322  }
3323 
3324 private:
3325  typedef typename IndexCheck<N_rows - MatrixLhsExpr::iNumRows>::CheckType check_iNumRowsLhs;
3326  typedef typename IndexCheck<N_cols - MatrixLhsExpr::iNumCols>::CheckType check_iNumColsLhs;
3327  typedef typename IndexCheck<N_cols - VectorRhsExpr::iNumRows>::CheckType check_iNumRowsRhs;
3330 
3331  const MatrixLhsExpr A;
3332  const VectorRhsExpr x;
3333  };
3334 
3335  template <typename MatrixLhsExpr, typename VectorRhsExpr>
3336  class MatrixVectorProduct<DYNAMIC_SIZE, DYNAMIC_SIZE, MatrixLhsExpr, VectorRhsExpr> {
3337  public:
3338  static const bool bAlias = MatrixLhsExpr::bAlias || VectorRhsExpr::bAlias;
3339  static const index_type iNumRows = MatrixLhsExpr::iNumRows;
3340  typedef typename MatrixLhsExpr::ScalarType MatrixLhsScalarExpr;
3341  typedef typename VectorRhsExpr::ScalarType VectorRhsScalarExpr;
3342 
3346 
3347  MatrixVectorProduct(const MatrixLhsExpr& A, const VectorRhsExpr& x)
3348  :A(A), x(x) {
3349 
3350  }
3351 
3353  MATVEC_ASSERT(i >= 1);
3354  MATVEC_ASSERT(i <= iGetNumRows());
3355  MATVEC_ASSERT(A.iGetNumCols() == x.iGetNumRows());
3356 
3357  ScalarType b_i(0);
3358 
3359  for (integer j = 1; j <= A.iGetNumCols(); ++j) {
3360  b_i += A(i, j) * x(j);
3361  }
3362 
3363  return b_i;
3364  }
3365 
3367  return A.iGetNumRows();
3368  }
3369 
3370  template <typename ScalarType2>
3371  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
3372  return A.bHaveReferenceTo(pFirst, pLast) || x.bHaveReferenceTo(pFirst, pLast);
3373  }
3374 
3375  private:
3376  const MatrixLhsExpr A;
3377  const VectorRhsExpr x;
3378 };
3379 
3380  template <index_type N_rows, typename MatrixLhsExpr, typename VectorRhsExpr>
3381  class MatrixVectorProduct<N_rows, DYNAMIC_SIZE, MatrixLhsExpr, VectorRhsExpr>: public MatrixVectorProduct<DYNAMIC_SIZE, DYNAMIC_SIZE, MatrixLhsExpr, VectorRhsExpr> {
3382  public:
3383  MatrixVectorProduct(const MatrixLhsExpr& A, const VectorRhsExpr& x)
3384  :MatrixVectorProduct<DYNAMIC_SIZE, DYNAMIC_SIZE, MatrixLhsExpr, VectorRhsExpr>(A, x) {
3385 
3386  }
3387  };
3388 
3389  template <index_type N_cols, typename MatrixLhsExpr, typename VectorRhsExpr>
3390  class MatrixVectorProduct<DYNAMIC_SIZE, N_cols, MatrixLhsExpr, VectorRhsExpr>: public MatrixVectorProduct<DYNAMIC_SIZE, DYNAMIC_SIZE, MatrixLhsExpr, VectorRhsExpr> {
3391  public:
3392  MatrixVectorProduct(const MatrixLhsExpr& A, const VectorRhsExpr& x)
3393  :MatrixVectorProduct<DYNAMIC_SIZE, DYNAMIC_SIZE, MatrixLhsExpr, VectorRhsExpr>(A, x) {
3394 
3395  }
3396  };
3397 
3398 template <typename MatrixLhsExpr, typename MatrixRhsExpr>
3400 public:
3401  static const bool bAlias = MatrixLhsExpr::bAlias || MatrixRhsExpr::bAlias;
3402  static const index_type iNumRows = MatrixLhsExpr::iNumRows;
3403  static const index_type iNumCols = MatrixRhsExpr::iNumCols;
3404  typedef typename MatrixLhsExpr::ScalarType MatrixLhsScalarExpr;
3405  typedef typename MatrixRhsExpr::ScalarType MatrixRhsScalarExpr;
3406  typedef typename MatrixLhsExpr::RowVectorType MatrixLhsRowVector;
3407  typedef typename MatrixRhsExpr::ColumnVectorType MatrixRhsColumnVector;
3413 
3414  MatrixMatrixProduct(const MatrixLhsExpr& A, const MatrixRhsExpr& B)
3415  :A(A), B(B) {
3416  MATVEC_ASSERT(A.iGetNumRows() == iNumRows);
3417  MATVEC_ASSERT(B.iGetNumCols() == iNumCols);
3418  MATVEC_ASSERT(A.iGetNumCols() == B.iGetNumRows());
3419  }
3420 
3422  MATVEC_ASSERT(i >= 1);
3423  MATVEC_ASSERT(i <= iGetNumRows());
3424  MATVEC_ASSERT(j >= 1);
3425  MATVEC_ASSERT(j <= iGetNumCols());
3426  return Dot(A.GetRow(i), B.GetCol(j));
3427  }
3428 
3430  MATVEC_ASSERT(iNumRows == A.iGetNumRows());
3431  return A.iGetNumRows();
3432  }
3433 
3435  MATVEC_ASSERT(iNumCols == B.iGetNumCols());
3436  return B.iGetNumCols();
3437  }
3438 
3440  MATVEC_ASSERT(iRow >= 1);
3441  MATVEC_ASSERT(iRow <= iNumRows);
3442  return RowVectorType(*this, iRow);
3443  }
3444 
3446  MATVEC_ASSERT(iCol >= 1);
3447  MATVEC_ASSERT(iCol <= iNumCols);
3448  return ColumnVectorType(*this, iCol);
3449  }
3450 
3451  template <typename ScalarType2>
3452  bool bHaveReferenceTo(const ScalarType2* pFirst, const ScalarType2* pLast) const {
3453  return A.bHaveReferenceTo(pFirst, pLast) || B.bHaveReferenceTo(pFirst, pLast);
3454  }
3455 
3456 private:
3457  const MatrixLhsExpr A;
3458  const MatrixRhsExpr B;
3459 };
3460 
3461 /****************************************************************************************************************
3462  * matrix vector product
3463  ****************************************************************************************************************/
3464 
3465 template <typename MatrixLhsExpr, index_type N_rows, index_type N_cols, typename VectorRhsExpr>
3467 operator* (const MatrixExpression<MatrixLhsExpr, N_rows, N_cols>& A, const VectorExpression<VectorRhsExpr, N_cols>& x) {
3468  return VectorExpression<MatrixVectorProduct<N_rows, N_cols, MatrixExpression<MatrixLhsExpr, N_rows, N_cols>, VectorExpression<VectorRhsExpr, N_cols> >, N_rows>(A, x);
3469 }
3470 
3471 template <typename T, index_type N_rows, index_type N_cols, typename VectorRhsExpr>
3472  inline VectorExpression<MatrixVectorProduct<N_rows, N_cols, MatrixExpression<MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, N_rows, N_cols>, VectorExpression<VectorRhsExpr, N_cols> >, N_rows>
3473 operator* (const Matrix<T, N_rows, N_cols>& A, const VectorExpression<VectorRhsExpr, N_cols>& x) {
3474  return VectorExpression<MatrixVectorProduct<N_rows, N_cols, MatrixExpression<MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, N_rows, N_cols>, VectorExpression<VectorRhsExpr, N_cols> >, N_rows>(Direct(A), x);
3475 }
3476 
3477 template <typename MatrixLhsExpr, index_type N_rows, index_type N_cols, typename T>
3478  inline VectorExpression<MatrixVectorProduct<N_rows, N_cols, MatrixExpression<MatrixLhsExpr, N_rows, N_cols>, VectorExpression<VectorDirectExpr<Vector<T, N_cols> >, N_cols> >, N_rows>
3480  return VectorExpression<MatrixVectorProduct<N_rows, N_cols, MatrixExpression<MatrixLhsExpr, N_rows, N_cols>, VectorExpression<VectorDirectExpr<Vector<T, N_cols> >, N_cols> >, N_rows>(A, Direct(x));
3481 }
3482 
3483 template <typename T1, typename T2, index_type N_rows, index_type N_cols>
3484  inline VectorExpression<MatrixVectorProduct<N_rows, N_cols, MatrixExpression<MatrixDirectExpr<Matrix<T1, N_rows, N_cols> >, N_rows, N_cols>, VectorExpression<VectorDirectExpr<Vector<T2, N_cols> >, N_cols> >, N_rows>
3486  return VectorExpression<MatrixVectorProduct<N_rows, N_cols, MatrixExpression<MatrixDirectExpr<Matrix<T1, N_rows, N_cols> >, N_rows, N_cols>, VectorExpression<VectorDirectExpr<Vector<T2, N_cols> >, N_cols> >, N_rows>(Direct(A), Direct(x));
3487 }
3488 
3489 template <typename T, index_type N_rows>
3490  inline VectorExpression<MatrixVectorProduct<N_rows, 3, MatrixExpression<MatrixDirectExpr<Matrix<T, N_rows, 3> >, N_rows, 3>, VectorExpression<Vec3DirectExpr, 3> >, N_rows>
3491 operator* (const Matrix<T, N_rows, 3>& A, const Vec3& x) {
3492  return VectorExpression<MatrixVectorProduct<N_rows, 3, MatrixExpression<MatrixDirectExpr<Matrix<T, N_rows, 3> >, N_rows, 3>, VectorExpression<Vec3DirectExpr, 3> >, N_rows>(Direct(A), Direct(x));
3493 }
3494 
3495 template <typename T>
3496  inline VectorExpression<MatrixVectorProduct<3, 3, MatrixExpression<Mat3x3DirectExpr, 3, 3>, VectorExpression<VectorDirectExpr<Vector<T, 3> >, 3> >, 3>
3497 operator* (const Mat3x3& A, const Vector<T, 3>& x) {
3498  return VectorExpression<MatrixVectorProduct<3, 3, MatrixExpression<Mat3x3DirectExpr, 3, 3>, VectorExpression<VectorDirectExpr<Vector<T, 3> >, 3> >, 3>(Direct(A), Direct(x));
3499  }
3500 
3501  template <typename T>
3504  return VectorExpression<MatrixVectorProduct<DYNAMIC_SIZE, DYNAMIC_SIZE, MatrixExpression<MatNxNDirectExpr, DYNAMIC_SIZE, DYNAMIC_SIZE>, VectorExpression<VectorDirectExpr<Vector<T, DYNAMIC_SIZE> >, DYNAMIC_SIZE> >, DYNAMIC_SIZE>(Direct(A), Direct(x));
3505  }
3506 
3507  template <typename T>
3508  inline VectorExpression<MatrixVectorProduct<3, DYNAMIC_SIZE, MatrixExpression<Mat3xNDirectExpr, 3, DYNAMIC_SIZE>, VectorExpression<VectorDirectExpr<Vector<T, DYNAMIC_SIZE> >, DYNAMIC_SIZE> >, 3> operator*(const Mat3xN& A, const Vector<T, DYNAMIC_SIZE>& x) {
3510  MATVEC_ASSERT(A.iGetNumRows() == 3);
3511  return VectorExpression<MatrixVectorProduct<3, DYNAMIC_SIZE, MatrixExpression<Mat3xNDirectExpr, 3, DYNAMIC_SIZE>, VectorExpression<VectorDirectExpr<Vector<T, DYNAMIC_SIZE> >, DYNAMIC_SIZE> >, 3> (Direct(A), Direct(x));
3512 }
3513 
3514 /****************************************************************************************************************
3515  * matrix matrix product
3516  ****************************************************************************************************************/
3517 
3518 template <typename MatrixLhsExpr, index_type N_rows_Lhs, index_type N_cols_Lhs, typename MatrixRhsExpr, index_type N_cols_Rhs>
3519 inline MatrixExpression<MatrixMatrixProduct<MatrixExpression<MatrixLhsExpr, N_rows_Lhs, N_cols_Lhs>, MatrixExpression<MatrixRhsExpr, N_cols_Lhs, N_cols_Rhs> >, N_rows_Lhs, N_cols_Rhs>
3520 operator* (const MatrixExpression<MatrixLhsExpr, N_rows_Lhs, N_cols_Lhs>& A, const MatrixExpression<MatrixRhsExpr, N_cols_Lhs, N_cols_Rhs>& B) {
3521  return MatrixExpression<MatrixMatrixProduct<MatrixExpression<MatrixLhsExpr, N_rows_Lhs, N_cols_Lhs>, MatrixExpression<MatrixRhsExpr, N_cols_Lhs, N_cols_Rhs> >, N_rows_Lhs, N_cols_Rhs>(A, B);
3522 }
3523 
3524 template <typename MatrixLhsExpr, index_type N_rows_Lhs, index_type N_cols_Lhs, typename T, index_type N_cols_Rhs>
3525 inline MatrixExpression<MatrixMatrixProduct<MatrixExpression<MatrixLhsExpr, N_rows_Lhs, N_cols_Lhs>, MatrixDirectExpr<Matrix<T, N_cols_Lhs, N_cols_Rhs> > >, N_rows_Lhs, N_cols_Rhs>
3527  return MatrixExpression<MatrixMatrixProduct<MatrixExpression<MatrixLhsExpr, N_rows_Lhs, N_cols_Lhs>, MatrixDirectExpr<Matrix<T, N_cols_Lhs, N_cols_Rhs> > >, N_rows_Lhs, N_cols_Rhs>(A, Direct(B));
3528 }
3529 
3530 template <typename T, index_type N_rows_Lhs, index_type N_cols_Lhs, typename MatrixRhsExpr, index_type N_cols_Rhs>
3531 inline MatrixExpression<MatrixMatrixProduct<MatrixDirectExpr<Matrix<T, N_rows_Lhs, N_cols_Lhs> >, MatrixExpression<MatrixRhsExpr, N_cols_Lhs, N_cols_Rhs> >, N_rows_Lhs, N_cols_Rhs>
3532 operator* (const Matrix<T, N_rows_Lhs, N_cols_Lhs>& A, const MatrixExpression<MatrixRhsExpr, N_cols_Lhs, N_cols_Rhs>& B) {
3533  return MatrixExpression<MatrixMatrixProduct<MatrixDirectExpr<Matrix<T, N_rows_Lhs, N_cols_Lhs> >, MatrixExpression<MatrixRhsExpr, N_cols_Lhs, N_cols_Rhs> >, N_rows_Lhs, N_cols_Rhs>(Direct(A), B);
3534 }
3535 
3536 template <typename T_Lhs, index_type N_rows_Lhs, index_type N_cols_Lhs, typename T_Rhs, index_type N_cols_Rhs>
3537 inline MatrixExpression<MatrixMatrixProduct<MatrixDirectExpr<Matrix<T_Lhs, N_rows_Lhs, N_cols_Lhs> >, MatrixDirectExpr<Matrix<T_Rhs, N_cols_Lhs, N_cols_Rhs> > >, N_rows_Lhs, N_cols_Rhs>
3539  return MatrixExpression<MatrixMatrixProduct<MatrixDirectExpr<Matrix<T_Lhs, N_rows_Lhs, N_cols_Lhs> >, MatrixDirectExpr<Matrix<T_Rhs, N_cols_Lhs, N_cols_Rhs> > >, N_rows_Lhs, N_cols_Rhs>(Direct(A), Direct(B));
3540 }
3541 
3542 template <typename T, index_type N_rows_Lhs>
3543 inline MatrixExpression<MatrixMatrixProduct<MatrixDirectExpr<Matrix<T, N_rows_Lhs, 3> >, Mat3x3DirectExpr>, N_rows_Lhs, 3>
3545  return MatrixExpression<MatrixMatrixProduct<MatrixDirectExpr<Matrix<T, N_rows_Lhs, 3> >, Mat3x3DirectExpr>, N_rows_Lhs, 3>(Direct(A), Direct(B));
3546 }
3547 
3548 template <typename T, index_type N_cols_Rhs>
3549 inline MatrixExpression<MatrixMatrixProduct<Mat3x3DirectExpr, MatrixDirectExpr<Matrix<T, 3, N_cols_Rhs> > >, 3, N_cols_Rhs>
3552 }
3553 
3554 #define VECTOR_VECTOR_DEFINE_BINARY_FUNCTION(ExpressionName, FunctionName, FunctionClass) \
3555  template <typename VectorLhsExpr, typename VectorRhsExpr, index_type N_rows> \
3556  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename VectorLhsExpr::ExpressionType, typename VectorRhsExpr::ExpressionType>, VectorLhsExpr, VectorRhsExpr>, N_rows> \
3557  FunctionName(const VectorExpression<VectorLhsExpr, N_rows>& u, const VectorExpression<VectorRhsExpr, N_rows>& v) { \
3558  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename VectorLhsExpr::ExpressionType, typename VectorRhsExpr::ExpressionType>, VectorLhsExpr, VectorRhsExpr>, N_rows>(u, v); \
3559  } \
3560  \
3561  template <typename T, index_type N_rows> \
3562  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, typename ScalarTypeTraits<T>::DirectExpressionType>, VectorDirectExpr<Vector<T, N_rows> >, VectorDirectExpr<Vector<T, N_rows> > >, N_rows> \
3563  FunctionName(const Vector<T, N_rows>& u, const Vector<T, N_rows>& v) { \
3564  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, typename ScalarTypeTraits<T>::DirectExpressionType>, VectorDirectExpr<Vector<T, N_rows> >, VectorDirectExpr<Vector<T, N_rows> > >, N_rows>(u, v); \
3565  } \
3566  \
3567  template <typename T1, typename T2, index_type N_rows> \
3568  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T1>::DirectExpressionType, typename ScalarTypeTraits<T2>::DirectExpressionType>, VectorDirectExpr<Vector<T1, N_rows> >, VectorDirectExpr<Vector<T2, N_rows> > >, N_rows> \
3569  FunctionName(const Vector<T1, N_rows>& u, const Vector<T2, N_rows>& v) { \
3570  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T1>::DirectExpressionType, typename ScalarTypeTraits<T2>::DirectExpressionType>, VectorDirectExpr<Vector<T1, N_rows> >, VectorDirectExpr<Vector<T2, N_rows> > >, N_rows>(u, v); \
3571  } \
3572  \
3573  template <typename T1, index_type N_rows> \
3574  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T1>::DirectExpressionType, typename ScalarTypeTraits<doublereal>::DirectExpressionType>, VectorDirectExpr<Vector<T1, N_rows> >, Vec3DirectExpr>, N_rows> \
3575  FunctionName(const Vector<T1, N_rows>& u, const Vec3& v) { \
3576  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T1>::DirectExpressionType, typename ScalarTypeTraits<doublereal>::DirectExpressionType>, VectorDirectExpr<Vector<T1, N_rows> >, Vec3DirectExpr>, N_rows>(u, v); \
3577  } \
3578  \
3579  template <index_type N_rows> \
3580  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, scalar_func_type, scalar_func_type>, VectorDirectExpr<Vector<scalar_func_type, N_rows> >, VectorDirectExpr<Vector<scalar_func_type, N_rows> > >, N_rows> \
3581  FunctionName(const Vector<scalar_func_type, N_rows>& u, const Vector<scalar_func_type, N_rows>& v) { \
3582  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, scalar_func_type, scalar_func_type>, VectorDirectExpr<Vector<scalar_func_type, N_rows> >, VectorDirectExpr<Vector<scalar_func_type, N_rows> > >, N_rows>(u, v); \
3583  } \
3584  \
3585  template <typename VectorLhsExpr, typename T, index_type N_rows> \
3586  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename VectorLhsExpr::ExpressionType, typename ScalarTypeTraits<T>::DirectExpressionType >, VectorLhsExpr, VectorDirectExpr<Vector<T, N_rows> > >, N_rows> \
3587  FunctionName(const VectorExpression<VectorLhsExpr, N_rows>& u, const Vector<T, N_rows>& v) { \
3588  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename VectorLhsExpr::ExpressionType, typename ScalarTypeTraits<T>::DirectExpressionType >, VectorLhsExpr, VectorDirectExpr<Vector<T, N_rows> > >, N_rows>(u, v); \
3589  } \
3590  \
3591  template <typename T, index_type N_rows, typename VectorRhsExpr> \
3592  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, typename VectorRhsExpr::ExpressionType>, VectorDirectExpr<Vector<T, N_rows> >, VectorRhsExpr>, N_rows> \
3593  FunctionName(const Vector<T, N_rows>& u, const VectorExpression<VectorRhsExpr, N_rows>& v) { \
3594  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, typename VectorRhsExpr::ExpressionType>, VectorDirectExpr<Vector<T, N_rows> >, VectorRhsExpr>, N_rows>(u, v); \
3595  }
3596 
3597 #define VECTOR_SCALAR_DEFINE_BINARY_FUNCTION(ExpressionName, FunctionName, FunctionClass) \
3598  template <typename VectorLhsExpr, typename ScalarRhsExpr, index_type N_rows> \
3599  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename VectorLhsExpr::ExpressionType, ScalarRhsExpr>, VectorLhsExpr, ScalarRhsExpr>, N_rows> \
3600  FunctionName(const VectorExpression<VectorLhsExpr, N_rows>& u, const GradientExpression<ScalarRhsExpr>& v) { \
3601  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename VectorLhsExpr::ExpressionType, ScalarRhsExpr>, VectorLhsExpr, ScalarRhsExpr>, N_rows>(u, v); \
3602  } \
3603  \
3604  template <typename VectorLhsExpr, index_type N_rows, index_type N_SIZE> \
3605  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename VectorLhsExpr::ExpressionType, DirectExpr<Gradient<N_SIZE> > >, VectorLhsExpr, DirectExpr<Gradient<N_SIZE> > >, N_rows> \
3606  FunctionName(const VectorExpression<VectorLhsExpr, N_rows>& u, const Gradient<N_SIZE>& v) { \
3607  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename VectorLhsExpr::ExpressionType, DirectExpr<Gradient<N_SIZE> > >, VectorLhsExpr, DirectExpr<Gradient<N_SIZE> > >, N_rows>(u, v); \
3608  } \
3609  \
3610  template <typename VectorLhsExpr, index_type N_rows> \
3611  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename VectorLhsExpr::ExpressionType, scalar_func_type>, VectorLhsExpr, scalar_func_type>, N_rows> \
3612  FunctionName(const VectorExpression<VectorLhsExpr, N_rows>& u, scalar_func_type v) { \
3613  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename VectorLhsExpr::ExpressionType, scalar_func_type >, VectorLhsExpr, scalar_func_type>, N_rows>(u, v); \
3614  } \
3615  \
3616  template <typename T, index_type N_rows> \
3617  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, scalar_func_type >, VectorDirectExpr<Vector<T, N_rows> >, scalar_func_type>, N_rows> \
3618  FunctionName(const Vector<T, N_rows>& u, scalar_func_type v) { \
3619  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, scalar_func_type >, VectorDirectExpr<Vector<T, N_rows> >, scalar_func_type>, N_rows>(u, v); \
3620  } \
3621  \
3622  template <typename T, index_type N_rows, index_type N_SIZE> \
3623  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, DirectExpr<Gradient<N_SIZE> > >, VectorDirectExpr<Vector<T, N_rows> >, DirectExpr<Gradient<N_SIZE> > >, N_rows> \
3624  FunctionName(const Vector<T, N_rows>& u, const Gradient<N_SIZE>& v) { \
3625  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, DirectExpr<Gradient<N_SIZE> > >, VectorDirectExpr<Vector<T, N_rows> >, DirectExpr<Gradient<N_SIZE> > >, N_rows>(u, v); \
3626  } \
3627  \
3628  template <typename T, index_type N_rows, typename RhsExpr> \
3629  inline VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, GradientExpression<RhsExpr> >, VectorDirectExpr<Vector<T, N_rows> >, GradientExpression<RhsExpr> >, N_rows> \
3630  FunctionName(const Vector<T, N_rows>& u, const GradientExpression<RhsExpr>& v) { \
3631  return VectorExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, GradientExpression<RhsExpr> >, VectorDirectExpr<Vector<T, N_rows> >, GradientExpression<RhsExpr> >, N_rows>(u, v); \
3632  }
3633 
3634 #define VECTOR_DEFINE_UNARY_FUNCTION(ExpressionName, FunctionName, FunctionClass) \
3635  template <typename VectorExpr, index_type N_rows> \
3636  inline VectorExpression<ExpressionName<ScalarUnaryOperation<FunctionClass, typename VectorExpr::ExpressionType>, VectorExpr>, N_rows> \
3637  FunctionName(const VectorExpression<VectorExpr, N_rows>& u) { \
3638  return VectorExpression<ExpressionName<ScalarUnaryOperation<FunctionClass, typename VectorExpr::ExpressionType>, VectorExpr>, N_rows>(u); \
3639  } \
3640  \
3641  template <typename T, index_type N_rows> \
3642  inline VectorExpression<ExpressionName<ScalarUnaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType>, VectorDirectExpr<Vector<T, N_rows> > >, N_rows> \
3643  FunctionName(const Vector<T, N_rows>& u) { \
3644  return VectorExpression<ExpressionName<ScalarUnaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType>, VectorDirectExpr<Vector<T, N_rows> > >, N_rows>(u); \
3645  }
3646 
3647 #define MATRIX_DEFINE_UNARY_FUNCTION(ExpressionName, FunctionName, FunctionClass) \
3648  template <typename MatrixExpr, index_type N_rows, index_type N_cols> \
3649  inline MatrixExpression<ExpressionName<ScalarUnaryOperation<FunctionClass, typename MatrixExpr::ExpressionType>, MatrixExpr>, N_rows, N_cols> \
3650  FunctionName(const MatrixExpression<MatrixExpr, N_rows, N_cols>& u) { \
3651  return MatrixExpression<ExpressionName<ScalarUnaryOperation<FunctionClass, typename MatrixExpr::ExpressionType>, MatrixExpr>, N_rows, N_cols>(u); \
3652  } \
3653  \
3654  template <typename T, index_type N_rows, index_type N_cols> \
3655  inline MatrixExpression<ExpressionName<ScalarUnaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType>, MatrixDirectExpr<Matrix<T, N_rows, N_cols> > >, N_rows, N_cols> \
3656  FunctionName(const Matrix<T, N_rows, N_cols>& u) { \
3657  return MatrixExpression<ExpressionName<ScalarUnaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType>, MatrixDirectExpr<Matrix<T, N_rows, N_cols> > >, N_rows, N_cols>(u); \
3658  }
3659 
3660 #define MATRIX_MATRIX_DEFINE_BINARY_FUNCTION(ExpressionName, FunctionName, FunctionClass) \
3661  template <typename MatrixLhsExpr, typename MatrixRhsExpr, index_type N_rows, index_type N_cols> \
3662  inline MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename MatrixLhsExpr::ExpressionType, typename MatrixRhsExpr::ExpressionType>, MatrixLhsExpr, MatrixRhsExpr>, N_rows, N_cols> \
3663  FunctionName(const MatrixExpression<MatrixLhsExpr, N_rows, N_cols>& u, const MatrixExpression<MatrixRhsExpr, N_rows, N_cols>& v) { \
3664  return MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename MatrixLhsExpr::ExpressionType, typename MatrixRhsExpr::ExpressionType>, MatrixLhsExpr, MatrixRhsExpr>, N_rows, N_cols>(u, v); \
3665  } \
3666  \
3667  template <typename T, index_type N_rows, index_type N_cols> \
3668  inline MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, typename ScalarTypeTraits<T>::DirectExpressionType>, MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, MatrixDirectExpr<Matrix<T, N_rows, N_cols> > >, N_rows, N_cols> \
3669  FunctionName(const Matrix<T, N_rows, N_cols>& u, const Matrix<T, N_rows, N_cols>& v) { \
3670  return MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, typename ScalarTypeTraits<T>::DirectExpressionType>, MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, MatrixDirectExpr<Matrix<T, N_rows, N_cols> > >, N_rows, N_cols>(u, v); \
3671  } \
3672  \
3673  template <typename MatrixLhsExpr, typename T, index_type N_rows, index_type N_cols> \
3674  inline MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename MatrixLhsExpr::ExpressionType, typename ScalarTypeTraits<T>::DirectExpressionType >, MatrixLhsExpr, MatrixDirectExpr<Matrix<T, N_rows, N_cols> > >, N_rows, N_cols> \
3675  FunctionName(const MatrixExpression<MatrixLhsExpr, N_rows, N_cols>& u, const Matrix<T, N_rows, N_cols>& v) { \
3676  return MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename MatrixLhsExpr::ExpressionType, typename ScalarTypeTraits<T>::DirectExpressionType >, MatrixLhsExpr, MatrixDirectExpr<Matrix<T, N_rows, N_cols> > >, N_rows, N_cols>(u, v); \
3677  } \
3678  \
3679  template <typename T, index_type N_rows, index_type N_cols, typename MatrixRhsExpr> \
3680  inline MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, typename MatrixRhsExpr::ExpressionType>, MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, MatrixRhsExpr>, N_rows, N_cols> \
3681  FunctionName(const Matrix<T, N_rows, N_cols>& u, const MatrixExpression<MatrixRhsExpr, N_rows, N_cols>& v) { \
3682  return MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, typename MatrixRhsExpr::ExpressionType>, MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, MatrixRhsExpr>, N_rows, N_cols>(u, v); \
3683  }
3684 
3685 #define MATRIX_SCALAR_DEFINE_BINARY_FUNCTION(ExpressionName, FunctionName, FunctionClass) \
3686  template <typename MatrixLhsExpr, typename ScalarRhsExpr, index_type N_rows, index_type N_cols> \
3687  inline MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename MatrixLhsExpr::ExpressionType, ScalarRhsExpr>, MatrixLhsExpr, ScalarRhsExpr>, N_rows, N_cols> \
3688  FunctionName(const MatrixExpression<MatrixLhsExpr, N_rows, N_cols>& u, const GradientExpression<ScalarRhsExpr>& v) { \
3689  return MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename MatrixLhsExpr::ExpressionType, ScalarRhsExpr>, MatrixLhsExpr, ScalarRhsExpr>, N_rows, N_cols>(u, v); \
3690  } \
3691  \
3692  template <typename MatrixLhsExpr, index_type N_rows, index_type N_cols, index_type N_SIZE> \
3693  inline MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename MatrixLhsExpr::ExpressionType, DirectExpr<Gradient<N_SIZE> > >, MatrixLhsExpr, DirectExpr<Gradient<N_SIZE> > >, N_rows, N_cols> \
3694  FunctionName(const MatrixExpression<MatrixLhsExpr, N_rows, N_cols>& u, const Gradient<N_SIZE>& v) { \
3695  return MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename MatrixLhsExpr::ExpressionType, DirectExpr<Gradient<N_SIZE> > >, MatrixLhsExpr, DirectExpr<Gradient<N_SIZE> > >, N_rows, N_cols>(u, v); \
3696  } \
3697  \
3698  template <typename MatrixLhsExpr, index_type N_rows, index_type N_cols> \
3699  inline MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename MatrixLhsExpr::ExpressionType, scalar_func_type>, MatrixLhsExpr, scalar_func_type>, N_rows, N_cols> \
3700  FunctionName(const MatrixExpression<MatrixLhsExpr, N_rows, N_cols>& u, scalar_func_type v) { \
3701  return MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename MatrixLhsExpr::ExpressionType, scalar_func_type >, MatrixLhsExpr, scalar_func_type>, N_rows, N_cols>(u, v); \
3702  } \
3703  \
3704  template <typename T, index_type N_rows, index_type N_cols> \
3705  inline MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, scalar_func_type >, MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, scalar_func_type>, N_rows, N_cols> \
3706  FunctionName(const Matrix<T, N_rows, N_cols>& u, scalar_func_type v) { \
3707  return MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, scalar_func_type >, MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, scalar_func_type>, N_rows, N_cols>(u, v); \
3708  } \
3709  \
3710  template <typename T, index_type N_rows, index_type N_cols, index_type N_SIZE> \
3711  inline MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, DirectExpr<Gradient<N_SIZE> > >, MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, DirectExpr<Gradient<N_SIZE> > >, N_rows, N_cols> \
3712  FunctionName(const Matrix<T, N_rows, N_cols>& u, const Gradient<N_SIZE>& v) { \
3713  return MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, DirectExpr<Gradient<N_SIZE> > >, MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, DirectExpr<Gradient<N_SIZE> > >, N_rows, N_cols>(u, v); \
3714  } \
3715  \
3716  template <typename T, index_type N_rows, index_type N_cols, typename RhsExpr> \
3717  inline MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, GradientExpression<RhsExpr> >, MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, GradientExpression<RhsExpr> >, N_rows, N_cols> \
3718  FunctionName(const Matrix<T, N_rows, N_cols>& u, const GradientExpression<RhsExpr>& v) { \
3719  return MatrixExpression<ExpressionName<ScalarBinaryOperation<FunctionClass, typename ScalarTypeTraits<T>::DirectExpressionType, GradientExpression<RhsExpr> >, MatrixDirectExpr<Matrix<T, N_rows, N_cols> >, GradientExpression<RhsExpr> >, N_rows, N_cols>(u, v); \
3720  }
3721 
3725 VECTOR_SCALAR_DEFINE_BINARY_FUNCTION(VectorScalarVectorBinaryExpr, operator/, FuncDiv)
3727 
3729 MATRIX_MATRIX_DEFINE_BINARY_FUNCTION(MatrixMatrixMatrixBinaryExpr, operator-, FuncMinus)
3731 MATRIX_SCALAR_DEFINE_BINARY_FUNCTION(MatrixScalarMatrixBinaryExpr, operator/, FuncDiv)
3733 
3734 #undef VECTOR_VECTOR_DEFINE_BINARY_FUNCTION
3735 #undef VECTOR_SCALAR_DEFINE_BINARY_FUNCTION
3736 #undef VECTOR_DEFINE_UNARY_FUNCTION
3737 
3738 #undef MATRIX_MATRIX_DEFINE_BINARY_FUNCTION
3739 #undef MATRIX_SCALAR_DEFINE_BINARY_FUNCTION
3740 #undef MATRIX_DEFINE_UNARY_FUNCTION
3741 }
3742 
3743 #endif
static const index_type iNumCols
Definition: matvec.h:1566
VectorExpr::ScalarType ScalarType
Definition: matvec.h:807
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:3234
RowVectorType GetRow(index_type iRow) const
Definition: matvec.h:3439
CommonScalarType< ScalarTypeLhs, ScalarTypeRhs >::ScalarType ScalarType
Definition: matvec.h:3092
static const index_type iNumRows
Definition: matvec.h:1413
static int angle(const MathParser::MathArgs &args)
Definition: modelns.cc:313
MatrixData< ScalarType, iNumRows, iNumCols > rgMat
Definition: matvec.h:2180
static const bool bAlias
Definition: matvec.h:732
Vector & operator+=(const VectorExpression< Expression, N_rows > &f)
Definition: matvec.h:2409
Expression::ScalarType ScalarType
Definition: matvec.h:640
void Resize(index_type iRows, index_type iCols)
Definition: matvec.h:1989
ExpressionType operator()(index_type i) const
Definition: matvec.h:3214
Vector(const VectorExpression< Expression, N_rows > &f)
Definition: matvec.h:2314
Mat3x3DirectExpr(const Mat3x3 &A)
Definition: matvec.h:1477
VectorVectorUnaryExpr(const VectorExpr &u)
Definition: matvec.h:777
VectorData(const VectorData &oVec)
Definition: matvec.h:1688
ColumnVectorType GetCol(index_type iCol) const
Definition: matvec.h:3445
ExpressionType operator()(index_type iCol) const
Definition: matvec.h:1086
ScalarUnaryFunc::ExpressionType ExpressionType
Definition: matvec.h:1360
IndexCheck< VectorRhsExpr::iNumRows-3 >::CheckType check_iNumRowsRhs
Definition: matvec.h:3195
static const index_type iNumRows
Definition: matvec.h:1565
const ScalarType *const pVec
Definition: matvec.h:930
Matrix(index_type iRows, index_type iCols)
Definition: matvec.h:1926
const T * pGetLastElem() const
Definition: matvec.h:1764
index_type iGetNumRows() const
Definition: matvec.h:2822
MatrixExpr::ScalarType ScalarType
Definition: matvec.h:1047
VectorExpression< ColumnVectorExpr< MatrixMatrixProduct >, iNumRows > ColumnVectorType
Definition: matvec.h:3409
index_type iGetStartIndexLocalVector() const
Definition: matvec.h:205
static const index_type iNumRows
Definition: matvec.h:1168
void Resize(index_type iRows, index_type iCols)
Definition: matvec.h:1755
void Resize(index_type iRows, index_type iCols)
Definition: matvec.h:1888
static const index_type iNumCols
Definition: matvec.h:1469
DotTraits< MatrixLhsRowVector, VectorRhsExpr, MatrixLhsExpr::iNumCols, MatrixLhsExpr::iNumCols >::ExpressionType ExpressionType
Definition: matvec.h:3301
void ApplyScalarFunc(const U &a)
Definition: matvec.h:2512
const ScalarType & operator()(index_type i, index_type j) const
Definition: matvec.h:1578
scalar_deriv_type dGetDerivativeLocal(index_type iLocalDof) const
Definition: matvec.h:189
index_type iGetNumRows() const
Definition: matvec.h:3314
Mat3xNDirectExpr(const Mat3xN &A)
Definition: matvec.h:1532
ExpressionType operator()(index_type i) const
Definition: matvec.h:815
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:925
T rgData[N_rows]
Definition: matvec.h:1678
MatrixExpression< TransposedMatrix< MatrixDirectExpr< Matrix< T, N_rows, N_cols > > >, N_cols, N_rows > Transpose(const Matrix< T, N_rows, N_cols > &A)
Definition: matvec.h:2206
ScalarTypeTraits< doublereal >::ScalarType ScalarType
Definition: matvec.h:874
#define MATRIX_SCALAR_DEFINE_BINARY_FUNCTION(ExpressionName, FunctionName, FunctionClass)
Definition: matvec.h:3685
GradientExpression< BinaryExpr< ScalarBinaryFunction, ScalarLhsExpr, ScalarRhsExpr > > ExpressionType
Definition: matvec.h:332
VectorExpr::ExpressionType ExpressionType
Definition: matvec.h:808
static ExpressionType Sum(const VectorExpression< VectorExpressionType, N_rows > &v)
Definition: matvec.h:3058
index_type iGetNumRows() const
Definition: matvec.h:709
SumTraits< VectorDirectExpr< Vector< T, N_rows > >, N_rows, N_rows >::ExpressionType Sum(const Vector< T, N_rows > &v)
Definition: matvec.h:3082
GenericBinaryExpression< ScalarBinaryFunction, T, ScalarLhsExpr, ScalarRhsExpr >::ScalarType ScalarType
Definition: matvec.h:289
static const index_type iNumRows
Definition: matvec.h:588
#define MBDYN_EXCEPT_ARGS
Definition: except.h:63
CommonScalarType< ScalarTypeLhs, ScalarTypeRhs >::ScalarType ScalarType
Definition: matvec.h:3118
Definition: matvec3.h:98
static const index_type DYNAMIC_SIZE
Definition: gradient.h:141
GradientExpression< UnaryExpr< FuncSin, Expr > > sin(const GradientExpression< Expr > &u)
Definition: gradient.h:2977
const ScalarType * pGetLastElem() const
Definition: matvec.h:2504
const ScalarType & operator()(index_type i) const
Definition: matvec.h:850
const Mat3x3 & A
Definition: matvec.h:1518
void ApplyMatrixFuncAlias(const MatrixExpression< Expression, N_rows, N_cols > &A, const Func &f)
Definition: matvec.h:2175
MatNxNDirectExpr(const MatNxN &A)
Definition: matvec.h:1573
void ApplyScalarFunc(const U &a)
Definition: matvec.h:2143
static const index_type iNumRows
Definition: matvec.h:873
integer iGetNumRows(void) const
Definition: matvec3n.h:562
VectorExpression< RowVectorExpr< MatrixMatrixProduct >, iNumCols > RowVectorType
Definition: matvec.h:3408
ScalarTypeTraits< T >::ScalarType ScalarType
Definition: matvec.h:976
index_type iGetNumCols() const
Definition: matvec.h:2639
RowVectorType GetRow(index_type iRow) const
Definition: matvec.h:1198
index_type iGetNumRows() const
Definition: matvec.h:748
ScalarTypeTraits< doublereal >::DirectExpressionType ExpressionType
Definition: matvec.h:1530
MatrixExpr::ExpressionType ExpressionType
Definition: matvec.h:1174
MaxSizeCheck< iNumRows!=DYNAMIC_SIZE >::CheckType check_iNumRows
Definition: matvec.h:968
index_type iGetNumCols() const
Definition: matvec.h:1194
GenericBinaryExpression(const ScalarLhsExpr &lhs, const ScalarRhsExpr &rhs)
Definition: matvec.h:112
Vector & operator=(const VectorExpression< Expression, N_rows > &f)
Definition: matvec.h:2367
index_type iGetNumRows() const
Definition: matvec.h:1378
const ScalarType * pGetFirstElem() const
Definition: matvec.h:2500
ExpressionType operator()(index_type i, index_type j) const
Definition: matvec.h:3421
bool bHaveReferenceTo(const void *p) const
Definition: matvec.h:152
static const index_type iNumCols
Definition: matvec.h:1358
MatrixLhsExpr::ScalarType MatrixLhsScalarExpr
Definition: matvec.h:3404
Vector(index_type N)
Definition: matvec.h:2267
ColumnVectorType GetCol(index_type iCol) const
Definition: matvec.h:1204
ScalarBinaryExpressionTraits< FuncMinus, ScalarType, ExprMult, ExprMult >::ExpressionType ExpressionType
Definition: matvec.h:3192
Matrix & operator+=(const Matrix< T_Rhs, N_rows, N_cols > &A)
Definition: matvec.h:2006
VectorRhsExpr::ScalarType ScalarTypeRhs
Definition: matvec.h:3180
static const bool bAlias
Definition: matvec.h:1076
GradientExpression< BinaryExpr< ScalarBinaryFunction, ScalarLhsExpr, DirectExpr< Gradient< N_SIZE > > > > ExpressionType
Definition: matvec.h:337
static void Eval(T &b, const U &a)
Definition: matvec.h:558
void Resize(index_type iNumRows)
Definition: matvec.h:1665
void Compute() const
Definition: matvec.h:225
BasicScalarType< Expr >::ScalarType ScalarType
Definition: matvec.h:269
ScalarTypeTraits< doublereal >::ScalarType ScalarType
Definition: matvec.h:1471
integer iGetNumRows(void) const
Definition: matvec3n.h:254
MatrixRhsExpr::ColumnVectorType MatrixRhsColumnVector
Definition: matvec.h:3407
VectorExpression(const Expression &e)
Definition: matvec.h:592
ScalarBinaryExpressionTraits< FuncPlus, ScalarType, typename DotTraits< VectorExprLhs, VectorExprRhs, N_rows, N_index-1 >::ExpressionType, MultExpressionType >::ExpressionType ExpressionType
Definition: matvec.h:3104
CommonScalarType< typename BasicScalarType< LhsExpr >::ScalarType, typename BasicScalarType< RhsExpr >::ScalarType >::ScalarType ScalarType
Definition: matvec.h:264
ScalarBinFunc::ScalarType ScalarType
Definition: matvec.h:1231
static const index_type iNumRows
Definition: matvec.h:1468
MatrixScalarMatrixBinaryExpr(const MatrixLhsExpr &u, const ScalarRhsExpr &v)
Definition: matvec.h:1304
Matrix(const Matrix &A)
Definition: matvec.h:1941
IndexCheck< iNumRows-VectorRhsExpr::iNumRows >::CheckType check_VectorRhsExpr
Definition: matvec.h:3243
const T * pGetFirstElem() const
Definition: matvec.h:1714
RowVectorType GetRow(index_type iRow) const
Definition: matvec.h:1498
IndexCheck< N_rows-MatrixLhsExpr::iNumRows >::CheckType check_iNumRowsLhs
Definition: matvec.h:3325
VectorExpression< VectorExpressionType, N_rows >::ScalarType ScalarType
Definition: matvec.h:3065
ScalarBinFunc::ScalarType ScalarType
Definition: matvec.h:695
const T * pGetFirstElem() const
Definition: matvec.h:1760
void Initialize(Vector< T, 3 > &unit) const
Definition: matvec.h:2888
index_type iGetNumRows() const
Definition: matvec.h:788
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:1595
VectorExpression< ColumnVectorExpr< MatrixMatrixUnaryExpr >, iNumRows > ColumnVectorType
Definition: matvec.h:1362
Matrix(const MatrixInit< InitClass, T, N_rows, N_cols > &func)
Definition: matvec.h:1948
RowVectorType GetRow(index_type iRow) const
Definition: matvec.h:1329
BasicScalarType< ScalarLhsExpr >::ScalarType ScalarTypeLhs
Definition: matvec.h:388
VectorData(index_type N, bool bZeroInit)
Definition: matvec.h:1632
const ScalarType & operator()(index_type iRow) const
Definition: matvec.h:915
MaxSizeCheck< N_cols!=DYNAMIC_SIZE >::CheckType check_iNumColsDynamic
Definition: matvec.h:3329
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:1274
MaxSizeCheck< iEndIndex<=VectorExpr::iNumRows >::CheckType check_iEndIndex;typedef typename MaxSizeCheck< iStartIndex >=1 >::CheckType check_iStartIndex
Definition: matvec.h:833
const ScalarType * pGetVec() const
Definition: matvec.h:2489
index_type iGetNumRows() const
Definition: matvec.h:1190
Vector(const GradientExpression< Expr1 > &v1, const GradientExpression< Expr2 > &v2)
Definition: matvec.h:2284
ExpressionType operator()(index_type iRow) const
Definition: matvec.h:1055
Vector & operator-=(const Vector< T_Rhs, N_rows > &v)
Definition: matvec.h:2399
index_type iGetNumRows() const
Definition: matvec.h:1545
const ScalarType & operator()(index_type i, index_type j) const
Definition: matvec.h:1121
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:1554
Matrix(const MatrixExpression< Expression, N_rows, N_cols > &A)
Definition: matvec.h:1954
VectorExpression(const LhsExpr &u, const RhsExpr &v)
Definition: matvec.h:608
#define MATVEC_ASSERT(expr)
Definition: matvec.h:64
index_type iGetNumRows() const
Definition: matvec.h:2118
ScalarBinaryExpressionTraits< FuncMult, ScalarType, ExpressionTypeLhs, ExpressionTypeRhs >::ExpressionType ExprMult
Definition: matvec.h:3187
Matrix & operator=(const MatrixExpression< Expression, N_rows, N_cols > &A)
Definition: matvec.h:1994
ScalarBinFunc::ExpressionType ExpressionType
Definition: matvec.h:1232
Matrix(const Matrix< T2, N_rows, N_cols > &A, LocalDofMap *pDofMap)
Definition: matvec.h:1962
BasicScalarType< Expression >::ScalarType ScalarType
Definition: matvec.h:274
index_type iGetStartIndexLocal() const
Definition: matvec.h:132
MatrixExpr::RowVectorType ColumnVectorType
Definition: matvec.h:1110
GenericUnaryExpression(const ScalarExpr &expr)
Definition: matvec.h:177
index_type iGetNumRows() const
Definition: matvec.h:3229
index_type iGetNumRows() const
Definition: matvec.h:1586
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:3452
SumTraitsN_minus_1::ExpressionType ExpressionTypeN_minus_1
Definition: matvec.h:3054
static const index_type iDimension
Definition: matvec.h:107
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:754
ScalarUnaryFunc::ExpressionType ExpressionType
Definition: matvec.h:775
static const index_type iDimension
Definition: matvec.h:172
MatRInit(const VectorExpression< VectorExpr, 3 > &g)
Definition: matvec.h:2651
TransposedMatrix(const MatrixExpr &A)
Definition: matvec.h:1115
Vector & operator-=(const VectorExpression< Expression, N_rows > &f)
Definition: matvec.h:2419
MatrixDirectExpr(const MatrixType &A)
Definition: matvec.h:1420
index_type iGetNumCols() const
Definition: matvec.h:1383
index_type iGetStartIndexLocalVector() const
Definition: matvec.h:140
Vec3 GetCol(unsigned short int i) const
Definition: matvec3.h:903
const MatrixExpr A
Definition: matvec.h:1100
static const index_type iNumRows
Definition: matvec.h:638
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:3371
LocalDofMap * pGetDofMap() const
Definition: matvec.h:213
ScalarBinaryExpressionTraits< ScalarBinaryFunction, ScalarType, ScalarLhsExpr, ScalarRhsExpr >::ExpressionType ExpressionType
Definition: matvec.h:397
ExpressionType operator()(index_type i, index_type j) const
Definition: matvec.h:1245
MatrixLhsExpr::RowVectorType MatrixLhsRowVector
Definition: matvec.h:3406
integer index_type
Definition: gradient.h:104
RowVectorType GetRow(index_type iRow) const
Definition: matvec.h:1444
const T * pGetFirstElem() const
Definition: matvec.h:1669
MatrixData(index_type iNumRows, index_type iNumCols, bool bZeroInit)
Definition: matvec.h:1867
GradientExpression< BinaryExpr< FuncMult, LhsExpr, RhsExpr > > operator*(const GradientExpression< LhsExpr > &u, const GradientExpression< RhsExpr > &v)
Definition: gradient.h:2959
void Initialize(Matrix< T, 3, 3 > &A) const
Definition: matvec.h:2622
SliceVector(const ScalarType *p, index_type iRows, index_type iOffset)
Definition: matvec.h:979
MaxSizeCheck< N_offset!=DYNAMIC_SIZE >::CheckType check_iOffset
Definition: matvec.h:933
GenericUnaryExpression< ScalarUnaryFunction, T, ScalarExpr >::ScalarType ScalarType
Definition: matvec.h:294
VectorCrossExpr(const VectorLhsExpr &u, const VectorRhsExpr &v)
Definition: matvec.h:3207
IndexCheck< N_cols-VectorRhsExpr::iNumRows >::CheckType check_iNumRowsRhs
Definition: matvec.h:3327
index_type iGetEndIndexLocalVector() const
Definition: matvec.h:209
index_type iGetNumRows() const
Definition: matvec.h:821
static const index_type iNumCols
Definition: matvec.h:1298
CrossTraits< VectorLhsExpr, VectorRhsExpr >::ScalarType ScalarType
Definition: matvec.h:3204
ScalarTypeTraits< T >::DirectExpressionType ExpressionType
Definition: matvec.h:1013
ScalarBinFunc::ScalarType ScalarType
Definition: matvec.h:734
MatrixVectorProduct(const MatrixLhsExpr &A, const VectorRhsExpr &x)
Definition: matvec.h:3303
static ExpressionType Dot(const VectorExpression< VectorExprLhs, N_rows > &u, const VectorExpression< VectorExprRhs, N_rows > &v)
Definition: matvec.h:3107
VectorInit(const InitArg &R)
Definition: matvec.h:1623
index_type iGetNumRows() const
Definition: matvec.h:1490
std::vector< T, GradientAllocator< T > > rgData
Definition: matvec.h:1723
index_type iGetNumCols() const
Definition: matvec.h:1549
CommonScalarType< typename BasicScalarType< MatrixLhsScalarExpr >::ScalarType, typename BasicScalarType< VectorRhsScalarExpr >::ScalarType >::ScalarType ScalarType
Definition: matvec.h:3300
#define MATRIX_MATRIX_DEFINE_BINARY_FUNCTION(ExpressionName, FunctionName, FunctionClass)
Definition: matvec.h:3660
static const bool bAlias
Definition: matvec.h:104
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:960
index_type iGetNumRows() const
Definition: matvec.h:3429
void ZeroInit< long double >(long double *first, long double *last)
Definition: matvec.h:465
const MatrixLhsExpr A
Definition: matvec.h:3331
static void Eval(T &b, const U &a)
Definition: matvec.h:579
const VectorLhsExpr oU
Definition: matvec.h:761
void Initialize(Matrix< T, 3, 3 > &G) const
Definition: matvec.h:2732
const ScalarType & operator()(index_type i) const
Definition: matvec.h:882
Vector(const GradientExpression< Expr1 > &v1, const GradientExpression< Expr2 > &v2, const GradientExpression< Expr3 > &v3)
Definition: matvec.h:2302
void ApplyMatrixFuncNoAlias(const VectorExpression< Expression, N_rows > &A, const Func &)
Definition: matvec.h:2526
GradientExpression< UnaryExpr< FuncFabs, Expr > > fabs(const GradientExpression< Expr > &u)
Definition: gradient.h:2973
index_type iGetNumCols() const
Definition: matvec.h:3434
const ScalarType & operator()(index_type iRow) const
Definition: matvec.h:985
index_type iGetNumRows() const
Definition: matvec.h:2676
T Det(const Matrix< T, 2, 2 > &A)
Definition: matvec.h:3277
index_type iGetNumCols() const
Definition: matvec.h:1134
void ApplyMatrixFunc(const MatrixExpression< Expression, N_rows, N_cols > &A)
Definition: matvec.h:2152
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:862
index_type iGetNumCols() const
Definition: matvec.h:2119
MatrixInit< MatGInit< T, VectorDirectExpr< Vector< T, 3 > > >, T, 3, 3 > MatGVec(const Vector< T, 3 > &g)
Definition: matvec.h:2791
Matrix & operator/=(const T_Rhs &a)
Definition: matvec.h:2060
static const index_type iNumRows
Definition: matvec.h:1108
const ScalarType * pGetMat() const
Definition: matvec.h:2120
void ApplyMatrixFuncAlias(const VectorExpression< Expression, N_rows > &A, const Func &f)
Definition: matvec.h:2537
Vector & operator*=(const GradientExpression< ScalarExpression > &f)
Definition: matvec.h:2447
void Resize(index_type iNumRows)
Definition: matvec.h:1708
void Print(std::ostream &os) const
Definition: matvec.h:2996
MatrixExpression< SubMatrixExpr< iRowStart, iRowEnd, iColStart, iColEnd, MatrixExpression< MatrixDirectExpr< Matrix< T, N_rows, N_cols > >, N_rows, N_cols > >, iRowEnd-iRowStart+1, iColEnd-iColStart+1 > SubMatrix(const Matrix< T, N_rows, N_cols > &A)
Definition: matvec.h:2601
MaxSizeCheck< N_rows!=DYNAMIC_SIZE >::CheckType check_iNumRowsDynamic
Definition: matvec.h:3328
static const index_type iNumRows
Definition: matvec.h:2258
MatrixInit< MatRInit< T, VectorDirectExpr< Vector< T, 3 > > >, T, 3, 3 > MatRVec(const Vector< T, 3 > &g)
Definition: matvec.h:2803
void Resize(index_type iRows, index_type iCols)
Definition: matvec.h:1797
static const bool bAlias
Definition: matvec.h:1564
const doublereal d
Definition: matvec.h:2642
RowVectorType GetRow(index_type iRow) const
Definition: matvec.h:1265
index_type iGetNumRows() const
Definition: matvec.h:2638
index_type iGetNumRows() const
Definition: matvec.h:1061
VectorData< T, DYNAMIC_SIZE > rgData
Definition: matvec.h:1814
static const index_type iMaxDerivatives
Definition: matvec.h:105
VectorExpression(const Expr &u)
Definition: matvec.h:600
const ScalarType & operator()(index_type i, index_type j) const
Definition: matvec.h:1537
const MatrixLhsExpr A
Definition: matvec.h:3457
index_type iGetNumRows() const
Definition: matvec.h:623
ColumnVectorType GetCol(index_type iCol) const
Definition: matvec.h:1450
SubVectorExpr(const VectorExpr &u)
Definition: matvec.h:810
static const bool bAlias
Definition: matvec.h:904
index_type iGetStartIndexLocal() const
Definition: matvec.h:197
index_type iGetNumRows() const
Definition: matvec.h:1129
static ExpressionType f(const ScalarLhsExpr &u, const ScalarRhsExpr &v)
Definition: matvec.h:399
const VectorExpr oU
Definition: matvec.h:799
VectorLhsExpr::ExpressionType ExpressionTypeLhs
Definition: matvec.h:3179
Expression::ScalarType ScalarType
Definition: matvec.h:589
MatrixExpression(const Expression &e)
Definition: matvec.h:643
GradientExpression< BinaryExpr< ScalarBinaryFunction, DirectExpr< Gradient< N_SIZE > >, ScalarRhsExpr > > ExpressionType
Definition: matvec.h:347
void func(const T &u, const T &v, const T &w, doublereal e, T &f)
MatrixMatrixUnaryExpr(const MatrixExpr &u)
Definition: matvec.h:1364
static void RotCo(const T &phi, T &cf)
Definition: matvec.h:2945
void Copy(const Matrix< T2, N_rows, N_cols > &A, LocalDofMap *pDofMap)
Definition: matvec.h:1968
void Reset()
Definition: matvec.h:2356
ScalarBinFunc::ExpressionType ExpressionType
Definition: matvec.h:1300
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:1513
ColumnVectorType GetCol(index_type iCol) const
Definition: matvec.h:2112
MatrixExpression< MatrixDirectExpr< Matrix< T, N_rows, N_cols > >, N_rows, N_cols > Direct(const Matrix< T, N_rows, N_cols > &A)
Definition: matvec.h:2185
bool bArrayOverlap(const ScalarType *pFirstArray1, const ScalarType *pLastArray1, const ScalarType *pFirstArray2, const ScalarType *pLastArray2)
Definition: matvec.h:430
static void ApplyMatrixFunc(MatrixType &A, const Expression &B, const Func &f)
Definition: matvec.h:544
const MatrixRhsExpr B
Definition: matvec.h:3458
const VectorExpression< VectorExpr, 3 > g
Definition: matvec.h:2680
MatrixVectorProduct(const MatrixLhsExpr &A, const VectorRhsExpr &x)
Definition: matvec.h:3392
T rgData[N_rows][N_cols]
Definition: matvec.h:1769
MaxSizeCheck< N_offset!=DYNAMIC_SIZE >::CheckType check_iOffset
Definition: matvec.h:1003
vector_deriv_type dGetDerivativeLocalVector(index_type iLocalVecDof) const
Definition: matvec.h:128
MatrixExpression(const Expr &u)
Definition: matvec.h:651
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:794
static const index_type iInitNumCols
Definition: matvec.h:1915
vector_deriv_type dGetDerivativeLocalVector(index_type iLocalVecDof) const
Definition: matvec.h:193
Definition: matvec6.h:37
index_type iGetEndIndexLocalVector() const
Definition: matvec.h:144
const MatNxN & A
Definition: matvec.h:1600
void Reset()
Definition: matvec.h:1981
Matrix & operator+=(const MatrixExpression< Expression, N_rows, N_cols > &A)
Definition: matvec.h:2029
static const index_type iNumRows
Definition: matvec.h:1912
static const index_type iNumRows
Definition: matvec.h:1524
VectorExpression< SliceVector< T, N_cols, 1 >, N_cols > RowVectorType
Definition: matvec.h:1919
Expression::ExpressionType ExpressionType
Definition: matvec.h:590
ScalarType & operator()(index_type iRow, index_type iCol)
Definition: matvec.h:2096
VectorExpression< VectorExprRhs, N_rows >::ScalarType ScalarTypeRhs
Definition: matvec.h:3090
VectorData(const VectorData &oData)
Definition: matvec.h:1640
void Copy(const Vector< T2, N_rows > &v, LocalDofMap *pDofMap)
Definition: matvec.h:2346
index_type iGetNumRows() const
Definition: matvec.h:1253
index_type iGetEndIndexLocal() const
Definition: matvec.h:201
static index_type iGetMaxDerivatives()
Definition: matvec.h:221
MaxSizeCheck< iRowEnd<=MatrixExpr::iNumRows >::CheckType check_iRowEnd;typedef typename MaxSizeCheck< iColStart >=1 >::CheckType check_iColStart
Definition: matvec.h:1218
static const bool bVectorize
Definition: matvec.h:171
bool bHaveReferenceTo(const void *p) const
Definition: matvec.h:217
MatrixExpr::ScalarType ScalarType
Definition: matvec.h:1173
VectorVectorVectorBinaryExpr(const VectorLhsExpr &u, const VectorRhsExpr &v)
Definition: matvec.h:698
ColumnVectorType GetCol(index_type iCol) const
Definition: matvec.h:1394
const VectorRhsExpr oV
Definition: matvec.h:722
VectorData(index_type N, bool)
Definition: matvec.h:1684
static void Eval(T &b, const U &a)
Definition: matvec.h:551
const MatrixExpr oU
Definition: matvec.h:1406
ScalarBinaryExpressionTraits< FuncMult, ScalarType, ScalarExprLhs, ScalarExprRhs >::ExpressionType MultExpressionType
Definition: matvec.h:3098
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:1457
Matrix< T, 2, 2 > Inv(const Matrix< T, 2, 2 > &A)
Definition: matvec.h:3282
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:827
static const bool bAlias
Definition: matvec.h:1523
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:995
const ExpressionType a
Definition: matvec.h:163
void Reset(scalar_func_type &d)
Definition: gradient.h:2836
VectorRhsExpr::ExpressionType ExpressionTypeRhs
Definition: matvec.h:3181
index_type iGetNumRows() const
Definition: matvec.h:856
MatrixExpr::ExpressionType ExpressionType
Definition: matvec.h:1079
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:2123
Definition: mbdyn.h:76
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:2492
void Initialize(Matrix< T, 3, 3 > &A) const
Definition: matvec.h:2694
static const index_type iNumRows
Definition: matvec.h:806
ScalarBinaryExpressionTraits< FuncMult, ScalarType, ScalarExprLhs, ScalarExprRhs >::ExpressionType ExpressionType
Definition: matvec.h:3123
VectorExpression< SubVectorExpr< iStartIndex, iEndIndex, VectorDirectExpr< Vector< T, N_rows > > >, iEndIndex-iStartIndex+1 > SubVector(const Vector< T, N_rows > &v)
Definition: matvec.h:2585
BasicScalarType< ScalarExpr >::ScalarType ScalarType
Definition: matvec.h:407
static const bool bAlias
Definition: matvec.h:841
VectorDirectExpr(const VectorType &u)
Definition: matvec.h:846
void Resize(index_type iRows, index_type iCols)
Definition: matvec.h:1842
const ScalarType & operator()(index_type iRow) const
Definition: matvec.h:950
VectorExpression< ColumnVectorExpr< MatrixScalarMatrixBinaryExpr >, iNumRows > ColumnVectorType
Definition: matvec.h:1302
ColumnVectorType GetCol(index_type iCol) const
Definition: matvec.h:1335
Vector & operator*=(const T_Rhs &g)
Definition: matvec.h:2429
MatrixType::RowVectorType RowVectorType
Definition: matvec.h:1417
static const bool bVectorize
Definition: matvec.h:106
const ScalarType & operator()(index_type iRow) const
Definition: matvec.h:2466
RowVectorType GetRow(index_type iRow) const
Definition: matvec.h:1388
IndexCheck< VectorLhsExpr::iNumRows-3 >::CheckType check_iNumRowsLhs
Definition: matvec.h:3194
MatrixExpr::ScalarType ScalarType
Definition: matvec.h:1078
void Initialize(Matrix< T, 3, 3 > &RDelta) const
Definition: matvec.h:2655
T & operator[](index_type i)
Definition: matvec.h:1649
static const bool bAlias
Definition: matvec.h:772
index_type iGetNumRows() const
Definition: matvec.h:1661
scalar_func_type dGetValue() const
Definition: matvec.h:120
ExpressionType operator()(index_type i, index_type j) const
Definition: matvec.h:1370
static const index_type iInitNumRows
Definition: matvec.h:1914
VectorData(const VectorData< U, N_rows > &oData)
Definition: matvec.h:1645
VectorExpression< VectorExprLhs, N_rows >::ExpressionType ScalarExprLhs
Definition: matvec.h:3115
ScalarTypeTraits< T >::DirectExpressionType ExpressionType
Definition: matvec.h:2261
GenericUnaryExpression< ScalarUnaryFunction, T, ScalarExpr > ExpressionType
Definition: matvec.h:372
VectorData< T, DYNAMIC_SIZE > rgData
Definition: matvec.h:1906
Definition: mbdyn.h:77
VectorExpression< VectorExpressionType, N_rows >::ExpressionType ExpressionType
Definition: matvec.h:3066
const MatrixExpr A
Definition: matvec.h:1069
index_type iGetNumRows() const
Definition: matvec.h:888
static void Eval(T &b, const U &a)
Definition: matvec.h:565
Vector(const T &v1, const T &v2, const T &v3)
Definition: matvec.h:2292
SumTraits< VectorExpressionType, N_rows, N_index-1 > SumTraitsN_minus_1
Definition: matvec.h:3053
MaxSizeCheck< iNumRows!=DYNAMIC_SIZE >::CheckType check_iNumRows
Definition: matvec.h:1039
integer iGetNumCols(void) const
Definition: matvec3n.h:570
ScalarTypeTraits< T >::DirectExpressionType ExpressionType
Definition: matvec.h:942
#define VECTOR_DEFINE_UNARY_FUNCTION(ExpressionName, FunctionName, FunctionClass)
Definition: matvec.h:3634
void Copy(scalar_func_type &d1, const scalar_func_type &d2, LocalDofMap *)
Definition: gradient.h:2827
VectorExpression(const ScalarType *pArray, index_type iRows, index_type iOffset)
Definition: matvec.h:616
MatrixLhsExpr::ScalarType MatrixLhsScalarExpr
Definition: matvec.h:3296
ScalarTypeTraits< doublereal >::ScalarType ScalarType
Definition: matvec.h:1570
static const index_type iNumCols
Definition: matvec.h:639
VectorExpression< RowVectorExpr< MatrixMatrixMatrixBinaryExpr >, iNumCols > RowVectorType
Definition: matvec.h:1233
VectorExpression< SliceVector< doublereal, iNumCols, iNumRows >, iNumCols > RowVectorType
Definition: matvec.h:1474
ExpressionType operator()(index_type i) const
Definition: matvec.h:742
BasicScalarType< ScalarRhsExpr >::ScalarType ScalarTypeRhs
Definition: matvec.h:389
scalar_deriv_type dGetDerivativeLocal(index_type iLocalDof) const
Definition: matvec.h:124
ScalarTypeTraits< doublereal >::DirectExpressionType ExpressionType
Definition: matvec.h:1472
MaxSizeCheck< iNumRows!=DYNAMIC_SIZE >::CheckType check_iNumRows
Definition: matvec.h:932
const VectorLhsExpr oU
Definition: matvec.h:721
index_type iGetNumCols() const
Definition: matvec.h:1439
MatrixData(index_type iNumRows, index_type iNumCols, bool bZeroInit)
Definition: matvec.h:1821
const ScalarType * pGetFirstElem() const
Definition: matvec.h:2131
const MatrixRhsExpr oV
Definition: matvec.h:1280
const ScalarType & operator()(index_type i, index_type j) const
Definition: matvec.h:1182
MatGInit(const VectorExpression< VectorExpr, 3 > &g)
Definition: matvec.h:2727
integer iGetNumCols(void) const
Definition: matvec3n.h:298
ScalarTypeTraits< T >::ScalarType ScalarType
Definition: matvec.h:1916
static const bool bAlias
Definition: matvec.h:1045
DotTraits< VectorExprLhs, VectorExprRhs, N_rows, N_rows >::ExpressionType Dot(const VectorExpression< VectorExprLhs, N_rows > &u, const VectorExpression< VectorExprRhs, N_rows > &v)
Definition: matvec.h:3133
CommonScalarType< typename BasicScalarType< MatrixLhsScalarExpr >::ScalarType, typename BasicScalarType< VectorRhsScalarExpr >::ScalarType >::ScalarType ScalarType
Definition: matvec.h:3344
static void Eval(T &b, const U &a)
Definition: matvec.h:572
const VectorRhsExpr x
Definition: matvec.h:3332
const VectorType & oU
Definition: matvec.h:867
static const index_type iNumCols
Definition: matvec.h:1913
const Vec3 & oU
Definition: matvec.h:898
const ScalarType * pGetLastElem() const
Definition: matvec.h:2135
index_type iGetNumRows() const
Definition: matvec.h:666
ColumnVectorType GetCol(index_type iCol) const
Definition: matvec.h:1145
index_type iGetEndIndexLocal() const
Definition: matvec.h:136
index_type iGetNumCols() const
Definition: matvec.h:2754
static const bool bAlias
Definition: matvec.h:1356
GradientExpression< UnaryExpr< FuncSqrt, Expr > > sqrt(const GradientExpression< Expr > &u)
Definition: gradient.h:2974
static const bool bAlias
Definition: matvec.h:1167
ScalarBinFunc::ExpressionType ExpressionType
Definition: matvec.h:735
ColumnVectorExpr(const MatrixExpr &A, index_type iCol)
Definition: matvec.h:1050
MatrixData(index_type iNumRows, index_type iNumCols, bool bZeroInit)
Definition: matvec.h:1731
MatrixRhsExpr::ScalarType MatrixRhsScalarExpr
Definition: matvec.h:3405
CrossTraits< VectorLhsExpr, VectorRhsExpr >::ExprMult ExprMult
Definition: matvec.h:3203
index_type iGetNumRows() const
Definition: matvec.h:1740
VectorExpression< VectorCrossExpr< VectorLhsExpr, VectorRhsExpr >, 3 > Cross(const VectorExpression< VectorLhsExpr, 3 > &u, const VectorExpression< VectorRhsExpr, 3 > &v)
Definition: matvec.h:3248
index_type iGetNumRows() const
Definition: matvec.h:1704
static const index_type iNumRows
Definition: matvec.h:773
MatrixExpr::ExpressionType ExpressionType
Definition: matvec.h:1113
MatrixVectorProduct(const MatrixLhsExpr &A, const VectorRhsExpr &x)
Definition: matvec.h:3383
MatrixInit(const InitArg1 &v1, const InitArg2 &a2)
Definition: matvec.h:1613
ScalarUnaryExpressionTraits< ScalarUnaryFunction, ScalarType, ScalarExpr >::ExpressionType ExpressionType
Definition: matvec.h:411
static const index_type iNumRows
Definition: matvec.h:842
index_type iGetNumRows() const
Definition: matvec.h:2753
const VectorRhsExpr oV
Definition: matvec.h:3240
static ExpressionType f(const ScalarExpr &u)
Definition: matvec.h:413
static const bool bAlias
Definition: matvec.h:1412
const ScalarType & operator()(index_type iRow, index_type iCol) const
Definition: matvec.h:2088
Vector(const T &v1, const T &v2)
Definition: matvec.h:2275
MatrixType::ExpressionType ExpressionType
Definition: matvec.h:1416
scalar_func_type dGetValue() const
Definition: matvec.h:185
IndexCheck< iNumRows-Expression::iNumRows >::CheckType check_iNumRows
Definition: matvec.h:631
ScalarType & operator()(index_type iRow)
Definition: matvec.h:2473
VectorExpression< SubVectorExpr< iColStart, iColEnd, typename MatrixExpr::RowVectorType >, iNumCols > RowVectorType
Definition: matvec.h:1171
VectorData< ScalarType, N_rows > rgVec
Definition: matvec.h:2542
IndexCheck< iNumRows-VectorLhsExpr::iNumRows >::CheckType check_VectorLhsExpr
Definition: matvec.h:3242
VectorInit< VecRotInit< T, MatrixDirectExpr< Matrix< T, 3, 3 > > >, T, 3 > VecRotMat(const Matrix< T, 3, 3 > &R)
Definition: matvec.h:2984
VectorExpression< VectorExprLhs, N_rows >::ScalarType ScalarTypeLhs
Definition: matvec.h:3114
#define VECTOR_SCALAR_DEFINE_BINARY_FUNCTION(ExpressionName, FunctionName, FunctionClass)
Definition: matvec.h:3597
static bool bHaveReferenceTo(const GradientExpression< Expression > &g, const PointerType *pFirst, const PointerType *pLast)
Definition: matvec.h:314
VectorExpression< RowVectorExpr< MatrixScalarMatrixBinaryExpr >, iNumCols > RowVectorType
Definition: matvec.h:1301
static const index_type iNumCols
Definition: matvec.h:1230
const ScalarRhsExpr oV
Definition: matvec.h:1350
index_type iGetNumRows() const
Definition: matvec.h:2480
MatCrossInit(const VectorExpression< VectorExpr, 3 > &v, doublereal d)
Definition: matvec.h:2617
VectorExpression< SliceVector< T, N_rows, N_cols >, N_rows > ColumnVectorType
Definition: matvec.h:1920
static const index_type iNumCols
Definition: matvec.h:1525
Matrix & operator/=(const GradientExpression< ScalarExpression > &a)
Definition: matvec.h:2078
ScalarBinFunc::ExpressionType ExpressionType
Definition: matvec.h:696
static const index_type iNumRows
Definition: matvec.h:1077
const index_type iCol
Definition: matvec.h:1070
static index_type iGetMaxDerivatives()
Definition: matvec.h:156
Vector(const Vector< T2, N_rows > &v)
Definition: matvec.h:2328
ScalarType * pGetVec()
Definition: matvec.h:2488
index_type iGetNumRows() const
Definition: matvec.h:1092
static const index_type iNumRows
Definition: matvec.h:1046
Vector(const VectorInit< InitClass, T, N_rows > &func)
Definition: matvec.h:2322
VectorExpression< VectorExpr, N_rows >::ScalarType Norm(const VectorExpression< VectorExpr, N_rows > &u)
Definition: matvec.h:3163
RowVectorType GetRow(index_type iRow) const
Definition: matvec.h:1139
CommonScalarType< ScalarTypeLhs, ScalarTypeRhs >::ScalarType ScalarType
Definition: matvec.h:392
Vector & operator+=(const Vec3 &v)
Definition: matvec.h:2387
const index_type iRow
Definition: matvec.h:1101
Matrix & operator*=(const GradientExpression< ScalarExpression > &a)
Definition: matvec.h:2069
const ScalarType & operator()(index_type i, index_type j) const
Definition: matvec.h:1482
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:3320
RowVectorExpr(const MatrixExpr &A, index_type iRow)
Definition: matvec.h:1081
ScalarBinaryExpressionTraits< FuncPlus, ScalarType, ExpressionTypeN_minus_1, ScalarExpressionType >::ExpressionType ExpressionType
Definition: matvec.h:3055
ScalarTypeTraits< T >::ScalarType ScalarType
Definition: matvec.h:906
RowVectorType GetRow(index_type iRow) const
Definition: matvec.h:2105
ScalarTypeTraits< T >::ScalarType ScalarType
Definition: matvec.h:1012
SubMatrixExpr(const MatrixExpr &A)
Definition: matvec.h:1176
GradientExpression< DirectExpr< Gradient< N_SIZE >, true > > Alias(const Gradient< N_SIZE > &g)
Definition: gradient.h:2866
index_type iGetNumCols() const
Definition: matvec.h:1590
ScalarTypeTraits< T >::DirectExpressionType ExpressionType
Definition: matvec.h:977
index_type iGetNumCols() const
Definition: matvec.h:2716
VectorExpression< VectorExprRhs, N_rows >::ExpressionType ScalarExprRhs
Definition: matvec.h:3117
static const index_type iNumCols
Definition: matvec.h:1414
const doublereal * pGetMat(void) const
Definition: matvec3.h:743
static ExpressionType Sum(const VectorExpression< VectorExpressionType, N_rows > &v)
Definition: matvec.h:3069
const ScalarType & operator()(index_type iRow) const
Definition: matvec.h:1020
index_type iGetNumRows() const
Definition: matvec.h:1319
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:893
IndexCheck< MatrixLhsExpr::iNumCols-MatrixRhsExpr::iNumCols >::CheckType check_iNumCols
Definition: matvec.h:1284
const T & operator[](integer i) const
Definition: matvec.h:1655
static const index_type iNumCols
Definition: matvec.h:1169
ScalarTypeTraits< T >::DirectExpressionType ExpressionType
Definition: matvec.h:907
VectorExpression< VectorExprRhs, N_rows >::ExpressionType ScalarExprRhs
Definition: matvec.h:3091
IndexCheck< iNumRows-Expression::iNumRows >::CheckType check_iNumRows
Definition: matvec.h:682
ScalarTypeTraits< T >::DirectExpressionType ExpressionType
Definition: matvec.h:1917
const Mat3xN & A
Definition: matvec.h:1559
SliceVector(const ScalarType *p, index_type iRows, index_type iOffset)
Definition: matvec.h:944
VectorExpression< VectorExprRhs, N_rows >::ScalarType ScalarTypeRhs
Definition: matvec.h:3116
VectorData< T, DYNAMIC_SIZE > rgData
Definition: matvec.h:1860
#define VECTOR_VECTOR_DEFINE_BINARY_FUNCTION(ExpressionName, FunctionName, FunctionClass)
Definition: matvec.h:3554
VectorExpression< VectorExprLhs, N_rows >::ExpressionType ScalarExprLhs
Definition: matvec.h:3089
static const bool bAlias
Definition: matvec.h:1467
index_type iGetNumRows() const
Definition: matvec.h:2715
const ScalarType & operator()(index_type i, index_type j) const
Definition: matvec.h:1426
MatrixLhsExpr::RowVectorType MatrixLhsRowVector
Definition: matvec.h:3298
const doublereal * pGetVec(void) const
Definition: matvec3.h:192
CommonScalarType< ScalarTypeLhs, ScalarTypeRhs >::ScalarType ScalarType
Definition: matvec.h:3182
index_type iGetNumCols() const
Definition: matvec.h:1324
ScalarUnaryFunc::ScalarType ScalarType
Definition: matvec.h:774
Expression::ExpressionType ExpressionType
Definition: matvec.h:641
const VectorExpression< VectorExpr, 3 > g
Definition: matvec.h:2756
ColumnVectorType GetCol(index_type iCol) const
Definition: matvec.h:1269
ScalarTypeTraits< T >::ScalarType ScalarType
Definition: matvec.h:941
VectorExpression< ColumnVectorExpr< MatrixMatrixMatrixBinaryExpr >, iNumRows > ColumnVectorType
Definition: matvec.h:1234
SliceVector(const ScalarType *p, index_type iRows, index_type iOffset)
Definition: matvec.h:1015
GradientExpression< UnaryExpr< FuncAcos, Expr > > acos(const GradientExpression< Expr > &u)
Definition: gradient.h:2984
TabularMatrixView< T, N_rows, N_cols > Tabular(const Matrix< T, N_rows, N_cols > &A, int iColWidth=10)
Definition: matvec.h:3012
ScalarTypeTraits< doublereal >::DirectExpressionType ExpressionType
Definition: matvec.h:1571
T & operator()(index_type i, index_type j)
Definition: matvec.h:1743
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:716
static const index_type iNumCols
Definition: matvec.h:1109
MaxSizeCheck< VectorExpr::iNumRows!=DYNAMIC_SIZE >::CheckType check_iStaticSize
Definition: matvec.h:834
static const bool bAlias
Definition: matvec.h:805
VectorExpression< RowVectorExpr< MatrixMatrixUnaryExpr >, iNumCols > RowVectorType
Definition: matvec.h:1361
const T & operator[](integer i) const
Definition: matvec.h:1698
GradientExpression< UnaryExpr< ScalarUnaryFunction, DirectExpr< Gradient< N_SIZE > > > > ExpressionType
Definition: matvec.h:382
MatrixExpr::ExpressionType ExpressionType
Definition: matvec.h:1048
static ExpressionType Dot(const VectorExpression< VectorExprLhs, N_rows > &u, const VectorExpression< VectorExprRhs, N_rows > &v)
Definition: matvec.h:3126
ExpressionType operator()(index_type i) const
Definition: matvec.h:3308
MatrixType::ColumnVectorType ColumnVectorType
Definition: matvec.h:1418
SumTraits< VectorExpressionType, N_rows, N_rows >::ExpressionType Sum(const VectorExpression< VectorExpressionType, N_rows > &v)
Definition: matvec.h:3076
index_type iGetNumRows() const
Definition: matvec.h:922
index_type iGetNumCols() const
Definition: matvec.h:1259
index_type iGetNumCols() const
Definition: matvec.h:1494
static const index_type iMaxDerivatives
Definition: matvec.h:170
index_type iGetNumCols() const
Definition: matvec.h:1741
VectorExpression< VectorExprLhs, N_rows >::ScalarType ScalarTypeLhs
Definition: matvec.h:3088
GradientExpression< UnaryExpr< FuncCos, Expr > > cos(const GradientExpression< Expr > &u)
Definition: gradient.h:2978
MatrixExpr::ColumnVectorType RowVectorType
Definition: matvec.h:1111
VecRotInit(const MatrixExpression< MatrixExpr, 3, 3 > &R)
Definition: matvec.h:2818
Vector & operator/=(const GradientExpression< ScalarExpression > &f)
Definition: matvec.h:2456
DotTraits< MatrixLhsRowVector, MatrixRhsColumnVector, MatrixLhsExpr::iNumCols, MatrixLhsExpr::iNumCols >::ExpressionType ExpressionType
Definition: matvec.h:3412
static const int N
Definition: matvec.h:73
static const doublereal a
Definition: hfluid_.h:289
static const bool bAlias
Definition: matvec.h:1107
static const index_type iNumRows
Definition: matvec.h:1357
const MatrixLhsExpr oU
Definition: matvec.h:1279
const ScalarRhsExpr oV
Definition: matvec.h:762
MatrixInit< MatCrossCrossInit< T, VectorDirectExpr< Vector< T, 3 > > >, T, 3, 3 > MatCrossCrossVec(const Vector< T, 3 > &v)
Definition: matvec.h:2779
VectorLhsExpr::ScalarType ScalarTypeLhs
Definition: matvec.h:3178
Matrix & operator-=(const MatrixExpression< Expression, N_rows, N_cols > &A)
Definition: matvec.h:2040
GradientExpression< BinaryExpr< ScalarBinaryFunction, DirectExpr< Gradient< N_SIZE > >, DirectExpr< Gradient< N_SIZE > > > > ExpressionType
Definition: matvec.h:357
static const bool bAlias
Definition: matvec.h:169
Vector(const Vector< T2, N_rows > &v, LocalDofMap *pDofMap)
Definition: matvec.h:2339
Matrix(const T &A11, const T &A21, const T &A12, const T &A22)
Definition: matvec.h:1930
doublereal scalar_func_type
Definition: gradient.h:346
const ExpressionType a
Definition: matvec.h:228
GradientExpression< UnaryExpr< ScalarUnaryFunction, ScalarExpr > > ExpressionType
Definition: matvec.h:377
IndexCheck< MatrixLhsExpr::iNumRows-MatrixRhsExpr::iNumRows >::CheckType check_iNumRows
Definition: matvec.h:1283
const MatrixExpr A
Definition: matvec.h:1157
void ApplyMatrixFunc(const VectorExpression< Expression, N_rows > &A)
Definition: matvec.h:2519
static bool bHaveReferenceTo(const Expression &, const PointerType *, const PointerType *)
Definition: matvec.h:303
T & operator[](index_type i)
Definition: matvec.h:1692
void Convert(scalar_func_type &d, const Gradient< N_SIZE > &g)
Definition: gradient.h:2855
GradientExpression< BinaryExpr< FuncAtan2, LhsExpr, RhsExpr > > atan2(const GradientExpression< LhsExpr > &u, const GradientExpression< RhsExpr > &v)
Definition: gradient.h:2962
GradientExpression< BinaryExpr< ScalarBinaryFunction, DirectExpr< Gradient< N_SIZE > >, ConstExpr< Gradient< N_SIZE > > > > ExpressionType
Definition: matvec.h:362
GradientExpression< BinaryExpr< ScalarBinaryFunction, ScalarLhsExpr, ConstExpr< Gradient< N_SIZE > > > > ExpressionType
Definition: matvec.h:342
double doublereal
Definition: colamd.c:52
static const bool bAlias
Definition: matvec.h:872
index_type iGetNumCols() const
Definition: matvec.h:669
const T * pGetLastElem() const
Definition: matvec.h:1718
ScalarTypeTraits< T >::ScalarType ScalarType
Definition: matvec.h:2260
MatrixExpression(const LhsExpr &u, const RhsExpr &v)
Definition: matvec.h:659
Matrix & operator*=(const T_Rhs &a)
Definition: matvec.h:2051
index_type iGetNumCols() const
Definition: matvec.h:2677
static const index_type iNumRows
Definition: matvec.h:1229
long int integer
Definition: colamd.c:51
const MatrixLhsExpr oU
Definition: matvec.h:1349
MatrixInit< MatCrossInit< T, VectorDirectExpr< Vector< T, 3 > > >, T, 3, 3 > MatCrossVec(const Vector< T, 3 > &v, doublereal d=0.)
Definition: matvec.h:2761
const MatrixType & A
Definition: matvec.h:1462
void ZeroInit(T *first, T *last)
Definition: matvec.h:450
LocalDofMap * pGetDofMap() const
Definition: matvec.h:148
const VectorExpression< VectorExpr, 3 > v
Definition: matvec.h:2718
#define MATRIX_DEFINE_UNARY_FUNCTION(ExpressionName, FunctionName, FunctionClass)
Definition: matvec.h:3647
VectorExpression< SliceVector< doublereal, iNumRows, 1 >, iNumRows > ColumnVectorType
Definition: matvec.h:1475
GradientExpression< BinaryExpr< ScalarBinaryFunction, ConstExpr< Gradient< N_SIZE > >, ScalarRhsExpr > > ExpressionType
Definition: matvec.h:352
ScalarTypeTraits< doublereal >::DirectExpressionType ExpressionType
Definition: matvec.h:875
CrossTraits< VectorLhsExpr, VectorRhsExpr >::ExpressionType ExpressionType
Definition: matvec.h:3205
TabularMatrixView(const Matrix< T, N_rows, N_cols > &A, int iColWidth)
Definition: matvec.h:2991
ColumnVectorType GetCol(index_type iCol) const
Definition: matvec.h:1505
ScalarTypeTraits< doublereal >::ScalarType ScalarType
Definition: matvec.h:1529
MatrixType::ScalarType ScalarType
Definition: matvec.h:1415
const Matrix< T, N_rows, N_cols > & A
Definition: matvec.h:3006
void Resize(index_type iNumRows)
Definition: matvec.h:2362
static void ApplyMatrixFunc(MatrixType &A, const Expression &B, const Func &f)
Definition: matvec.h:536
void ApplyMatrixFuncNoAlias(const MatrixExpression< Expression, N_rows, N_cols > &A, const Func &)
Definition: matvec.h:2159
VectorType::ExpressionType ExpressionType
Definition: matvec.h:844
const VectorLhsExpr oU
Definition: matvec.h:3239
ExpressionType operator()(index_type i) const
Definition: matvec.h:782
ScalarUnaryFunc::ScalarType ScalarType
Definition: matvec.h:1359
static const index_type iNumRows
Definition: matvec.h:733
static const index_type iNumRows
Definition: matvec.h:1297
IndexCheck< iNumCols-Expression::iNumCols >::CheckType check_iNumCols
Definition: matvec.h:683
CommonScalarType< typename BasicScalarType< MatrixLhsScalarExpr >::ScalarType, typename BasicScalarType< MatrixRhsScalarExpr >::ScalarType >::ScalarType ScalarType
Definition: matvec.h:3411
VectorRhsExpr::ScalarType VectorRhsScalarExpr
Definition: matvec.h:3297
MatrixMatrixProduct(const MatrixLhsExpr &A, const MatrixRhsExpr &B)
Definition: matvec.h:3414
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:1095
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:1342
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:1401
MatrixData(index_type iNumRows, index_type iNumCols, bool bZeroInit)
Definition: matvec.h:1776
VectorExpression< SubVectorExpr< iRowStart, iRowEnd, typename MatrixExpr::ColumnVectorType >, iNumRows > ColumnVectorType
Definition: matvec.h:1172
void array_fill(T *first, T *const last, const T &val)
Definition: gradient.h:309
VectorType::ScalarType ScalarType
Definition: matvec.h:843
void ZeroInit< double >(double *first, double *last)
Definition: matvec.h:460
IndexCheck< N_cols-MatrixLhsExpr::iNumCols >::CheckType check_iNumColsLhs
Definition: matvec.h:3326
SliceVector(const ScalarType *p, index_type iRows, index_type iOffset)
Definition: matvec.h:909
index_type iGetNumRows() const
Definition: matvec.h:1434
Mat3x3 R
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:1152
ExpressionType operator()(index_type i) const
Definition: matvec.h:703
GradientExpression< BinaryExpr< ScalarBinaryFunction, ConstExpr< Gradient< N_SIZE > >, DirectExpr< Gradient< N_SIZE > > > > ExpressionType
Definition: matvec.h:367
Vector & operator/=(const T_Rhs &g)
Definition: matvec.h:2438
Vector(const Vector &oVec)
Definition: matvec.h:2271
MatrixMatrixMatrixBinaryExpr(const MatrixLhsExpr &u, const MatrixRhsExpr &v)
Definition: matvec.h:1236
static const bool bAlias
Definition: matvec.h:637
Matrix & operator-=(const Matrix< T_Rhs, N_rows, N_cols > &A)
Definition: matvec.h:2017
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:1030
static const index_type iNumRows
Definition: matvec.h:694
MatCrossCrossInit(const VectorExpression< VectorExpr, 3 > &v)
Definition: matvec.h:2689
ExpressionType operator()(index_type i, index_type j) const
Definition: matvec.h:1311
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:1064
GenericBinaryExpression< ScalarBinaryFunction, T, ScalarLhsExpr, ScalarRhsExpr > ExpressionType
Definition: matvec.h:327
VectorExpression< VectorExpressionType, N_rows >::ExpressionType ScalarExpressionType
Definition: matvec.h:3052
bool bHaveReferenceTo(const ScalarType2 *pFirst, const ScalarType2 *pLast) const
Definition: matvec.h:1211
void ZeroInit< float >(float *first, float *last)
Definition: matvec.h:455
GradientExpression< DirectExpr< Gradient< N_SIZE > > > DirectExpressionType
Definition: matvec.h:311
VectorExpression< VectorExpressionType, N_rows >::ScalarType ScalarType
Definition: matvec.h:3051
const VectorExpr oU
Definition: matvec.h:835
VectorScalarVectorBinaryExpr(const VectorLhsExpr &u, const ScalarRhsExpr &v)
Definition: matvec.h:737
MatrixInit(const InitArg &v)
Definition: matvec.h:1607
Vector & operator+=(const Vector< T_Rhs, N_rows > &v)
Definition: matvec.h:2378
static const bool bAlias
Definition: matvec.h:693
MatrixExpr::ScalarType ScalarType
Definition: matvec.h:1112
const T * pGetLastElem() const
Definition: matvec.h:1673
Vec3DirectExpr(const Vec3 &u)
Definition: matvec.h:877
static const index_type iNumRows
Definition: matvec.h:905
const VectorExpression< VectorExpr, 3 > v
Definition: matvec.h:2641
ScalarBinFunc::ScalarType ScalarType
Definition: matvec.h:1299