From: David Wells Date: Sun, 16 Jul 2017 19:50:21 +0000 (-0400) Subject: Reorganize (some of) the internal grid functions. X-Git-Tag: v9.0.0-rc1~1415^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=992d445834f4d05f019e10708e6525e03454ba11;p=dealii.git Reorganize (some of) the internal grid functions. This is similar to 6310b2dee3 and 52bef21ac2: these functions should be in specific namespaces to avoid clashes if we want to support unity builds. --- diff --git a/source/grid/grid_refinement.cc b/source/grid/grid_refinement.cc index da97212ba4..1fe51cb306 100644 --- a/source/grid/grid_refinement.cc +++ b/source/grid/grid_refinement.cc @@ -35,102 +35,101 @@ DEAL_II_NAMESPACE_OPEN -namespace +namespace internal { - namespace internal + namespace GridRefinement { - template - inline - number - max_element (const Vector &criteria) + namespace { - return *std::max_element(criteria.begin(), criteria.end()); - } + template + inline + number + max_element (const dealii::Vector &criteria) + { + return *std::max_element(criteria.begin(), criteria.end()); + } - template - inline - number - min_element (const Vector &criteria) - { - return *std::min_element(criteria.begin(), criteria.end()); - } + template + inline + number + min_element (const dealii::Vector &criteria) + { + return *std::min_element(criteria.begin(), criteria.end()); + } - // Silence a (bogus) warning in clang-3.6 about the following four - // functions being unused: - DEAL_II_DISABLE_EXTRA_DIAGNOSTICS + // Silence a (bogus) warning in clang-3.6 about the following four + // functions being unused: + DEAL_II_DISABLE_EXTRA_DIAGNOSTICS #ifdef DEAL_II_WITH_TRILINOS - inline - TrilinosScalar - max_element (const TrilinosWrappers::MPI::Vector &criteria) - { - TrilinosScalar m = 0; - criteria.trilinos_vector().MaxValue(&m); - return m; - } + inline + TrilinosScalar + max_element (const dealii::TrilinosWrappers::MPI::Vector &criteria) + { + TrilinosScalar m = 0; + criteria.trilinos_vector().MaxValue(&m); + return m; + } - inline - TrilinosScalar - min_element (const TrilinosWrappers::MPI::Vector &criteria) - { - TrilinosScalar m = 0; - criteria.trilinos_vector().MinValue(&m); - return m; - } + inline + TrilinosScalar + min_element (const dealii::TrilinosWrappers::MPI::Vector &criteria) + { + TrilinosScalar m = 0; + criteria.trilinos_vector().MinValue(&m); + return m; + } #endif - DEAL_II_ENABLE_EXTRA_DIAGNOSTICS + DEAL_II_ENABLE_EXTRA_DIAGNOSTICS - } /* namespace internal */ - - - template - typename constraint_and_return_value::value, - typename VectorType::value_type>::type - min_element (const VectorType &criteria) - { - return internal::min_element (criteria); - } - - - template - typename constraint_and_return_value::value, - typename VectorType::value_type>::type - max_element (const VectorType &criteria) - { - return internal::max_element (criteria); - } + template + typename constraint_and_return_value::value, + typename VectorType::value_type>::type + min_element (const VectorType &criteria) + { + return min_element (criteria); + } - template - typename constraint_and_return_value::value, - typename VectorType::value_type>::type - min_element (const VectorType &criteria) - { - typename VectorType::value_type t = internal::min_element(criteria.block(0)); - for (unsigned int b=1; b + typename constraint_and_return_value::value, + typename VectorType::value_type>::type + max_element (const VectorType &criteria) + { + return max_element (criteria); + } - return t; - } + template + typename constraint_and_return_value::value, + typename VectorType::value_type>::type + min_element (const VectorType &criteria) + { + typename VectorType::value_type t = min_element(criteria.block(0)); + for (unsigned int b=1; b - typename constraint_and_return_value::value, - typename VectorType::value_type>::type - max_element (const VectorType &criteria) - { - typename VectorType::value_type t = internal::max_element(criteria.block(0)); - for (unsigned int b=1; b + typename constraint_and_return_value::value, + typename VectorType::value_type>::type + max_element (const VectorType &criteria) + { + typename VectorType::value_type t = max_element(criteria.block(0)); + for (unsigned int b=1; b @@ -439,7 +438,7 @@ GridRefinement::refine_and_coarsen_fixed_fraction (Triangulation & // threshold if it equals the // largest indicator and the // top_fraction!=1 - if ((top_threshold == max_element(criteria)) && + if ((top_threshold == internal::GridRefinement::max_element(criteria)) && (top_fraction != 1)) top_threshold *= 0.999; @@ -447,10 +446,10 @@ GridRefinement::refine_and_coarsen_fixed_fraction (Triangulation & bottom_threshold = 0.999*top_threshold; // actually flag cells - if (top_threshold < max_element(criteria)) + if (top_threshold < internal::GridRefinement::max_element(criteria)) refine (tria, criteria, top_threshold, pp - tmp.begin()); - if (bottom_threshold > min_element(criteria)) + if (bottom_threshold > internal::GridRefinement::min_element(criteria)) coarsen (tria, criteria, bottom_threshold); } diff --git a/source/grid/manifold_lib.cc b/source/grid/manifold_lib.cc index 6a806d0cda..6b1ee08145 100644 --- a/source/grid/manifold_lib.cc +++ b/source/grid/manifold_lib.cc @@ -300,37 +300,43 @@ get_new_point (const std::vector > &vertices, // CylindricalManifold // ============================================================ -namespace +namespace internal { - // helper function to compute a vector orthogonal to a given one. - template - Point - compute_normal(const Tensor<1,spacedim> &vector) + namespace CylindricalManifold { - AssertThrow(vector.norm() != 0., - ExcMessage("The direction parameter must not be zero!")); - Point<3> normal; - if (std::abs(vector[0]) >= std::abs(vector[1]) - && std::abs(vector[0]) >= std::abs(vector[2])) - { - normal[1]=-1.; - normal[2]=-1.; - normal[0]=(vector[1]+vector[2])/vector[0]; - } - else if (std::abs(vector[1]) >= std::abs(vector[0]) - && std::abs(vector[1]) >= std::abs(vector[2])) - { - normal[0]=-1.; - normal[2]=-1.; - normal[1]=(vector[0]+vector[2])/vector[1]; - } - else + namespace + { + // helper function to compute a vector orthogonal to a given one. + template + Point + compute_normal(const Tensor<1,spacedim> &vector) { - normal[0]=-1.; - normal[1]=-1.; - normal[2]=(vector[0]+vector[1])/vector[2]; + AssertThrow(vector.norm() != 0., + ExcMessage("The direction parameter must not be zero!")); + Point<3> normal; + if (std::abs(vector[0]) >= std::abs(vector[1]) + && std::abs(vector[0]) >= std::abs(vector[2])) + { + normal[1]=-1.; + normal[2]=-1.; + normal[0]=(vector[1]+vector[2])/vector[0]; + } + else if (std::abs(vector[1]) >= std::abs(vector[0]) + && std::abs(vector[1]) >= std::abs(vector[2])) + { + normal[0]=-1.; + normal[2]=-1.; + normal[1]=(vector[0]+vector[2])/vector[1]; + } + else + { + normal[0]=-1.; + normal[1]=-1.; + normal[2]=(vector[0]+vector[1])/vector[2]; + } + return normal; } - return normal; + } } } @@ -351,7 +357,7 @@ CylindricalManifold::CylindricalManifold(const Point &d const Point &point_on_axis_, const double tolerance) : ChartManifold(Tensor<1,3>({0,2.*numbers::PI,0})), - normal_direction(compute_normal(direction_)), + normal_direction(internal::CylindricalManifold::compute_normal(direction_)), direction (direction_/direction_.norm()), point_on_axis (point_on_axis_), tolerance(tolerance)