type_vector.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_TYPE_VECTOR_H
21 #define LIBMESH_TYPE_VECTOR_H
22 
23 // Local includes
24 #include "libmesh/libmesh_common.h"
25 #include "libmesh/compare_types.h"
26 #include "libmesh/tensor_tools.h"
27 
28 // C++ includes
29 #include <cstdlib> // *must* precede <cmath> for proper std:abs() on PGI, Sun Studio CC
30 #include <cmath>
31 #include <complex>
32 
33 namespace libMesh
34 {
35 
36 // Forward declarations
37 template <typename T> class TypeTensor;
38 template <typename T> class VectorValue;
39 template <typename T> class TensorValue;
40 
52 template <typename T>
53 class TypeVector
54 {
55  template <typename T2>
56  friend class TypeVector;
57 
58  friend class TypeTensor<T>;
59 
60 protected:
61 
66  TypeVector ();
67 
68 
73  TypeVector (const T x,
74  const T y=0,
75  const T z=0);
76 
77 
82  template <typename Scalar1, typename Scalar2, typename Scalar3>
83  TypeVector (typename
85  const Scalar1>::type x,
86  typename
88  const Scalar2>::type y=0,
89  typename
91  const Scalar3>::type z=0);
92 
93 
100  template <typename Scalar>
101  TypeVector (const Scalar x,
102  typename
104  const Scalar>::type * sfinae = libmesh_nullptr);
105 
106 public:
107 
111  typedef T value_type;
112 
116  template <typename T2>
117  TypeVector (const TypeVector<T2> & p);
118 
122  ~TypeVector ();
123 
127  template <typename T2>
128  void assign (const TypeVector<T2> &);
129 
133  template <typename Scalar>
134  typename boostcopy::enable_if_c<
136  TypeVector &>::type
137  operator = (const Scalar & libmesh_dbg_var(p))
138  { libmesh_assert_equal_to (p, Scalar(0)); this->zero(); return *this; }
139 
143  const T & operator () (const unsigned int i) const;
144 
145  const T & slice (const unsigned int i) const { return (*this)(i); }
146 
150  T & operator () (const unsigned int i);
151 
152  T & slice (const unsigned int i) { return (*this)(i); }
153 
157  template <typename T2>
159  operator + (const TypeVector<T2> &) const;
160 
164  template <typename T2>
165  const TypeVector<T> & operator += (const TypeVector<T2> &);
166 
170  template <typename T2>
171  void add (const TypeVector<T2> &);
172 
177  template <typename T2>
178  void add_scaled (const TypeVector<T2> &, const T);
179 
183  template <typename T2>
185  operator - (const TypeVector<T2> &) const;
186 
190  template <typename T2>
191  const TypeVector<T> & operator -= (const TypeVector<T2> &);
192 
196  template <typename T2>
197  void subtract (const TypeVector<T2> &);
198 
203  template <typename T2>
204  void subtract_scaled (const TypeVector<T2> &, const T);
205 
209  TypeVector<T> operator - () const;
210 
214  template <typename Scalar>
215  typename boostcopy::enable_if_c<
216  ScalarTraits<Scalar>::value,
218  operator * (const Scalar) const;
219 
223  const TypeVector<T> & operator *= (const T);
224 
228  template <typename Scalar>
229  typename boostcopy::enable_if_c<
230  ScalarTraits<Scalar>::value,
232  operator / (const Scalar) const;
233 
237  const TypeVector<T> & operator /= (const T);
238 
243  template <typename T2>
245  operator * (const TypeVector<T2> &) const;
246 
251  template <typename T2>
253  contract (const TypeVector<T2> &) const;
254 
258  template <typename T2>
260  cross(const TypeVector<T2> &) const;
261 
266  TypeVector<T> unit() const;
267 
273  Real size() const;
274 
279  Real norm() const;
280 
286  Real size_sq() const;
287 
292  Real norm_sq() const;
293 
297  void zero();
298 
303  bool relative_fuzzy_equals(const TypeVector<T> & rhs,
304  Real tol = TOLERANCE) const;
305 
310  bool absolute_fuzzy_equals(const TypeVector<T> & rhs,
311  Real tol = TOLERANCE) const;
312 
318  bool operator == (const TypeVector<T> & rhs) const;
319 
324  bool operator != (const TypeVector<T> & rhs) const;
325 
332  bool operator < (const TypeVector<T> & rhs) const;
333 
340  bool operator <= (const TypeVector<T> & rhs) const;
341 
348  bool operator > (const TypeVector<T> & rhs) const;
349 
356  bool operator >= (const TypeVector<T> & rhs) const;
357 
361  void print(std::ostream & os = libMesh::out) const;
362 
368  friend std::ostream & operator << (std::ostream & os, const TypeVector<T> & t)
369  {
370  t.print(os);
371  return os;
372  }
373 
379  void write_unformatted (std::ostream & out, const bool newline = true) const;
380 
381 protected:
382 
386  T _coords[LIBMESH_DIM];
387 };
388 
389 
390 
391 
392 
393 //------------------------------------------------------
394 // Inline functions
395 
396 template <typename T>
397 inline
399 {
400  _coords[0] = 0;
401 
402 #if LIBMESH_DIM > 1
403  _coords[1] = 0;
404 #endif
405 
406 #if LIBMESH_DIM > 2
407  _coords[2] = 0;
408 #endif
409 }
410 
411 
412 
413 template <typename T>
414 inline
416  const T y,
417  const T z)
418 {
419  _coords[0] = x;
420 
421 #if LIBMESH_DIM > 1
422  _coords[1] = y;
423 #else
424  libmesh_assert_equal_to (y, 0);
425 #endif
426 
427 #if LIBMESH_DIM > 2
428  _coords[2] = z;
429 #else
430  libmesh_assert_equal_to (z, 0);
431 #endif
432 }
433 
434 
435 template <typename T>
436 template <typename Scalar1, typename Scalar2, typename Scalar3>
437 inline
440  const Scalar1>::type x,
441  typename
443  const Scalar2>::type y,
444  typename
446  const Scalar3>::type z)
447 {
448  _coords[0] = x;
449 
450 #if LIBMESH_DIM > 1
451  _coords[1] = y;
452 #else
453  libmesh_assert_equal_to (y, 0);
454 #endif
455 
456 #if LIBMESH_DIM > 2
457  _coords[2] = z;
458 #else
459  libmesh_assert_equal_to (z, 0);
460 #endif
461 }
462 
463 
464 
465 template <typename T>
466 template <typename Scalar>
467 inline
469  typename
470  boostcopy::enable_if_c<ScalarTraits<Scalar>::value,
471  const Scalar>::type * /*sfinae*/)
472 {
473  _coords[0] = x;
474 
475 #if LIBMESH_DIM > 1
476  _coords[1] = 0;
477 #endif
478 
479 #if LIBMESH_DIM > 2
480  _coords[2] = 0;
481 #endif
482 }
483 
484 
485 
486 template <typename T>
487 template <typename T2>
488 inline
490 {
491  // copy the nodes from vector p to me
492  for (unsigned int i=0; i<LIBMESH_DIM; i++)
493  _coords[i] = p._coords[i];
494 }
495 
496 
497 
498 template <typename T>
499 inline
501 {
502 }
503 
504 
505 
506 template <typename T>
507 template <typename T2>
508 inline
510 {
511  for (unsigned int i=0; i<LIBMESH_DIM; i++)
512  _coords[i] = p._coords[i];
513 }
514 
515 
516 
517 template <typename T>
518 inline
519 const T & TypeVector<T>::operator () (const unsigned int i) const
520 {
521  libmesh_assert_less (i, LIBMESH_DIM);
522 
523  return _coords[i];
524 }
525 
526 
527 
528 template <typename T>
529 inline
530 T & TypeVector<T>::operator () (const unsigned int i)
531 {
532  libmesh_assert_less (i, LIBMESH_DIM);
533 
534  return _coords[i];
535 }
536 
537 
538 
539 template <typename T>
540 template <typename T2>
541 inline
544 {
545  typedef typename CompareTypes<T, T2>::supertype TS;
546 #if LIBMESH_DIM == 1
547  return TypeVector<TS> (_coords[0] + p._coords[0]);
548 #endif
549 
550 #if LIBMESH_DIM == 2
551  return TypeVector<TS> (_coords[0] + p._coords[0],
552  _coords[1] + p._coords[1]);
553 #endif
554 
555 #if LIBMESH_DIM == 3
556  return TypeVector<TS> (_coords[0] + p._coords[0],
557  _coords[1] + p._coords[1],
558  _coords[2] + p._coords[2]);
559 #endif
560 
561 }
562 
563 
564 
565 template <typename T>
566 template <typename T2>
567 inline
569 {
570  this->add (p);
571 
572  return *this;
573 }
574 
575 
576 
577 template <typename T>
578 template <typename T2>
579 inline
581 {
582 #if LIBMESH_DIM == 1
583  _coords[0] += p._coords[0];
584 #endif
585 
586 #if LIBMESH_DIM == 2
587  _coords[0] += p._coords[0];
588  _coords[1] += p._coords[1];
589 #endif
590 
591 #if LIBMESH_DIM == 3
592  _coords[0] += p._coords[0];
593  _coords[1] += p._coords[1];
594  _coords[2] += p._coords[2];
595 #endif
596 
597 }
598 
599 
600 
601 template <typename T>
602 template <typename T2>
603 inline
604 void TypeVector<T>::add_scaled (const TypeVector<T2> & p, const T factor)
605 {
606 #if LIBMESH_DIM == 1
607  _coords[0] += factor*p(0);
608 #endif
609 
610 #if LIBMESH_DIM == 2
611  _coords[0] += factor*p(0);
612  _coords[1] += factor*p(1);
613 #endif
614 
615 #if LIBMESH_DIM == 3
616  _coords[0] += factor*p(0);
617  _coords[1] += factor*p(1);
618  _coords[2] += factor*p(2);
619 #endif
620 
621 }
622 
623 
624 
625 template <typename T>
626 template <typename T2>
627 inline
630 {
631  typedef typename CompareTypes<T, T2>::supertype TS;
632 
633 #if LIBMESH_DIM == 1
634  return TypeVector<TS>(_coords[0] - p._coords[0]);
635 #endif
636 
637 #if LIBMESH_DIM == 2
638  return TypeVector<TS>(_coords[0] - p._coords[0],
639  _coords[1] - p._coords[1]);
640 #endif
641 
642 #if LIBMESH_DIM == 3
643  return TypeVector<TS>(_coords[0] - p._coords[0],
644  _coords[1] - p._coords[1],
645  _coords[2] - p._coords[2]);
646 #endif
647 
648 }
649 
650 
651 
652 template <typename T>
653 template <typename T2>
654 inline
656 {
657  this->subtract (p);
658 
659  return *this;
660 }
661 
662 
663 
664 template <typename T>
665 template <typename T2>
666 inline
668 {
669  for (unsigned int i=0; i<LIBMESH_DIM; i++)
670  _coords[i] -= p._coords[i];
671 }
672 
673 
674 
675 template <typename T>
676 template <typename T2>
677 inline
678 void TypeVector<T>::subtract_scaled (const TypeVector<T2> & p, const T factor)
679 {
680  for (unsigned int i=0; i<LIBMESH_DIM; i++)
681  _coords[i] -= factor*p(i);
682 }
683 
684 
685 
686 template <typename T>
687 inline
689 {
690 
691 #if LIBMESH_DIM == 1
692  return TypeVector(-_coords[0]);
693 #endif
694 
695 #if LIBMESH_DIM == 2
696  return TypeVector(-_coords[0],
697  -_coords[1]);
698 #endif
699 
700 #if LIBMESH_DIM == 3
701  return TypeVector(-_coords[0],
702  -_coords[1],
703  -_coords[2]);
704 #endif
705 
706 }
707 
708 
709 
710 template <typename T>
711 template <typename Scalar>
712 inline
713 typename boostcopy::enable_if_c<
714  ScalarTraits<Scalar>::value,
716 TypeVector<T>::operator * (const Scalar factor) const
717 {
718  typedef typename CompareTypes<T, Scalar>::supertype SuperType;
719 
720 #if LIBMESH_DIM == 1
721  return TypeVector<SuperType>(_coords[0]*factor);
722 #endif
723 
724 #if LIBMESH_DIM == 2
725  return TypeVector<SuperType>(_coords[0]*factor,
726  _coords[1]*factor);
727 #endif
728 
729 #if LIBMESH_DIM == 3
730  return TypeVector<SuperType>(_coords[0]*factor,
731  _coords[1]*factor,
732  _coords[2]*factor);
733 #endif
734 }
735 
736 
737 
738 template <typename T, typename Scalar>
739 inline
740 typename boostcopy::enable_if_c<
741  ScalarTraits<Scalar>::value,
743 operator * (const Scalar factor,
744  const TypeVector<T> & v)
745 {
746  return v * factor;
747 }
748 
749 
750 
751 template <typename T>
752 inline
754 {
755 #if LIBMESH_DIM == 1
756  _coords[0] *= factor;
757 #endif
758 
759 #if LIBMESH_DIM == 2
760  _coords[0] *= factor;
761  _coords[1] *= factor;
762 #endif
763 
764 #if LIBMESH_DIM == 3
765  _coords[0] *= factor;
766  _coords[1] *= factor;
767  _coords[2] *= factor;
768 #endif
769 
770  return *this;
771 }
772 
773 
774 
775 template <typename T>
776 template <typename Scalar>
777 inline
778 typename boostcopy::enable_if_c<
779  ScalarTraits<Scalar>::value,
781 TypeVector<T>::operator / (const Scalar factor) const
782 {
783  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
784 
785  typedef typename CompareTypes<T, Scalar>::supertype TS;
786 
787 #if LIBMESH_DIM == 1
788  return TypeVector<TS>(_coords[0]/factor);
789 #endif
790 
791 #if LIBMESH_DIM == 2
792  return TypeVector<TS>(_coords[0]/factor,
793  _coords[1]/factor);
794 #endif
795 
796 #if LIBMESH_DIM == 3
797  return TypeVector<TS>(_coords[0]/factor,
798  _coords[1]/factor,
799  _coords[2]/factor);
800 #endif
801 
802 }
803 
804 
805 
806 
807 template <typename T>
808 inline
809 const TypeVector<T> &
811 {
812  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
813 
814  for (unsigned int i=0; i<LIBMESH_DIM; i++)
815  _coords[i] /= factor;
816 
817  return *this;
818 }
819 
820 
821 
822 
823 template <typename T>
824 template <typename T2>
825 inline
828 {
829 #if LIBMESH_DIM == 1
830  return _coords[0]*p._coords[0];
831 #endif
832 
833 #if LIBMESH_DIM == 2
834  return (_coords[0]*p._coords[0] +
835  _coords[1]*p._coords[1]);
836 #endif
837 
838 #if LIBMESH_DIM == 3
839  return (_coords[0]*p(0) +
840  _coords[1]*p(1) +
841  _coords[2]*p(2));
842 #endif
843 }
844 
845 template <typename T>
846 template <typename T2>
847 inline
850 {
851  return (*this)*(p);
852 }
853 
854 
855 
856 template <typename T>
857 template <typename T2>
860 {
861  typedef typename CompareTypes<T, T2>::supertype TS;
862  libmesh_assert_equal_to (LIBMESH_DIM, 3);
863 
864  // | i j k |
865  // |(*this)(0) (*this)(1) (*this)(2)|
866  // | p(0) p(1) p(2) |
867 
868  return TypeVector<TS>( _coords[1]*p._coords[2] - _coords[2]*p._coords[1],
869  -_coords[0]*p._coords[2] + _coords[2]*p._coords[0],
870  _coords[0]*p._coords[1] - _coords[1]*p._coords[0]);
871 }
872 
873 
874 
875 template <typename T>
876 inline
878 {
879  libmesh_deprecated();
880  return this->norm();
881 }
882 
883 
884 
885 template <typename T>
886 inline
888 {
889  return std::sqrt(this->norm_sq());
890 }
891 
892 
893 
894 template <typename T>
895 inline
897 {
898  for (unsigned int i=0; i<LIBMESH_DIM; i++)
899  _coords[i] = 0.;
900 }
901 
902 
903 
904 template <typename T>
905 inline
907 {
908  libmesh_deprecated();
909  return this->norm_sq();
910 }
911 
912 
913 
914 template <typename T>
915 inline
917 {
918 #if LIBMESH_DIM == 1
919  return (TensorTools::norm_sq(_coords[0]));
920 #endif
921 
922 #if LIBMESH_DIM == 2
923  return (TensorTools::norm_sq(_coords[0]) +
925 #endif
926 
927 #if LIBMESH_DIM == 3
928  return (TensorTools::norm_sq(_coords[0]) +
931 #endif
932 }
933 
934 
935 
936 template <typename T>
937 inline
939 {
940 #if LIBMESH_DIM == 1
941  return (std::abs(_coords[0] - rhs._coords[0])
942  <= tol);
943 #endif
944 
945 #if LIBMESH_DIM == 2
946  return (std::abs(_coords[0] - rhs._coords[0]) +
947  std::abs(_coords[1] - rhs._coords[1])
948  <= tol);
949 #endif
950 
951 #if LIBMESH_DIM == 3
952  return (std::abs(_coords[0] - rhs._coords[0]) +
953  std::abs(_coords[1] - rhs._coords[1]) +
954  std::abs(_coords[2] - rhs._coords[2])
955  <= tol);
956 #endif
957 }
958 
959 
960 
961 template <typename T>
962 inline
964 {
965 #if LIBMESH_DIM == 1
966  return this->absolute_fuzzy_equals(rhs, tol *
967  (std::abs(_coords[0]) + std::abs(rhs._coords[0])));
968 #endif
969 
970 #if LIBMESH_DIM == 2
971  return this->absolute_fuzzy_equals(rhs, tol *
972  (std::abs(_coords[0]) + std::abs(rhs._coords[0]) +
973  std::abs(_coords[1]) + std::abs(rhs._coords[1])));
974 #endif
975 
976 #if LIBMESH_DIM == 3
977  return this->absolute_fuzzy_equals(rhs, tol *
978  (std::abs(_coords[0]) + std::abs(rhs._coords[0]) +
979  std::abs(_coords[1]) + std::abs(rhs._coords[1]) +
980  std::abs(_coords[2]) + std::abs(rhs._coords[2])));
981 #endif
982 }
983 
984 
985 
986 template <typename T>
987 inline
989 {
990 #if LIBMESH_DIM == 1
991  return (_coords[0] == rhs._coords[0]);
992 #endif
993 
994 #if LIBMESH_DIM == 2
995  return (_coords[0] == rhs._coords[0] &&
996  _coords[1] == rhs._coords[1]);
997 #endif
998 
999 #if LIBMESH_DIM == 3
1000  return (_coords[0] == rhs._coords[0] &&
1001  _coords[1] == rhs._coords[1] &&
1002  _coords[2] == rhs._coords[2]);
1003 #endif
1004 }
1005 
1006 
1007 
1008 template <typename T>
1009 inline
1011 {
1012  return (!(*this == rhs));
1013 }
1014 
1015 
1016 //------------------------------------------------------
1017 // Non-member functions on TypeVectors
1018 
1019 // Compute a * (b.cross(c)) without creating a temporary
1020 // for the cross product. Equivalent to the determinant
1021 // of the 3x3 tensor:
1022 // [a0, a1, a2]
1023 // [b0, b1, b2]
1024 // [c0, c1, c2]
1025 template <typename T>
1026 inline
1028  const TypeVector<T> & b,
1029  const TypeVector<T> & c)
1030 {
1031  // We only support cross products when LIBMESH_DIM==3, same goes for this.
1032  libmesh_assert_equal_to (LIBMESH_DIM, 3);
1033 
1034  return
1035  a(0)*(b(1)*c(2) - b(2)*c(1)) -
1036  a(1)*(b(0)*c(2) - b(2)*c(0)) +
1037  a(2)*(b(0)*c(1) - b(1)*c(0));
1038 }
1039 
1040 
1041 
1046 template <typename T>
1047 inline
1049  const TypeVector<T> & c)
1050 {
1051  // We only support cross products when LIBMESH_DIM==3, same goes for this.
1052  libmesh_assert_equal_to (LIBMESH_DIM, 3);
1053 
1054  T x = b(1)*c(2) - b(2)*c(1),
1055  y = b(0)*c(2) - b(2)*c(0),
1056  z = b(0)*c(1) - b(1)*c(0);
1057 
1058  return x*x + y*y + z*z;
1059 }
1060 
1061 
1062 
1066 template <typename T>
1067 inline
1069  const TypeVector<T> & c)
1070 {
1071  // We only support cross products when LIBMESH_DIM==3, same goes for this.
1072  libmesh_assert_equal_to (LIBMESH_DIM, 3);
1073 
1074  return std::sqrt(cross_norm_sq(b,c));
1075 }
1076 
1077 } // namespace libMesh
1078 
1079 #endif // LIBMESH_TYPE_VECTOR_H
bool operator>=(const TypeVector< T > &rhs) const
Definition: type_vector.C:152
T _coords[LIBMESH_DIM]
Definition: type_vector.h:386
double abs(double a)
bool operator==(const TypeVector< T > &rhs) const
Definition: type_vector.h:988
const T & slice(const unsigned int i) const
Definition: type_vector.h:145
bool operator>(const TypeVector< T > &rhs) const
Definition: type_vector.C:138
CompareTypes< T, T2 >::supertype contract(const TypeVector< T2 > &) const
Definition: type_vector.h:849
void add_scaled(const TypeVector< T2 > &, const T)
Definition: type_vector.h:604
T cross_norm(const TypeVector< T > &b, const TypeVector< T > &c)
Definition: type_vector.h:1068
void write_unformatted(std::ostream &out, const bool newline=true) const
Definition: type_vector.C:92
Real norm() const
Definition: type_vector.h:887
const class libmesh_nullptr_t libmesh_nullptr
static const Real TOLERANCE
TypeVector< T > operator-() const
Definition: type_vector.h:688
TypeVector< typename CompareTypes< T, T2 >::supertype > operator+(const TypeVector< T2 > &) const
Definition: type_vector.h:543
const TypeVector< T > & operator+=(const TypeVector< T2 > &)
Definition: type_vector.h:568
TypeVector< typename CompareTypes< T, T2 >::supertype > cross(const TypeVector< T2 > &) const
Definition: type_vector.h:859
void add(const TypeVector< T2 > &)
Definition: type_vector.h:580
const TypeVector< T > & operator/=(const T)
Definition: type_vector.h:810
T cross_norm_sq(const TypeVector< T > &b, const TypeVector< T > &c)
Definition: type_vector.h:1048
T triple_product(const TypeVector< T > &a, const TypeVector< T > &b, const TypeVector< T > &c)
Definition: type_vector.h:1027
PetscErrorCode Vec x
void subtract(const TypeVector< T2 > &)
Definition: type_vector.h:667
void print(std::ostream &os=libMesh::out) const
Definition: type_vector.C:64
bool operator!=(const TypeVector< T > &rhs) const
Definition: type_vector.h:1010
PetscErrorCode Vec Mat libmesh_dbg_var(j)
Real norm_sq() const
Definition: type_vector.h:916
void subtract_scaled(const TypeVector< T2 > &, const T)
Definition: type_vector.h:678
TypeVector< T > unit() const
Definition: type_vector.C:37
const TypeVector< T > & operator-=(const TypeVector< T2 > &)
Definition: type_vector.h:655
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool absolute_fuzzy_equals(const TypeVector< T > &rhs, Real tol=TOLERANCE) const
Definition: type_vector.h:938
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeVector< typename CompareTypes< T, Scalar >::supertype > >::type operator/(const Scalar) const
Definition: type_vector.h:781
void assign(const TypeVector< T2 > &)
Definition: type_vector.h:509
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeVector< typename CompareTypes< T, Scalar >::supertype > >::type operator*(const Scalar) const
Definition: type_vector.h:716
const T & operator()(const unsigned int i) const
Definition: type_vector.h:519
T & slice(const unsigned int i)
Definition: type_vector.h:152
Real size_sq() const
Definition: type_vector.h:906
Real size() const
Definition: type_vector.h:877
const TypeVector< T > & operator*=(const T)
Definition: type_vector.h:753
OStreamProxy out(std::cout)
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeVector & >::type operator=(const Scalar &libmesh_dbg_var(p))
Definition: type_vector.h:137
bool relative_fuzzy_equals(const TypeVector< T > &rhs, Real tol=TOLERANCE) const
Definition: type_vector.h:963