From: Rene Gassmoeller Date: Mon, 6 Nov 2017 17:19:13 +0000 (-0700) Subject: Address some comments X-Git-Tag: v9.0.0-rc1~811^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=306c971489cf075355c0cd896e58ddb3af8454ce;p=dealii.git Address some comments --- diff --git a/source/grid/manifold.cc b/source/grid/manifold.cc index 54168ea40f..69e70fd527 100644 --- a/source/grid/manifold.cc +++ b/source/grid/manifold.cc @@ -156,16 +156,22 @@ normal_vector (const Triangulation<3, 3>::face_iterator &face, const Point<3> &p) const { const int spacedim=3; + + // value of i used to compute t1 + unsigned int i_t1 = numbers::invalid_unsigned_int; Tensor<1,spacedim> t1,t2; Tensor<1,spacedim> normal; + double normal_direction = numbers::signaling_nan(); // Look for a combination of tangent vectors that // are of approximately equal length and not linearly dependent - for (unsigned int i=0,j=1 ; i<4 && j<4; ++j) + for (unsigned int i=0,j=1 ; + i::vertices_per_face && j::vertices_per_face; + ++j) { // if p is too close to vertex i try again with different i and j if ((p - face->vertex(i)).norm_square() < - std::numeric_limits::epsilon() * (p - face->vertex(j)).norm_square()) + 1e4 * std::numeric_limits::epsilon() * (p - face->vertex(j)).norm_square()) { ++i; continue; @@ -173,23 +179,29 @@ normal_vector (const Triangulation<3, 3>::face_iterator &face, // if p is too close to vertex j try again with different j if ((p - face->vertex(j)).norm_square() < - std::numeric_limits::epsilon() * (p - face->vertex(i)).norm_square()) + 1e4 * std::numeric_limits::epsilon() * (p - face->vertex(i)).norm_square()) continue; - t1 = get_tangent_vector(p, face->vertex(i)); - t2 = get_tangent_vector(p, face->vertex(j)); + // i might not have changed, and get_tangent_vector is potentially expensive. + // Only compute a new t1 if necessary + if (i != i_t1) + { + i_t1 = i; + t1 = get_tangent_vector(p, face->vertex(i)); + } + t2 = get_tangent_vector(p, face->vertex(j)); normal = cross_product_3d(t1,t2); // if t1 and t2 are (nearly) linearly dependent try again with different j / t2 - if (normal.norm_square() < std::numeric_limits::epsilon() * + if (normal.norm_square() < 1e4 * std::numeric_limits::epsilon() * t1.norm_square() * t2.norm_square()) continue; break; } - Assert(normal.norm_square() >= std::numeric_limits::epsilon() * + Assert(normal.norm_square() >= 1e4 * std::numeric_limits::epsilon() * t1.norm_square() * t2.norm_square(), ExcMessage("Manifold::normal_vector was unable to find a suitable combination " "of vertices to compute a normal on this face. Check for distorted " diff --git a/tests/manifold/flat_manifold_08.cc b/tests/manifold/flat_manifold_08.cc index 34a750e7fb..5072807204 100644 --- a/tests/manifold/flat_manifold_08.cc +++ b/tests/manifold/flat_manifold_08.cc @@ -1,15 +1,16 @@ -//---------------------------- manifold_id_01.cc --------------------------- -// Copyright (C) 2011 - 2015 by the mathLab team. +//------------------------------------------------------------------- +// Copyright (C) 2017 by the deal.II authors. // // This file is subject to LGPL and may not be distributed // without copyright and license information. Please refer // to the file deal.II/doc/license.html for the text and // further information on this license. // -//---------------------------- flat_manifold_01.cc --------------------------- +//------------------------------------------------------------------- -// Test that the flat manifold does what it should +// Test that Manifold::get_normal_vector also works for a special case that +// used to crash with another algorithm. #include "../tests.h"