FVAdvectedMinmodWeightBased.C

Go to the documentation of this file.

//* This file is part of the MOOSE framework
//* https://mooseframework.inl.gov
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html
 
#include "FVAdvectedMinmodWeightBased.h"
 
#include "MathFVUtils.h"
 
#include <algorithm>
#include <limits>
 
registerMooseObject("MooseApp", FVAdvectedMinmodWeightBased);
 
InputParameters
FVAdvectedMinmodWeightBased::validParams()
{
  InputParameters params = FVInterpolationMethod::validParams();
  params.addClassDescription(
      "Minmod interpolation for advected quantities implemented as limited blending weights "
      "(no MUSCL reconstruction, no deferred correction).");
  params.addParam<bool>(
      "limit_to_linear",
      true,
      "Whether to limit the scheme to be no more downwind-biased than linear interpolation.");
  params.addRangeCheckedParam<Real>(
      "blending_factor",
      1.0,
      "blending_factor>=0 & blending_factor<=1",
      "Scales the high-order blending strength; 0 gives pure upwind, 1 gives the full limited "
      "blending. Values < 1 can improve linear solver robustness for fully implicit assembly.");
  return params;
}
 
FVAdvectedMinmodWeightBased::FVAdvectedMinmodWeightBased(const InputParameters & params)
  : FVInterpolationMethod(params),
    _limit_to_linear(getParam<bool>("limit_to_linear")),
    _blending_factor(getParam<Real>("blending_factor"))
{
}
 
FVAdvectedInterpolationMethod::AdvectedSystemContribution
FVAdvectedMinmodWeightBased::advectedInterpolate(const FaceInfo & face,
                                                 Real elem_value,
                                                 Real neighbor_value,
                                                 const VectorValue<Real> * elem_grad,
                                                 const VectorValue<Real> * neighbor_grad,
                                                 Real mass_flux) const
{
  mooseAssert(elem_grad && neighbor_grad,
              "Minmod advected interpolation requires both element and neighbor gradients.");
 
  const bool upwind_is_elem = mass_flux >= 0.0;
 
  const Real phi_upwind = upwind_is_elem ? elem_value : neighbor_value;
  const Real phi_downwind = upwind_is_elem ? neighbor_value : elem_value;
 
  const VectorValue<Real> grad_upwind = upwind_is_elem ? *elem_grad : *neighbor_grad;
  const Point upwind_to_downwind = upwind_is_elem ? face.dCN() : Point(-face.dCN());
 
  const auto r_f = Moose::FV::rF(phi_upwind, phi_downwind, grad_upwind, upwind_to_downwind);
  const Real beta = std::max(Real(0.0), std::min(Real(1.0), r_f));
 
  // Geometric weight associated with the upwind cell for this face.
  const Real w_f = upwind_is_elem ? face.gC() : (1.0 - face.gC());
 
  // Use a Greenshields-style blending:
  //   phi_f = (1-g)*phi_upwind + g*phi_downwind, where g = beta*(1-w_f)
  const Real g_unclamped = _blending_factor * beta * (1.0 - w_f);
  const Real g_clamped = std::min(std::max(g_unclamped, 0.0), 1.0 - w_f);
  // Clamp to [0, 1-w_f] so the weights do not become more downwind-biased than linear.
  const Real g = _limit_to_linear ? g_clamped : g_unclamped;
 
  const Real w_upwind = 1.0 - g;
  const Real w_downwind = g;
 
  // Map (upwind, downwind) weights back to (elem, neighbor) ordering.
  const Real w_elem = upwind_is_elem ? w_upwind : w_downwind;
  const Real w_neighbor = upwind_is_elem ? w_downwind : w_upwind;
 
  AdvectedSystemContribution result;
  result.weights_matrix = std::make_pair(w_elem, w_neighbor);
 
  return result;
}
 
Real
FVAdvectedMinmodWeightBased::advectedInterpolateValue(const FaceInfo & face,
                                                      Real elem_value,
                                                      Real neighbor_value,
                                                      const VectorValue<Real> * elem_grad,
                                                      const VectorValue<Real> * neighbor_grad,
                                                      Real mass_flux) const
{
  const auto result =
      advectedInterpolate(face, elem_value, neighbor_value, elem_grad, neighbor_grad, mass_flux);
  return result.weights_matrix.first * elem_value + result.weights_matrix.second * neighbor_value;
}