]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use DEAL_II_IF_CONSTEXPR as attribute
authorDaniel Arndt <arndtd@ornl.gov>
Wed, 25 Jan 2023 21:04:52 +0000 (16:04 -0500)
committerDaniel Arndt <arndtd@ornl.gov>
Wed, 25 Jan 2023 21:28:53 +0000 (16:28 -0500)
.clang-format
include/deal.II/base/config.h.in
include/deal.II/base/partitioner.templates.h
include/deal.II/base/utilities.h
include/deal.II/lac/sparse_matrix.templates.h
include/deal.II/lac/vector_operations_internal.h
source/base/polynomials_rt_bubbles.cc

index 9b627cecdff422b8fd0f85a7d2c048e76e7dc497..2975f20f50074d4e06b633b385b715d0841d9ece 100644 (file)
@@ -65,9 +65,6 @@ DerivePointerAlignment: false
 
 FixNamespaceComments: true
 
-# Should be IfMacros in clang-format 13 and higher
-ForEachMacros: ['DEAL_II_IF_CONSTEXPR']
-
 IncludeBlocks: Regroup
 IncludeCategories:
 # config.h must always be first:
index dc29536bcc46e9375d443de6c3c620203f553d08..8cc111bedd5edb0b24a1808d0b16c73331a69b08 100644 (file)
@@ -548,13 +548,14 @@ _Pragma("GCC diagnostic pop")
 /***********************************************************************
  * Some places in deal.II require or benefit from using "if constexpr".
  * However, this feature is only available in C++17 and above.
- * The macro defined here expands to "if constexpr(...)" when C++17 is supported
- * and to "if (...)" otherwise.
+ * The macro defined here can be used so that
+ * "if DEAL_II_CONSTEXPR_IN_CONDITIONAL (...)" expands to "if constexpr(...)"
+ * when C++17 is supported and to "if (...)" otherwise.
  */
 #ifdef DEAL_II_HAVE_CXX17
-  #define DEAL_II_IF_CONSTEXPR(x) if constexpr x
+  #define DEAL_II_CONSTEXPR_IN_CONDITIONAL constexpr
 #else
-  #define DEAL_II_IF_CONSTEXPR(x) if x
+  #define DEAL_II_CONSTEXPR_IN_CONDITIONAL
 #endif
 
 
index 7cbb1f02d61ddeb407acce9dcdfd67313ddc9a33..6381d6481c60ab95e61ab2748d458ae4394cc7d0 100644 (file)
@@ -533,18 +533,15 @@ namespace Utilities
           else
 #      endif
             {
-              DEAL_II_IF_CONSTEXPR ((std::is_trivial<Number>::value))
-                {
-                  std::memset(ghost_array.data(),
-                              0,
-                              sizeof(Number) * ghost_array.size());
-                }
+              if DEAL_II_CONSTEXPR_IN_CONDITIONAL (std::is_trivial<
+                                                     Number>::value)
+                std::memset(ghost_array.data(),
+                            0,
+                            sizeof(Number) * ghost_array.size());
               else
-              {
                 std::fill(ghost_array.data(),
                           ghost_array.data() + ghost_array.size(),
                           0);
-              }
             }
           return;
         }
@@ -736,18 +733,17 @@ namespace Utilities
           else
 #    endif
             {
-              DEAL_II_IF_CONSTEXPR ((std::is_trivial<Number>::value))
+              if DEAL_II_CONSTEXPR_IN_CONDITIONAL (std::is_trivial<
+                                                     Number>::value)
                 {
                   std::memset(ghost_array.data(),
                               0,
                               sizeof(Number) * n_ghost_indices());
                 }
               else
-              {
                 std::fill(ghost_array.data(),
                           ghost_array.data() + n_ghost_indices(),
                           0);
-              }
             }
         }
 
index 9c8e2de96911530804ad9ee9f45ecf36e18dfa6f..ef8f22e2cbf3043ac7a85c863293206e10afb3f1 100644 (file)
@@ -1341,7 +1341,8 @@ namespace Utilities
     // see if the object is small and copyable via memcpy. if so, use
     // this fast path. otherwise, we have to go through the BOOST
     // serialization machinery
-    DEAL_II_IF_CONSTEXPR ((std::is_trivially_copyable<T>() && sizeof(T) < 256))
+    if DEAL_II_CONSTEXPR_IN_CONDITIONAL (std::is_trivially_copyable<T>() &&
+                                         sizeof(T) < 256)
       {
         // Determine the size. There are places where we would like to use a
         // truly empty type, for which we use std::tuple<> (i.e., a tuple
@@ -1365,40 +1366,40 @@ namespace Utilities
     // the data.
     else if (internal::IsVectorOfTriviallyCopyable<T>::value &&
              (allow_compression == false))
-    {
-      const std::size_t previous_size = dest_buffer.size();
+      {
+        const std::size_t previous_size = dest_buffer.size();
 
-      // When we have DEAL_II_HAVE_CXX17 set by default, we can just
-      // inline the code of the following function here and make the 'if'
-      // above a 'if constexpr'. Without the 'constexpr', we need to keep
-      // the general template of the function that throws an exception.
-      internal::append_vector_of_trivially_copyable_to_buffer(object,
-                                                              dest_buffer);
+        // When we have DEAL_II_HAVE_CXX17 set by default, we can just
+        // inline the code of the following function here and make the 'if'
+        // above a 'if constexpr'. Without the 'constexpr', we need to keep
+        // the general template of the function that throws an exception.
+        internal::append_vector_of_trivially_copyable_to_buffer(object,
+                                                                dest_buffer);
 
-      size = dest_buffer.size() - previous_size;
-    }
+        size = dest_buffer.size() - previous_size;
+      }
     else
-    {
-      // use buffer as the target of a compressing
-      // stream into which we serialize the current object
-      const std::size_t previous_size = dest_buffer.size();
       {
-        boost::iostreams::filtering_ostreambuf fosb;
+        // use buffer as the target of a compressing
+        // stream into which we serialize the current object
+        const std::size_t previous_size = dest_buffer.size();
+        {
+          boost::iostreams::filtering_ostreambuf fosb;
 #ifdef DEAL_II_WITH_ZLIB
-        if (allow_compression)
-          fosb.push(boost::iostreams::gzip_compressor());
+          if (allow_compression)
+            fosb.push(boost::iostreams::gzip_compressor());
 #else
-        (void)allow_compression;
+          (void)allow_compression;
 #endif
-        fosb.push(boost::iostreams::back_inserter(dest_buffer));
+          fosb.push(boost::iostreams::back_inserter(dest_buffer));
 
-        boost::archive::binary_oarchive boa(fosb);
-        boa << object;
-        // the stream object has to be destroyed before the return statement
-        // to ensure that all data has been written in the buffer
+          boost::archive::binary_oarchive boa(fosb);
+          boa << object;
+          // the stream object has to be destroyed before the return statement
+          // to ensure that all data has been written in the buffer
+        }
+        size = dest_buffer.size() - previous_size;
       }
-      size = dest_buffer.size() - previous_size;
-    }
 
     return size;
   }
@@ -1424,7 +1425,8 @@ namespace Utilities
     // see if the object is small and copyable via memcpy. if so, use
     // this fast path. otherwise, we have to go through the BOOST
     // serialization machinery
-    DEAL_II_IF_CONSTEXPR ((std::is_trivially_copyable<T>() && sizeof(T) < 256))
+    if DEAL_II_CONSTEXPR_IN_CONDITIONAL (std::is_trivially_copyable<T>() &&
+                                         sizeof(T) < 256)
       {
         // Determine the size. There are places where we would like to use a
         // truly empty type, for which we use std::tuple<> (i.e., a tuple
@@ -1452,35 +1454,35 @@ namespace Utilities
     // are not asked to compress the data.
     else if (internal::IsVectorOfTriviallyCopyable<T>::value &&
              (allow_compression == false))
-    {
-      // When we have DEAL_II_HAVE_CXX17 set by default, we can just
-      // inline the code of the following function here and make the 'if'
-      // above a 'if constexpr'. Without the 'constexpr', we need to keep
-      // the general template of the function that throws an exception.
-      T object;
-      internal::create_vector_of_trivially_copyable_from_buffer(cbegin,
-                                                                cend,
-                                                                object);
-      return object;
-    }
+      {
+        // When we have DEAL_II_HAVE_CXX17 set by default, we can just
+        // inline the code of the following function here and make the 'if'
+        // above a 'if constexpr'. Without the 'constexpr', we need to keep
+        // the general template of the function that throws an exception.
+        T object;
+        internal::create_vector_of_trivially_copyable_from_buffer(cbegin,
+                                                                  cend,
+                                                                  object);
+        return object;
+      }
     else
-    {
-      // decompress the buffer section into the object
-      boost::iostreams::filtering_istreambuf fisb;
+      {
+        // decompress the buffer section into the object
+        boost::iostreams::filtering_istreambuf fisb;
 #ifdef DEAL_II_WITH_ZLIB
-      if (allow_compression)
-        fisb.push(boost::iostreams::gzip_decompressor());
+        if (allow_compression)
+          fisb.push(boost::iostreams::gzip_decompressor());
 #else
-      (void)allow_compression;
+        (void)allow_compression;
 #endif
-      fisb.push(boost::iostreams::array_source(&*cbegin, cend - cbegin));
+        fisb.push(boost::iostreams::array_source(&*cbegin, cend - cbegin));
 
-      boost::archive::binary_iarchive bia(fisb);
+        boost::archive::binary_iarchive bia(fisb);
 
-      T object;
-      bia >> object;
-      return object;
-    }
+        T object;
+        bia >> object;
+        return object;
+      }
 
     return T();
   }
@@ -1504,28 +1506,28 @@ namespace Utilities
     // see if the object is small and copyable via memcpy. if so, use
     // this fast path. otherwise, we have to go through the BOOST
     // serialization machinery
-    DEAL_II_IF_CONSTEXPR (
-      (std::is_trivially_copyable<T>() && sizeof(T) * N < 256))
+    if DEAL_II_CONSTEXPR_IN_CONDITIONAL (std::is_trivially_copyable<T>() &&
+                                         sizeof(T) * N < 256)
       {
         Assert(std::distance(cbegin, cend) == sizeof(T) * N,
                ExcInternalError());
         std::memcpy(unpacked_object, &*cbegin, sizeof(T) * N);
       }
     else
-    {
-      // decompress the buffer section into the object
-      boost::iostreams::filtering_istreambuf fisb;
+      {
+        // decompress the buffer section into the object
+        boost::iostreams::filtering_istreambuf fisb;
 #ifdef DEAL_II_WITH_ZLIB
-      if (allow_compression)
-        fisb.push(boost::iostreams::gzip_decompressor());
+        if (allow_compression)
+          fisb.push(boost::iostreams::gzip_decompressor());
 #else
-      (void)allow_compression;
+        (void)allow_compression;
 #endif
-      fisb.push(boost::iostreams::array_source(&*cbegin, cend - cbegin));
+        fisb.push(boost::iostreams::array_source(&*cbegin, cend - cbegin));
 
-      boost::archive::binary_iarchive bia(fisb);
-      bia >> unpacked_object;
-    }
+        boost::archive::binary_iarchive bia(fisb);
+        bia >> unpacked_object;
+      }
   }
 
 
index 06986c0b19105938fce124fc87e1f25769de76ff..bfc5c39534b7ec9558f431ff5e8de371dbfe8a7b 100644 (file)
@@ -219,14 +219,10 @@ SparseMatrix<number>::operator=(const double d)
       grain_size);
   else if (matrix_size > 0)
     {
-      DEAL_II_IF_CONSTEXPR ((std::is_trivial<number>::value))
-        {
-          std::memset(val.get(), 0, matrix_size * sizeof(number));
-        }
+      if DEAL_II_CONSTEXPR_IN_CONDITIONAL (std::is_trivial<number>::value)
+        std::memset(val.get(), 0, matrix_size * sizeof(number));
       else
-      {
         std::fill(val.get(), val.get() + matrix_size, 0);
-      }
     }
 
   return *this;
index 3ffb97788153ada19fa9afd4f7c22747a7034b83..f2a67f618c36fbbbb68ac7d7e1b934fb92d70c63 100644 (file)
@@ -214,7 +214,7 @@ namespace internal
 
         if (value == Number())
           {
-            DEAL_II_IF_CONSTEXPR ((std::is_trivial<Number>::value))
+            if DEAL_II_CONSTEXPR_IN_CONDITIONAL (std::is_trivial<Number>::value)
               {
                 std::memset(dst + begin, 0, sizeof(Number) * (end - begin));
                 return;
@@ -243,15 +243,17 @@ namespace internal
       {
         Assert(end >= begin, ExcInternalError());
 
-        DEAL_II_IF_CONSTEXPR ((std::is_trivially_copyable<Number>() &&
-                               std::is_same<Number, OtherNumber>::value))
+        if DEAL_II_CONSTEXPR_IN_CONDITIONAL (std::is_trivially_copyable<
+                                               Number>() &&
+                                             std::is_same<Number,
+                                                          OtherNumber>::value)
           std::memcpy(dst + begin, src + begin, (end - begin) * sizeof(Number));
         else
-        {
-          DEAL_II_OPENMP_SIMD_PRAGMA
-          for (size_type i = begin; i < end; ++i)
-            dst[i] = src[i];
-        }
+          {
+            DEAL_II_OPENMP_SIMD_PRAGMA
+            for (size_type i = begin; i < end; ++i)
+              dst[i] = src[i];
+          }
       }
 
       const OtherNumber *const src;
index 95270e6ea1728f077bf40a7ca2311fbd71885190..6b0b188aea6f334273ddd1ece2ff37b35f510b74 100644 (file)
@@ -110,14 +110,14 @@ PolynomialsRT_Bubbles<dim>::evaluate(
   double             monoval_i[dim][n_derivatives + 1];
 
 
-  DEAL_II_IF_CONSTEXPR ((dim <= 1))
+  if DEAL_II_CONSTEXPR_IN_CONDITIONAL (dim <= 1)
     {
       (void)monoval_plus;
       (void)monoval_i;
     }
 
   unsigned int start = n_sub;
-  DEAL_II_IF_CONSTEXPR ((dim == 2))
+  if DEAL_II_CONSTEXPR_IN_CONDITIONAL (dim == 2)
     {
       // In 2d the curl part of the space is spanned by the vectors
       // of two types. The first one is
@@ -200,639 +200,641 @@ PolynomialsRT_Bubbles<dim>::evaluate(
         }
       Assert(start == this->n() - my_degree - 1, ExcInternalError());
     }
-  else DEAL_II_IF_CONSTEXPR ((dim == 3))
-  {
-    double monoval[dim][n_derivatives + 1];
-    double monoval_j[dim][n_derivatives + 1];
-    double monoval_jplus[dim][n_derivatives + 1];
-
-    // In 3d the first type of basis vector is
-    //  [ x^i * y^j * z^k * (j+k+2) ]
-    //  [  -[x^i]' * y^(j+1) * z^k  ]
-    //  [  -[x^i]' * y^j * z^(k+1)  ],
-    // For the second type of basis vector y and z
-    // are swapped. Then for each of these,
-    // two more are obtained by the cyclic rotation
-    // of the coordinates.
-    //  monoval = x^k,   monoval_plus = x^(k+1)
-    //  monoval_* = x^*, monoval_jplus = x^(j+1)
-    for (unsigned int d = 0; d < dim; ++d)
-      {
-        monomials[my_degree + 1].value(unit_point(d),
-                                       n_derivatives,
-                                       monoval_plus[d]);
-        monomials[my_degree].value(unit_point(d), n_derivatives, monoval[d]);
-      }
-
-    const unsigned int n_curls = (my_degree + 1) * (2 * my_degree + 1);
-    // Span of $\tilde{B}$
-    for (unsigned int i = 0; i <= my_degree; ++i)
-      {
-        for (unsigned int d = 0; d < dim; ++d)
-          monomials[i].value(unit_point(d), n_derivatives, monoval_i[d]);
-
-        for (unsigned int j = 0; j <= my_degree; ++j)
-          {
-            for (unsigned int d = 0; d < dim; ++d)
-              {
-                monomials[j].value(unit_point(d), n_derivatives, monoval_j[d]);
-                monomials[j + 1].value(unit_point(d),
-                                       n_derivatives,
-                                       monoval_jplus[d]);
-              }
-
-            if (values.size() != 0)
-              {
-                values[start][0] = monoval_i[0][0] * monoval_j[1][0] *
-                                   monoval[2][0] *
-                                   static_cast<double>(j + my_degree + 2);
-                values[start][1] =
-                  -monoval_i[0][1] * monoval_jplus[1][0] * monoval[2][0];
-                values[start][2] =
-                  -monoval_i[0][1] * monoval_j[1][0] * monoval_plus[2][0];
-
-                values[start + n_curls][0] =
-                  -monoval_jplus[0][0] * monoval_i[1][1] * monoval[2][0];
-                values[start + n_curls][1] =
-                  monoval_j[0][0] * monoval_i[1][0] * monoval[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                values[start + n_curls][2] =
-                  -monoval_j[0][0] * monoval_i[1][1] * monoval_plus[2][0];
-
-                values[start + 2 * n_curls][0] =
-                  -monoval_jplus[0][0] * monoval[1][0] * monoval_i[2][1];
-                values[start + 2 * n_curls][1] =
-                  -monoval_j[0][0] * monoval_plus[1][0] * monoval_i[2][1];
-                values[start + 2 * n_curls][2] =
-                  monoval_j[0][0] * monoval[1][0] * monoval_i[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-
-                // Only unique triples of powers (i j k)
-                // and (i k j) are allowed, 0 <= i,j <= k
-                if (j != my_degree)
-                  {
-                    values[start + 1][0] =
-                      monoval_i[0][0] * monoval[1][0] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    values[start + 1][1] =
-                      -monoval_i[0][1] * monoval_plus[1][0] * monoval_j[2][0];
-                    values[start + 1][2] =
-                      -monoval_i[0][1] * monoval[1][0] * monoval_jplus[2][0];
-
-                    values[start + n_curls + 1][0] =
-                      -monoval_plus[0][0] * monoval_i[1][1] * monoval_j[2][0];
-                    values[start + n_curls + 1][1] =
-                      monoval[0][0] * monoval_i[1][0] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    values[start + n_curls + 1][2] =
-                      -monoval[0][0] * monoval_i[1][1] * monoval_jplus[2][0];
-
-                    values[start + 2 * n_curls + 1][0] =
-                      -monoval_plus[0][0] * monoval_j[1][0] * monoval_i[2][1];
-                    values[start + 2 * n_curls + 1][1] =
-                      -monoval[0][0] * monoval_jplus[1][0] * monoval_i[2][1];
-                    values[start + 2 * n_curls + 1][2] =
-                      monoval[0][0] * monoval_j[1][0] * monoval_i[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                  }
-              }
-
-            if (grads.size() != 0)
-              {
-                grads[start][0][0] = monoval_i[0][1] * monoval_j[1][0] *
-                                     monoval[2][0] *
-                                     static_cast<double>(j + my_degree + 2);
-                grads[start][0][1] = monoval_i[0][0] * monoval_j[1][1] *
+  else if DEAL_II_CONSTEXPR_IN_CONDITIONAL (dim == 3)
+    {
+      double monoval[dim][n_derivatives + 1];
+      double monoval_j[dim][n_derivatives + 1];
+      double monoval_jplus[dim][n_derivatives + 1];
+
+      // In 3d the first type of basis vector is
+      //  [ x^i * y^j * z^k * (j+k+2) ]
+      //  [  -[x^i]' * y^(j+1) * z^k  ]
+      //  [  -[x^i]' * y^j * z^(k+1)  ],
+      // For the second type of basis vector y and z
+      // are swapped. Then for each of these,
+      // two more are obtained by the cyclic rotation
+      // of the coordinates.
+      //  monoval = x^k,   monoval_plus = x^(k+1)
+      //  monoval_* = x^*, monoval_jplus = x^(j+1)
+      for (unsigned int d = 0; d < dim; ++d)
+        {
+          monomials[my_degree + 1].value(unit_point(d),
+                                         n_derivatives,
+                                         monoval_plus[d]);
+          monomials[my_degree].value(unit_point(d), n_derivatives, monoval[d]);
+        }
+
+      const unsigned int n_curls = (my_degree + 1) * (2 * my_degree + 1);
+      // Span of $\tilde{B}$
+      for (unsigned int i = 0; i <= my_degree; ++i)
+        {
+          for (unsigned int d = 0; d < dim; ++d)
+            monomials[i].value(unit_point(d), n_derivatives, monoval_i[d]);
+
+          for (unsigned int j = 0; j <= my_degree; ++j)
+            {
+              for (unsigned int d = 0; d < dim; ++d)
+                {
+                  monomials[j].value(unit_point(d),
+                                     n_derivatives,
+                                     monoval_j[d]);
+                  monomials[j + 1].value(unit_point(d),
+                                         n_derivatives,
+                                         monoval_jplus[d]);
+                }
+
+              if (values.size() != 0)
+                {
+                  values[start][0] = monoval_i[0][0] * monoval_j[1][0] *
                                      monoval[2][0] *
                                      static_cast<double>(j + my_degree + 2);
-                grads[start][0][2] = monoval_i[0][0] * monoval_j[1][0] *
-                                     monoval[2][1] *
-                                     static_cast<double>(j + my_degree + 2);
-                grads[start][1][0] =
-                  -monoval_i[0][2] * monoval_jplus[1][0] * monoval[2][0];
-                grads[start][1][1] =
-                  -monoval_i[0][1] * monoval_jplus[1][1] * monoval[2][0];
-                grads[start][1][2] =
-                  -monoval_i[0][1] * monoval_jplus[1][0] * monoval[2][1];
-                grads[start][2][0] =
-                  -monoval_i[0][2] * monoval_j[1][0] * monoval_plus[2][0];
-                grads[start][2][1] =
-                  -monoval_i[0][1] * monoval_j[1][1] * monoval_plus[2][0];
-                grads[start][2][2] =
-                  -monoval_i[0][1] * monoval_j[1][0] * monoval_plus[2][1];
-
-                grads[start + n_curls][0][0] =
-                  -monoval_jplus[0][1] * monoval_i[1][1] * monoval[2][0];
-                grads[start + n_curls][0][1] =
-                  -monoval_jplus[0][0] * monoval_i[1][2] * monoval[2][0];
-                grads[start + n_curls][0][2] =
-                  -monoval_jplus[0][0] * monoval_i[1][1] * monoval[2][1];
-                grads[start + n_curls][1][0] =
-                  monoval_j[0][1] * monoval_i[1][0] * monoval[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grads[start + n_curls][1][1] =
-                  monoval_j[0][0] * monoval_i[1][1] * monoval[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grads[start + n_curls][1][2] =
-                  monoval_j[0][0] * monoval_i[1][0] * monoval[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grads[start + n_curls][2][0] =
-                  -monoval_j[0][1] * monoval_i[1][1] * monoval_plus[2][0];
-                grads[start + n_curls][2][1] =
-                  -monoval_j[0][0] * monoval_i[1][2] * monoval_plus[2][0];
-                grads[start + n_curls][2][2] =
-                  -monoval_j[0][0] * monoval_i[1][1] * monoval_plus[2][1];
-
-                grads[start + 2 * n_curls][0][0] =
-                  -monoval_jplus[0][1] * monoval[1][0] * monoval_i[2][1];
-                grads[start + 2 * n_curls][0][1] =
-                  -monoval_jplus[0][0] * monoval[1][1] * monoval_i[2][1];
-                grads[start + 2 * n_curls][0][2] =
-                  -monoval_jplus[0][0] * monoval[1][0] * monoval_i[2][2];
-                grads[start + 2 * n_curls][1][0] =
-                  -monoval_j[0][1] * monoval_plus[1][0] * monoval_i[2][1];
-                grads[start + 2 * n_curls][1][1] =
-                  -monoval_j[0][0] * monoval_plus[1][1] * monoval_i[2][1];
-                grads[start + 2 * n_curls][1][2] =
-                  -monoval_j[0][0] * monoval_plus[1][0] * monoval_i[2][2];
-                grads[start + 2 * n_curls][2][0] =
-                  monoval_j[0][1] * monoval[1][0] * monoval_i[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grads[start + 2 * n_curls][2][1] =
-                  monoval_j[0][0] * monoval[1][1] * monoval_i[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grads[start + 2 * n_curls][2][2] =
-                  monoval_j[0][0] * monoval[1][0] * monoval_i[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-
-                if (j != my_degree)
-                  {
-                    grads[start + 1][0][0] =
-                      monoval_i[0][1] * monoval[1][0] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grads[start + 1][0][1] =
-                      monoval_i[0][0] * monoval[1][1] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grads[start + 1][0][2] =
-                      monoval_i[0][0] * monoval[1][0] * monoval_j[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grads[start + 1][1][0] =
-                      -monoval_i[0][2] * monoval_plus[1][0] * monoval_j[2][0];
-                    grads[start + 1][1][1] =
-                      -monoval_i[0][1] * monoval_plus[1][1] * monoval_j[2][0];
-                    grads[start + 1][1][2] =
-                      -monoval_i[0][1] * monoval_plus[1][0] * monoval_j[2][1];
-                    grads[start + 1][2][0] =
-                      -monoval_i[0][2] * monoval[1][0] * monoval_jplus[2][0];
-                    grads[start + 1][2][1] =
-                      -monoval_i[0][1] * monoval[1][1] * monoval_jplus[2][0];
-                    grads[start + 1][2][2] =
-                      -monoval_i[0][1] * monoval[1][0] * monoval_jplus[2][1];
-
-                    grads[start + n_curls + 1][0][0] =
-                      -monoval_plus[0][1] * monoval_i[1][1] * monoval_j[2][0];
-                    grads[start + n_curls + 1][0][1] =
-                      -monoval_plus[0][0] * monoval_i[1][2] * monoval_j[2][0];
-                    grads[start + n_curls + 1][0][2] =
-                      -monoval_plus[0][0] * monoval_i[1][1] * monoval_j[2][1];
-                    grads[start + n_curls + 1][1][0] =
-                      monoval[0][1] * monoval_i[1][0] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grads[start + n_curls + 1][1][1] =
-                      monoval[0][0] * monoval_i[1][1] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grads[start + n_curls + 1][1][2] =
-                      monoval[0][0] * monoval_i[1][0] * monoval_j[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grads[start + n_curls + 1][2][0] =
-                      -monoval[0][1] * monoval_i[1][1] * monoval_jplus[2][0];
-                    grads[start + n_curls + 1][2][1] =
-                      -monoval[0][0] * monoval_i[1][2] * monoval_jplus[2][0];
-                    grads[start + n_curls + 1][2][2] =
-                      -monoval[0][0] * monoval_i[1][1] * monoval_jplus[2][1];
-
-                    grads[start + 2 * n_curls + 1][0][0] =
-                      -monoval_plus[0][1] * monoval_j[1][0] * monoval_i[2][1];
-                    grads[start + 2 * n_curls + 1][0][1] =
-                      -monoval_plus[0][0] * monoval_j[1][1] * monoval_i[2][1];
-                    grads[start + 2 * n_curls + 1][0][2] =
-                      -monoval_plus[0][0] * monoval_j[1][0] * monoval_i[2][2];
-                    grads[start + 2 * n_curls + 1][1][0] =
-                      -monoval[0][1] * monoval_jplus[1][0] * monoval_i[2][1];
-                    grads[start + 2 * n_curls + 1][1][1] =
-                      -monoval[0][0] * monoval_jplus[1][1] * monoval_i[2][1];
-                    grads[start + 2 * n_curls + 1][1][2] =
-                      -monoval[0][0] * monoval_jplus[1][0] * monoval_i[2][2];
-                    grads[start + 2 * n_curls + 1][2][0] =
-                      monoval[0][1] * monoval_j[1][0] * monoval_i[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grads[start + 2 * n_curls + 1][2][1] =
-                      monoval[0][0] * monoval_j[1][1] * monoval_i[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grads[start + 2 * n_curls + 1][2][2] =
-                      monoval[0][0] * monoval_j[1][0] * monoval_i[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                  }
-              }
-
-            if (grad_grads.size() != 0)
-              {
-                grad_grads[start][0][0][0] =
-                  monoval_i[0][2] * monoval_j[1][0] * monoval[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start][0][0][1] =
-                  monoval_i[0][1] * monoval_j[1][1] * monoval[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start][0][0][2] =
-                  monoval_i[0][1] * monoval_j[1][0] * monoval[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start][0][1][0] =
-                  monoval_i[0][1] * monoval_j[1][1] * monoval[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start][0][1][1] =
-                  monoval_i[0][0] * monoval_j[1][2] * monoval[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start][0][1][2] =
-                  monoval_i[0][0] * monoval_j[1][1] * monoval[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start][0][2][0] =
-                  monoval_i[0][1] * monoval_j[1][0] * monoval[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start][0][2][1] =
-                  monoval_i[0][0] * monoval_j[1][1] * monoval[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start][0][2][2] =
-                  monoval_i[0][0] * monoval_j[1][0] * monoval[2][2] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start][1][0][0] =
-                  -monoval_i[0][3] * monoval_jplus[1][0] * monoval[2][0];
-                grad_grads[start][1][0][1] =
-                  -monoval_i[0][2] * monoval_jplus[1][1] * monoval[2][0];
-                grad_grads[start][1][0][2] =
-                  -monoval_i[0][2] * monoval_jplus[1][0] * monoval[2][1];
-                grad_grads[start][1][1][0] =
-                  -monoval_i[0][2] * monoval_jplus[1][1] * monoval[2][0];
-                grad_grads[start][1][1][1] =
-                  -monoval_i[0][1] * monoval_jplus[1][2] * monoval[2][0];
-                grad_grads[start][1][1][2] =
-                  -monoval_i[0][1] * monoval_jplus[1][1] * monoval[2][1];
-                grad_grads[start][1][2][0] =
-                  -monoval_i[0][2] * monoval_jplus[1][0] * monoval[2][1];
-                grad_grads[start][1][2][1] =
-                  -monoval_i[0][1] * monoval_jplus[1][1] * monoval[2][1];
-                grad_grads[start][1][2][2] =
-                  -monoval_i[0][1] * monoval_jplus[1][0] * monoval[2][2];
-                grad_grads[start][2][0][0] =
-                  -monoval_i[0][3] * monoval_j[1][0] * monoval_plus[2][0];
-                grad_grads[start][2][0][1] =
-                  -monoval_i[0][2] * monoval_j[1][1] * monoval_plus[2][0];
-                grad_grads[start][2][0][2] =
-                  -monoval_i[0][2] * monoval_j[1][0] * monoval_plus[2][1];
-                grad_grads[start][2][1][0] =
-                  -monoval_i[0][2] * monoval_j[1][1] * monoval_plus[2][0];
-                grad_grads[start][2][1][1] =
-                  -monoval_i[0][1] * monoval_j[1][2] * monoval_plus[2][0];
-                grad_grads[start][2][1][2] =
-                  -monoval_i[0][1] * monoval_j[1][1] * monoval_plus[2][1];
-                grad_grads[start][2][2][0] =
-                  -monoval_i[0][2] * monoval_j[1][0] * monoval_plus[2][1];
-                grad_grads[start][2][2][1] =
-                  -monoval_i[0][1] * monoval_j[1][1] * monoval_plus[2][1];
-                grad_grads[start][2][2][2] =
-                  -monoval_i[0][1] * monoval_j[1][0] * monoval_plus[2][2];
-
-                grad_grads[start + n_curls][0][0][0] =
-                  -monoval_jplus[0][2] * monoval_i[1][1] * monoval[2][0];
-                grad_grads[start + n_curls][0][0][1] =
-                  -monoval_jplus[0][1] * monoval_i[1][2] * monoval[2][0];
-                grad_grads[start + n_curls][0][0][2] =
-                  -monoval_jplus[0][1] * monoval_i[1][1] * monoval[2][1];
-                grad_grads[start + n_curls][0][1][0] =
-                  -monoval_jplus[0][1] * monoval_i[1][2] * monoval[2][0];
-                grad_grads[start + n_curls][0][1][1] =
-                  -monoval_jplus[0][0] * monoval_i[1][3] * monoval[2][0];
-                grad_grads[start + n_curls][0][1][2] =
-                  -monoval_jplus[0][0] * monoval_i[1][2] * monoval[2][1];
-                grad_grads[start + n_curls][0][2][0] =
-                  -monoval_jplus[0][1] * monoval_i[1][1] * monoval[2][1];
-                grad_grads[start + n_curls][0][2][1] =
-                  -monoval_jplus[0][0] * monoval_i[1][2] * monoval[2][1];
-                grad_grads[start + n_curls][0][2][2] =
-                  -monoval_jplus[0][0] * monoval_i[1][1] * monoval[2][2];
-                grad_grads[start + n_curls][1][0][0] =
-                  monoval_j[0][2] * monoval_i[1][0] * monoval[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + n_curls][1][0][1] =
-                  monoval_j[0][1] * monoval_i[1][1] * monoval[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + n_curls][1][0][2] =
-                  monoval_j[0][1] * monoval_i[1][0] * monoval[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + n_curls][1][1][0] =
-                  monoval_j[0][1] * monoval_i[1][1] * monoval[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + n_curls][1][1][1] =
-                  monoval_j[0][0] * monoval_i[1][2] * monoval[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + n_curls][1][1][2] =
-                  monoval_j[0][0] * monoval_i[1][1] * monoval[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + n_curls][1][2][0] =
-                  monoval_j[0][1] * monoval_i[1][0] * monoval[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + n_curls][1][2][1] =
-                  monoval_j[0][0] * monoval_i[1][1] * monoval[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + n_curls][1][2][2] =
-                  monoval_j[0][0] * monoval_i[1][0] * monoval[2][2] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + n_curls][2][0][0] =
-                  -monoval_j[0][2] * monoval_i[1][1] * monoval_plus[2][0];
-                grad_grads[start + n_curls][2][0][1] =
-                  -monoval_j[0][1] * monoval_i[1][2] * monoval_plus[2][0];
-                grad_grads[start + n_curls][2][0][2] =
-                  -monoval_j[0][1] * monoval_i[1][1] * monoval_plus[2][1];
-                grad_grads[start + n_curls][2][1][0] =
-                  -monoval_j[0][1] * monoval_i[1][2] * monoval_plus[2][0];
-                grad_grads[start + n_curls][2][1][1] =
-                  -monoval_j[0][0] * monoval_i[1][3] * monoval_plus[2][0];
-                grad_grads[start + n_curls][2][1][2] =
-                  -monoval_j[0][0] * monoval_i[1][2] * monoval_plus[2][1];
-                grad_grads[start + n_curls][2][2][0] =
-                  -monoval_j[0][1] * monoval_i[1][1] * monoval_plus[2][1];
-                grad_grads[start + n_curls][2][2][1] =
-                  -monoval_j[0][0] * monoval_i[1][2] * monoval_plus[2][1];
-                grad_grads[start + n_curls][2][2][2] =
-                  -monoval_j[0][0] * monoval_i[1][1] * monoval_plus[2][2];
-
-                grad_grads[start + 2 * n_curls][0][0][0] =
-                  -monoval_jplus[0][2] * monoval[1][0] * monoval_i[2][1];
-                grad_grads[start + 2 * n_curls][0][0][1] =
-                  -monoval_jplus[0][1] * monoval[1][1] * monoval_i[2][1];
-                grad_grads[start + 2 * n_curls][0][0][2] =
-                  -monoval_jplus[0][1] * monoval[1][0] * monoval_i[2][2];
-                grad_grads[start + 2 * n_curls][0][1][0] =
-                  -monoval_jplus[0][1] * monoval[1][1] * monoval_i[2][1];
-                grad_grads[start + 2 * n_curls][0][1][1] =
-                  -monoval_jplus[0][0] * monoval[1][2] * monoval_i[2][1];
-                grad_grads[start + 2 * n_curls][0][1][2] =
-                  -monoval_jplus[0][0] * monoval[1][1] * monoval_i[2][2];
-                grad_grads[start + 2 * n_curls][0][2][0] =
-                  -monoval_jplus[0][1] * monoval[1][0] * monoval_i[2][2];
-                grad_grads[start + 2 * n_curls][0][2][1] =
-                  -monoval_jplus[0][0] * monoval[1][1] * monoval_i[2][2];
-                grad_grads[start + 2 * n_curls][0][2][2] =
-                  -monoval_jplus[0][0] * monoval[1][0] * monoval_i[2][3];
-                grad_grads[start + 2 * n_curls][1][0][0] =
-                  -monoval_j[0][2] * monoval_plus[1][0] * monoval_i[2][1];
-                grad_grads[start + 2 * n_curls][1][0][1] =
-                  -monoval_j[0][1] * monoval_plus[1][1] * monoval_i[2][1];
-                grad_grads[start + 2 * n_curls][1][0][2] =
-                  -monoval_j[0][1] * monoval_plus[1][0] * monoval_i[2][2];
-                grad_grads[start + 2 * n_curls][1][1][0] =
-                  -monoval_j[0][1] * monoval_plus[1][1] * monoval_i[2][1];
-                grad_grads[start + 2 * n_curls][1][1][1] =
-                  -monoval_j[0][0] * monoval_plus[1][2] * monoval_i[2][1];
-                grad_grads[start + 2 * n_curls][1][1][2] =
-                  -monoval_j[0][0] * monoval_plus[1][1] * monoval_i[2][2];
-                grad_grads[start + 2 * n_curls][1][2][0] =
-                  -monoval_j[0][1] * monoval_plus[1][0] * monoval_i[2][2];
-                grad_grads[start + 2 * n_curls][1][2][1] =
-                  -monoval_j[0][0] * monoval_plus[1][1] * monoval_i[2][2];
-                grad_grads[start + 2 * n_curls][1][2][2] =
-                  -monoval_j[0][0] * monoval_plus[1][0] * monoval_i[2][3];
-                grad_grads[start + 2 * n_curls][2][0][0] =
-                  monoval_j[0][2] * monoval[1][0] * monoval_i[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + 2 * n_curls][2][0][1] =
-                  monoval_j[0][1] * monoval[1][1] * monoval_i[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + 2 * n_curls][2][0][2] =
-                  monoval_j[0][1] * monoval[1][0] * monoval_i[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + 2 * n_curls][2][1][0] =
-                  monoval_j[0][1] * monoval[1][1] * monoval_i[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + 2 * n_curls][2][1][1] =
-                  monoval_j[0][0] * monoval[1][2] * monoval_i[2][0] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + 2 * n_curls][2][1][2] =
-                  monoval_j[0][0] * monoval[1][1] * monoval_i[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + 2 * n_curls][2][2][0] =
-                  monoval_j[0][1] * monoval[1][0] * monoval_i[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + 2 * n_curls][2][2][1] =
-                  monoval_j[0][0] * monoval[1][1] * monoval_i[2][1] *
-                  static_cast<double>(j + my_degree + 2);
-                grad_grads[start + 2 * n_curls][2][2][2] =
-                  monoval_j[0][0] * monoval[1][0] * monoval_i[2][2] *
-                  static_cast<double>(j + my_degree + 2);
-
-                if (j != my_degree)
-                  {
-                    grad_grads[start + 1][0][0][0] =
-                      monoval_i[0][2] * monoval[1][0] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 1][0][0][1] =
-                      monoval_i[0][1] * monoval[1][1] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 1][0][0][2] =
-                      monoval_i[0][1] * monoval[1][0] * monoval_j[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 1][0][1][0] =
-                      monoval_i[0][1] * monoval[1][1] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 1][0][1][1] =
-                      monoval_i[0][0] * monoval[1][2] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 1][0][1][2] =
-                      monoval_i[0][0] * monoval[1][1] * monoval_j[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 1][0][2][0] =
-                      monoval_i[0][1] * monoval[1][0] * monoval_j[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 1][0][2][1] =
-                      monoval_i[0][0] * monoval[1][1] * monoval_j[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 1][0][2][2] =
-                      monoval_i[0][0] * monoval[1][0] * monoval_j[2][2] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 1][1][0][0] =
-                      -monoval_i[0][3] * monoval_plus[1][0] * monoval_j[2][0];
-                    grad_grads[start + 1][1][0][1] =
-                      -monoval_i[0][2] * monoval_plus[1][1] * monoval_j[2][0];
-                    grad_grads[start + 1][1][0][2] =
-                      -monoval_i[0][2] * monoval_plus[1][0] * monoval_j[2][1];
-                    grad_grads[start + 1][1][1][0] =
-                      -monoval_i[0][2] * monoval_plus[1][1] * monoval_j[2][0];
-                    grad_grads[start + 1][1][1][1] =
-                      -monoval_i[0][1] * monoval_plus[1][2] * monoval_j[2][0];
-                    grad_grads[start + 1][1][1][2] =
-                      -monoval_i[0][1] * monoval_plus[1][1] * monoval_j[2][1];
-                    grad_grads[start + 1][1][2][0] =
-                      -monoval_i[0][2] * monoval_plus[1][0] * monoval_j[2][1];
-                    grad_grads[start + 1][1][2][1] =
-                      -monoval_i[0][1] * monoval_plus[1][1] * monoval_j[2][1];
-                    grad_grads[start + 1][1][2][2] =
-                      -monoval_i[0][1] * monoval_plus[1][0] * monoval_j[2][2];
-                    grad_grads[start + 1][2][0][0] =
-                      -monoval_i[0][3] * monoval[1][0] * monoval_jplus[2][0];
-                    grad_grads[start + 1][2][0][1] =
-                      -monoval_i[0][2] * monoval[1][1] * monoval_jplus[2][0];
-                    grad_grads[start + 1][2][0][2] =
-                      -monoval_i[0][2] * monoval[1][0] * monoval_jplus[2][1];
-                    grad_grads[start + 1][2][1][0] =
-                      -monoval_i[0][2] * monoval[1][1] * monoval_jplus[2][0];
-                    grad_grads[start + 1][2][1][1] =
-                      -monoval_i[0][1] * monoval[1][2] * monoval_jplus[2][0];
-                    grad_grads[start + 1][2][1][2] =
-                      -monoval_i[0][1] * monoval[1][1] * monoval_jplus[2][1];
-                    grad_grads[start + 1][2][2][0] =
-                      -monoval_i[0][2] * monoval[1][0] * monoval_jplus[2][1];
-                    grad_grads[start + 1][2][2][1] =
-                      -monoval_i[0][1] * monoval[1][1] * monoval_jplus[2][1];
-                    grad_grads[start + 1][2][2][2] =
-                      -monoval_i[0][1] * monoval[1][0] * monoval_jplus[2][2];
-
-                    grad_grads[start + n_curls + 1][0][0][0] =
-                      -monoval_plus[0][2] * monoval_i[1][1] * monoval_j[2][0];
-                    grad_grads[start + n_curls + 1][0][0][1] =
-                      -monoval_plus[0][1] * monoval_i[1][2] * monoval_j[2][0];
-                    grad_grads[start + n_curls + 1][0][0][2] =
-                      -monoval_plus[0][1] * monoval_i[1][1] * monoval_j[2][1];
-                    grad_grads[start + n_curls + 1][0][1][0] =
-                      -monoval_plus[0][1] * monoval_i[1][2] * monoval_j[2][0];
-                    grad_grads[start + n_curls + 1][0][1][1] =
-                      -monoval_plus[0][0] * monoval_i[1][3] * monoval_j[2][0];
-                    grad_grads[start + n_curls + 1][0][1][2] =
-                      -monoval_plus[0][0] * monoval_i[1][2] * monoval_j[2][1];
-                    grad_grads[start + n_curls + 1][0][2][0] =
-                      -monoval_plus[0][1] * monoval_i[1][1] * monoval_j[2][1];
-                    grad_grads[start + n_curls + 1][0][2][1] =
-                      -monoval_plus[0][0] * monoval_i[1][2] * monoval_j[2][1];
-                    grad_grads[start + n_curls + 1][0][2][2] =
-                      -monoval_plus[0][0] * monoval_i[1][1] * monoval_j[2][2];
-                    grad_grads[start + n_curls + 1][1][0][0] =
-                      monoval[0][2] * monoval_i[1][0] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + n_curls + 1][1][0][1] =
-                      monoval[0][1] * monoval_i[1][1] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + n_curls + 1][1][0][2] =
-                      monoval[0][1] * monoval_i[1][0] * monoval_j[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + n_curls + 1][1][1][0] =
-                      monoval[0][1] * monoval_i[1][1] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + n_curls + 1][1][1][1] =
-                      monoval[0][0] * monoval_i[1][2] * monoval_j[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + n_curls + 1][1][1][2] =
-                      monoval[0][0] * monoval_i[1][1] * monoval_j[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + n_curls + 1][1][2][0] =
-                      monoval[0][1] * monoval_i[1][0] * monoval_j[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + n_curls + 1][1][2][1] =
-                      monoval[0][0] * monoval_i[1][1] * monoval_j[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + n_curls + 1][1][2][2] =
-                      monoval[0][0] * monoval_i[1][0] * monoval_j[2][2] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + n_curls + 1][2][0][0] =
-                      -monoval[0][2] * monoval_i[1][1] * monoval_jplus[2][0];
-                    grad_grads[start + n_curls + 1][2][0][1] =
-                      -monoval[0][1] * monoval_i[1][2] * monoval_jplus[2][0];
-                    grad_grads[start + n_curls + 1][2][0][2] =
-                      -monoval[0][1] * monoval_i[1][1] * monoval_jplus[2][1];
-                    grad_grads[start + n_curls + 1][2][1][0] =
-                      -monoval[0][1] * monoval_i[1][2] * monoval_jplus[2][0];
-                    grad_grads[start + n_curls + 1][2][1][1] =
-                      -monoval[0][0] * monoval_i[1][3] * monoval_jplus[2][0];
-                    grad_grads[start + n_curls + 1][2][1][2] =
-                      -monoval[0][0] * monoval_i[1][2] * monoval_jplus[2][1];
-                    grad_grads[start + n_curls + 1][2][2][0] =
-                      -monoval[0][1] * monoval_i[1][1] * monoval_jplus[2][1];
-                    grad_grads[start + n_curls + 1][2][2][1] =
-                      -monoval[0][0] * monoval_i[1][2] * monoval_jplus[2][1];
-                    grad_grads[start + n_curls + 1][2][2][2] =
-                      -monoval[0][0] * monoval_i[1][1] * monoval_jplus[2][2];
-
-                    grad_grads[start + 2 * n_curls + 1][0][0][0] =
-                      -monoval_plus[0][2] * monoval_j[1][0] * monoval_i[2][1];
-                    grad_grads[start + 2 * n_curls + 1][0][0][1] =
-                      -monoval_plus[0][1] * monoval_j[1][1] * monoval_i[2][1];
-                    grad_grads[start + 2 * n_curls + 1][0][0][2] =
-                      -monoval_plus[0][1] * monoval_j[1][0] * monoval_i[2][2];
-                    grad_grads[start + 2 * n_curls + 1][0][1][0] =
-                      -monoval_plus[0][1] * monoval_j[1][1] * monoval_i[2][1];
-                    grad_grads[start + 2 * n_curls + 1][0][1][1] =
-                      -monoval_plus[0][0] * monoval_j[1][2] * monoval_i[2][1];
-                    grad_grads[start + 2 * n_curls + 1][0][1][2] =
-                      -monoval_plus[0][0] * monoval_j[1][1] * monoval_i[2][2];
-                    grad_grads[start + 2 * n_curls + 1][0][2][0] =
-                      -monoval_plus[0][1] * monoval_j[1][0] * monoval_i[2][2];
-                    grad_grads[start + 2 * n_curls + 1][0][2][1] =
-                      -monoval_plus[0][0] * monoval_j[1][1] * monoval_i[2][2];
-                    grad_grads[start + 2 * n_curls + 1][0][2][2] =
-                      -monoval_plus[0][0] * monoval_j[1][0] * monoval_i[2][3];
-                    grad_grads[start + 2 * n_curls + 1][1][0][0] =
-                      -monoval[0][2] * monoval_jplus[1][0] * monoval_i[2][1];
-                    grad_grads[start + 2 * n_curls + 1][1][0][1] =
-                      -monoval[0][1] * monoval_jplus[1][1] * monoval_i[2][1];
-                    grad_grads[start + 2 * n_curls + 1][1][0][2] =
-                      -monoval[0][1] * monoval_jplus[1][0] * monoval_i[2][2];
-                    grad_grads[start + 2 * n_curls + 1][1][1][0] =
-                      -monoval[0][1] * monoval_jplus[1][1] * monoval_i[2][1];
-                    grad_grads[start + 2 * n_curls + 1][1][1][1] =
-                      -monoval[0][0] * monoval_jplus[1][2] * monoval_i[2][1];
-                    grad_grads[start + 2 * n_curls + 1][1][1][2] =
-                      -monoval[0][0] * monoval_jplus[1][1] * monoval_i[2][2];
-                    grad_grads[start + 2 * n_curls + 1][1][2][0] =
-                      -monoval[0][1] * monoval_jplus[1][0] * monoval_i[2][2];
-                    grad_grads[start + 2 * n_curls + 1][1][2][1] =
-                      -monoval[0][0] * monoval_jplus[1][1] * monoval_i[2][2];
-                    grad_grads[start + 2 * n_curls + 1][1][2][2] =
-                      -monoval[0][0] * monoval_jplus[1][0] * monoval_i[2][3];
-                    grad_grads[start + 2 * n_curls + 1][2][0][0] =
-                      monoval[0][2] * monoval_j[1][0] * monoval_i[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 2 * n_curls + 1][2][0][1] =
-                      monoval[0][1] * monoval_j[1][1] * monoval_i[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 2 * n_curls + 1][2][0][2] =
-                      monoval[0][1] * monoval_j[1][0] * monoval_i[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 2 * n_curls + 1][2][1][0] =
-                      monoval[0][1] * monoval_j[1][1] * monoval_i[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 2 * n_curls + 1][2][1][1] =
-                      monoval[0][0] * monoval_j[1][2] * monoval_i[2][0] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 2 * n_curls + 1][2][1][2] =
-                      monoval[0][0] * monoval_j[1][1] * monoval_i[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 2 * n_curls + 1][2][2][0] =
-                      monoval[0][1] * monoval_j[1][0] * monoval_i[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 2 * n_curls + 1][2][2][1] =
-                      monoval[0][0] * monoval_j[1][1] * monoval_i[2][1] *
-                      static_cast<double>(j + my_degree + 2);
-                    grad_grads[start + 2 * n_curls + 1][2][2][2] =
-                      monoval[0][0] * monoval_j[1][0] * monoval_i[2][2] *
-                      static_cast<double>(j + my_degree + 2);
-                  }
-              }
-
-            if (j == my_degree)
-              start += 1;
-            else
-              start += 2;
-          }
-      }
-    Assert(start == this->n() - 2 * n_curls, ExcInternalError());
-  }
+                  values[start][1] =
+                    -monoval_i[0][1] * monoval_jplus[1][0] * monoval[2][0];
+                  values[start][2] =
+                    -monoval_i[0][1] * monoval_j[1][0] * monoval_plus[2][0];
+
+                  values[start + n_curls][0] =
+                    -monoval_jplus[0][0] * monoval_i[1][1] * monoval[2][0];
+                  values[start + n_curls][1] =
+                    monoval_j[0][0] * monoval_i[1][0] * monoval[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  values[start + n_curls][2] =
+                    -monoval_j[0][0] * monoval_i[1][1] * monoval_plus[2][0];
+
+                  values[start + 2 * n_curls][0] =
+                    -monoval_jplus[0][0] * monoval[1][0] * monoval_i[2][1];
+                  values[start + 2 * n_curls][1] =
+                    -monoval_j[0][0] * monoval_plus[1][0] * monoval_i[2][1];
+                  values[start + 2 * n_curls][2] =
+                    monoval_j[0][0] * monoval[1][0] * monoval_i[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+
+                  // Only unique triples of powers (i j k)
+                  // and (i k j) are allowed, 0 <= i,j <= k
+                  if (j != my_degree)
+                    {
+                      values[start + 1][0] =
+                        monoval_i[0][0] * monoval[1][0] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      values[start + 1][1] =
+                        -monoval_i[0][1] * monoval_plus[1][0] * monoval_j[2][0];
+                      values[start + 1][2] =
+                        -monoval_i[0][1] * monoval[1][0] * monoval_jplus[2][0];
+
+                      values[start + n_curls + 1][0] =
+                        -monoval_plus[0][0] * monoval_i[1][1] * monoval_j[2][0];
+                      values[start + n_curls + 1][1] =
+                        monoval[0][0] * monoval_i[1][0] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      values[start + n_curls + 1][2] =
+                        -monoval[0][0] * monoval_i[1][1] * monoval_jplus[2][0];
+
+                      values[start + 2 * n_curls + 1][0] =
+                        -monoval_plus[0][0] * monoval_j[1][0] * monoval_i[2][1];
+                      values[start + 2 * n_curls + 1][1] =
+                        -monoval[0][0] * monoval_jplus[1][0] * monoval_i[2][1];
+                      values[start + 2 * n_curls + 1][2] =
+                        monoval[0][0] * monoval_j[1][0] * monoval_i[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                    }
+                }
+
+              if (grads.size() != 0)
+                {
+                  grads[start][0][0] = monoval_i[0][1] * monoval_j[1][0] *
+                                       monoval[2][0] *
+                                       static_cast<double>(j + my_degree + 2);
+                  grads[start][0][1] = monoval_i[0][0] * monoval_j[1][1] *
+                                       monoval[2][0] *
+                                       static_cast<double>(j + my_degree + 2);
+                  grads[start][0][2] = monoval_i[0][0] * monoval_j[1][0] *
+                                       monoval[2][1] *
+                                       static_cast<double>(j + my_degree + 2);
+                  grads[start][1][0] =
+                    -monoval_i[0][2] * monoval_jplus[1][0] * monoval[2][0];
+                  grads[start][1][1] =
+                    -monoval_i[0][1] * monoval_jplus[1][1] * monoval[2][0];
+                  grads[start][1][2] =
+                    -monoval_i[0][1] * monoval_jplus[1][0] * monoval[2][1];
+                  grads[start][2][0] =
+                    -monoval_i[0][2] * monoval_j[1][0] * monoval_plus[2][0];
+                  grads[start][2][1] =
+                    -monoval_i[0][1] * monoval_j[1][1] * monoval_plus[2][0];
+                  grads[start][2][2] =
+                    -monoval_i[0][1] * monoval_j[1][0] * monoval_plus[2][1];
+
+                  grads[start + n_curls][0][0] =
+                    -monoval_jplus[0][1] * monoval_i[1][1] * monoval[2][0];
+                  grads[start + n_curls][0][1] =
+                    -monoval_jplus[0][0] * monoval_i[1][2] * monoval[2][0];
+                  grads[start + n_curls][0][2] =
+                    -monoval_jplus[0][0] * monoval_i[1][1] * monoval[2][1];
+                  grads[start + n_curls][1][0] =
+                    monoval_j[0][1] * monoval_i[1][0] * monoval[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grads[start + n_curls][1][1] =
+                    monoval_j[0][0] * monoval_i[1][1] * monoval[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grads[start + n_curls][1][2] =
+                    monoval_j[0][0] * monoval_i[1][0] * monoval[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grads[start + n_curls][2][0] =
+                    -monoval_j[0][1] * monoval_i[1][1] * monoval_plus[2][0];
+                  grads[start + n_curls][2][1] =
+                    -monoval_j[0][0] * monoval_i[1][2] * monoval_plus[2][0];
+                  grads[start + n_curls][2][2] =
+                    -monoval_j[0][0] * monoval_i[1][1] * monoval_plus[2][1];
+
+                  grads[start + 2 * n_curls][0][0] =
+                    -monoval_jplus[0][1] * monoval[1][0] * monoval_i[2][1];
+                  grads[start + 2 * n_curls][0][1] =
+                    -monoval_jplus[0][0] * monoval[1][1] * monoval_i[2][1];
+                  grads[start + 2 * n_curls][0][2] =
+                    -monoval_jplus[0][0] * monoval[1][0] * monoval_i[2][2];
+                  grads[start + 2 * n_curls][1][0] =
+                    -monoval_j[0][1] * monoval_plus[1][0] * monoval_i[2][1];
+                  grads[start + 2 * n_curls][1][1] =
+                    -monoval_j[0][0] * monoval_plus[1][1] * monoval_i[2][1];
+                  grads[start + 2 * n_curls][1][2] =
+                    -monoval_j[0][0] * monoval_plus[1][0] * monoval_i[2][2];
+                  grads[start + 2 * n_curls][2][0] =
+                    monoval_j[0][1] * monoval[1][0] * monoval_i[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grads[start + 2 * n_curls][2][1] =
+                    monoval_j[0][0] * monoval[1][1] * monoval_i[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grads[start + 2 * n_curls][2][2] =
+                    monoval_j[0][0] * monoval[1][0] * monoval_i[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+
+                  if (j != my_degree)
+                    {
+                      grads[start + 1][0][0] =
+                        monoval_i[0][1] * monoval[1][0] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grads[start + 1][0][1] =
+                        monoval_i[0][0] * monoval[1][1] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grads[start + 1][0][2] =
+                        monoval_i[0][0] * monoval[1][0] * monoval_j[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grads[start + 1][1][0] =
+                        -monoval_i[0][2] * monoval_plus[1][0] * monoval_j[2][0];
+                      grads[start + 1][1][1] =
+                        -monoval_i[0][1] * monoval_plus[1][1] * monoval_j[2][0];
+                      grads[start + 1][1][2] =
+                        -monoval_i[0][1] * monoval_plus[1][0] * monoval_j[2][1];
+                      grads[start + 1][2][0] =
+                        -monoval_i[0][2] * monoval[1][0] * monoval_jplus[2][0];
+                      grads[start + 1][2][1] =
+                        -monoval_i[0][1] * monoval[1][1] * monoval_jplus[2][0];
+                      grads[start + 1][2][2] =
+                        -monoval_i[0][1] * monoval[1][0] * monoval_jplus[2][1];
+
+                      grads[start + n_curls + 1][0][0] =
+                        -monoval_plus[0][1] * monoval_i[1][1] * monoval_j[2][0];
+                      grads[start + n_curls + 1][0][1] =
+                        -monoval_plus[0][0] * monoval_i[1][2] * monoval_j[2][0];
+                      grads[start + n_curls + 1][0][2] =
+                        -monoval_plus[0][0] * monoval_i[1][1] * monoval_j[2][1];
+                      grads[start + n_curls + 1][1][0] =
+                        monoval[0][1] * monoval_i[1][0] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grads[start + n_curls + 1][1][1] =
+                        monoval[0][0] * monoval_i[1][1] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grads[start + n_curls + 1][1][2] =
+                        monoval[0][0] * monoval_i[1][0] * monoval_j[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grads[start + n_curls + 1][2][0] =
+                        -monoval[0][1] * monoval_i[1][1] * monoval_jplus[2][0];
+                      grads[start + n_curls + 1][2][1] =
+                        -monoval[0][0] * monoval_i[1][2] * monoval_jplus[2][0];
+                      grads[start + n_curls + 1][2][2] =
+                        -monoval[0][0] * monoval_i[1][1] * monoval_jplus[2][1];
+
+                      grads[start + 2 * n_curls + 1][0][0] =
+                        -monoval_plus[0][1] * monoval_j[1][0] * monoval_i[2][1];
+                      grads[start + 2 * n_curls + 1][0][1] =
+                        -monoval_plus[0][0] * monoval_j[1][1] * monoval_i[2][1];
+                      grads[start + 2 * n_curls + 1][0][2] =
+                        -monoval_plus[0][0] * monoval_j[1][0] * monoval_i[2][2];
+                      grads[start + 2 * n_curls + 1][1][0] =
+                        -monoval[0][1] * monoval_jplus[1][0] * monoval_i[2][1];
+                      grads[start + 2 * n_curls + 1][1][1] =
+                        -monoval[0][0] * monoval_jplus[1][1] * monoval_i[2][1];
+                      grads[start + 2 * n_curls + 1][1][2] =
+                        -monoval[0][0] * monoval_jplus[1][0] * monoval_i[2][2];
+                      grads[start + 2 * n_curls + 1][2][0] =
+                        monoval[0][1] * monoval_j[1][0] * monoval_i[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grads[start + 2 * n_curls + 1][2][1] =
+                        monoval[0][0] * monoval_j[1][1] * monoval_i[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grads[start + 2 * n_curls + 1][2][2] =
+                        monoval[0][0] * monoval_j[1][0] * monoval_i[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                    }
+                }
+
+              if (grad_grads.size() != 0)
+                {
+                  grad_grads[start][0][0][0] =
+                    monoval_i[0][2] * monoval_j[1][0] * monoval[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start][0][0][1] =
+                    monoval_i[0][1] * monoval_j[1][1] * monoval[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start][0][0][2] =
+                    monoval_i[0][1] * monoval_j[1][0] * monoval[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start][0][1][0] =
+                    monoval_i[0][1] * monoval_j[1][1] * monoval[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start][0][1][1] =
+                    monoval_i[0][0] * monoval_j[1][2] * monoval[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start][0][1][2] =
+                    monoval_i[0][0] * monoval_j[1][1] * monoval[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start][0][2][0] =
+                    monoval_i[0][1] * monoval_j[1][0] * monoval[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start][0][2][1] =
+                    monoval_i[0][0] * monoval_j[1][1] * monoval[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start][0][2][2] =
+                    monoval_i[0][0] * monoval_j[1][0] * monoval[2][2] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start][1][0][0] =
+                    -monoval_i[0][3] * monoval_jplus[1][0] * monoval[2][0];
+                  grad_grads[start][1][0][1] =
+                    -monoval_i[0][2] * monoval_jplus[1][1] * monoval[2][0];
+                  grad_grads[start][1][0][2] =
+                    -monoval_i[0][2] * monoval_jplus[1][0] * monoval[2][1];
+                  grad_grads[start][1][1][0] =
+                    -monoval_i[0][2] * monoval_jplus[1][1] * monoval[2][0];
+                  grad_grads[start][1][1][1] =
+                    -monoval_i[0][1] * monoval_jplus[1][2] * monoval[2][0];
+                  grad_grads[start][1][1][2] =
+                    -monoval_i[0][1] * monoval_jplus[1][1] * monoval[2][1];
+                  grad_grads[start][1][2][0] =
+                    -monoval_i[0][2] * monoval_jplus[1][0] * monoval[2][1];
+                  grad_grads[start][1][2][1] =
+                    -monoval_i[0][1] * monoval_jplus[1][1] * monoval[2][1];
+                  grad_grads[start][1][2][2] =
+                    -monoval_i[0][1] * monoval_jplus[1][0] * monoval[2][2];
+                  grad_grads[start][2][0][0] =
+                    -monoval_i[0][3] * monoval_j[1][0] * monoval_plus[2][0];
+                  grad_grads[start][2][0][1] =
+                    -monoval_i[0][2] * monoval_j[1][1] * monoval_plus[2][0];
+                  grad_grads[start][2][0][2] =
+                    -monoval_i[0][2] * monoval_j[1][0] * monoval_plus[2][1];
+                  grad_grads[start][2][1][0] =
+                    -monoval_i[0][2] * monoval_j[1][1] * monoval_plus[2][0];
+                  grad_grads[start][2][1][1] =
+                    -monoval_i[0][1] * monoval_j[1][2] * monoval_plus[2][0];
+                  grad_grads[start][2][1][2] =
+                    -monoval_i[0][1] * monoval_j[1][1] * monoval_plus[2][1];
+                  grad_grads[start][2][2][0] =
+                    -monoval_i[0][2] * monoval_j[1][0] * monoval_plus[2][1];
+                  grad_grads[start][2][2][1] =
+                    -monoval_i[0][1] * monoval_j[1][1] * monoval_plus[2][1];
+                  grad_grads[start][2][2][2] =
+                    -monoval_i[0][1] * monoval_j[1][0] * monoval_plus[2][2];
+
+                  grad_grads[start + n_curls][0][0][0] =
+                    -monoval_jplus[0][2] * monoval_i[1][1] * monoval[2][0];
+                  grad_grads[start + n_curls][0][0][1] =
+                    -monoval_jplus[0][1] * monoval_i[1][2] * monoval[2][0];
+                  grad_grads[start + n_curls][0][0][2] =
+                    -monoval_jplus[0][1] * monoval_i[1][1] * monoval[2][1];
+                  grad_grads[start + n_curls][0][1][0] =
+                    -monoval_jplus[0][1] * monoval_i[1][2] * monoval[2][0];
+                  grad_grads[start + n_curls][0][1][1] =
+                    -monoval_jplus[0][0] * monoval_i[1][3] * monoval[2][0];
+                  grad_grads[start + n_curls][0][1][2] =
+                    -monoval_jplus[0][0] * monoval_i[1][2] * monoval[2][1];
+                  grad_grads[start + n_curls][0][2][0] =
+                    -monoval_jplus[0][1] * monoval_i[1][1] * monoval[2][1];
+                  grad_grads[start + n_curls][0][2][1] =
+                    -monoval_jplus[0][0] * monoval_i[1][2] * monoval[2][1];
+                  grad_grads[start + n_curls][0][2][2] =
+                    -monoval_jplus[0][0] * monoval_i[1][1] * monoval[2][2];
+                  grad_grads[start + n_curls][1][0][0] =
+                    monoval_j[0][2] * monoval_i[1][0] * monoval[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + n_curls][1][0][1] =
+                    monoval_j[0][1] * monoval_i[1][1] * monoval[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + n_curls][1][0][2] =
+                    monoval_j[0][1] * monoval_i[1][0] * monoval[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + n_curls][1][1][0] =
+                    monoval_j[0][1] * monoval_i[1][1] * monoval[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + n_curls][1][1][1] =
+                    monoval_j[0][0] * monoval_i[1][2] * monoval[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + n_curls][1][1][2] =
+                    monoval_j[0][0] * monoval_i[1][1] * monoval[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + n_curls][1][2][0] =
+                    monoval_j[0][1] * monoval_i[1][0] * monoval[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + n_curls][1][2][1] =
+                    monoval_j[0][0] * monoval_i[1][1] * monoval[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + n_curls][1][2][2] =
+                    monoval_j[0][0] * monoval_i[1][0] * monoval[2][2] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + n_curls][2][0][0] =
+                    -monoval_j[0][2] * monoval_i[1][1] * monoval_plus[2][0];
+                  grad_grads[start + n_curls][2][0][1] =
+                    -monoval_j[0][1] * monoval_i[1][2] * monoval_plus[2][0];
+                  grad_grads[start + n_curls][2][0][2] =
+                    -monoval_j[0][1] * monoval_i[1][1] * monoval_plus[2][1];
+                  grad_grads[start + n_curls][2][1][0] =
+                    -monoval_j[0][1] * monoval_i[1][2] * monoval_plus[2][0];
+                  grad_grads[start + n_curls][2][1][1] =
+                    -monoval_j[0][0] * monoval_i[1][3] * monoval_plus[2][0];
+                  grad_grads[start + n_curls][2][1][2] =
+                    -monoval_j[0][0] * monoval_i[1][2] * monoval_plus[2][1];
+                  grad_grads[start + n_curls][2][2][0] =
+                    -monoval_j[0][1] * monoval_i[1][1] * monoval_plus[2][1];
+                  grad_grads[start + n_curls][2][2][1] =
+                    -monoval_j[0][0] * monoval_i[1][2] * monoval_plus[2][1];
+                  grad_grads[start + n_curls][2][2][2] =
+                    -monoval_j[0][0] * monoval_i[1][1] * monoval_plus[2][2];
+
+                  grad_grads[start + 2 * n_curls][0][0][0] =
+                    -monoval_jplus[0][2] * monoval[1][0] * monoval_i[2][1];
+                  grad_grads[start + 2 * n_curls][0][0][1] =
+                    -monoval_jplus[0][1] * monoval[1][1] * monoval_i[2][1];
+                  grad_grads[start + 2 * n_curls][0][0][2] =
+                    -monoval_jplus[0][1] * monoval[1][0] * monoval_i[2][2];
+                  grad_grads[start + 2 * n_curls][0][1][0] =
+                    -monoval_jplus[0][1] * monoval[1][1] * monoval_i[2][1];
+                  grad_grads[start + 2 * n_curls][0][1][1] =
+                    -monoval_jplus[0][0] * monoval[1][2] * monoval_i[2][1];
+                  grad_grads[start + 2 * n_curls][0][1][2] =
+                    -monoval_jplus[0][0] * monoval[1][1] * monoval_i[2][2];
+                  grad_grads[start + 2 * n_curls][0][2][0] =
+                    -monoval_jplus[0][1] * monoval[1][0] * monoval_i[2][2];
+                  grad_grads[start + 2 * n_curls][0][2][1] =
+                    -monoval_jplus[0][0] * monoval[1][1] * monoval_i[2][2];
+                  grad_grads[start + 2 * n_curls][0][2][2] =
+                    -monoval_jplus[0][0] * monoval[1][0] * monoval_i[2][3];
+                  grad_grads[start + 2 * n_curls][1][0][0] =
+                    -monoval_j[0][2] * monoval_plus[1][0] * monoval_i[2][1];
+                  grad_grads[start + 2 * n_curls][1][0][1] =
+                    -monoval_j[0][1] * monoval_plus[1][1] * monoval_i[2][1];
+                  grad_grads[start + 2 * n_curls][1][0][2] =
+                    -monoval_j[0][1] * monoval_plus[1][0] * monoval_i[2][2];
+                  grad_grads[start + 2 * n_curls][1][1][0] =
+                    -monoval_j[0][1] * monoval_plus[1][1] * monoval_i[2][1];
+                  grad_grads[start + 2 * n_curls][1][1][1] =
+                    -monoval_j[0][0] * monoval_plus[1][2] * monoval_i[2][1];
+                  grad_grads[start + 2 * n_curls][1][1][2] =
+                    -monoval_j[0][0] * monoval_plus[1][1] * monoval_i[2][2];
+                  grad_grads[start + 2 * n_curls][1][2][0] =
+                    -monoval_j[0][1] * monoval_plus[1][0] * monoval_i[2][2];
+                  grad_grads[start + 2 * n_curls][1][2][1] =
+                    -monoval_j[0][0] * monoval_plus[1][1] * monoval_i[2][2];
+                  grad_grads[start + 2 * n_curls][1][2][2] =
+                    -monoval_j[0][0] * monoval_plus[1][0] * monoval_i[2][3];
+                  grad_grads[start + 2 * n_curls][2][0][0] =
+                    monoval_j[0][2] * monoval[1][0] * monoval_i[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + 2 * n_curls][2][0][1] =
+                    monoval_j[0][1] * monoval[1][1] * monoval_i[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + 2 * n_curls][2][0][2] =
+                    monoval_j[0][1] * monoval[1][0] * monoval_i[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + 2 * n_curls][2][1][0] =
+                    monoval_j[0][1] * monoval[1][1] * monoval_i[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + 2 * n_curls][2][1][1] =
+                    monoval_j[0][0] * monoval[1][2] * monoval_i[2][0] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + 2 * n_curls][2][1][2] =
+                    monoval_j[0][0] * monoval[1][1] * monoval_i[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + 2 * n_curls][2][2][0] =
+                    monoval_j[0][1] * monoval[1][0] * monoval_i[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + 2 * n_curls][2][2][1] =
+                    monoval_j[0][0] * monoval[1][1] * monoval_i[2][1] *
+                    static_cast<double>(j + my_degree + 2);
+                  grad_grads[start + 2 * n_curls][2][2][2] =
+                    monoval_j[0][0] * monoval[1][0] * monoval_i[2][2] *
+                    static_cast<double>(j + my_degree + 2);
+
+                  if (j != my_degree)
+                    {
+                      grad_grads[start + 1][0][0][0] =
+                        monoval_i[0][2] * monoval[1][0] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 1][0][0][1] =
+                        monoval_i[0][1] * monoval[1][1] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 1][0][0][2] =
+                        monoval_i[0][1] * monoval[1][0] * monoval_j[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 1][0][1][0] =
+                        monoval_i[0][1] * monoval[1][1] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 1][0][1][1] =
+                        monoval_i[0][0] * monoval[1][2] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 1][0][1][2] =
+                        monoval_i[0][0] * monoval[1][1] * monoval_j[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 1][0][2][0] =
+                        monoval_i[0][1] * monoval[1][0] * monoval_j[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 1][0][2][1] =
+                        monoval_i[0][0] * monoval[1][1] * monoval_j[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 1][0][2][2] =
+                        monoval_i[0][0] * monoval[1][0] * monoval_j[2][2] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 1][1][0][0] =
+                        -monoval_i[0][3] * monoval_plus[1][0] * monoval_j[2][0];
+                      grad_grads[start + 1][1][0][1] =
+                        -monoval_i[0][2] * monoval_plus[1][1] * monoval_j[2][0];
+                      grad_grads[start + 1][1][0][2] =
+                        -monoval_i[0][2] * monoval_plus[1][0] * monoval_j[2][1];
+                      grad_grads[start + 1][1][1][0] =
+                        -monoval_i[0][2] * monoval_plus[1][1] * monoval_j[2][0];
+                      grad_grads[start + 1][1][1][1] =
+                        -monoval_i[0][1] * monoval_plus[1][2] * monoval_j[2][0];
+                      grad_grads[start + 1][1][1][2] =
+                        -monoval_i[0][1] * monoval_plus[1][1] * monoval_j[2][1];
+                      grad_grads[start + 1][1][2][0] =
+                        -monoval_i[0][2] * monoval_plus[1][0] * monoval_j[2][1];
+                      grad_grads[start + 1][1][2][1] =
+                        -monoval_i[0][1] * monoval_plus[1][1] * monoval_j[2][1];
+                      grad_grads[start + 1][1][2][2] =
+                        -monoval_i[0][1] * monoval_plus[1][0] * monoval_j[2][2];
+                      grad_grads[start + 1][2][0][0] =
+                        -monoval_i[0][3] * monoval[1][0] * monoval_jplus[2][0];
+                      grad_grads[start + 1][2][0][1] =
+                        -monoval_i[0][2] * monoval[1][1] * monoval_jplus[2][0];
+                      grad_grads[start + 1][2][0][2] =
+                        -monoval_i[0][2] * monoval[1][0] * monoval_jplus[2][1];
+                      grad_grads[start + 1][2][1][0] =
+                        -monoval_i[0][2] * monoval[1][1] * monoval_jplus[2][0];
+                      grad_grads[start + 1][2][1][1] =
+                        -monoval_i[0][1] * monoval[1][2] * monoval_jplus[2][0];
+                      grad_grads[start + 1][2][1][2] =
+                        -monoval_i[0][1] * monoval[1][1] * monoval_jplus[2][1];
+                      grad_grads[start + 1][2][2][0] =
+                        -monoval_i[0][2] * monoval[1][0] * monoval_jplus[2][1];
+                      grad_grads[start + 1][2][2][1] =
+                        -monoval_i[0][1] * monoval[1][1] * monoval_jplus[2][1];
+                      grad_grads[start + 1][2][2][2] =
+                        -monoval_i[0][1] * monoval[1][0] * monoval_jplus[2][2];
+
+                      grad_grads[start + n_curls + 1][0][0][0] =
+                        -monoval_plus[0][2] * monoval_i[1][1] * monoval_j[2][0];
+                      grad_grads[start + n_curls + 1][0][0][1] =
+                        -monoval_plus[0][1] * monoval_i[1][2] * monoval_j[2][0];
+                      grad_grads[start + n_curls + 1][0][0][2] =
+                        -monoval_plus[0][1] * monoval_i[1][1] * monoval_j[2][1];
+                      grad_grads[start + n_curls + 1][0][1][0] =
+                        -monoval_plus[0][1] * monoval_i[1][2] * monoval_j[2][0];
+                      grad_grads[start + n_curls + 1][0][1][1] =
+                        -monoval_plus[0][0] * monoval_i[1][3] * monoval_j[2][0];
+                      grad_grads[start + n_curls + 1][0][1][2] =
+                        -monoval_plus[0][0] * monoval_i[1][2] * monoval_j[2][1];
+                      grad_grads[start + n_curls + 1][0][2][0] =
+                        -monoval_plus[0][1] * monoval_i[1][1] * monoval_j[2][1];
+                      grad_grads[start + n_curls + 1][0][2][1] =
+                        -monoval_plus[0][0] * monoval_i[1][2] * monoval_j[2][1];
+                      grad_grads[start + n_curls + 1][0][2][2] =
+                        -monoval_plus[0][0] * monoval_i[1][1] * monoval_j[2][2];
+                      grad_grads[start + n_curls + 1][1][0][0] =
+                        monoval[0][2] * monoval_i[1][0] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + n_curls + 1][1][0][1] =
+                        monoval[0][1] * monoval_i[1][1] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + n_curls + 1][1][0][2] =
+                        monoval[0][1] * monoval_i[1][0] * monoval_j[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + n_curls + 1][1][1][0] =
+                        monoval[0][1] * monoval_i[1][1] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + n_curls + 1][1][1][1] =
+                        monoval[0][0] * monoval_i[1][2] * monoval_j[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + n_curls + 1][1][1][2] =
+                        monoval[0][0] * monoval_i[1][1] * monoval_j[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + n_curls + 1][1][2][0] =
+                        monoval[0][1] * monoval_i[1][0] * monoval_j[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + n_curls + 1][1][2][1] =
+                        monoval[0][0] * monoval_i[1][1] * monoval_j[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + n_curls + 1][1][2][2] =
+                        monoval[0][0] * monoval_i[1][0] * monoval_j[2][2] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + n_curls + 1][2][0][0] =
+                        -monoval[0][2] * monoval_i[1][1] * monoval_jplus[2][0];
+                      grad_grads[start + n_curls + 1][2][0][1] =
+                        -monoval[0][1] * monoval_i[1][2] * monoval_jplus[2][0];
+                      grad_grads[start + n_curls + 1][2][0][2] =
+                        -monoval[0][1] * monoval_i[1][1] * monoval_jplus[2][1];
+                      grad_grads[start + n_curls + 1][2][1][0] =
+                        -monoval[0][1] * monoval_i[1][2] * monoval_jplus[2][0];
+                      grad_grads[start + n_curls + 1][2][1][1] =
+                        -monoval[0][0] * monoval_i[1][3] * monoval_jplus[2][0];
+                      grad_grads[start + n_curls + 1][2][1][2] =
+                        -monoval[0][0] * monoval_i[1][2] * monoval_jplus[2][1];
+                      grad_grads[start + n_curls + 1][2][2][0] =
+                        -monoval[0][1] * monoval_i[1][1] * monoval_jplus[2][1];
+                      grad_grads[start + n_curls + 1][2][2][1] =
+                        -monoval[0][0] * monoval_i[1][2] * monoval_jplus[2][1];
+                      grad_grads[start + n_curls + 1][2][2][2] =
+                        -monoval[0][0] * monoval_i[1][1] * monoval_jplus[2][2];
+
+                      grad_grads[start + 2 * n_curls + 1][0][0][0] =
+                        -monoval_plus[0][2] * monoval_j[1][0] * monoval_i[2][1];
+                      grad_grads[start + 2 * n_curls + 1][0][0][1] =
+                        -monoval_plus[0][1] * monoval_j[1][1] * monoval_i[2][1];
+                      grad_grads[start + 2 * n_curls + 1][0][0][2] =
+                        -monoval_plus[0][1] * monoval_j[1][0] * monoval_i[2][2];
+                      grad_grads[start + 2 * n_curls + 1][0][1][0] =
+                        -monoval_plus[0][1] * monoval_j[1][1] * monoval_i[2][1];
+                      grad_grads[start + 2 * n_curls + 1][0][1][1] =
+                        -monoval_plus[0][0] * monoval_j[1][2] * monoval_i[2][1];
+                      grad_grads[start + 2 * n_curls + 1][0][1][2] =
+                        -monoval_plus[0][0] * monoval_j[1][1] * monoval_i[2][2];
+                      grad_grads[start + 2 * n_curls + 1][0][2][0] =
+                        -monoval_plus[0][1] * monoval_j[1][0] * monoval_i[2][2];
+                      grad_grads[start + 2 * n_curls + 1][0][2][1] =
+                        -monoval_plus[0][0] * monoval_j[1][1] * monoval_i[2][2];
+                      grad_grads[start + 2 * n_curls + 1][0][2][2] =
+                        -monoval_plus[0][0] * monoval_j[1][0] * monoval_i[2][3];
+                      grad_grads[start + 2 * n_curls + 1][1][0][0] =
+                        -monoval[0][2] * monoval_jplus[1][0] * monoval_i[2][1];
+                      grad_grads[start + 2 * n_curls + 1][1][0][1] =
+                        -monoval[0][1] * monoval_jplus[1][1] * monoval_i[2][1];
+                      grad_grads[start + 2 * n_curls + 1][1][0][2] =
+                        -monoval[0][1] * monoval_jplus[1][0] * monoval_i[2][2];
+                      grad_grads[start + 2 * n_curls + 1][1][1][0] =
+                        -monoval[0][1] * monoval_jplus[1][1] * monoval_i[2][1];
+                      grad_grads[start + 2 * n_curls + 1][1][1][1] =
+                        -monoval[0][0] * monoval_jplus[1][2] * monoval_i[2][1];
+                      grad_grads[start + 2 * n_curls + 1][1][1][2] =
+                        -monoval[0][0] * monoval_jplus[1][1] * monoval_i[2][2];
+                      grad_grads[start + 2 * n_curls + 1][1][2][0] =
+                        -monoval[0][1] * monoval_jplus[1][0] * monoval_i[2][2];
+                      grad_grads[start + 2 * n_curls + 1][1][2][1] =
+                        -monoval[0][0] * monoval_jplus[1][1] * monoval_i[2][2];
+                      grad_grads[start + 2 * n_curls + 1][1][2][2] =
+                        -monoval[0][0] * monoval_jplus[1][0] * monoval_i[2][3];
+                      grad_grads[start + 2 * n_curls + 1][2][0][0] =
+                        monoval[0][2] * monoval_j[1][0] * monoval_i[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 2 * n_curls + 1][2][0][1] =
+                        monoval[0][1] * monoval_j[1][1] * monoval_i[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 2 * n_curls + 1][2][0][2] =
+                        monoval[0][1] * monoval_j[1][0] * monoval_i[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 2 * n_curls + 1][2][1][0] =
+                        monoval[0][1] * monoval_j[1][1] * monoval_i[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 2 * n_curls + 1][2][1][1] =
+                        monoval[0][0] * monoval_j[1][2] * monoval_i[2][0] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 2 * n_curls + 1][2][1][2] =
+                        monoval[0][0] * monoval_j[1][1] * monoval_i[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 2 * n_curls + 1][2][2][0] =
+                        monoval[0][1] * monoval_j[1][0] * monoval_i[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 2 * n_curls + 1][2][2][1] =
+                        monoval[0][0] * monoval_j[1][1] * monoval_i[2][1] *
+                        static_cast<double>(j + my_degree + 2);
+                      grad_grads[start + 2 * n_curls + 1][2][2][2] =
+                        monoval[0][0] * monoval_j[1][0] * monoval_i[2][2] *
+                        static_cast<double>(j + my_degree + 2);
+                    }
+                }
+
+              if (j == my_degree)
+                start += 1;
+              else
+                start += 2;
+            }
+        }
+      Assert(start == this->n() - 2 * n_curls, ExcInternalError());
+    }
 }
 
 

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.