Line data    Source code

       1             : //* This file is part of the MOOSE framework
       2             : //* https://mooseframework.inl.gov
       3             : //*
       4             : //* All rights reserved, see COPYRIGHT for full restrictions
       5             : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
       6             : //*
       7             : //* Licensed under LGPL 2.1, please see LICENSE for details
       8             : //* https://www.gnu.org/licenses/lgpl-2.1.html
       9             : 
      10             : #include "FVAdvectedMinmodWeightBased.h"
      11             : 
      12             : #include "MathFVUtils.h"
      13             : 
      14             : #include <algorithm>
      15             : #include <limits>
      16             : 
      17             : registerMooseObject("MooseApp", FVAdvectedMinmodWeightBased);
      18             : 
      19             : InputParameters
      20        3759 : FVAdvectedMinmodWeightBased::validParams()
      21             : {
      22        3759 :   InputParameters params = FVInterpolationMethod::validParams();
      23        7518 :   params.addClassDescription(
      24             :       "Minmod interpolation for advected quantities implemented as limited blending weights "
      25             :       "(no MUSCL reconstruction, no deferred correction).");
      26       11277 :   params.addParam<bool>(
      27             :       "limit_to_linear",
      28        7518 :       true,
      29             :       "Whether to limit the scheme to be no more downwind-biased than linear interpolation.");
      30       15036 :   params.addRangeCheckedParam<Real>(
      31             :       "blending_factor",
      32        7518 :       1.0,
      33             :       "blending_factor>=0 & blending_factor<=1",
      34             :       "Scales the high-order blending strength; 0 gives pure upwind, 1 gives the full limited "
      35             :       "blending. Values < 1 can improve linear solver robustness for fully implicit assembly.");
      36        3759 :   return params;
      37           0 : }
      38             : 
      39         349 : FVAdvectedMinmodWeightBased::FVAdvectedMinmodWeightBased(const InputParameters & params)
      40             :   : FVInterpolationMethod(params),
      41         349 :     _limit_to_linear(getParam<bool>("limit_to_linear")),
      42        1047 :     _blending_factor(getParam<Real>("blending_factor"))
      43             : {
      44         349 : }
      45             : 
      46             : FVAdvectedInterpolationMethod::AdvectedSystemContribution
      47    10724228 : FVAdvectedMinmodWeightBased::advectedInterpolate(const FaceInfo & face,
      48             :                                                  Real elem_value,
      49             :                                                  Real neighbor_value,
      50             :                                                  const VectorValue<Real> * elem_grad,
      51             :                                                  const VectorValue<Real> * neighbor_grad,
      52             :                                                  Real mass_flux) const
      53             : {
      54             :   mooseAssert(elem_grad && neighbor_grad,
      55             :               "Minmod advected interpolation requires both element and neighbor gradients.");
      56             : 
      57    10724228 :   const bool upwind_is_elem = mass_flux >= 0.0;
      58             : 
      59    10724228 :   const Real phi_upwind = upwind_is_elem ? elem_value : neighbor_value;
      60    10724228 :   const Real phi_downwind = upwind_is_elem ? neighbor_value : elem_value;
      61             : 
      62    10724228 :   const VectorValue<Real> grad_upwind = upwind_is_elem ? *elem_grad : *neighbor_grad;
      63    10724228 :   const Point upwind_to_downwind = upwind_is_elem ? face.dCN() : Point(-face.dCN());
      64             : 
      65    10724228 :   const auto r_f = Moose::FV::rF(phi_upwind, phi_downwind, grad_upwind, upwind_to_downwind);
      66    10724228 :   const Real beta = std::max(Real(0.0), std::min(Real(1.0), r_f));
      67             : 
      68             :   // Geometric weight associated with the upwind cell for this face.
      69    10724228 :   const Real w_f = upwind_is_elem ? face.gC() : (1.0 - face.gC());
      70             : 
      71             :   // Use a Greenshields-style blending:
      72             :   //   phi_f = (1-g)*phi_upwind + g*phi_downwind, where g = beta*(1-w_f)
      73    10724228 :   const Real g_unclamped = _blending_factor * beta * (1.0 - w_f);
      74    10724228 :   const Real g_clamped = std::min(std::max(g_unclamped, 0.0), 1.0 - w_f);
      75             :   // Clamp to [0, 1-w_f] so the weights do not become more downwind-biased than linear.
      76    10724228 :   const Real g = _limit_to_linear ? g_clamped : g_unclamped;
      77             : 
      78    10724228 :   const Real w_upwind = 1.0 - g;
      79    10724228 :   const Real w_downwind = g;
      80             : 
      81             :   // Map (upwind, downwind) weights back to (elem, neighbor) ordering.
      82    10724228 :   const Real w_elem = upwind_is_elem ? w_upwind : w_downwind;
      83    10724228 :   const Real w_neighbor = upwind_is_elem ? w_downwind : w_upwind;
      84             : 
      85    10724228 :   AdvectedSystemContribution result;
      86    10724228 :   result.weights_matrix = std::make_pair(w_elem, w_neighbor);
      87             : 
      88    21448456 :   return result;
      89             : }
      90             : 
      91             : Real
      92           0 : FVAdvectedMinmodWeightBased::advectedInterpolateValue(const FaceInfo & face,
      93             :                                                       Real elem_value,
      94             :                                                       Real neighbor_value,
      95             :                                                       const VectorValue<Real> * elem_grad,
      96             :                                                       const VectorValue<Real> * neighbor_grad,
      97             :                                                       Real mass_flux) const
      98             : {
      99             :   const auto result =
     100           0 :       advectedInterpolate(face, elem_value, neighbor_value, elem_grad, neighbor_grad, mass_flux);
     101           0 :   return result.weights_matrix.first * elem_value + result.weights_matrix.second * neighbor_value;
     102             : }