type_vector.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 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 Scalar>
83  TypeVector (const Scalar x,
84  const Scalar y=0,
85  typename
87  const Scalar>::type z=0);
88 
89 public:
90 
94  template <typename T2>
95  TypeVector (const TypeVector<T2>& p);
96 
100  ~TypeVector ();
101 
105  template <typename T2>
106  void assign (const TypeVector<T2> &);
107 
111  template <typename Scalar>
112  typename boostcopy::enable_if_c<
114  TypeVector&>::type
115  operator = (const Scalar& p)
116  { libmesh_assert_equal_to (p, Scalar(0)); this->zero(); return *this; }
117 
121  const T & operator () (const unsigned int i) const;
122 
123  const T & slice (const unsigned int i) const { return (*this)(i); }
124 
128  T & operator () (const unsigned int i);
129 
130  T & slice (const unsigned int i) { return (*this)(i); }
131 
135  template <typename T2>
137  operator + (const TypeVector<T2> &) const;
138 
142  template <typename T2>
143  const TypeVector<T> & operator += (const TypeVector<T2> &);
144 
148  template <typename T2>
149  void add (const TypeVector<T2> &);
150 
155  template <typename T2>
156  void add_scaled (const TypeVector<T2> &, const T);
157 
161  template <typename T2>
163  operator - (const TypeVector<T2> &) const;
164 
168  template <typename T2>
169  const TypeVector<T> & operator -= (const TypeVector<T2> &);
170 
174  template <typename T2>
175  void subtract (const TypeVector<T2> &);
176 
181  template <typename T2>
182  void subtract_scaled (const TypeVector<T2> &, const T);
183 
187  TypeVector<T> operator - () const;
188 
192  template <typename Scalar>
193  typename boostcopy::enable_if_c<
196  operator * (const Scalar) const;
197 
201  const TypeVector<T> & operator *= (const T);
202 
206  template <typename Scalar>
207  typename boostcopy::enable_if_c<
210  operator / (const Scalar) const;
211 
215  const TypeVector<T> & operator /= (const T);
216 
221  template <typename T2>
223  operator * (const TypeVector<T2> &) const;
224 
229  template <typename T2>
231  contract (const TypeVector<T2> &) const;
232 
236  template <typename T2>
238  cross(const TypeVector<T2> &) const;
239 
244  TypeVector<T> unit() const;
245 
250  Real size() const;
251 
256  Real size_sq() const;
257 
261  void zero();
262 
267  bool relative_fuzzy_equals(const TypeVector<T>& rhs, Real tol = TOLERANCE) const;
268 
273  bool absolute_fuzzy_equals(const TypeVector<T>& rhs, Real tol = TOLERANCE) const;
274 
279  bool operator == (const TypeVector<T>& rhs) const;
280 
285  bool operator != (const TypeVector<T>& rhs) const;
286 
293  bool operator < (const TypeVector<T>& rhs) const;
294 
301  bool operator <= (const TypeVector<T>& rhs) const;
302 
309  bool operator > (const TypeVector<T>& rhs) const;
310 
317  bool operator >= (const TypeVector<T>& rhs) const;
318 
322  void print(std::ostream& os = libMesh::out) const;
323 
329  friend std::ostream& operator << (std::ostream& os, const TypeVector<T>& t)
330  {
331  t.print(os);
332  return os;
333  }
334 
340  void write_unformatted (std::ostream &out, const bool newline = true) const;
341 
342  protected:
343 
347  T _coords[LIBMESH_DIM];
348 };
349 
350 
351 
352 
353 
354 //------------------------------------------------------
355 // Inline functions
356 
357 template <typename T>
358 inline
360 {
361  _coords[0] = 0;
362 
363 #if LIBMESH_DIM > 1
364  _coords[1] = 0;
365 #endif
366 
367 #if LIBMESH_DIM > 2
368  _coords[2] = 0;
369 #endif
370 }
371 
372 
373 
374 template <typename T>
375 inline
377  const T y,
378  const T z)
379 {
380  _coords[0] = x;
381 
382 #if LIBMESH_DIM > 1
383  _coords[1] = y;
384 #else
385  libmesh_assert_equal_to (y, 0);
386 #endif
387 
388 #if LIBMESH_DIM > 2
389  _coords[2] = z;
390 #else
391  libmesh_assert_equal_to (z, 0);
392 #endif
393 }
394 
395 
396 template <typename T>
397 template <typename Scalar>
398 inline
400  const Scalar y,
401  typename
403  const Scalar>::type z)
404 {
405  _coords[0] = x;
406 
407 #if LIBMESH_DIM > 1
408  _coords[1] = y;
409 #else
410  libmesh_assert_equal_to (y, 0);
411 #endif
412 
413 #if LIBMESH_DIM > 2
414  _coords[2] = z;
415 #else
416  libmesh_assert_equal_to (z, 0);
417 #endif
418 }
419 
420 
421 
422 template <typename T>
423 template <typename T2>
424 inline
426 {
427  // copy the nodes from vector p to me
428  for (unsigned int i=0; i<LIBMESH_DIM; i++)
429  _coords[i] = p._coords[i];
430 }
431 
432 
433 
434 template <typename T>
435 inline
437 {
438 }
439 
440 
441 
442 template <typename T>
443 template <typename T2>
444 inline
446 {
447  for (unsigned int i=0; i<LIBMESH_DIM; i++)
448  _coords[i] = p._coords[i];
449 }
450 
451 
452 
453 template <typename T>
454 inline
455 const T & TypeVector<T>::operator () (const unsigned int i) const
456 {
457  libmesh_assert_less (i, LIBMESH_DIM);
458 
459  return _coords[i];
460 }
461 
462 
463 
464 template <typename T>
465 inline
466 T & TypeVector<T>::operator () (const unsigned int i)
467 {
468  libmesh_assert_less (i, LIBMESH_DIM);
469 
470  return _coords[i];
471 }
472 
473 
474 
475 template <typename T>
476 template <typename T2>
477 inline
480 {
481  typedef typename CompareTypes<T, T2>::supertype TS;
482 #if LIBMESH_DIM == 1
483  return TypeVector<TS> (_coords[0] + p._coords[0]);
484 #endif
485 
486 #if LIBMESH_DIM == 2
487  return TypeVector<TS> (_coords[0] + p._coords[0],
488  _coords[1] + p._coords[1]);
489 #endif
490 
491 #if LIBMESH_DIM == 3
492  return TypeVector<TS> (_coords[0] + p._coords[0],
493  _coords[1] + p._coords[1],
494  _coords[2] + p._coords[2]);
495 #endif
496 
497 }
498 
499 
500 
501 template <typename T>
502 template <typename T2>
503 inline
505 {
506  this->add (p);
507 
508  return *this;
509 }
510 
511 
512 
513 template <typename T>
514 template <typename T2>
515 inline
517 {
518 #if LIBMESH_DIM == 1
519  _coords[0] += p._coords[0];
520 #endif
521 
522 #if LIBMESH_DIM == 2
523  _coords[0] += p._coords[0];
524  _coords[1] += p._coords[1];
525 #endif
526 
527 #if LIBMESH_DIM == 3
528  _coords[0] += p._coords[0];
529  _coords[1] += p._coords[1];
530  _coords[2] += p._coords[2];
531 #endif
532 
533 }
534 
535 
536 
537 template <typename T>
538 template <typename T2>
539 inline
540 void TypeVector<T>::add_scaled (const TypeVector<T2> &p, const T factor)
541 {
542 #if LIBMESH_DIM == 1
543  _coords[0] += factor*p(0);
544 #endif
545 
546 #if LIBMESH_DIM == 2
547  _coords[0] += factor*p(0);
548  _coords[1] += factor*p(1);
549 #endif
550 
551 #if LIBMESH_DIM == 3
552  _coords[0] += factor*p(0);
553  _coords[1] += factor*p(1);
554  _coords[2] += factor*p(2);
555 #endif
556 
557 }
558 
559 
560 
561 template <typename T>
562 template <typename T2>
563 inline
566 {
567  typedef typename CompareTypes<T, T2>::supertype TS;
568 
569 #if LIBMESH_DIM == 1
570  return TypeVector<TS>(_coords[0] - p._coords[0]);
571 #endif
572 
573 #if LIBMESH_DIM == 2
574  return TypeVector<TS>(_coords[0] - p._coords[0],
575  _coords[1] - p._coords[1]);
576 #endif
577 
578 #if LIBMESH_DIM == 3
579  return TypeVector<TS>(_coords[0] - p._coords[0],
580  _coords[1] - p._coords[1],
581  _coords[2] - p._coords[2]);
582 #endif
583 
584 }
585 
586 
587 
588 template <typename T>
589 template <typename T2>
590 inline
592 {
593  this->subtract (p);
594 
595  return *this;
596 }
597 
598 
599 
600 template <typename T>
601 template <typename T2>
602 inline
604 {
605  for (unsigned int i=0; i<LIBMESH_DIM; i++)
606  _coords[i] -= p._coords[i];
607 }
608 
609 
610 
611 template <typename T>
612 template <typename T2>
613 inline
614 void TypeVector<T>::subtract_scaled (const TypeVector<T2> &p, const T factor)
615 {
616  for (unsigned int i=0; i<LIBMESH_DIM; i++)
617  _coords[i] -= factor*p(i);
618 }
619 
620 
621 
622 template <typename T>
623 inline
625 {
626 
627 #if LIBMESH_DIM == 1
628  return TypeVector(-_coords[0]);
629 #endif
630 
631 #if LIBMESH_DIM == 2
632  return TypeVector(-_coords[0],
633  -_coords[1]);
634 #endif
635 
636 #if LIBMESH_DIM == 3
637  return TypeVector(-_coords[0],
638  -_coords[1],
639  -_coords[2]);
640 #endif
641 
642 }
643 
644 
645 
646 template <typename T>
647 template <typename Scalar>
648 inline
649 typename boostcopy::enable_if_c<
652 TypeVector<T>::operator * (const Scalar factor) const
653 {
654  typedef typename CompareTypes<T, Scalar>::supertype SuperType;
655 
656 #if LIBMESH_DIM == 1
657  return TypeVector<SuperType>(_coords[0]*factor);
658 #endif
659 
660 #if LIBMESH_DIM == 2
661  return TypeVector<SuperType>(_coords[0]*factor,
662  _coords[1]*factor);
663 #endif
664 
665 #if LIBMESH_DIM == 3
666  return TypeVector<SuperType>(_coords[0]*factor,
667  _coords[1]*factor,
668  _coords[2]*factor);
669 #endif
670 }
671 
672 
673 
674 template <typename T, typename Scalar>
675 inline
676 typename boostcopy::enable_if_c<
679 operator * (const Scalar factor,
680  const TypeVector<T> &v)
681 {
682  return v * factor;
683 }
684 
685 
686 
687 template <typename T>
688 inline
690 {
691 #if LIBMESH_DIM == 1
692  _coords[0] *= factor;
693 #endif
694 
695 #if LIBMESH_DIM == 2
696  _coords[0] *= factor;
697  _coords[1] *= factor;
698 #endif
699 
700 #if LIBMESH_DIM == 3
701  _coords[0] *= factor;
702  _coords[1] *= factor;
703  _coords[2] *= factor;
704 #endif
705 
706  return *this;
707 }
708 
709 
710 
711 template <typename T>
712 template <typename Scalar>
713 inline
714 typename boostcopy::enable_if_c<
717 TypeVector<T>::operator / (const Scalar factor) const
718 {
719  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
720 
721  typedef typename CompareTypes<T, Scalar>::supertype TS;
722 
723 #if LIBMESH_DIM == 1
724  return TypeVector<TS>(_coords[0]/factor);
725 #endif
726 
727 #if LIBMESH_DIM == 2
728  return TypeVector<TS>(_coords[0]/factor,
729  _coords[1]/factor);
730 #endif
731 
732 #if LIBMESH_DIM == 3
733  return TypeVector<TS>(_coords[0]/factor,
734  _coords[1]/factor,
735  _coords[2]/factor);
736 #endif
737 
738 }
739 
740 
741 
742 
743 template <typename T>
744 inline
745 const TypeVector<T> &
747 {
748  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
749 
750  for (unsigned int i=0; i<LIBMESH_DIM; i++)
751  _coords[i] /= factor;
752 
753  return *this;
754 }
755 
756 
757 
758 
759 template <typename T>
760 template <typename T2>
761 inline
764 {
765 #if LIBMESH_DIM == 1
766  return _coords[0]*p._coords[0];
767 #endif
768 
769 #if LIBMESH_DIM == 2
770  return (_coords[0]*p._coords[0] +
771  _coords[1]*p._coords[1]);
772 #endif
773 
774 #if LIBMESH_DIM == 3
775  return (_coords[0]*p(0) +
776  _coords[1]*p(1) +
777  _coords[2]*p(2));
778 #endif
779 }
780 
781 template <typename T>
782 template <typename T2>
783 inline
786 {
787  return (*this)*(p);
788 }
789 
790 
791 
792 template <typename T>
793 template <typename T2>
796 {
797  typedef typename CompareTypes<T, T2>::supertype TS;
798  libmesh_assert_equal_to (LIBMESH_DIM, 3);
799 
800  // | i j k |
801  // |(*this)(0) (*this)(1) (*this)(2)|
802  // | p(0) p(1) p(2) |
803 
804  return TypeVector<TS>( _coords[1]*p._coords[2] - _coords[2]*p._coords[1],
805  -_coords[0]*p._coords[2] + _coords[2]*p._coords[0],
806  _coords[0]*p._coords[1] - _coords[1]*p._coords[0]);
807 }
808 
809 
810 
811 template <typename T>
812 inline
814 {
815  return std::sqrt(this->size_sq());
816 }
817 
818 
819 
820 template <typename T>
821 inline
823 {
824  for (unsigned int i=0; i<LIBMESH_DIM; i++)
825  _coords[i] = 0.;
826 }
827 
828 
829 
830 template <typename T>
831 inline
833 {
834 #if LIBMESH_DIM == 1
835  return (TensorTools::norm_sq(_coords[0]));
836 #endif
837 
838 #if LIBMESH_DIM == 2
839  return (TensorTools::norm_sq(_coords[0]) +
840  TensorTools::norm_sq(_coords[1]));
841 #endif
842 
843 #if LIBMESH_DIM == 3
844  return (TensorTools::norm_sq(_coords[0]) +
845  TensorTools::norm_sq(_coords[1]) +
846  TensorTools::norm_sq(_coords[2]));
847 #endif
848 }
849 
850 
851 
852 template <typename T>
853 inline
855 {
856 #if LIBMESH_DIM == 1
857  return (std::abs(_coords[0] - rhs._coords[0])
858  <= tol);
859 #endif
860 
861 #if LIBMESH_DIM == 2
862  return (std::abs(_coords[0] - rhs._coords[0]) +
863  std::abs(_coords[1] - rhs._coords[1])
864  <= tol);
865 #endif
866 
867 #if LIBMESH_DIM == 3
868  return (std::abs(_coords[0] - rhs._coords[0]) +
869  std::abs(_coords[1] - rhs._coords[1]) +
870  std::abs(_coords[2] - rhs._coords[2])
871  <= tol);
872 #endif
873 }
874 
875 
876 
877 template <typename T>
878 inline
880 {
881 #if LIBMESH_DIM == 1
882  return this->absolute_fuzzy_equals(rhs, tol *
883  (std::abs(_coords[0]) + std::abs(rhs._coords[0])));
884 #endif
885 
886 #if LIBMESH_DIM == 2
887  return this->absolute_fuzzy_equals(rhs, tol *
888  (std::abs(_coords[0]) + std::abs(rhs._coords[0]) +
889  std::abs(_coords[1]) + std::abs(rhs._coords[1])));
890 #endif
891 
892 #if LIBMESH_DIM == 3
893  return this->absolute_fuzzy_equals(rhs, tol *
894  (std::abs(_coords[0]) + std::abs(rhs._coords[0]) +
895  std::abs(_coords[1]) + std::abs(rhs._coords[1]) +
896  std::abs(_coords[2]) + std::abs(rhs._coords[2])));
897 #endif
898 }
899 
900 
901 
902 template <typename T>
903 inline
905 {
906 #if LIBMESH_DIM == 1
907  return (_coords[0] == rhs._coords[0]);
908 #endif
909 
910 #if LIBMESH_DIM == 2
911  return (_coords[0] == rhs._coords[0] &&
912  _coords[1] == rhs._coords[1]);
913 #endif
914 
915 #if LIBMESH_DIM == 3
916  return (_coords[0] == rhs._coords[0] &&
917  _coords[1] == rhs._coords[1] &&
918  _coords[2] == rhs._coords[2]);
919 #endif
920 }
921 
922 
923 
924 template <>
925 inline
927 {
928  return (!(*this == rhs));
929 }
930 
931 } // namespace libMesh
932 
933 #endif // LIBMESH_TYPE_VECTOR_H

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:07 UTC

Hosted By:
SourceForge.net Logo