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 "FVPorousFlowFluidMass.h"
      11             : 
      12             : registerMooseObject("PorousFlowApp", FVPorousFlowFluidMass);
      13             : 
      14             : InputParameters
      15         927 : FVPorousFlowFluidMass::validParams()
      16             : {
      17         927 :   InputParameters params = ElementIntegralPostprocessor::validParams();
      18        1854 :   params.addParam<unsigned int>(
      19             :       "fluid_component",
      20        1854 :       0,
      21             :       "The index corresponding to the fluid component that this Postprocessor acts on");
      22        1854 :   params.addRequiredParam<UserObjectName>(
      23             :       "PorousFlowDictator", "The UserObject that holds the list of PorousFlow variable names.");
      24        1854 :   params.addParam<std::vector<unsigned int>>("phase",
      25             :                                              {},
      26             :                                              "The index of the fluid phase that this "
      27             :                                              "Postprocessor is restricted to.  Multiple "
      28             :                                              "indices can be entered");
      29        2781 :   params.addRangeCheckedParam<Real>("saturation_threshold",
      30        1854 :                                     1.0,
      31             :                                     "saturation_threshold >= 0 & saturation_threshold <= 1",
      32             :                                     "The saturation threshold below which the mass is calculated "
      33             :                                     "for a specific phase. Default is 1.0. Note: only one "
      34             :                                     "phase_index can be entered");
      35        1854 :   params.addParam<std::string>(
      36             :       "base_name",
      37             :       "For non-mechanically-coupled systems with no TensorMechanics strain calculators, base_name "
      38             :       "need not be set.  For mechanically-coupled systems, base_name should be the same base_name "
      39             :       "as given to the TensorMechanics object that computes strain, so that this Postprocessor can "
      40             :       "correctly account for changes in mesh volume.  For non-mechanically-coupled systems, "
      41             :       "base_name should not be the base_name of any TensorMechanics strain calculators.");
      42         927 :   params.set<bool>("use_displaced_mesh") = false;
      43         927 :   params.suppressParameter<bool>("use_displaced_mesh");
      44         927 :   params.addClassDescription("Calculates the mass of a fluid component in a region");
      45         927 :   return params;
      46           0 : }
      47             : 
      48         498 : FVPorousFlowFluidMass::FVPorousFlowFluidMass(const InputParameters & parameters)
      49             :   : ElementIntegralPostprocessor(parameters),
      50         498 :     _dictator(getUserObject<PorousFlowDictator>("PorousFlowDictator")),
      51         996 :     _fluid_component(getParam<unsigned int>("fluid_component")),
      52         996 :     _phase_index(getParam<std::vector<unsigned int>>("phase")),
      53        1000 :     _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : ""),
      54         498 :     _has_total_strain(hasMaterialProperty<RankTwoTensor>(_base_name + "total_strain")),
      55         996 :     _total_strain(_has_total_strain
      56         498 :                       ? &getMaterialProperty<RankTwoTensor>(_base_name + "total_strain")
      57             :                       : nullptr),
      58         996 :     _porosity(getADMaterialProperty<Real>("PorousFlow_porosity_qp")),
      59         996 :     _fluid_density(getADMaterialProperty<std::vector<Real>>("PorousFlow_fluid_phase_density_qp")),
      60         996 :     _fluid_saturation(getADMaterialProperty<std::vector<Real>>("PorousFlow_saturation_qp")),
      61         498 :     _mass_fraction(
      62         498 :         getADMaterialProperty<std::vector<std::vector<Real>>>("PorousFlow_mass_frac_qp")),
      63        1494 :     _saturation_threshold(getParam<Real>("saturation_threshold"))
      64             : {
      65         498 :   const unsigned int num_phases = _dictator.numPhases();
      66         498 :   const unsigned int num_components = _dictator.numComponents();
      67             : 
      68             :   // Check that the number of components entered is not greater than the total number of components
      69         498 :   if (_fluid_component >= num_components)
      70           0 :     paramError(
      71             :         "fluid_component",
      72             :         "The Dictator proclaims that the number of components in this simulation is ",
      73             :         num_components,
      74             :         " whereas you have used a component index of ",
      75             :         _fluid_component,
      76             :         ". Remember that indexing starts at 0. The Dictator does not take such mistakes lightly.");
      77             : 
      78             :   // Check that the number of phases entered is not more than the total possible phases
      79         498 :   if (_phase_index.size() > num_phases)
      80           0 :     paramError("phase",
      81             :                "The Dictator decrees that the number of phases in this simulation is ",
      82             :                num_phases,
      83             :                " but you have entered ",
      84             :                _phase_index.size(),
      85             :                " phases.");
      86             : 
      87             :   // Also check that the phase indices entered are not greater than the number of phases
      88             :   // to avoid a segfault. Note that the input parser takes care of negative inputs so we
      89             :   // don't need to guard against them
      90         498 :   if (!_phase_index.empty())
      91             :   {
      92         484 :     unsigned int max_phase_num = *std::max_element(_phase_index.begin(), _phase_index.end());
      93         484 :     if (max_phase_num > num_phases - 1)
      94           0 :       paramError("phase",
      95             :                  "The Dictator proclaims that the phase index ",
      96             :                  max_phase_num,
      97             :                  " is greater than the largest phase index possible, which is ",
      98             :                  num_phases - 1);
      99             :   }
     100             : 
     101             :   // Using saturation_threshold only makes sense for a specific phase_index
     102         498 :   if (_saturation_threshold < 1.0 && _phase_index.size() != 1)
     103           0 :     paramError("saturation_threshold",
     104             :                "A single phase_index must be entered when prescribing a saturation_threshold");
     105             : 
     106             :   // If _phase_index is empty, create vector of all phase numbers to calculate mass over all phases
     107         498 :   if (_phase_index.empty())
     108          28 :     for (unsigned int i = 0; i < num_phases; ++i)
     109          14 :       _phase_index.push_back(i);
     110             : 
     111             :   // Error if a strain base_name is provided but doesn't exist
     112         996 :   if (parameters.isParamSetByUser("base_name") && !_has_total_strain)
     113           2 :     paramError("base_name", "A strain base_name ", _base_name, " does not exist");
     114         496 : }
     115             : 
     116             : Real
     117       88683 : FVPorousFlowFluidMass::computeQpIntegral()
     118             : {
     119             :   // The fluid mass for the finite volume case is much simpler
     120             :   Real mass = 0.0;
     121       88683 :   const Real strain = (_has_total_strain ? (*_total_strain)[_qp].trace() : 0.0);
     122             : 
     123      221542 :   for (auto ph : _phase_index)
     124             :   {
     125      132859 :     if (MetaPhysicL::raw_value(_fluid_saturation[_qp][ph]) <= _saturation_threshold)
     126      265718 :       mass += MetaPhysicL::raw_value(_fluid_density[_qp][ph] * _fluid_saturation[_qp][ph] *
     127      132859 :                                      _mass_fraction[_qp][ph][_fluid_component]);
     128             :   }
     129             : 
     130       88683 :   return MetaPhysicL::raw_value(_porosity[_qp]) * (1.0 + strain) * mass;
     131             : }