From 256ce343566352467221f7bbf26963870e829d3c Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Sat, 31 Oct 2020 19:10:06 +0100 Subject: [PATCH] Transfinite interpolation: Fix refresh of Jacobian --- source/grid/manifold_lib.cc | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) 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); -- 2.39.5