From b24f0fdd723a717446409c1aed6c69767da4936d Mon Sep 17 00:00:00 2001 From: David Wells Date: Thu, 18 Jan 2018 16:42:00 -0500 Subject: [PATCH] Make the TFI chart point calculation more robust. If the guess computed by the shortcut is not accurate enough then we should try to use the affine approximation instead. --- source/grid/manifold_lib.cc | 16 +++++++++++++++- 1 file changed, 15 insertions(+), 1 deletion(-) diff --git a/source/grid/manifold_lib.cc b/source/grid/manifold_lib.cc index 43989a311a..7a0bf39e92 100644 --- a/source/grid/manifold_lib.cc +++ b/source/grid/manifold_lib.cc @@ -28,6 +28,12 @@ DEAL_II_NAMESPACE_OPEN namespace internal { + // The pull_back function fails regularly in the compute_chart_points + // method, and, instead of throwing an exception, returns a point outside + // the unit cell. The individual coordinates of that point are given by the + // value below. + static constexpr double invalid_pull_back_coordinate = 20.0; + // Rotate a given unit vector u around the axis dir // where the angle is given by the length of dir. // This is the exponential map for a sphere. @@ -1600,7 +1606,7 @@ TransfiniteInterpolationManifold { Point outside; for (unsigned int d=0; d chart_point = GeometryInfo::project_to_unit_cell(initial_guess); @@ -1887,7 +1893,15 @@ TransfiniteInterpolationManifold Assert(false, ExcInternalError()); } + // This guess should be pretty good, but if the pull_back fails + // then try again with affine approximation (which is more + // expensive) chart_points[i] = pull_back(cell, surrounding_points[i], guess); + if (chart_points[i][0] == internal::invalid_pull_back_coordinate) + { + chart_points[i] = pull_back(cell, surrounding_points[i], + cell->real_to_unit_cell_affine_approximation(surrounding_points[i])); + } } else chart_points[i] = pull_back(cell, surrounding_points[i], -- 2.39.5