From: Martin Kronbichler Date: Sat, 31 Oct 2020 18:10:06 +0000 (+0100) Subject: Transfinite interpolation: Fix refresh of Jacobian X-Git-Tag: v9.3.0-rc1~946^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F11126%2Fhead;p=dealii.git Transfinite interpolation: Fix refresh of Jacobian --- diff --git a/source/grid/manifold_lib.cc b/source/grid/manifold_lib.cc index a2110f226b..d1f0813f34 100644 --- a/source/grid/manifold_lib.cc +++ b/source/grid/manifold_lib.cc @@ -2154,13 +2154,14 @@ TransfiniteInterpolationManifold::pull_back( const Tensor<1, spacedim> old_residual = residual; while (alpha > 1e-4) { - Point guess = chart_point + alpha * update; - residual = + Point guess = chart_point + alpha * update; + const Tensor<1, spacedim> residual_guess = point - compute_transfinite_interpolation( *cell, guess, coarse_cell_is_flat[cell->index()]); - const double residual_norm_new = residual.norm_square(); + const double residual_norm_new = residual_guess.norm_square(); if (residual_norm_new < residual_norm_square) { + residual = residual_guess; residual_norm_square = residual_norm_new; chart_point += alpha * update; break; @@ -2198,7 +2199,7 @@ TransfiniteInterpolationManifold::pull_back( // prevent division by zero. This number should be scale-invariant // because Jinv_deltaf carries no units and x is in reference // coordinates. - if (std::abs(delta_x * Jinv_deltaf) > 1e-12) + if (std::abs(delta_x * Jinv_deltaf) > 1e-12 && !must_recompute_jacobian) { const Tensor<1, dim> factor = (delta_x - Jinv_deltaf) / (delta_x * Jinv_deltaf);