type_tensor.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_TENSOR_H
21 #define LIBMESH_TYPE_TENSOR_H
22 
23 // Local includes
24 #include "libmesh/libmesh_common.h"
25 #include "libmesh/type_vector.h"
26 
27 // C++ includes
28 #include <cstdlib> // *must* precede <cmath> for proper std:abs() on PGI, Sun Studio CC
29 #include <cmath>
30 
31 namespace libMesh
32 {
33 
34 // Forward declarations
35 template <typename T> class TypeTensorColumn;
36 template <typename T> class ConstTypeTensorColumn;
37 template <unsigned int N, typename T> class TypeNTensor;
38 
39 
40 
50 template <typename T>
51 class TypeTensor
52 {
53 template <typename T2>
54 friend class TypeTensor;
55 
56 protected:
57 
62  TypeTensor ();
63 
70  explicit TypeTensor (const T xx,
71  const T xy=0,
72  const T xz=0,
73  const T yx=0,
74  const T yy=0,
75  const T yz=0,
76  const T zx=0,
77  const T zy=0,
78  const T zz=0);
79 
80 
84  template <typename Scalar>
85  explicit TypeTensor (const Scalar xx,
86  const Scalar xy=0,
87  const Scalar xz=0,
88  const Scalar yx=0,
89  const Scalar yy=0,
90  const Scalar yz=0,
91  const Scalar zx=0,
92  const Scalar zy=0,
93  typename
95  const Scalar>::type zz=0);
96 
102  template<typename T2>
103  TypeTensor(const TypeVector<T2>& vx);
104 
105  template<typename T2>
106  TypeTensor(const TypeVector<T2>& vx, const TypeVector<T2> &vy);
107 
108  template<typename T2>
109  TypeTensor(const TypeVector<T2>& vx, const TypeVector<T2> &vy, const TypeVector<T2> &vz);
110 
111 public:
112 
116  template<typename T2>
117  TypeTensor(const TypeTensor<T2>& p);
118 
122  ~TypeTensor();
123 
127  template<typename T2>
128  void assign (const TypeTensor<T2> &);
129 
133  template <typename Scalar>
134  typename boostcopy::enable_if_c<
136  TypeTensor&>::type
137  operator = (const Scalar& p)
138  { libmesh_assert_equal_to (p, Scalar(0)); this->zero(); return *this; }
139 
143  const T & operator () (const unsigned int i, const unsigned int j) const;
144 
149  T & operator () (const unsigned int i, const unsigned int j);
150 
154  ConstTypeTensorColumn<T> slice (const unsigned int i) const;
155 
159  TypeTensorColumn<T> slice (const unsigned int i);
160 
164  TypeVector<T> row(const unsigned int r);
165 
169  template<typename T2>
171  operator + (const TypeTensor<T2> &) const;
172 
176  template<typename T2>
177  const TypeTensor<T> & operator += (const TypeTensor<T2> &);
178 
182  template<typename T2>
183  void add (const TypeTensor<T2> &);
184 
189  template <typename T2>
190  void add_scaled (const TypeTensor<T2> &, const T);
191 
195  template<typename T2>
197  operator - (const TypeTensor<T2> &) const;
198 
202  template<typename T2>
203  const TypeTensor<T> & operator -= (const TypeTensor<T2> &);
204 
208  template<typename T2>
209  void subtract (const TypeTensor<T2> &);
210 
215  template <typename T2>
216  void subtract_scaled (const TypeTensor<T2> &, const T);
217 
221  TypeTensor<T> operator - () const;
222 
226  template <typename Scalar>
227  typename boostcopy::enable_if_c<
230  operator * (const Scalar) const;
231 
235  template <typename Scalar>
236  const TypeTensor<T> & operator *= (const Scalar);
237 
241  template <typename Scalar>
242  typename boostcopy::enable_if_c<
245  operator / (const Scalar) const;
246 
250  const TypeTensor<T> & operator /= (const T);
251 
256  template <typename T2>
258 
264  template <typename T2>
266  contract (const TypeTensor<T2> &) const;
267 
272  template <typename T2>
274  operator * (const TypeVector<T2> &) const;
275 
279  TypeTensor<T> transpose() const;
280 
285  Real size() const;
286 
291  Real size_sq() const;
292 
298  T det() const;
299 
303  T tr() const;
304 
308  void zero();
309 
313  bool operator == (const TypeTensor<T>& rhs) const;
314 
319  bool operator < (const TypeTensor<T>& rhs) const;
320 
325  bool operator > (const TypeTensor<T>& rhs) const;
326 
330  void print(std::ostream& os = libMesh::out) const;
331 
336  friend std::ostream& operator << (std::ostream& os, const TypeTensor<T>& t)
337  {
338  t.print(os);
339  return os;
340  }
341 
346  void write_unformatted (std::ostream &out, const bool newline = true) const;
347 
348  // protected:
349 
353  T _coords[LIBMESH_DIM*LIBMESH_DIM];
354 };
355 
356 
357 template <typename T>
358 class TypeTensorColumn
359 {
360 public:
361 
363  unsigned int j) :
364  _tensor(&tensor), _j(j) {}
365 
370  T & operator () (const unsigned int i)
371  { return (*_tensor)(i,_j); }
372 
373  T & slice (const unsigned int i)
374  { return (*_tensor)(i,_j); }
375 
380  {
381  for (unsigned int i=0; i != LIBMESH_DIM; ++i)
382  (*this)(i) = rhs(i);
383  return *this;
384  }
385 
386 private:
388  const unsigned int _j;
389 };
390 
391 
392 template <typename T>
394 {
395 public:
396 
398  unsigned int j) :
399  _tensor(&tensor), _j(j) {}
400 
404  const T & operator () (const unsigned int i) const
405  { return (*_tensor)(i,_j); }
406 
407  const T & slice (const unsigned int i) const
408  { return (*_tensor)(i,_j); }
409 
410 private:
412  const unsigned int _j;
413 };
414 
415 //------------------------------------------------------
416 // Inline functions
417 template <typename T>
418 inline
420 {
421  _coords[0] = 0;
422 
423 #if LIBMESH_DIM > 1
424  _coords[1] = 0;
425  _coords[2] = 0;
426  _coords[3] = 0;
427 #endif
428 
429 #if LIBMESH_DIM > 2
430  _coords[4] = 0;
431  _coords[5] = 0;
432  _coords[6] = 0;
433  _coords[7] = 0;
434  _coords[8] = 0;
435 #endif
436 }
437 
438 
439 
440 template <typename T>
441 inline
443  (const T xx,
444  const T xy,
445  const T xz,
446  const T yx,
447  const T yy,
448  const T yz,
449  const T zx,
450  const T zy,
451  const T zz)
452 {
453  _coords[0] = xx;
454 
455 #if LIBMESH_DIM == 2
456  _coords[1] = xy;
457  _coords[2] = yx;
458  _coords[3] = yy;
459 #endif
460 
461 #if LIBMESH_DIM == 3
462  _coords[1] = xy;
463  _coords[2] = xz;
464  _coords[3] = yx;
465  _coords[4] = yy;
466  _coords[5] = yz;
467  _coords[6] = zx;
468  _coords[7] = zy;
469  _coords[8] = zz;
470 #endif
471 }
472 
473 
474 template <typename T>
475 template <typename Scalar>
476 inline
478  (const Scalar xx,
479  const Scalar xy,
480  const Scalar xz,
481  const Scalar yx,
482  const Scalar yy,
483  const Scalar yz,
484  const Scalar zx,
485  const Scalar zy,
486  typename
488  const Scalar>::type zz)
489 {
490  _coords[0] = xx;
491 
492 #if LIBMESH_DIM == 2
493  _coords[1] = xy;
494  _coords[2] = yx;
495  _coords[3] = yy;
496 #endif
497 
498 #if LIBMESH_DIM == 3
499  _coords[1] = xy;
500  _coords[2] = xz;
501  _coords[3] = yx;
502  _coords[4] = yy;
503  _coords[5] = yz;
504  _coords[6] = zx;
505  _coords[7] = zy;
506  _coords[8] = zz;
507 #endif
508 }
509 
510 
511 
512 template <typename T>
513 template<typename T2>
514 inline
516 {
517  // copy the nodes from vector p to me
518  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
519  _coords[i] = p._coords[i];
520 }
521 
522 
523 template <typename T>
524 template <typename T2>
526 {
527  libmesh_assert_equal_to (LIBMESH_DIM, 1);
528  _coords[0] = vx(0);
529 }
530 
531 template <typename T>
532 template <typename T2>
534 {
535  libmesh_assert_equal_to (LIBMESH_DIM, 2);
536  _coords[0] = vx(0);
537  _coords[1] = vx(1);
538  _coords[2] = vy(0);
539  _coords[3] = vy(1);
540 }
541 
542 template <typename T>
543 template <typename T2>
545 {
546  libmesh_assert_equal_to (LIBMESH_DIM, 3);
547  _coords[0] = vx(0);
548  _coords[1] = vx(1);
549  _coords[2] = vx(2);
550  _coords[3] = vy(0);
551  _coords[4] = vy(1);
552  _coords[5] = vy(2);
553  _coords[6] = vz(0);
554  _coords[7] = vz(1);
555  _coords[8] = vz(2);
556 }
557 
558 
559 
560 
561 template <typename T>
562 inline
564 {
565 }
566 
567 
568 
569 template <typename T>
570 template<typename T2>
571 inline
573 {
574  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
575  _coords[i] = p._coords[i];
576 }
577 
578 
579 
580 template <typename T>
581 inline
582 const T & TypeTensor<T>::operator () (const unsigned int i,
583  const unsigned int j) const
584 {
585  libmesh_assert_less (i, 3);
586  libmesh_assert_less (j, 3);
587 
588 #if LIBMESH_DIM < 3
589  const static T my_zero = 0;
590  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
591  return my_zero;
592 #endif
593 
594  return _coords[i*LIBMESH_DIM+j];
595 }
596 
597 
598 
599 template <typename T>
600 inline
601 T & TypeTensor<T>::operator () (const unsigned int i,
602  const unsigned int j)
603 {
604 #if LIBMESH_DIM < 3
605 
606  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
607  {
608 // libMesh::err << "ERROR: You are assigning to a tensor component" << std::endl
609 // << "that is out of range for the compiled LIBMESH_DIM!" << std::endl
610 // << " LIBMESH_DIM=" << LIBMESH_DIM << " , i=" << i << " , j=" << j << std::endl;
611  libmesh_error();
612  }
613 
614 #endif
615 
616  libmesh_assert_less (i, LIBMESH_DIM);
617  libmesh_assert_less (j, LIBMESH_DIM);
618 
619  return _coords[i*LIBMESH_DIM+j];
620 }
621 
622 
623 template <typename T>
624 inline
626 TypeTensor<T>::slice (const unsigned int i) const
627 {
628  libmesh_assert_less (i, LIBMESH_DIM);
629  return ConstTypeTensorColumn<T>(*this, i);
630 }
631 
632 
633 template <typename T>
634 inline
636 TypeTensor<T>::slice (const unsigned int i)
637 {
638  libmesh_assert_less (i, LIBMESH_DIM);
639  return TypeTensorColumn<T>(*this, i);
640 }
641 
642 
643 template <typename T>
644 inline
646 TypeTensor<T>::row(const unsigned int r)
647 {
648  TypeVector<T> return_vector;
649 
650  for(unsigned int j=0; j<LIBMESH_DIM; j++)
651  return_vector._coords[j] = _coords[r*LIBMESH_DIM + j];
652 
653  return return_vector;
654 }
655 
656 template <typename T>
657 template<typename T2>
658 inline
661 {
662 
663 #if LIBMESH_DIM == 1
664  return TypeTensor(_coords[0] + p._coords[0]);
665 #endif
666 
667 #if LIBMESH_DIM == 2
668  return TypeTensor(_coords[0] + p._coords[0],
669  _coords[1] + p._coords[1],
670  0.,
671  _coords[2] + p._coords[2],
672  _coords[3] + p._coords[3]);
673 #endif
674 
675 #if LIBMESH_DIM == 3
676  return TypeTensor(_coords[0] + p._coords[0],
677  _coords[1] + p._coords[1],
678  _coords[2] + p._coords[2],
679  _coords[3] + p._coords[3],
680  _coords[4] + p._coords[4],
681  _coords[5] + p._coords[5],
682  _coords[6] + p._coords[6],
683  _coords[7] + p._coords[7],
684  _coords[8] + p._coords[8]);
685 #endif
686 
687 }
688 
689 
690 
691 template <typename T>
692 template<typename T2>
693 inline
695 {
696  this->add (p);
697 
698  return *this;
699 }
700 
701 
702 
703 template <typename T>
704 template<typename T2>
705 inline
707 {
708  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
709  _coords[i] += p._coords[i];
710 }
711 
712 
713 
714 template <typename T>
715 template <typename T2>
716 inline
717 void TypeTensor<T>::add_scaled (const TypeTensor<T2> &p, const T factor)
718 {
719  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
720  _coords[i] += factor*p._coords[i];
721 
722 }
723 
724 
725 
726 template <typename T>
727 template<typename T2>
728 inline
731 {
732 
733 #if LIBMESH_DIM == 1
734  return TypeTensor(_coords[0] - p._coords[0]);
735 #endif
736 
737 #if LIBMESH_DIM == 2
738  return TypeTensor(_coords[0] - p._coords[0],
739  _coords[1] - p._coords[1],
740  0.,
741  _coords[2] - p._coords[2],
742  _coords[3] - p._coords[3]);
743 #endif
744 
745 #if LIBMESH_DIM == 3
746  return TypeTensor(_coords[0] - p._coords[0],
747  _coords[1] - p._coords[1],
748  _coords[2] - p._coords[2],
749  _coords[3] - p._coords[3],
750  _coords[4] - p._coords[4],
751  _coords[5] - p._coords[5],
752  _coords[6] - p._coords[6],
753  _coords[7] - p._coords[7],
754  _coords[8] - p._coords[8]);
755 #endif
756 
757 }
758 
759 
760 
761 template <typename T>
762 template <typename T2>
763 inline
765 {
766  this->subtract (p);
767 
768  return *this;
769 }
770 
771 
772 
773 template <typename T>
774 template <typename T2>
775 inline
777 {
778  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
779  _coords[i] -= p._coords[i];
780 }
781 
782 
783 
784 template <typename T>
785 template <typename T2>
786 inline
787 void TypeTensor<T>::subtract_scaled (const TypeTensor<T2> &p, const T factor)
788 {
789  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
790  _coords[i] -= factor*p._coords[i];
791 }
792 
793 
794 
795 template <typename T>
796 inline
798 {
799 
800 #if LIBMESH_DIM == 1
801  return TypeTensor(-_coords[0]);
802 #endif
803 
804 #if LIBMESH_DIM == 2
805  return TypeTensor(-_coords[0],
806  -_coords[1],
807  -_coords[2],
808  -_coords[3]);
809 #endif
810 
811 #if LIBMESH_DIM == 3
812  return TypeTensor(-_coords[0],
813  -_coords[1],
814  -_coords[2],
815  -_coords[3],
816  -_coords[4],
817  -_coords[5],
818  -_coords[6],
819  -_coords[7],
820  -_coords[8]);
821 #endif
822 
823 }
824 
825 
826 
827 template <typename T>
828 template <typename Scalar>
829 inline
830 typename boostcopy::enable_if_c<
833 TypeTensor<T>::operator * (const Scalar factor) const
834 {
835  typedef typename CompareTypes<T, Scalar>::supertype TS;
836 
837 
838 #if LIBMESH_DIM == 1
839  return TypeTensor<TS>(_coords[0]*factor);
840 #endif
841 
842 #if LIBMESH_DIM == 2
843  return TypeTensor<TS>(_coords[0]*factor,
844  _coords[1]*factor,
845  _coords[2]*factor,
846  _coords[3]*factor);
847 #endif
848 
849 #if LIBMESH_DIM == 3
850  return TypeTensor<TS>(_coords[0]*factor,
851  _coords[1]*factor,
852  _coords[2]*factor,
853  _coords[3]*factor,
854  _coords[4]*factor,
855  _coords[5]*factor,
856  _coords[6]*factor,
857  _coords[7]*factor,
858  _coords[8]*factor);
859 #endif
860 }
861 
862 
863 
864 template <typename T, typename Scalar>
865 inline
866 typename boostcopy::enable_if_c<
869 operator * (const Scalar factor,
870  const TypeTensor<T> &t)
871 {
872  return t * factor;
873 }
874 
875 
876 
877 template <typename T>
878 template <typename Scalar>
879 inline
880 const TypeTensor<T> & TypeTensor<T>::operator *= (const Scalar factor)
881 {
882  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
883  _coords[i] *= factor;
884 
885  return *this;
886 }
887 
888 
889 
890 
891 template <typename T>
892 template <typename Scalar>
893 inline
894 typename boostcopy::enable_if_c<
897 TypeTensor<T>::operator / (const Scalar factor) const
898 {
899  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
900 
901  typedef typename CompareTypes<T, Scalar>::supertype TS;
902 
903 #if LIBMESH_DIM == 1
904  return TypeTensor<TS>(_coords[0]/factor);
905 #endif
906 
907 #if LIBMESH_DIM == 2
908  return TypeTensor<TS>(_coords[0]/factor,
909  _coords[1]/factor,
910  _coords[2]/factor,
911  _coords[3]/factor);
912 #endif
913 
914 #if LIBMESH_DIM == 3
915  return TypeTensor<TS>(_coords[0]/factor,
916  _coords[1]/factor,
917  _coords[2]/factor,
918  _coords[3]/factor,
919  _coords[4]/factor,
920  _coords[5]/factor,
921  _coords[6]/factor,
922  _coords[7]/factor,
923  _coords[8]/factor);
924 #endif
925 
926 }
927 
928 
929 
930 template <typename T>
931 inline
933 {
934 #if LIBMESH_DIM == 1
935  return TypeTensor(_coords[0]);
936 #endif
937 
938 #if LIBMESH_DIM == 2
939  return TypeTensor(_coords[0],
940  _coords[2],
941  _coords[1],
942  _coords[3]);
943 #endif
944 
945 #if LIBMESH_DIM == 3
946  return TypeTensor(_coords[0],
947  _coords[3],
948  _coords[6],
949  _coords[1],
950  _coords[4],
951  _coords[7],
952  _coords[2],
953  _coords[5],
954  _coords[8]);
955 #endif
956 }
957 
958 
959 
960 template <typename T>
961 inline
963 {
964  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
965 
966  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
967  _coords[i] /= factor;
968 
969  return *this;
970 }
971 
972 
973 
974 
975 template <typename T>
976 template <typename T2>
977 inline
980 {
982  for (unsigned int i=0; i<LIBMESH_DIM; i++)
983  for (unsigned int j=0; j<LIBMESH_DIM; j++)
984  returnval(i) += (*this)(i,j)*p(j);
985 
986  return returnval;
987 }
988 
989 
990 
991 template <typename T>
992 template <typename T2>
993 inline
995 {
996  TypeTensor<T> returnval;
997  for (unsigned int i=0; i<LIBMESH_DIM; i++)
998  for (unsigned int j=0; j<LIBMESH_DIM; j++)
999  for (unsigned int k=0; k<LIBMESH_DIM; k++)
1000  returnval(i,j) += (*this)(i,k)*p(k,j);
1001 
1002  return returnval;
1003 }
1004 
1005 
1006 
1011 template <typename T>
1012 template <typename T2>
1013 inline
1016 {
1017  typename CompareTypes<T,T2>::supertype sum = 0.;
1018  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1019  sum += _coords[i]*t._coords[i];
1020  return sum;
1021 }
1022 
1023 
1024 
1025 template <typename T>
1026 inline
1028 {
1029  return std::sqrt(this->size_sq());
1030 }
1031 
1032 
1033 
1034 template <typename T>
1035 inline
1037 {
1038 #if LIBMESH_DIM == 1
1039  return _coords[0];
1040 #endif
1041 
1042 #if LIBMESH_DIM == 2
1043  return (_coords[0] * _coords[3]
1044  - _coords[1] * _coords[2]);
1045 #endif
1046 
1047 #if LIBMESH_DIM == 3
1048  return (_coords[0] * _coords[4] * _coords[8]
1049  + _coords[1] * _coords[5] * _coords[6]
1050  + _coords[2] * _coords[3] * _coords[7]
1051  - _coords[0] * _coords[5] * _coords[7]
1052  - _coords[1] * _coords[3] * _coords[8]
1053  - _coords[2] * _coords[4] * _coords[6]);
1054 #endif
1055 }
1056 
1057 template <typename T>
1058 inline
1060 {
1061 #if LIBMESH_DIM == 1
1062  return _coords[0];
1063 #endif
1064 
1065 #if LIBMESH_DIM == 2
1066  return _coords[0] + _coords[3];
1067 #endif
1068 
1069 #if LIBMESH_DIM == 3
1070  return _coords[0] + _coords[4] + _coords[8];
1071 #endif
1072 }
1073 
1074 template <typename T>
1075 inline
1077 {
1078  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1079  _coords[i] = 0.;
1080 }
1081 
1082 
1083 
1084 template <typename T>
1085 inline
1087 {
1088  Real sum = 0.;
1089  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1090  sum += TensorTools::norm_sq(_coords[i]);
1091  return sum;
1092 }
1093 
1094 
1095 
1096 template <typename T>
1097 inline
1099 {
1100 #if LIBMESH_DIM == 1
1101  return (std::abs(_coords[0] - rhs._coords[0])
1102  < TOLERANCE);
1103 #endif
1104 
1105 #if LIBMESH_DIM == 2
1106  return ((std::abs(_coords[0] - rhs._coords[0]) +
1107  std::abs(_coords[1] - rhs._coords[1]) +
1108  std::abs(_coords[2] - rhs._coords[2]) +
1109  std::abs(_coords[3] - rhs._coords[3]))
1110  < 4.*TOLERANCE);
1111 #endif
1112 
1113 #if LIBMESH_DIM == 3
1114  return ((std::abs(_coords[0] - rhs._coords[0]) +
1115  std::abs(_coords[1] - rhs._coords[1]) +
1116  std::abs(_coords[2] - rhs._coords[2]) +
1117  std::abs(_coords[3] - rhs._coords[3]) +
1118  std::abs(_coords[4] - rhs._coords[4]) +
1119  std::abs(_coords[5] - rhs._coords[5]) +
1120  std::abs(_coords[6] - rhs._coords[6]) +
1121  std::abs(_coords[7] - rhs._coords[7]) +
1122  std::abs(_coords[8] - rhs._coords[8]))
1123  < 9.*TOLERANCE);
1124 #endif
1125 
1126 }
1127 
1128 } // namespace libMesh
1129 
1130 #endif // LIBMESH_TYPE_TENSOR_H

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

Hosted By:
SourceForge.net Logo