]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Cleanup use of DEAL_II_CONSTEXPR
authorDaniel Arndt <arndtd@ornl.gov>
Tue, 6 Aug 2019 22:54:24 +0000 (16:54 -0600)
committerDaniel Arndt <arndtd@ornl.gov>
Tue, 6 Aug 2019 23:07:19 +0000 (17:07 -0600)
include/deal.II/base/numbers.h
include/deal.II/base/symmetric_tensor.h
include/deal.II/base/tensor.h
tests/tensors/constexpr_symmetric_tensor.cc
tests/tensors/constexpr_tensor.cc
tests/tensors/constexpr_tensor.output

index e264687af3930580852476fdd7175911cc993b62..332d0d7ef45140c951491915b2cf76207111d10c 100644 (file)
@@ -336,7 +336,7 @@ namespace numbers
    * of @p value_1.
    */
   template <typename Number1, typename Number2>
-  bool
+  constexpr bool
   values_are_equal(const Number1 &value_1, const Number2 &value_2);
 
   /**
@@ -361,7 +361,7 @@ namespace numbers
    * by the input arguments.
    */
   template <typename Number>
-  bool
+  constexpr bool
   value_is_zero(const Number &value);
 
   /**
@@ -937,7 +937,7 @@ namespace numbers
 
 
   template <typename Number1, typename Number2>
-  inline bool
+  constexpr bool
   values_are_equal(const Number1 &value_1, const Number2 &value_2)
   {
     return (value_1 == internal::NumberType<Number1>::value(value_2));
@@ -953,7 +953,7 @@ namespace numbers
 
 
   template <typename Number>
-  inline bool
+  constexpr bool
   value_is_zero(const Number &value)
   {
     return values_are_equal(value, 0.0);
index b9a7adc70386363824315af8aba119b348cc0438..a552189bb38aa1533629d7eff0c04d9e1733a318 100644 (file)
@@ -593,10 +593,12 @@ public:
    * is symmetric only up to round-off, then you may want to call the
    * <tt>symmetrize</tt> function first. If you aren't sure, it is good
    * practice to check before calling <tt>symmetrize</tt>.
+   *
+   * If you want to use this function in a `constexpr` context,
+   * use symmetrize() instead.
    */
   template <typename OtherNumber>
-  DEAL_II_CONSTEXPR explicit SymmetricTensor(
-    const Tensor<2, dim, OtherNumber> &t);
+  explicit SymmetricTensor(const Tensor<2, dim, OtherNumber> &t);
 
   /**
    * A constructor that creates a symmetric tensor from an array holding its
@@ -628,26 +630,26 @@ public:
   /**
    * Return a pointer to the first element of the underlying storage.
    */
-  DEAL_II_CONSTEXPR Number *
-                    begin_raw();
+  Number *
+  begin_raw();
 
   /**
    * Return a const pointer to the first element of the underlying storage.
    */
-  constexpr const Number *
+  const Number *
   begin_raw() const;
 
   /**
    * Return a pointer to the element past the end of the underlying storage.
    */
-  DEAL_II_CONSTEXPR Number *
-                    end_raw();
+  Number *
+  end_raw();
 
   /**
    * Return a const pointer to the element past the end of the underlying
    * storage.
    */
-  constexpr const Number *
+  const Number *
   end_raw() const;
 
   /**
@@ -1022,7 +1024,7 @@ namespace internal
 
 template <int rank_, int dim, typename Number>
 template <typename OtherNumber>
-DEAL_II_CONSTEXPR inline DEAL_II_ALWAYS_INLINE
+inline DEAL_II_ALWAYS_INLINE
 SymmetricTensor<rank_, dim, Number>::SymmetricTensor(
   const Tensor<2, dim, OtherNumber> &t)
 {
@@ -2070,7 +2072,7 @@ DEAL_II_CONSTEXPR inline DEAL_II_ALWAYS_INLINE Number &
 
 
 template <int rank_, int dim, typename Number>
-DEAL_II_CONSTEXPR inline Number *
+inline Number *
 SymmetricTensor<rank_, dim, Number>::begin_raw()
 {
   return std::addressof(this->access_raw_entry(0));
@@ -2079,7 +2081,7 @@ SymmetricTensor<rank_, dim, Number>::begin_raw()
 
 
 template <int rank_, int dim, typename Number>
-constexpr const Number *
+inline const Number *
 SymmetricTensor<rank_, dim, Number>::begin_raw() const
 {
   return std::addressof(this->access_raw_entry(0));
@@ -2088,7 +2090,7 @@ SymmetricTensor<rank_, dim, Number>::begin_raw() const
 
 
 template <int rank_, int dim, typename Number>
-DEAL_II_CONSTEXPR inline Number *
+inline Number *
 SymmetricTensor<rank_, dim, Number>::end_raw()
 {
   return begin_raw() + n_independent_components;
@@ -2097,7 +2099,7 @@ SymmetricTensor<rank_, dim, Number>::end_raw()
 
 
 template <int rank_, int dim, typename Number>
-constexpr const Number *
+inline const Number *
 SymmetricTensor<rank_, dim, Number>::end_raw() const
 {
   return begin_raw() + n_independent_components;
index 90fc129269d9603ba8aed3dd7be4c0afc3d96ef7..8d96f7d4a278079822429cce6029761ae25be3f0 100644 (file)
@@ -167,26 +167,26 @@ public:
   /**
    * Return a pointer to the first element of the underlying storage.
    */
-  DEAL_II_CONSTEXPR Number *
-                    begin_raw();
+  Number *
+  begin_raw();
 
   /**
    * Return a const pointer to the first element of the underlying storage.
    */
-  constexpr const Number *
+  const Number *
   begin_raw() const;
 
   /**
    * Return a pointer to the element past the end of the underlying storage.
    */
-  DEAL_II_CONSTEXPR Number *
-                    end_raw();
+  Number *
+  end_raw();
 
   /**
    * Return a const pointer to the element past the end of the underlying
    * storage.
    */
-  constexpr const Number *
+  const Number *
   end_raw() const;
 
   /**
@@ -322,8 +322,8 @@ public:
    * the absolute squares of all entries. For the present case of rank-1
    * tensors, this equals the usual <tt>l<sub>2</sub></tt> norm of the vector.
    */
-  DEAL_II_CONSTEXPR real_type
-                    norm() const;
+  real_type
+  norm() const;
 
   /**
    * Return the square of the Frobenius-norm of a tensor, i.e. the sum of the
@@ -331,8 +331,8 @@ public:
    *
    * @note This function can also be used in CUDA device code.
    */
-  DEAL_II_CUDA_HOST_DEV real_type
-                        norm_square() const;
+  DEAL_II_CONSTEXPR DEAL_II_CUDA_HOST_DEV real_type
+                                          norm_square() const;
 
   /**
    * Read or write the data of this object to or from a stream for the purpose
@@ -521,25 +521,25 @@ public:
   /**
    * Return a pointer to the first element of the underlying storage.
    */
-  DEAL_II_CONSTEXPR Number *
-                    begin_raw();
+  Number *
+  begin_raw();
 
   /**
    * Return a const pointer to the first element of the underlying storage.
    */
-  constexpr const Number *
+  const Number *
   begin_raw() const;
 
   /**
    * Return a pointer to the element past the end of the underlying storage.
    */
-  DEAL_II_CONSTEXPR Number *
-                    end_raw();
+  Number *
+  end_raw();
 
   /**
    * Return a pointer to the element past the end of the underlying storage.
    */
-  constexpr const Number *
+  const Number *
   end_raw() const;
 
   /**
@@ -643,10 +643,9 @@ public:
    *
    * @note This function can also be used in CUDA device code.
    */
-
-  DEAL_II_CONSTEXPR DEAL_II_CUDA_HOST_DEV
-    typename numbers::NumberTraits<Number>::real_type
-    norm() const;
+  DEAL_II_CUDA_HOST_DEV
+  typename numbers::NumberTraits<Number>::real_type
+  norm() const;
 
   /**
    * Return the square of the Frobenius-norm of a tensor, i.e. the sum of the
@@ -666,7 +665,7 @@ public:
    * fastest.
    */
   template <typename OtherNumber>
-  DEAL_II_CONSTEXPR void
+  void
   unroll(Vector<OtherNumber> &result) const;
 
   /**
@@ -716,7 +715,7 @@ private:
    * Internal helper function for unroll.
    */
   template <typename OtherNumber>
-  DEAL_II_CONSTEXPR void
+  void
   unroll_recursion(Vector<OtherNumber> &result,
                    unsigned int &       start_index) const;
 
@@ -798,7 +797,7 @@ constexpr DEAL_II_ALWAYS_INLINE DEAL_II_CUDA_HOST_DEV
 
 
 template <int dim, typename Number>
-DEAL_II_CONSTEXPR inline Number *
+inline Number *
 Tensor<0, dim, Number>::begin_raw()
 {
   return std::addressof(value);
@@ -807,7 +806,7 @@ Tensor<0, dim, Number>::begin_raw()
 
 
 template <int dim, typename Number>
-constexpr const Number *
+inline const Number *
 Tensor<0, dim, Number>::begin_raw() const
 {
   return std::addressof(value);
@@ -816,7 +815,7 @@ Tensor<0, dim, Number>::begin_raw() const
 
 
 template <int dim, typename Number>
-DEAL_II_CONSTEXPR inline Number *
+inline Number *
 Tensor<0, dim, Number>::end_raw()
 {
   return begin_raw() + n_independent_components;
@@ -825,7 +824,7 @@ Tensor<0, dim, Number>::end_raw()
 
 
 template <int dim, typename Number>
-constexpr const Number *
+const Number *
 Tensor<0, dim, Number>::end_raw() const
 {
   return begin_raw() + n_independent_components;
@@ -997,9 +996,8 @@ Tensor<0, dim, Number>::operator-() const
 
 
 template <int dim, typename Number>
-DEAL_II_CONSTEXPR inline DEAL_II_ALWAYS_INLINE
-  typename Tensor<0, dim, Number>::real_type
-  Tensor<0, dim, Number>::norm() const
+inline DEAL_II_ALWAYS_INLINE typename Tensor<0, dim, Number>::real_type
+Tensor<0, dim, Number>::norm() const
 {
   Assert(dim != 0,
          ExcMessage("Cannot access an object of type Tensor<0,0,Number>"));
@@ -1008,7 +1006,7 @@ DEAL_II_CONSTEXPR inline DEAL_II_ALWAYS_INLINE
 
 
 template <int dim, typename Number>
-DEAL_II_CUDA_HOST_DEV inline DEAL_II_ALWAYS_INLINE
+DEAL_II_CONSTEXPR DEAL_II_CUDA_HOST_DEV inline
   typename Tensor<0, dim, Number>::real_type
   Tensor<0, dim, Number>::norm_square() const
 {
@@ -1023,7 +1021,7 @@ DEAL_II_CUDA_HOST_DEV inline DEAL_II_ALWAYS_INLINE
 
 template <int dim, typename Number>
 template <typename OtherNumber>
-DEAL_II_CONSTEXPR inline void
+inline DEAL_II_CONSTEXPR void
 Tensor<0, dim, Number>::unroll_recursion(Vector<OtherNumber> &result,
                                          unsigned int &       index) const
 {
@@ -1193,7 +1191,7 @@ DEAL_II_CONSTEXPR inline DEAL_II_ALWAYS_INLINE Number &
 
 
 template <int rank_, int dim, typename Number>
-DEAL_II_CONSTEXPR inline Number *
+inline Number *
 Tensor<rank_, dim, Number>::begin_raw()
 {
   return std::addressof(
@@ -1203,7 +1201,7 @@ Tensor<rank_, dim, Number>::begin_raw()
 
 
 template <int rank_, int dim, typename Number>
-constexpr const Number *
+inline const Number *
 Tensor<rank_, dim, Number>::begin_raw() const
 {
   return std::addressof(
@@ -1213,7 +1211,7 @@ Tensor<rank_, dim, Number>::begin_raw() const
 
 
 template <int rank_, int dim, typename Number>
-DEAL_II_CONSTEXPR inline Number *
+inline Number *
 Tensor<rank_, dim, Number>::end_raw()
 {
   return begin_raw() + n_independent_components;
@@ -1222,7 +1220,7 @@ Tensor<rank_, dim, Number>::end_raw()
 
 
 template <int rank_, int dim, typename Number>
-constexpr const Number *
+inline const Number *
 Tensor<rank_, dim, Number>::end_raw() const
 {
   return begin_raw() + n_independent_components;
@@ -1235,8 +1233,9 @@ template <typename OtherNumber>
 DEAL_II_CONSTEXPR inline DEAL_II_ALWAYS_INLINE Tensor<rank_, dim, Number> &
 Tensor<rank_, dim, Number>::operator=(const Tensor<rank_, dim, OtherNumber> &t)
 {
-  if (dim > 0)
-    std::copy(&t.values[0], &t.values[0] + dim, &values[0]);
+  // std::copy is constexpr only from C++20 on.
+  for (unsigned int i = 0; i < dim; ++i)
+    values[i] = t.values[i];
   return *this;
 }
 
@@ -1357,7 +1356,7 @@ DEAL_II_CONSTEXPR inline DEAL_II_ALWAYS_INLINE
 
 
 template <int rank_, int dim, typename Number>
-DEAL_II_CONSTEXPR inline typename numbers::NumberTraits<Number>::real_type
+inline typename numbers::NumberTraits<Number>::real_type
 Tensor<rank_, dim, Number>::norm() const
 {
   return std::sqrt(norm_square());
@@ -1380,7 +1379,7 @@ DEAL_II_CONSTEXPR inline DEAL_II_ALWAYS_INLINE DEAL_II_CUDA_HOST_DEV
 
 template <int rank_, int dim, typename Number>
 template <typename OtherNumber>
-DEAL_II_CONSTEXPR inline void
+inline void
 Tensor<rank_, dim, Number>::unroll(Vector<OtherNumber> &result) const
 {
   AssertDimension(result.size(),
@@ -1393,7 +1392,7 @@ Tensor<rank_, dim, Number>::unroll(Vector<OtherNumber> &result) const
 
 template <int rank_, int dim, typename Number>
 template <typename OtherNumber>
-DEAL_II_CONSTEXPR inline void
+inline void
 Tensor<rank_, dim, Number>::unroll_recursion(Vector<OtherNumber> &result,
                                              unsigned int &       index) const
 {
@@ -1422,7 +1421,7 @@ namespace internal
   // and rank=2. Make sure we don't have compiler warnings.
 
   template <int dim>
-  inline unsigned int
+  inline DEAL_II_CONSTEXPR unsigned int
   mod(const unsigned int x)
   {
     return x % dim;
@@ -1432,12 +1431,12 @@ namespace internal
   inline unsigned int
   mod<0>(const unsigned int x)
   {
-    Assert(false, ExcInternalError());
+    // Assert(false, ExcInternalError());
     return x;
   }
 
   template <int dim>
-  inline unsigned int
+  inline DEAL_II_CONSTEXPR unsigned int
   div(const unsigned int x)
   {
     return x / dim;
@@ -1799,10 +1798,10 @@ DEAL_II_CONSTEXPR DEAL_II_CUDA_HOST_DEV inline DEAL_II_ALWAYS_INLINE
  * @relatesalso Tensor
  */
 template <int dim, typename Number, typename OtherNumber>
-inline DEAL_II_ALWAYS_INLINE
-  Tensor<0, dim, typename ProductType<Number, OtherNumber>::type>
-  schur_product(const Tensor<0, dim, Number> &     src1,
-                const Tensor<0, dim, OtherNumber> &src2)
+inline DEAL_II_CONSTEXPR DEAL_II_ALWAYS_INLINE
+                         Tensor<0, dim, typename ProductType<Number, OtherNumber>::type>
+                         schur_product(const Tensor<0, dim, Number> &     src1,
+                                       const Tensor<0, dim, OtherNumber> &src2)
 {
   Tensor<0, dim, typename ProductType<Number, OtherNumber>::type> tmp(src1);
 
@@ -1817,10 +1816,10 @@ inline DEAL_II_ALWAYS_INLINE
  * @relatesalso Tensor
  */
 template <int dim, typename Number, typename OtherNumber>
-inline DEAL_II_ALWAYS_INLINE
-  Tensor<1, dim, typename ProductType<Number, OtherNumber>::type>
-  schur_product(const Tensor<1, dim, Number> &     src1,
-                const Tensor<1, dim, OtherNumber> &src2)
+inline DEAL_II_CONSTEXPR DEAL_II_ALWAYS_INLINE
+                         Tensor<1, dim, typename ProductType<Number, OtherNumber>::type>
+                         schur_product(const Tensor<1, dim, Number> &     src1,
+                                       const Tensor<1, dim, OtherNumber> &src2)
 {
   Tensor<1, dim, typename ProductType<Number, OtherNumber>::type> tmp(src1);
 
@@ -1847,10 +1846,10 @@ inline DEAL_II_ALWAYS_INLINE
  * @relatesalso Tensor
  */
 template <int rank, int dim, typename Number, typename OtherNumber>
-inline DEAL_II_ALWAYS_INLINE
-  Tensor<rank, dim, typename ProductType<Number, OtherNumber>::type>
-  schur_product(const Tensor<rank, dim, Number> &     src1,
-                const Tensor<rank, dim, OtherNumber> &src2)
+inline DEAL_II_CONSTEXPR DEAL_II_ALWAYS_INLINE
+                         Tensor<rank, dim, typename ProductType<Number, OtherNumber>::type>
+                         schur_product(const Tensor<rank, dim, Number> &     src1,
+                                       const Tensor<rank, dim, OtherNumber> &src2)
 {
   Tensor<rank, dim, typename ProductType<Number, OtherNumber>::type> tmp(src1);
 
@@ -2088,7 +2087,7 @@ DEAL_II_CONSTEXPR inline
   typename Tensor<rank_1 + rank_2 - 4,
                   dim,
                   typename ProductType<Number, OtherNumber>::type>::tensor_type
-    result;
+    result{};
   TensorAccessors::contract<2, rank_1, rank_2, dim>(result, reord_3, reord_4);
   return result;
 }
@@ -2113,7 +2112,7 @@ DEAL_II_CONSTEXPR inline DEAL_II_ALWAYS_INLINE
   scalar_product(const Tensor<rank, dim, Number> &     left,
                  const Tensor<rank, dim, OtherNumber> &right)
 {
-  typename ProductType<Number, OtherNumber>::type result;
+  typename ProductType<Number, OtherNumber>::type result{};
   TensorAccessors::contract<rank, rank, rank, dim>(result, left, right);
   return result;
 }
@@ -2491,8 +2490,8 @@ cofactor(const Tensor<2, dim, Number> &t)
  * @author Wolfgang Bangerth, 2012
  */
 template <int dim, typename Number>
-inline Number
-l1_norm(const Tensor<2, dim, Number> &t)
+inline DEAL_II_CONSTEXPR Number
+                         l1_norm(const Tensor<2, dim, Number> &t)
 {
   Number max = internal::NumberType<Number>::value(0.0);
   for (unsigned int j = 0; j < dim; ++j)
@@ -2517,8 +2516,8 @@ l1_norm(const Tensor<2, dim, Number> &t)
  * @author Wolfgang Bangerth, 2012
  */
 template <int dim, typename Number>
-inline Number
-linfty_norm(const Tensor<2, dim, Number> &t)
+inline DEAL_II_CONSTEXPR Number
+                         linfty_norm(const Tensor<2, dim, Number> &t)
 {
   Number max = internal::NumberType<Number>::value(0.0);
   for (unsigned int i = 0; i < dim; ++i)
index 6d5fca5ac7009171abe4af3e585f3ee3ba6a274d..bf4ef812a82108c689d794639db35f06b1a2797e 100644 (file)
@@ -111,6 +111,37 @@ main()
   DEAL_II_CONSTEXPR const auto B = get_tensor_4();
   deallog << "SymmetricTensor<4,2> = " << B << std::endl;
 
+  {
+    DEAL_II_CONSTEXPR double a_init[3][3] = {{1., 0., 0.},
+                                             {2., 1., 0.},
+                                             {3., 2., 1.}};
+    DEAL_II_CONSTEXPR Tensor<2, 3> dummy_a{a_init};
+    DEAL_II_CONSTEXPR SymmetricTensor<2, 3> a = symmetrize(dummy_a);
+    std::cout << a << std::endl;
+    DEAL_II_CONSTEXPR auto inverted = invert(a);
+    std::cout << inverted << std::endl;
+    DEAL_II_CONSTEXPR double ref_init[3][3] = {{0., -2., 2.},
+                                               {-2., 5., -2.},
+                                               {2., -2., 0.}};
+    DEAL_II_CONSTEXPR Tensor<2, 3> dummy_ref{ref_init};
+    DEAL_II_CONSTEXPR SymmetricTensor<2, 3> ref = symmetrize(dummy_ref);
+    std::cout << ref << std::endl;
+    Assert(inverted == ref, ExcInternalError());
+  }
+  {
+    DEAL_II_CONSTEXPR double a_init[3][3] = {{1., 2., 3.},
+                                             {2., 1., 2.},
+                                             {3., 2., 1.}};
+    DEAL_II_CONSTEXPR Tensor<2, 3> dummy_a{a_init};
+    DEAL_II_CONSTEXPR SymmetricTensor<2, 3> a          = symmetrize(dummy_a);
+    DEAL_II_CONSTEXPR auto                  transposed = transpose(a);
+    Assert(transposed == a, ExcInternalError());
+    constexpr auto dummy   = scalar_product(a, a);
+    constexpr auto dummy_5 = a * a;
+    constexpr auto middle  = outer_product(a, a);
+    constexpr auto dummy_6 = contract3(a, middle, a);
+  }
+
   deallog << "OK" << std::endl;
   return 0;
 }
index 1afa257a92aa73f1beeb0b58055306c085090c15..a5cc8696f6bd74df42385e2227ba692228771a3d 100644 (file)
 
 template <int rank, int dim, typename Number>
 void
-test_consexpr_tensor_constructors()
+test_constexpr_tensor_constructors()
 {
-  constexpr dealii::Tensor<rank, dim, Number> A;
-  constexpr dealii::Tensor<rank, dim, Number> B(A);
-  constexpr dealii::Tensor<rank, dim, Number> C = A;
+  constexpr dealii::Tensor<rank, dim, Number> a;
+  constexpr dealii::Tensor<rank, dim, Number> b(a);
+  constexpr dealii::Tensor<rank, dim, Number> c = a;
   deallog << " Tensor<" << rank << "," << dim << ">" << std::endl;
-  deallog << A << std::endl;
-  deallog << B << std::endl;
-  deallog << C << std::endl;
+  deallog << a << std::endl;
+  deallog << b << std::endl;
+  deallog << c << std::endl;
+
+  constexpr auto d = a * 2.;
+  constexpr auto e = a / 2.;
+  constexpr auto f = d + e;
+  constexpr auto g = d - e;
+  deallog << d << std::endl;
+  deallog << e << std::endl;
+  deallog << f << std::endl;
+  deallog << g << std::endl;
 }
 
 DEAL_II_CONSTEXPR double
@@ -42,11 +51,29 @@ tensor_rank0_constexpr(const double value)
 DEAL_II_CONSTEXPR double
 tensor_rank2_constexpr()
 {
-  constexpr double       a[3][3] = {{4., 2., 0.}, {2., 3., 0.}, {0., 0., 5.}};
-  constexpr Tensor<2, 3> A(a);
-  DEAL_II_CONSTEXPR const auto det = determinant(A);
-  DEAL_II_CONSTEXPR const auto tr  = trace(transpose(A) * A);
-  return det + tr;
+  constexpr double a_init[3][3] = {{4., 2., 0.}, {2., 3., 0.}, {0., 0., 5.}};
+  constexpr Tensor<2, 3>       a(a_init);
+  DEAL_II_CONSTEXPR const auto det   = determinant(a);
+  DEAL_II_CONSTEXPR const auto tr    = trace(transpose(a) * a);
+  DEAL_II_CONSTEXPR const auto dummy = symmetrize(a);
+
+  auto d = a;
+  d /= 2.;
+  const auto test_1 = d[0][0];
+  d *= 2.;
+  const auto test_2 = d[0][0];
+  d += a;
+  const auto test_3 = d[0][0];
+  d -= a;
+  const auto test_4 = d[0][0];
+  d.clear();
+  const auto test_5 = d[0][0];
+  d                 = 0.;
+  const auto test_6 = d[0][0];
+  const auto test_7 = d.norm_square();
+
+  return det + tr + test_1 + test_2 + test_3 + test_4 + test_5 + test_6 +
+         test_7;
 }
 
 int
@@ -54,32 +81,31 @@ main()
 {
   initlog();
 
-  deallog << "Cheching constexpr default constructor of Tensor<rank,dim,Number>"
+  deallog << "Checking constexpr default constructor of Tensor<rank,dim,Number>"
           << std::endl;
   {
     dealii::LogStream::Prefix p("float");
-    test_consexpr_tensor_constructors<0, 1, float>();
-    test_consexpr_tensor_constructors<0, 2, float>();
-    test_consexpr_tensor_constructors<0, 3, float>();
-    test_consexpr_tensor_constructors<1, 1, float>();
-    test_consexpr_tensor_constructors<1, 2, float>();
-    test_consexpr_tensor_constructors<1, 3, float>();
-    test_consexpr_tensor_constructors<2, 1, float>();
-    test_consexpr_tensor_constructors<2, 2, float>();
-    test_consexpr_tensor_constructors<2, 3, float>();
+    test_constexpr_tensor_constructors<0, 1, float>();
+    test_constexpr_tensor_constructors<0, 2, float>();
+    test_constexpr_tensor_constructors<0, 3, float>();
+    test_constexpr_tensor_constructors<1, 1, float>();
+    test_constexpr_tensor_constructors<1, 2, float>();
+    test_constexpr_tensor_constructors<1, 3, float>();
+    test_constexpr_tensor_constructors<2, 1, float>();
+    test_constexpr_tensor_constructors<2, 2, float>();
+    test_constexpr_tensor_constructors<2, 3, float>();
   }
-
   {
     dealii::LogStream::Prefix p("double");
-    test_consexpr_tensor_constructors<0, 1, double>();
-    test_consexpr_tensor_constructors<0, 2, double>();
-    test_consexpr_tensor_constructors<0, 3, double>();
-    test_consexpr_tensor_constructors<1, 1, double>();
-    test_consexpr_tensor_constructors<1, 2, double>();
-    test_consexpr_tensor_constructors<1, 3, double>();
-    test_consexpr_tensor_constructors<2, 1, double>();
-    test_consexpr_tensor_constructors<2, 2, double>();
-    test_consexpr_tensor_constructors<2, 3, double>();
+    test_constexpr_tensor_constructors<0, 1, double>();
+    test_constexpr_tensor_constructors<0, 2, double>();
+    test_constexpr_tensor_constructors<0, 3, double>();
+    test_constexpr_tensor_constructors<1, 1, double>();
+    test_constexpr_tensor_constructors<1, 2, double>();
+    test_constexpr_tensor_constructors<1, 3, double>();
+    test_constexpr_tensor_constructors<2, 1, double>();
+    test_constexpr_tensor_constructors<2, 2, double>();
+    test_constexpr_tensor_constructors<2, 3, double>();
   }
 
   deallog << "Using Tensor within constexpr functions" << std::endl;
@@ -90,6 +116,75 @@ main()
   DEAL_II_CONSTEXPR const auto tensor_rank2_result = tensor_rank2_constexpr();
   deallog << tensor_rank2_result << std::endl;
 
+  {
+    DEAL_II_CONSTEXPR double initializer[2] = {1., -1.};
+    DEAL_II_CONSTEXPR Tensor<1, 2> c{initializer};
+    DEAL_II_CONSTEXPR Tensor<1, 2> c_cross            = cross_product_2d(c);
+    DEAL_II_CONSTEXPR double       initializer_ref[2] = {-1., -1.};
+    DEAL_II_CONSTEXPR Tensor<1, 2> ref{initializer_ref};
+    DEAL_II_CONSTEXPR bool         is_same     = (c_cross == ref);
+    DEAL_II_CONSTEXPR bool         is_not_same = (c_cross != ref);
+    Assert(is_same && !is_not_same, ExcInternalError());
+  }
+  {
+    DEAL_II_CONSTEXPR double initializer_1[3] = {1., 0., 0.};
+    DEAL_II_CONSTEXPR Tensor<1, 3> c_1{initializer_1};
+    DEAL_II_CONSTEXPR double       initializer_2[3] = {0., 1., 0.};
+    DEAL_II_CONSTEXPR Tensor<1, 3> c_2{initializer_2};
+
+    DEAL_II_CONSTEXPR auto   c_cross            = cross_product_3d(c_1, c_2);
+    DEAL_II_CONSTEXPR double initializer_ref[3] = {0., 0., 1.};
+    DEAL_II_CONSTEXPR Tensor<1, 3> ref{initializer_ref};
+    DEAL_II_CONSTEXPR bool         is_same     = (c_cross == ref);
+    DEAL_II_CONSTEXPR bool         is_not_same = (c_cross != ref);
+    Assert(is_same && !is_not_same, ExcInternalError());
+  }
+
+  DEAL_II_CONSTEXPR auto table_indices =
+    Tensor<2, 3>::unrolled_to_component_indices(0);
+  DEAL_II_CONSTEXPR auto index =
+    Tensor<2, 3>::component_to_unrolled_index(TableIndices<2>{});
+  Assert(index == 0, ExcInternalError());
+
+  DEAL_II_CONSTEXPR auto used_memory = Tensor<2, 3>::memory_consumption();
+  deallog << "Used memory: " << used_memory << std::endl;
+
+  {
+    DEAL_II_CONSTEXPR double a_init[3][3] = {{1., 0., 0.},
+                                             {2., 1., 0.},
+                                             {3., 2., 1.}};
+    DEAL_II_CONSTEXPR Tensor<2, 3> a{a_init};
+    DEAL_II_CONSTEXPR auto         inverted       = invert(a);
+    DEAL_II_CONSTEXPR double       ref_init[3][3] = {{1., 0., 0.},
+                                               {-2., 1., 0.},
+                                               {1., -2., 1}};
+    DEAL_II_CONSTEXPR Tensor<2, 3> ref{ref_init};
+    Assert(inverted == ref, ExcInternalError());
+  }
+  {
+    DEAL_II_CONSTEXPR double a_init[3][3] = {{1., 0., 0.},
+                                             {2., 1., 0.},
+                                             {3., 2., 1.}};
+    DEAL_II_CONSTEXPR Tensor<2, 3> a{a_init};
+    DEAL_II_CONSTEXPR auto         transposed     = transpose(a);
+    DEAL_II_CONSTEXPR double       ref_init[3][3] = {{1., 2., 3.},
+                                               {0., 1., 2.},
+                                               {0., 0., 1}};
+    DEAL_II_CONSTEXPR Tensor<2, 3> ref{ref_init};
+    Assert(transposed == ref, ExcInternalError());
+    constexpr auto dummy    = scalar_product(a, ref);
+    constexpr auto dummy_2  = contract<0, 0>(a, ref);
+    constexpr auto dummy_3  = double_contract<0, 0, 1, 1>(a, ref);
+    constexpr auto dummy_4  = schur_product(a, ref);
+    constexpr auto dummy_5  = a * ref;
+    constexpr auto middle   = outer_product(a, a);
+    constexpr auto dummy_6  = contract3(a, middle, a);
+    constexpr auto dummy_7  = adjugate(a);
+    constexpr auto dummy_8  = cofactor(a);
+    constexpr auto dummy_9  = l1_norm(a);
+    constexpr auto dummy_10 = linfty_norm(a);
+  }
+
   deallog << "OK" << std::endl;
   return 0;
 }
index 7c10bec9a8ffe0242f3a40d83b2fded00ef00085..218e060df40a10f4370bbe8aae793ff886708073 100644 (file)
 
-DEAL::Cheching constexpr default constructor of Tensor<rank,dim,Number>
+DEAL::Checking constexpr default constructor of Tensor<rank,dim,Number>
 DEAL:float:: Tensor<0,1>
 DEAL:float::0.00000
 DEAL:float::0.00000
 DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
 DEAL:float:: Tensor<0,2>
 DEAL:float::0.00000
 DEAL:float::0.00000
 DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
 DEAL:float:: Tensor<0,3>
 DEAL:float::0.00000
 DEAL:float::0.00000
 DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
 DEAL:float:: Tensor<1,1>
 DEAL:float::0.00000
 DEAL:float::0.00000
 DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
 DEAL:float:: Tensor<1,2>
 DEAL:float::0.00000 0.00000
 DEAL:float::0.00000 0.00000
 DEAL:float::0.00000 0.00000
+DEAL:float::0.00000 0.00000
+DEAL:float::0.00000 0.00000
+DEAL:float::0.00000 0.00000
+DEAL:float::0.00000 0.00000
 DEAL:float:: Tensor<1,3>
 DEAL:float::0.00000 0.00000 0.00000
 DEAL:float::0.00000 0.00000 0.00000
 DEAL:float::0.00000 0.00000 0.00000
+DEAL:float::0.00000 0.00000 0.00000
+DEAL:float::0.00000 0.00000 0.00000
+DEAL:float::0.00000 0.00000 0.00000
+DEAL:float::0.00000 0.00000 0.00000
 DEAL:float:: Tensor<2,1>
 DEAL:float::0.00000
 DEAL:float::0.00000
 DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
+DEAL:float::0.00000
 DEAL:float:: Tensor<2,2>
 DEAL:float::0.00000 0.00000 0.00000 0.00000
 DEAL:float::0.00000 0.00000 0.00000 0.00000
 DEAL:float::0.00000 0.00000 0.00000 0.00000
+DEAL:float::0.00000 0.00000 0.00000 0.00000
+DEAL:float::0.00000 0.00000 0.00000 0.00000
+DEAL:float::0.00000 0.00000 0.00000 0.00000
+DEAL:float::0.00000 0.00000 0.00000 0.00000
 DEAL:float:: Tensor<2,3>
 DEAL:float::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
 DEAL:float::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
 DEAL:float::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:float::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:float::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:float::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:float::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
 DEAL:double:: Tensor<0,1>
 DEAL:double::0.00000
 DEAL:double::0.00000
 DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
 DEAL:double:: Tensor<0,2>
 DEAL:double::0.00000
 DEAL:double::0.00000
 DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
 DEAL:double:: Tensor<0,3>
 DEAL:double::0.00000
 DEAL:double::0.00000
 DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
 DEAL:double:: Tensor<1,1>
 DEAL:double::0.00000
 DEAL:double::0.00000
 DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
 DEAL:double:: Tensor<1,2>
 DEAL:double::0.00000 0.00000
 DEAL:double::0.00000 0.00000
 DEAL:double::0.00000 0.00000
+DEAL:double::0.00000 0.00000
+DEAL:double::0.00000 0.00000
+DEAL:double::0.00000 0.00000
+DEAL:double::0.00000 0.00000
 DEAL:double:: Tensor<1,3>
 DEAL:double::0.00000 0.00000 0.00000
 DEAL:double::0.00000 0.00000 0.00000
 DEAL:double::0.00000 0.00000 0.00000
+DEAL:double::0.00000 0.00000 0.00000
+DEAL:double::0.00000 0.00000 0.00000
+DEAL:double::0.00000 0.00000 0.00000
+DEAL:double::0.00000 0.00000 0.00000
 DEAL:double:: Tensor<2,1>
 DEAL:double::0.00000
 DEAL:double::0.00000
 DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
+DEAL:double::0.00000
 DEAL:double:: Tensor<2,2>
 DEAL:double::0.00000 0.00000 0.00000 0.00000
 DEAL:double::0.00000 0.00000 0.00000 0.00000
 DEAL:double::0.00000 0.00000 0.00000 0.00000
+DEAL:double::0.00000 0.00000 0.00000 0.00000
+DEAL:double::0.00000 0.00000 0.00000 0.00000
+DEAL:double::0.00000 0.00000 0.00000 0.00000
+DEAL:double::0.00000 0.00000 0.00000 0.00000
 DEAL:double:: Tensor<2,3>
 DEAL:double::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
 DEAL:double::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
 DEAL:double::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:double::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:double::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:double::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
+DEAL:double::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
 DEAL::Using Tensor within constexpr functions
 DEAL::15.2000
-DEAL::98.0000
+DEAL::116.000
+DEAL::Used memory: 72
 DEAL::OK

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.