]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Address comments
authorDaniel Arndt <arndtd@ornl.gov>
Wed, 7 Aug 2019 15:43:00 +0000 (09:43 -0600)
committerDaniel Arndt <arndtd@ornl.gov>
Wed, 7 Aug 2019 15:44:46 +0000 (09:44 -0600)
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 a552189bb38aa1533629d7eff0c04d9e1733a318..83901d0d48463410fc4c8ee3fc07bd036f5e6022 100644 (file)
@@ -596,6 +596,9 @@ public:
    *
    * If you want to use this function in a `constexpr` context,
    * use symmetrize() instead.
+   *
+   * Because we check for symmetry via a non-constexpr function call, you will
+   * have to use the symmetrize() function in constexpr contexts instead.
    */
   template <typename OtherNumber>
   explicit SymmetricTensor(const Tensor<2, dim, OtherNumber> &t);
index 8d96f7d4a278079822429cce6029761ae25be3f0..49725bfcee55e85cfc8d2f5326cb82afa6812a2f 100644 (file)
@@ -358,7 +358,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;
 
@@ -1006,7 +1006,7 @@ Tensor<0, dim, Number>::norm() const
 
 
 template <int dim, typename Number>
-DEAL_II_CONSTEXPR DEAL_II_CUDA_HOST_DEV inline
+DEAL_II_CONSTEXPR DEAL_II_CUDA_HOST_DEV inline DEAL_II_ALWAYS_INLINE
   typename Tensor<0, dim, Number>::real_type
   Tensor<0, dim, Number>::norm_square() const
 {
@@ -1021,7 +1021,7 @@ DEAL_II_CONSTEXPR DEAL_II_CUDA_HOST_DEV inline
 
 template <int dim, typename Number>
 template <typename OtherNumber>
-inline DEAL_II_CONSTEXPR void
+inline void
 Tensor<0, dim, Number>::unroll_recursion(Vector<OtherNumber> &result,
                                          unsigned int &       index) const
 {
@@ -1233,7 +1233,8 @@ 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)
 {
-  // std::copy is constexpr only from C++20 on.
+  // The following loop could be written more concisely using std::copy, but
+  // that function is only constexpr from C++20 on.
   for (unsigned int i = 0; i < dim; ++i)
     values[i] = t.values[i];
   return *this;
@@ -1431,7 +1432,7 @@ namespace internal
   inline unsigned int
   mod<0>(const unsigned int x)
   {
-    // Assert(false, ExcInternalError());
+    Assert(false, ExcInternalError());
     return x;
   }
 
@@ -1563,7 +1564,7 @@ operator<<(std::ostream &out, const Tensor<0, dim, Number> &p)
  * @relatesalso Tensor<0,dim,Number>
  */
 template <int dim, typename Number, typename Other>
-constexpr DEAL_II_CUDA_HOST_DEV DEAL_II_ALWAYS_INLINE
+DEAL_II_CONSTEXPR DEAL_II_CUDA_HOST_DEV DEAL_II_ALWAYS_INLINE
   typename ProductType<Other, Number>::type
   operator*(const Other &object, const Tensor<0, dim, Number> &t)
 {
@@ -1583,7 +1584,7 @@ constexpr DEAL_II_CUDA_HOST_DEV DEAL_II_ALWAYS_INLINE
  * @relatesalso Tensor<0,dim,Number>
  */
 template <int dim, typename Number, typename Other>
-constexpr DEAL_II_CUDA_HOST_DEV DEAL_II_ALWAYS_INLINE
+DEAL_II_CONSTEXPR DEAL_II_CUDA_HOST_DEV DEAL_II_ALWAYS_INLINE
   typename ProductType<Number, Other>::type
   operator*(const Tensor<0, dim, Number> &t, const Other &object)
 {
@@ -1707,12 +1708,12 @@ DEAL_II_CONSTEXPR DEAL_II_CUDA_HOST_DEV inline DEAL_II_ALWAYS_INLINE
  * @relatesalso Tensor
  */
 template <int rank, int dim, typename Number, typename OtherNumber>
-DEAL_II_CUDA_HOST_DEV constexpr DEAL_II_ALWAYS_INLINE
-  Tensor<rank,
+DEAL_II_CUDA_HOST_DEV DEAL_II_CONSTEXPR DEAL_II_ALWAYS_INLINE
+                                        Tensor<rank,
          dim,
          typename ProductType<typename EnableIfScalar<Number>::type,
                               OtherNumber>::type>
-  operator*(const Number &factor, const Tensor<rank, dim, OtherNumber> &t)
+                                        operator*(const Number &factor, const Tensor<rank, dim, OtherNumber> &t)
 {
   // simply forward to the operator above
   return t * factor;
index bf4ef812a82108c689d794639db35f06b1a2797e..c700deaf06b82345b5e1ec3964bd58d3d7202964 100644 (file)
@@ -112,34 +112,27 @@ main()
   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;
+    constexpr double a_init[3][3] = {{1., 0., 0.}, {2., 1., 0.}, {3., 2., 1.}};
+    constexpr Tensor<2, 3>  dummy_a{a_init};
+    DEAL_II_CONSTEXPR const SymmetricTensor<2, 3> a = symmetrize(dummy_a);
+    DEAL_II_CONSTEXPR const auto                  inverted = invert(a);
+    constexpr double        ref_init[3][3]                 = {{0., -2., 2.},
+                                       {-2., 5., -2.},
+                                       {2., -2., 0.}};
+    constexpr Tensor<2, 3>  dummy_ref{ref_init};
+    DEAL_II_CONSTEXPR const SymmetricTensor<2, 3> ref = symmetrize(dummy_ref);
     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);
+    constexpr double a_init[3][3] = {{1., 2., 3.}, {2., 1., 2.}, {3., 2., 1.}};
+    constexpr Tensor<2, 3>  dummy_a{a_init};
+    DEAL_II_CONSTEXPR const SymmetricTensor<2, 3> a = symmetrize(dummy_a);
+    DEAL_II_CONSTEXPR const 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);
+    DEAL_II_CONSTEXPR const auto dummy   = scalar_product(a, a);
+    DEAL_II_CONSTEXPR const auto dummy_5 = a * a;
+    DEAL_II_CONSTEXPR const auto middle  = outer_product(a, a);
+    DEAL_II_CONSTEXPR const auto dummy_6 = contract3(a, middle, a);
   }
 
   deallog << "OK" << std::endl;
index a5cc8696f6bd74df42385e2227ba692228771a3d..35c96f433c6fc9bbe92f83711487e49f2d2a3e51 100644 (file)
@@ -31,14 +31,16 @@ test_constexpr_tensor_constructors()
   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;
+  DEAL_II_CONSTEXPR const auto d = a * 2.;
+  DEAL_II_CONSTEXPR const auto e = 2. * a;
+  DEAL_II_CONSTEXPR const auto f = a / 2.;
+  DEAL_II_CONSTEXPR const auto g = d + e;
+  DEAL_II_CONSTEXPR const auto h = d - e;
   deallog << d << std::endl;
   deallog << e << std::endl;
   deallog << f << std::endl;
   deallog << g << std::endl;
+  deallog << h << std::endl;
 }
 
 DEAL_II_CONSTEXPR double
@@ -117,72 +119,72 @@ main()
   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);
+    constexpr double        initializer[2] = {1., -1.};
+    constexpr Tensor<1, 2>  c{initializer};
+    DEAL_II_CONSTEXPR const Tensor<1, 2> c_cross = cross_product_2d(c);
+    constexpr double                     initializer_ref[2] = {-1., -1.};
+    constexpr Tensor<1, 2>               ref{initializer_ref};
+    DEAL_II_CONSTEXPR const bool         is_same     = (c_cross == ref);
+    DEAL_II_CONSTEXPR const 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);
+    DEAL_II_CONSTEXPR const double initializer_1[3] = {1., 0., 0.};
+    DEAL_II_CONSTEXPR const Tensor<1, 3> c_1{initializer_1};
+    DEAL_II_CONSTEXPR const double       initializer_2[3] = {0., 1., 0.};
+    DEAL_II_CONSTEXPR const Tensor<1, 3> c_2{initializer_2};
+
+    DEAL_II_CONSTEXPR const auto   c_cross = cross_product_3d(c_1, c_2);
+    DEAL_II_CONSTEXPR const double initializer_ref[3] = {0., 0., 1.};
+    DEAL_II_CONSTEXPR const Tensor<1, 3> ref{initializer_ref};
+    DEAL_II_CONSTEXPR const bool         is_same     = (c_cross == ref);
+    DEAL_II_CONSTEXPR const bool         is_not_same = (c_cross != ref);
     Assert(is_same && !is_not_same, ExcInternalError());
   }
 
-  DEAL_II_CONSTEXPR auto table_indices =
+  DEAL_II_CONSTEXPR const auto table_indices =
     Tensor<2, 3>::unrolled_to_component_indices(0);
-  DEAL_II_CONSTEXPR auto index =
+  DEAL_II_CONSTEXPR const 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();
+  DEAL_II_CONSTEXPR const 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};
+    DEAL_II_CONSTEXPR const double a_init[3][3] = {{1., 0., 0.},
+                                                   {2., 1., 0.},
+                                                   {3., 2., 1.}};
+    DEAL_II_CONSTEXPR const Tensor<2, 3> a{a_init};
+    DEAL_II_CONSTEXPR const auto         inverted       = invert(a);
+    DEAL_II_CONSTEXPR const double       ref_init[3][3] = {{1., 0., 0.},
+                                                     {-2., 1., 0.},
+                                                     {1., -2., 1}};
+    DEAL_II_CONSTEXPR const 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};
+    DEAL_II_CONSTEXPR const double a_init[3][3] = {{1., 0., 0.},
+                                                   {2., 1., 0.},
+                                                   {3., 2., 1.}};
+    DEAL_II_CONSTEXPR const Tensor<2, 3> a{a_init};
+    DEAL_II_CONSTEXPR const auto         transposed     = transpose(a);
+    DEAL_II_CONSTEXPR const double       ref_init[3][3] = {{1., 2., 3.},
+                                                     {0., 1., 2.},
+                                                     {0., 0., 1}};
+    DEAL_II_CONSTEXPR const 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);
+    DEAL_II_CONSTEXPR const auto dummy    = scalar_product(a, ref);
+    DEAL_II_CONSTEXPR const auto dummy_2  = contract<0, 0>(a, ref);
+    DEAL_II_CONSTEXPR const auto dummy_3  = double_contract<0, 0, 1, 1>(a, ref);
+    DEAL_II_CONSTEXPR const auto dummy_4  = schur_product(a, ref);
+    DEAL_II_CONSTEXPR const auto dummy_5  = a * ref;
+    DEAL_II_CONSTEXPR const auto middle   = outer_product(a, a);
+    DEAL_II_CONSTEXPR const auto dummy_6  = contract3(a, middle, a);
+    DEAL_II_CONSTEXPR const auto dummy_7  = adjugate(a);
+    DEAL_II_CONSTEXPR const auto dummy_8  = cofactor(a);
+    DEAL_II_CONSTEXPR const auto dummy_9  = l1_norm(a);
+    DEAL_II_CONSTEXPR const auto dummy_10 = linfty_norm(a);
   }
 
   deallog << "OK" << std::endl;
index 218e060df40a10f4370bbe8aae793ff886708073..5b07862e373ac4b6df30f9973496dbcdc1f8485d 100644 (file)
@@ -8,6 +8,7 @@ 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
@@ -16,6 +17,7 @@ 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
@@ -24,6 +26,7 @@ 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
@@ -32,6 +35,7 @@ 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
@@ -40,6 +44,7 @@ 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
@@ -48,6 +53,7 @@ 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
@@ -56,6 +62,7 @@ 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
@@ -64,6 +71,7 @@ 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
@@ -72,6 +80,7 @@ DEAL:float::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00
 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
@@ -80,6 +89,7 @@ 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
@@ -88,6 +98,7 @@ 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
@@ -96,6 +107,7 @@ 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
@@ -104,6 +116,7 @@ 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
@@ -112,6 +125,7 @@ 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
@@ -120,6 +134,7 @@ 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
@@ -128,6 +143,7 @@ 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
@@ -136,6 +152,7 @@ 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
@@ -144,6 +161,7 @@ DEAL:double::0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.0
 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::116.000

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.