From b6941fe52e07f568ac4191dcc13deb4ee48b9da2 Mon Sep 17 00:00:00 2001 From: Bruno Turcksin Date: Thu, 6 Feb 2025 15:10:20 +0000 Subject: [PATCH] Fix compiler error when using nvcc with C++20 --- .../algorithms/theta_timestepping.templates.h | 6 ++- include/deal.II/base/tensor.h | 2 +- include/deal.II/grid/reference_cell.h | 7 +++- source/fe/fe_system.cc | 39 ++++++++++++++----- source/fe/mapping_q_cache.cc | 6 ++- 5 files changed, 46 insertions(+), 14 deletions(-) diff --git a/include/deal.II/algorithms/theta_timestepping.templates.h b/include/deal.II/algorithms/theta_timestepping.templates.h index 2ec44e5208..1f55692a64 100644 --- a/include/deal.II/algorithms/theta_timestepping.templates.h +++ b/include/deal.II/algorithms/theta_timestepping.templates.h @@ -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(d_explicit.time) < control.final(); + ++count) { const bool step_change = control.advance(); d_implicit.time = control.now(); diff --git a/include/deal.II/base/tensor.h b/include/deal.II/base/tensor.h index f3888740ce..28c772abac 100644 --- a/include/deal.II/base/tensor.h +++ b/include/deal.II/base/tensor.h @@ -1253,7 +1253,7 @@ Tensor::Tensor(const ArrayLike &initializer, } -# ifdef DEAL_II_HAVE_CXX20 +# if defined(DEAL_II_HAVE_CXX20) && !defined(__NVCC__) template constexpr DEAL_II_HOST_DEVICE_ALWAYS_INLINE diff --git a/include/deal.II/grid/reference_cell.h b/include/deal.II/grid/reference_cell.h index f179143b9e..4b0a7a357b 100644 --- a/include/deal.II/grid/reference_cell.h +++ b/include/deal.II/grid/reference_cell.h @@ -3540,7 +3540,12 @@ ReferenceCell::get_combined_orientation( return o; } - Assert(false, (internal::NoPermutation(*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(ref_cell, vertices_0, vertices_1))); return std::numeric_limits::max(); }; diff --git a/source/fe/fe_system.cc b/source/fe/fe_system.cc index 5c29ee4dcd..04d982d245 100644 --- a/source/fe/fe_system.cc +++ b/source/fe/fe_system.cc @@ -1799,7 +1799,12 @@ FESystem::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::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::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::initialize( index += temp2.size(); } } - Assert(index == this->n_dofs_per_line(), ExcInternalError()); + Assert(index == n_dofs_per_line, ExcInternalError()); }); diff --git a/source/fe/mapping_q_cache.cc b/source/fe/mapping_q_cache.cc index 16dd292b50..b6ed5e1271 100644 --- a/source/fe/mapping_q_cache.cc +++ b/source/fe/mapping_q_cache.cc @@ -159,9 +159,13 @@ MappingQCache::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(), /* scratch_data */ nullptr, -- 2.39.5