]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Extended tensor operations contract3. Fixed a bug I put in one contract3 (embarressed)
authorToby D. Young <tyoung@ippt.pan.pl>
Thu, 13 Jan 2011 12:44:41 +0000 (12:44 +0000)
committerToby D. Young <tyoung@ippt.pan.pl>
Thu, 13 Jan 2011 12:44:41 +0000 (12:44 +0000)
git-svn-id: https://svn.dealii.org/trunk@23182 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/base/tensor.h

index 4af007f0692e9ba43943319b7729d6db203c4c6e..e7f9a3d09bccea851b7a772ed106d5ca4cf23af9 100644 (file)
@@ -333,7 +333,6 @@ Tensor<rank_,dim>::operator = (const Tensor<rank_,dim> &t)
 }
 
 
-
 template <int rank_, int dim>
 inline
 Tensor<rank_,dim> &
@@ -347,7 +346,6 @@ Tensor<rank_,dim>::operator = (const double d)
 }
 
 
-
 template <int rank_, int dim>
 inline
 bool
@@ -454,7 +452,6 @@ Tensor<rank_,dim>::operator - () const
 }
 
 
-
 template <int rank_, int dim>
 inline
 double Tensor<rank_,dim>::norm () const
@@ -463,7 +460,6 @@ double Tensor<rank_,dim>::norm () const
 }
 
 
-
 template <int rank_, int dim>
 inline
 double Tensor<rank_,dim>::norm_square () const
@@ -476,7 +472,6 @@ double Tensor<rank_,dim>::norm_square () const
 }
 
 
-
 template <int rank_, int dim>
 inline
 void Tensor<rank_,dim>::clear ()
@@ -486,7 +481,6 @@ void Tensor<rank_,dim>::clear ()
 }
 
 
-
 template <int rank_, int dim>
 inline
 unsigned int
@@ -495,6 +489,7 @@ Tensor<rank_,dim>::memory_consumption ()
   return sizeof(Tensor<rank_,dim>);
 }
 
+
 template <int rank_, int dim>
 template <class Archive>
 inline
@@ -564,7 +559,6 @@ double contract (const Tensor<1,dim> &src1,
 }
 
 
-
 /**
  * Multiplication operator performing a contraction of the last index
  * of the first argument and the first index of the second
@@ -591,7 +585,6 @@ operator * (const Tensor<1,dim> &src1,
 }
 
 
-
 /**
  * Contract a tensor of rank 2 with a tensor of rank 1. The result is
  * <tt>dest[i] = sum_j src1[i][j] src2[j]</tt>.
@@ -612,7 +605,6 @@ void contract (Tensor<1,dim>       &dest,
 }
 
 
-
 /**
  * Multiplication operator performing a contraction of the last index
  * of the first argument and the first index of the second
@@ -642,7 +634,6 @@ operator * (const Tensor<2,dim> &src1,
 }
 
 
-
 /**
  * Contract a tensor of rank 1 with a tensor of rank 2. The result is
  * <tt>dest[i] = sum_j src1[j] src2[j][i]</tt>.
@@ -663,7 +654,6 @@ void contract (Tensor<1,dim>       &dest,
 }
 
 
-
 /**
  * Multiplication operator performing a contraction of the last index
  * of the first argument and the first index of the second
@@ -694,7 +684,6 @@ operator * (const Tensor<1,dim> &src1,
 }
 
 
-
 /**
  * Contract a tensor of rank 2 with a tensor of rank 2. The result is
  * <tt>dest[i][k] = sum_j src1[i][j] src2[j][k]</tt>.
@@ -716,7 +705,6 @@ void contract (Tensor<2,dim>       &dest,
 }
 
 
-
 /**
  * Multiplication operator performing a contraction of the last index
  * of the first argument and the first index of the second
@@ -748,7 +736,6 @@ operator * (const Tensor<2,dim> &src1,
 }
 
 
-
 /**
  * Contract a tensor of rank 2 with a tensor of rank 2. The
  * contraction is performed over index <tt>index1</tt> of the first tensor,
@@ -822,7 +809,6 @@ void contract (Tensor<2,dim>       &dest,
 }
 
 
-
 /**
  * Contract a tensor of rank 3 with a tensor of rank 1. The
  * contraction is performed over index <tt>index1</tt> of the first
@@ -872,7 +858,6 @@ void contract (Tensor<2,dim>       &dest,
 }
 
 
-
 /**
  * Contract a tensor of rank 3 with a tensor of rank 2. The result is
  * <tt>dest[i][j][l] = sum_k src1[i][j][k] src2[k][l]</tt>.
@@ -895,7 +880,6 @@ void contract (Tensor<3,dim>       &dest,
 }
 
 
-
 /**
  * Contract a tensor of rank 3 with a tensor of rank 2. The
  * contraction is performed over index <tt>index1</tt> of the first tensor,
@@ -994,7 +978,6 @@ void contract (Tensor<3,dim>       &dest,
     }
 }
 
-
    
 /**
  * Multiplication operator performing a contraction of the last index
@@ -1028,7 +1011,6 @@ operator * (const Tensor<3,dim> &src1,
 }
 
 
-
 /**
  * Contract a tensor of rank 2 with a tensor of rank 3. The result is
  * <tt>dest[i][j][l] = sum_k src1[i][k] src2[k][j][l]</tt>.
@@ -1051,7 +1033,6 @@ void contract (Tensor<3,dim>       &dest,
 }
 
 
-
 /**
  * Multiplication operator performing a contraction of the last index
  * of the first argument and the first index of the second
@@ -1131,7 +1112,6 @@ void double_contract (Tensor<2,dim>       &dest,
 }
 
 
-
 /**
  * Contract three tensors, corresponding to the matrix vector product
  * <i>u<sup>T</sup> A v</i>.
@@ -1155,39 +1135,63 @@ double contract3 (const Tensor<1,dim>& u,
 
 
 /**
- * Compute the contraction of three tensors $s=\sum_{i,j,k,l}
- * a_{ij}b_{ijkl}c_{kl}$.
+ * Compute the contraction of three tensors $s=\sum_{i,j,k}
+ * a_{i}b_{ijk}c_{jk}$.
  *
  * @relates Tensor
+ * @author Toby D. Young, 2011
  */
 template <int dim>
 inline
 double
-contract3 (const Tensor<2,dim> &t1,
-          const Tensor<4,dim> &t2,
+contract3 (const Tensor<1,dim> &t1,
+          const Tensor<3,dim> &t2,
           const Tensor<2,dim> &t3)
 {
   double s = 0;
   for (unsigned int i=0; i<dim; ++i)
     for (unsigned int j=0; j<dim; ++j)
       for (unsigned int k=0; k<dim; ++k)
-       for (unsigned int l=0; l<dim; ++l)
-         s += t1[i][j] * t2[i][j][k][l] * t3[k][l];
+         s += t1[i] * t2[i][j][k] * t3[j][k];
   return s;
 }
 
 
 /**
- * Compute the contraction of three tensors $s=\sum_{i,j,k,l}
- * a_{i}b_{ijk}c_{kl}$.
+ * Compute the contraction of three tensors $s=\sum_{i,j,k}
+ * a_{ij}b_{ijk}c_{k}$.
  *
  * @relates Tensor
+ * @author Toby D. Young, 2011
  */
 template <int dim>
 inline
 double
-contract3 (const Tensor<1,dim> &t1,
+contract3 (const Tensor<2,dim> &t1,
           const Tensor<3,dim> &t2,
+          const Tensor<1,dim> &t3)
+{
+  double s = 0;
+  for (unsigned int i=0; i<dim; ++i)
+    for (unsigned int j=0; j<dim; ++j)
+      for (unsigned int k=0; k<dim; ++k)
+         s += t1[i][j] * t2[i][j][k] * t3[k];
+  return s;
+}
+
+
+/**
+ * Compute the contraction of three tensors $s=\sum_{i,j,k,l}
+ * a_{ij}b_{ijkl}c_{kl}$.
+ *
+ * @relates Tensor
+ * @author Toby D. Young, 2011
+ */
+template <int dim>
+inline
+double
+contract3 (const Tensor<2,dim> &t1,
+          const Tensor<4,dim> &t2,
           const Tensor<2,dim> &t3)
 {
   double s = 0;
@@ -1195,7 +1199,7 @@ contract3 (const Tensor<1,dim> &t1,
     for (unsigned int j=0; j<dim; ++j)
       for (unsigned int k=0; k<dim; ++k)
        for (unsigned int l=0; l<dim; ++l)
-         s += t1[i] * t2[i][j][k] * t3[k][l];
+         s += t1[i][j] * t2[i][j][k][l] * t3[k][l];
   return s;
 }
 
@@ -1218,7 +1222,6 @@ void outer_product (Tensor<2,dim>       &dst,
 }
 
 
-
 /**
  * Form the outer product of two tensors of rank 1 and 2, i.e.
  * <tt>dst[i][j][k] = src1[i] * src2[j][k]</tt>.
@@ -1238,7 +1241,6 @@ void outer_product (Tensor<3,dim>       &dst,
 }
 
 
-
 /**
  * Form the outer product of two tensors of rank 2 and 1, i.e.
  * <tt>dst[i][j][k] = src1[i][j] * src2[k]</tt>.
@@ -1258,7 +1260,6 @@ void outer_product (Tensor<3,dim>       &dst,
 }
 
 
-
 /**
  * Form the outer product of two tensors of rank 0 and 1, i.e.
  * <tt>dst[i] = src1 * src2[i]</tt>. Of course, this is only a scaling of
@@ -1302,7 +1303,6 @@ void outer_product (Tensor<1,dim>       &dst,
 }
 
 
-
 /**
  * Cross-product in 2d. This is just a rotation by 90 degrees
  * clockwise to compute the outer normal from a tangential
@@ -1328,7 +1328,6 @@ cross_product (Tensor<1,dim>       &dst,
 }
 
 
-
 /**
  * Cross-product of 2 vectors in 3d. This function is defined for all
  * space dimensions to allow for dimension independent programming
@@ -1354,7 +1353,6 @@ cross_product (Tensor<1,dim>       &dst,
 }
 
 
-
 /**
  * Compute the scalar product $a:b=\sum_{i,j} a_{ij}b_{ij}$ between two
  * tensors $a,b$ of rank 2. We don't use <code>operator*</code> for this
@@ -1411,7 +1409,6 @@ double determinant (const Tensor<1,1> &t)
 }
 
 
-
 /**
  * Compute the determinant of a tensor of rank two and dimension
  * one. Since this is a number, the return value is, of course, the
@@ -1442,8 +1439,6 @@ double determinant (const Tensor<2,2> &t)
 }
 
 
-
-
 /**
  * Compute the determinant of a tensor or rank 2, here for <tt>dim==3</tt>.
  *
@@ -1468,7 +1463,6 @@ double determinant (const Tensor<2,3> &t)
 }
 
 
-
 /**
  * Compute the determinant of a tensor or rank 2, here for all dimensions for
  * which no explicit specialization is available above.

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.