]> https://gitweb.dealii.org/ - dealii.git/commitdiff
tensor.h: Deprecate "contract" function
authorMatthias Maier <tamiko@43-1.org>
Tue, 15 Sep 2015 06:18:28 +0000 (01:18 -0500)
committerMatthias Maier <tamiko@43-1.org>
Tue, 15 Sep 2015 22:19:06 +0000 (17:19 -0500)
include/deal.II/base/tensor.h
include/deal.II/base/tensor_deprecated.h

index 060f47019a406e6c5469eac315b5879fe0979308..cc53f565c8d6a63597db54ab6b75bc24bcbca5ef 100644 (file)
@@ -1537,27 +1537,6 @@ outer_product(const Tensor<rank_1, dim, Number> &src1,
 //@{
 
 
-/**
- * Contract a tensor of rank 1 with a tensor of rank 1. The result is
- * <tt>sum_j src1[j] src2[j]</tt>.
- *
- * @relates Tensor
- */
-template <int dim, typename Number, typename OtherNumber>
-inline
-typename ProductType<Number,OtherNumber>::type
-contract (const Tensor<1,dim,Number> &src1,
-          const Tensor<1,dim,OtherNumber> &src2)
-{
-  typename ProductType<Number,OtherNumber>::type res
-    = typename ProductType<Number,OtherNumber>::type();
-  for (unsigned int i=0; i<dim; ++i)
-    res += src1[i] * src2[i];
-
-  return res;
-}
-
-
 /**
  * Double contract two tensors of rank 2, thus computing the Frobenius inner
  * product <tt> sum<sub>i,j</sub> src1[i][j]*src2[i][j]</tt>.
@@ -1578,73 +1557,6 @@ Number double_contract (const Tensor<2, dim, Number> &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>.
- *
- * @relates Tensor
- * @author Wolfgang Bangerth, 1998
- */
-template <int dim, typename Number>
-inline
-void contract (Tensor<1,dim,Number>       &dest,
-               const Tensor<2,dim,Number> &src1,
-               const Tensor<1,dim,Number> &src2)
-{
-  for (unsigned int i=0; i<dim; ++i)
-    {
-      dest[i] = src1[i][0] * src2[0];
-      for (unsigned int j=1; j<dim; ++j)
-        dest[i] += src1[i][j] * src2[j];
-    }
-}
-
-
-/**
- * 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>.
- *
- * @relates Tensor
- * @author Guido Kanschat, 2001
- */
-template <int dim, typename Number>
-inline
-void contract (Tensor<1,dim,Number>       &dest,
-               const Tensor<1,dim,Number> &src1,
-               const Tensor<2,dim,Number> &src2)
-{
-  for (unsigned int i=0; i<dim; ++i)
-    {
-      dest[i] = src1[0] * src2[0][i];
-      for (unsigned int j=1; j<dim; ++j)
-        dest[i] += src1[j] * src2[j][i];
-    }
-}
-
-
-/**
- * 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>.
- *
- * @relates Tensor
- * @author Wolfgang Bangerth, 1998
- */
-template <int dim, typename Number>
-inline
-void contract (Tensor<2,dim,Number>       &dest,
-               const Tensor<2,dim,Number> &src1,
-               const Tensor<2,dim,Number> &src2)
-{
-  for (unsigned int i=0; i<dim; ++i)
-    for (unsigned int j=0; j<dim; ++j)
-      {
-        dest[i][j] = src1[i][0] * src2[0][j];
-        for (unsigned int k=1; k<dim; ++k)
-          dest[i][j] += src1[i][k] * src2[k][j];
-      }
-}
-
-
 /**
  * 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, and
@@ -1766,28 +1678,6 @@ void contract (Tensor<2,dim,Number>       &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>.
- *
- * @relates Tensor
- * @author Wolfgang Bangerth, 1998
- */
-template <int dim, typename Number>
-inline
-void contract (Tensor<3,dim,Number>       &dest,
-               const Tensor<3,dim,Number> &src1,
-               const Tensor<2,dim,Number> &src2)
-{
-  dest.clear ();
-  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)
-          dest[i][j][k] += src1[i][j][l] * src2[l][k];
-}
-
-
 /**
  * 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, and
@@ -1887,28 +1777,6 @@ void contract (Tensor<3,dim,Number>       &dest,
 }
 
 
-/**
- * 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>.
- *
- * @relates Tensor
- * @author Wolfgang Bangerth, 1998
- */
-template <int dim, typename Number>
-inline
-void contract (Tensor<3,dim,Number>       &dest,
-               const Tensor<2,dim,Number> &src1,
-               const Tensor<3,dim,Number> &src2)
-{
-  dest.clear ();
-  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)
-          dest[i][j][k] += src1[i][l] * src2[l][j][k];
-}
-
-
 /**
  * Contract the last two indices of <tt>src1</tt> with the two indices
  * <tt>src2</tt>, creating a rank-2 tensor. This is the matrix-vector product
index f8eab6f503b8087ef7bc152e2248e9ae490b8981..2aaf1b08e6f9435d0a9513c87d876c20f63c765f 100644 (file)
@@ -28,6 +28,33 @@ DEAL_II_NAMESPACE_OPEN
  */
 //@{
 
+/**
+ * Single contraction for tensors: contract the last index of a tensor @p
+ * src1 of rank @p rank_1 with the first index of a tensor @p src2 of rank
+ * @p rank_2.
+ *
+ * @deprecated Use operator* instead. It denotes a single contraction.
+ * @relates Tensor
+ */
+template <int rank_1, int rank_2, int dim, typename Number>
+inline
+void contract (Tensor<rank_1 + rank_2 - 2, dim, Number> &dest,
+               const Tensor<rank_1 ,dim, Number>        &src1,
+               const Tensor<rank_2 ,dim, Number>        &src2) DEAL_II_DEPRECATED;
+
+/**
+ * Contract a tensor of rank 1 with a tensor of rank 1 and return the
+ * result.
+ *
+ * @deprecated Use operator* instead. It denotes a single contraction.
+ * @relates Tensor
+ */
+template <int dim, typename Number, typename OtherNumber>
+inline
+typename ProductType<Number,OtherNumber>::type
+contract (const Tensor<1,dim,Number> &src1,
+          const Tensor<1,dim,OtherNumber> &src2) DEAL_II_DEPRECATED;
+
 /**
  * The cross-product of 2 vectors in 3d.
  *
@@ -96,6 +123,31 @@ Number determinant (const Tensor<1,1,Number> &t) DEAL_II_DEPRECATED;
 
 #ifndef DOXYGEN
 
+template <int rank_1, int rank_2, int dim, typename Number>
+inline
+void contract (Tensor<rank_1 + rank_2 - 2, dim, Number> &dest,
+               const Tensor<rank_1 ,dim, Number>        &src1,
+               const Tensor<rank_2 ,dim, Number>        &src2)
+{
+  TensorAccessors::internal::ReorderedIndexView<0, rank_2, const Tensor<rank_2, dim, Number> >
+  reordered = TensorAccessors::reordered_index_view<0, rank_2>(src2);
+  TensorAccessors::contract<1, rank_1, rank_2, dim>(dest, src1, reordered);
+}
+
+template <int dim, typename Number, typename OtherNumber>
+inline
+typename ProductType<Number,OtherNumber>::type
+contract (const Tensor<1,dim,Number> &src1,
+          const Tensor<1,dim,OtherNumber> &src2)
+{
+  typename ProductType<Number,OtherNumber>::type res
+    = typename ProductType<Number,OtherNumber>::type();
+  for (unsigned int i=0; i<dim; ++i)
+    res += src1[i] * src2[i];
+
+  return res;
+}
+
 template <int dim, typename Number>
 inline
 void

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.