From 9dca30e4abd8049196c18e59fe25db1fd0153f24 Mon Sep 17 00:00:00 2001 From: Nicolas Aunai Date: Sun, 20 Oct 2024 12:53:45 +0200 Subject: [PATCH] safer nu --- src/core/numerics/ohm/ohm.hpp | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/src/core/numerics/ohm/ohm.hpp b/src/core/numerics/ohm/ohm.hpp index 070364320..27d265e6d 100644 --- a/src/core/numerics/ohm/ohm.hpp +++ b/src/core/numerics/ohm/ohm.hpp @@ -360,7 +360,8 @@ class Ohm : public LayoutHolder auto spatial_hyperresistive_(VecField const& J, VecField const& B, Field const& n, MeshIndex index) const { // TODO : https://github.com/PHAREHUB/PHARE/issues/3 - auto const lvlCoeff = 1. / std::pow(4, layout_->levelNumber()); + auto const lvlCoeff = 1. / std::pow(4, layout_->levelNumber()); + auto constexpr min_density = 0.1; auto computeHR = [&](auto BxProj, auto ByProj, auto BzProj, auto nProj) { auto const BxOnE = GridLayout::project(B(Component::X), index, BxProj); @@ -368,7 +369,8 @@ class Ohm : public LayoutHolder auto const BzOnE = GridLayout::project(B(Component::Z), index, BzProj); auto const nOnE = GridLayout::project(n, index, nProj); auto b = std::sqrt(BxOnE * BxOnE + ByOnE * ByOnE + BzOnE * BzOnE); - return -nu_ * b / nOnE * lvlCoeff * layout_->laplacian(J(component), index); + return -nu_ * (b / (nOnE + min_density) + 1) * lvlCoeff + * layout_->laplacian(J(component), index); }; if constexpr (component == Component::X)