]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix compiler error when using nvcc with C++20
authorBruno Turcksin <bruno.turcksin@gmail.com>
Thu, 6 Feb 2025 15:10:20 +0000 (15:10 +0000)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Thu, 6 Feb 2025 15:10:20 +0000 (15:10 +0000)
include/deal.II/algorithms/theta_timestepping.templates.h
include/deal.II/base/tensor.h
include/deal.II/grid/reference_cell.h
source/fe/fe_system.cc
source/fe/mapping_q_cache.cc

index 2ec44e5208b67425f8e176d82a79a6e41a1b8025..1f55692a643a1163ad065073a6f662dcb864fca0 100644 (file)
@@ -118,7 +118,11 @@ namespace Algorithms
     if (output != nullptr)
       (*output) << 0U << out;
 
-    for (unsigned int count = 1; d_explicit.time < control.final(); ++count)
+    // When using nvcc 12.6 with C++20, the compiler has trouble determining the
+    // type of d_explicit.time. We need to use a static_cast to help it.
+    for (unsigned int count = 1;
+         static_cast<double>(d_explicit.time) < control.final();
+         ++count)
       {
         const bool step_change = control.advance();
         d_implicit.time        = control.now();
index f3888740cef6bed00040da8eeaa15013f08b3e34..28c772abacf7746514926f8ffc34cecd3d35941f 100644 (file)
@@ -1253,7 +1253,7 @@ Tensor<rank_, dim, Number>::Tensor(const ArrayLike &initializer,
 }
 
 
-#  ifdef DEAL_II_HAVE_CXX20
+#  if defined(DEAL_II_HAVE_CXX20) && !defined(__NVCC__)
 
 template <int rank_, int dim, typename Number>
 constexpr DEAL_II_HOST_DEVICE_ALWAYS_INLINE
index f179143b9eae15e9b77032d8931c5a869bf04199..4b0a7a357bb6f599de394bbf19fe269b6d5ce63e 100644 (file)
@@ -3540,7 +3540,12 @@ ReferenceCell::get_combined_orientation(
           return o;
       }
 
-    Assert(false, (internal::NoPermutation<T>(*this, vertices_0, vertices_1)));
+    // Do not use this in Assert because nvcc when using C++20 assumes that this
+    // is an integer and we get the following error: invalid type argument of
+    // unary '*' (have 'int')
+    [[maybe_unused]] const auto &ref_cell = *this;
+    Assert(false,
+           (internal::NoPermutation<T>(ref_cell, vertices_0, vertices_1)));
     return std::numeric_limits<types::geometric_orientation>::max();
   };
 
index 5c29ee4dcdf50efaf66b608dfc02c74ff16632f1..04d982d2455f915e2b291414448b9b7f4b50aa1e 100644 (file)
@@ -1799,7 +1799,12 @@ FESystem<dim, spacedim>::initialize(
       {
         const unsigned int base = this->system_to_base_table[i].first.first,
                            base_index = this->system_to_base_table[i].second;
-        Assert(base < this->n_base_elements(), ExcInternalError());
+        // Do not use this in Assert because nvcc when using C++20 assumes that
+        // this is an integer and we get the following error: base operand of
+        // '->' is not a pointer
+        [[maybe_unused]] const unsigned int n_base_elements =
+          this->n_base_elements();
+        Assert(base < n_base_elements, ExcInternalError());
         Assert(base_index < base_element(base).unit_support_points.size(),
                ExcInternalError());
         this->unit_support_points[i] =
@@ -1957,10 +1962,18 @@ FESystem<dim, spacedim>::initialize(
         {
           // the array into which we want to write should have the correct size
           // already.
-          Assert(this->adjust_quad_dof_index_for_face_orientation_table[face_no]
-                     .n_elements() ==
-                   this->reference_cell().n_face_orientations(face_no) *
-                     this->n_dofs_per_quad(face_no),
+          // Do not use this in Assert because nvcc when using C++20 assumes
+          // that this is an integer and we get the following error: base
+          // operand of
+          // '->' is not a pointer
+          [[maybe_unused]] const unsigned int n_elements =
+            this->adjust_quad_dof_index_for_face_orientation_table[face_no]
+              .n_elements();
+          [[maybe_unused]] const unsigned int n_face_orientations =
+            this->reference_cell().n_face_orientations(face_no);
+          [[maybe_unused]] const unsigned int n_dofs_per_quad =
+            this->n_dofs_per_quad(face_no);
+          Assert(n_elements == n_face_orientations * n_dofs_per_quad,
                  ExcInternalError());
 
           // to obtain the shifts for this composed element, copy the shift
@@ -1983,13 +1996,19 @@ FESystem<dim, spacedim>::initialize(
                   index += temp.size(0);
                 }
             }
-          Assert(index == this->n_dofs_per_quad(face_no), ExcInternalError());
+          // error: base operand of '->' is not a pointer
+          Assert(index == n_dofs_per_quad, ExcInternalError());
         }
 
       // additionally compose the permutation information for lines
-      Assert(this->adjust_line_dof_index_for_line_orientation_table.size() ==
-               this->n_dofs_per_line(),
-             ExcInternalError());
+      // Do not use this in Assert because nvcc when using C++20 assumes that
+      // this is an integer and we get the following error: base operand of
+      // '->' is not a pointer
+      [[maybe_unused]] const unsigned int table_size =
+        this->adjust_line_dof_index_for_line_orientation_table.size();
+      [[maybe_unused]] const unsigned int n_dofs_per_line =
+        this->n_dofs_per_line();
+      Assert(table_size == n_dofs_per_line, ExcInternalError());
       unsigned int index = 0;
       for (unsigned int b = 0; b < this->n_base_elements(); ++b)
         {
@@ -2006,7 +2025,7 @@ FESystem<dim, spacedim>::initialize(
               index += temp2.size();
             }
         }
-      Assert(index == this->n_dofs_per_line(), ExcInternalError());
+      Assert(index == n_dofs_per_line, ExcInternalError());
     });
 
 
index 16dd292b50c6de2001fc5a72ec17e24b9cb0a4f8..b6ed5e1271bc75b2ab11f4becc3f4f87bba19d80 100644 (file)
@@ -159,9 +159,13 @@ MappingQCache<dim, spacedim>::initialize(
         void *) {
       (*support_point_cache)[cell->level()][cell->index()] =
         compute_points_on_cell(cell);
+      // Do not use this in Assert because nvcc when using C++20 assumes that
+      // this is an integer and we get the following error: invalid type
+      // argument of unary '*' (have 'int')
+      [[maybe_unused]] const unsigned int d = this->get_degree() + 1;
       AssertDimension(
         (*support_point_cache)[cell->level()][cell->index()].size(),
-        Utilities::pow(this->get_degree() + 1, dim));
+        Utilities::pow(d, dim));
     },
     /* copier */ std::function<void(void *)>(),
     /* scratch_data */ nullptr,

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.