tensor(2rheolef)						    rheolef-6.1 						  tensor(2rheolef)

NAME

       tensor - a N*N tensor, N=1,2,3

SYNOPSYS

       The tensor class defines a 3*3 tensor, as the value of a tensorial valued field. Basic algebra with scalars, vectors of R^3 (i.e. the point
       class) and tensor objects are supported.

IMPLEMENTATION

       template<class T>
       class tensor_basic {
	   public:

	       typedef size_t size_type;
	       typedef T  element_type;

       // allocators:

	       tensor_basic (const T& init_val = 0);
	       tensor_basic (T x[3][3]);
	       tensor_basic (const tensor_basic<T>& a);

       // affectation:

	       tensor_basic<T>& operator = (const tensor_basic<T>& a);
	       tensor_basic<T>& operator = (const T& val);

       // modifiers:

	       void fill (const T& init_val);
	       void reset ();
	       void set_row    (const point_basic<T>& r, size_t i, size_t d = 3);
	       void set_column (const point_basic<T>& c, size_t j, size_t d = 3);

       // accessors:

	       T& operator()(size_type i, size_type j);
	       T  operator()(size_type i, size_type j) const;
	       point_basic<T>  row(size_type i) const;
	       point_basic<T>  col(size_type i) const;
	       size_t nrow() const; // = 3, for template matrix compatibility
	       size_t ncol() const;

       // inputs/outputs:

	       std::ostream& put (std::ostream& s, size_type d = 3) const;
	       std::istream& get (std::istream&);

       // algebra:

	       bool operator== (const tensor_basic<T>&) const;
	       bool operator!= (const tensor_basic<T>& b) const { return ! operator== (b); }
	       template <class U>
	       friend tensor_basic<U> operator- (const tensor_basic<U>&);
	       template <class U>
	       friend tensor_basic<U> operator+ (const tensor_basic<U>&, const tensor_basic<U>&);
	       template <class U>
	       friend tensor_basic<U> operator- (const tensor_basic<U>&, const tensor_basic<U>&);
	       template <class U>
	       friend tensor_basic<U> operator* (int k, const tensor_basic<U>& a);
	       template <class U>
	       friend tensor_basic<U> operator* (const U& k, const tensor_basic<U>& a);
	       template <class U>
	       friend tensor_basic<U> operator* (const tensor_basic<U>& a, int k);
	       template <class U>
	       friend tensor_basic<U> operator* (const tensor_basic<U>& a, const U& k);
	       template <class U>
	       friend tensor_basic<U> operator/ (const tensor_basic<U>& a, int k);
	       template <class U>
	       friend tensor_basic<U> operator/ (const tensor_basic<U>& a, const U& k);
	       template <class U>
	       friend point_basic<U>  operator* (const tensor_basic<U>&, const point_basic<U>&);
	       template <class U>
	       friend point_basic<U>  operator* (const point_basic<U>& yt, const tensor_basic<U>& a);
		      point_basic<T>  trans_mult (const point_basic<T>& x) const;
	       template <class U>
	       friend tensor_basic<U> trans	 (const tensor_basic<U>& a, size_t d = 3);
	       template <class U>
	       friend tensor_basic<U> operator* (const tensor_basic<U>& a, const tensor_basic<U>& b);
	       template <class U>
	       friend void prod (const tensor_basic<U>& a, const tensor_basic<U>& b, tensor_basic<U>& result,
		       size_t di=3, size_t dj=3, size_t dk=3);
	       template <class U>
	       friend tensor_basic<U> inv	 (const tensor_basic<U>& a, size_t d = 3);
	       template <class U>
	       friend tensor_basic<U> diag (const point_basic<U>& d);
	       template <class U>
	       friend tensor_basic<U> identity (size_t d=3);
	       template <class U>
	       friend tensor_basic<U> dyadic (const point_basic<U>& u, const point_basic<U>& v, size_t d=3);

       // metric and geometric transformations:

	       template <class U>
	       friend U dotdot (const tensor_basic<U>&, const tensor_basic<U>&);
	       template <class U>
	       friend U norm2 (const tensor_basic<U>& a) { return dotdot(a,a); }
	       template <class U>
	       friend U dist2 (const tensor_basic<U>& a, const tensor_basic<U>& b) { return norm2(a-b); }
	       template <class U>
	       friend U norm  (const tensor_basic<U>& a) { return ::sqrt(norm2(a)); }
	       template <class U>
	       friend U dist  (const tensor_basic<U>& a, const tensor_basic<U>& b) { return norm(a-b); }
	       T determinant (size_type d = 3) const;
	       template <class U>
	       friend U determinant (const tensor_basic<U>& A, size_t d = 3);
	       template <class U>
	       friend bool invert_3x3 (const tensor_basic<U>& A, tensor_basic<U>& result);

       // spectral:
	       // eigenvalues & eigenvectors:
	       // a = q*d*q^T
	       // a may be symmetric
	       // where q=(q1,q2,q3) are eigenvectors in rows (othonormal matrix)
	       // and	d=(d1,d2,d3) are eigenvalues, sorted in decreasing order d1 >= d2 >= d3
	       // return d
	       point_basic<T> eig (tensor_basic<T>& q, size_t dim = 3) const;
	       point_basic<T> eig (size_t dim = 3) const;

	       // singular value decomposition:
	       // a = u*s*v^T
	       // a can be unsymmetric
	       // where u=(u1,u2,u3) are left pseudo-eigenvectors in rows (othonormal matrix)
	       //	v=(v1,v2,v3) are right pseudo-eigenvectors in rows (othonormal matrix)
	       // and	s=(s1,s2,s3) are eigenvalues, sorted in decreasing order s1 >= s2 >= s3
	       // return s
	       point_basic<T> svd (tensor_basic<T>& u, tensor_basic<T>& v, size_t dim = 3) const;

       // data:
	       T _x[3][3];
       };
       typedef tensor_basic<Float> tensor;

       // inputs/outputs:
       template<class T>
       inline
       std::istream& operator>> (std::istream& in, tensor_basic<T>& a)
       {
	   return a.get (in);
       }
       template<class T>
       inline
       std::ostream& operator<< (std::ostream& out, const tensor_basic<T>& a)
       {
	   return a.put (out);
       }
       // t = a otimes b
       template<class T>
       tensor_basic<T> otimes (const point_basic<T>& a, const point_basic<T>& b, size_t na = 3);
       // t += a otimes b
       template<class T>
       void cumul_otimes (tensor_basic<T>& t, const point_basic<T>& a, const point_basic<T>& b, size_t na = 3);
       template<class T>
       void cumul_otimes (tensor_basic<T>& t, const point_basic<T>& a, const point_basic<T>& b, size_t na, size_t nb);

rheolef-6.1							    rheolef-6.1 						  tensor(2rheolef)