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 "InternalSideIndicatorBase.h"
      11             : 
      12             : // MOOSE includes
      13             : #include "Assembly.h"
      14             : #include "MooseTypes.h"
      15             : #include "MooseVariableFE.h"
      16             : #include "Problem.h"
      17             : #include "SubProblem.h"
      18             : #include "SystemBase.h"
      19             : 
      20             : #include "libmesh/dof_map.h"
      21             : #include "libmesh/dense_vector.h"
      22             : #include "libmesh/numeric_vector.h"
      23             : #include "libmesh/dense_subvector.h"
      24             : #include "libmesh/libmesh_common.h"
      25             : #include "libmesh/quadrature.h"
      26             : 
      27             : const BoundaryID InternalSideIndicatorBase::InternalBndId = 12345;
      28             : 
      29             : InputParameters
      30       58197 : InternalSideIndicatorBase::validParams()
      31             : {
      32       58197 :   InputParameters params = Indicator::validParams();
      33       58197 :   params.addRequiredParam<VariableName>(
      34             :       "variable", "The name of the variable that this side indicator applies to");
      35      174591 :   params.addParam<bool>("scale_by_flux_faces",
      36      116394 :                         false,
      37             :                         "Whether or not to scale the error values by "
      38             :                         "the number of flux faces.  This attempts to "
      39             :                         "not penalize elements on boundaries for "
      40             :                         "having less neighbors.");
      41             : 
      42       58197 :   params.addPrivateParam<BoundaryID>("_boundary_id", InternalSideIndicatorBase::InternalBndId);
      43       58197 :   return params;
      44           0 : }
      45             : 
      46         595 : InternalSideIndicatorBase::InternalSideIndicatorBase(const InputParameters & parameters)
      47             :   : Indicator(parameters),
      48             :     NeighborCoupleable(this, false, false),
      49             :     ScalarCoupleable(this),
      50         595 :     _field_var(_subproblem.getStandardVariable(_tid, name())),
      51         595 :     _current_elem(_assembly.elem()),
      52         595 :     _neighbor_elem(_assembly.neighbor()),
      53         595 :     _current_side(_assembly.side()),
      54         595 :     _current_side_elem(_assembly.sideElem()),
      55         595 :     _coord_sys(_assembly.coordSystem()),
      56         595 :     _qrule(_assembly.qRuleFace()),
      57         595 :     _q_point(_assembly.qPointsFace()),
      58         595 :     _JxW(_assembly.JxWFace()),
      59         595 :     _coord(_assembly.coordTransformation()),
      60         595 :     _boundary_id(parameters.get<BoundaryID>("_boundary_id")),
      61         595 :     _scale_by_flux_faces(parameters.get<bool>("scale_by_flux_faces")),
      62        1190 :     _normals(_assembly.normals())
      63             : {
      64         595 :   const std::vector<MooseVariableFieldBase *> & coupled_vars = getCoupledMooseVars();
      65         595 :   for (const auto & var : coupled_vars)
      66           0 :     addMooseVariableDependency(var);
      67             : 
      68             :   // Not supported with the base linear lagrange case
      69         595 :   if (_use_displaced_mesh)
      70           0 :     paramError("use_displaced_mesh",
      71             :                "Internal side indicators do not support using the displaced mesh at this time. "
      72             :                "They can be used on the undisplaced mesh in a Problem with displaced mesh");
      73             :   // Access into the solution vector assumes constant monomial
      74         595 :   if (_field_var.feType() != libMesh::FEType(CONSTANT, MONOMIAL))
      75           0 :     mooseError("Only constant monomial variables for the internal side indicator are supported");
      76         595 : }
      77             : 
      78             : void
      79     1570220 : InternalSideIndicatorBase::computeIndicator()
      80             : {
      81             :   // Derived class decides the contribution at each qp to the indicator value
      82     1570220 :   Real sum = 0;
      83     4438666 :   for (_qp = 0; _qp < _qrule->n_points(); _qp++)
      84     2868446 :     sum += _JxW[_qp] * _coord[_qp] * computeQpIntegral();
      85             : 
      86             :   // Contribution is added to the two elements
      87             :   {
      88     1570220 :     const auto sys_num = _field_var.sys().number();
      89     1570220 :     const auto var_num = _field_var.number();
      90     1570220 :     _solution.add(_current_elem->dof_number(sys_num, var_num, 0), sum * _current_elem->hmax());
      91     1570220 :     if (_field_var.hasBlocks(_neighbor_elem->subdomain_id()))
      92     1569861 :       _solution.add(_neighbor_elem->dof_number(sys_num, var_num, 0), sum * _neighbor_elem->hmax());
      93             :   }
      94     1570220 : }
      95             : 
      96             : void
      97      783545 : InternalSideIndicatorBase::finalize()
      98             : {
      99      783545 :   unsigned int n_flux_faces = 0;
     100             : 
     101      783545 :   if (_scale_by_flux_faces)
     102             :   {
     103        5448 :     if (isVarFV())
     104           0 :       paramError("scale_by_flux_faces", "Unsupported at this time for finite volume variables");
     105             : 
     106             :     // Figure out the total number of sides contributing to the error.
     107             :     // We'll scale by this so boundary elements are less penalized
     108       28536 :     for (unsigned int side = 0; side < _current_elem->n_sides(); side++)
     109       23088 :       if (_current_elem->neighbor_ptr(side) != nullptr)
     110       20328 :         n_flux_faces++;
     111             :   }
     112             :   else
     113      778097 :     n_flux_faces = 1;
     114             : 
     115             :   // The 0 is because CONSTANT MONOMIALS only have one coefficient per element
     116      783545 :   Real value = _field_var.dofValues()[0];
     117             : 
     118             :   {
     119      783545 :     const auto sys_num = _field_var.sys().number();
     120      783545 :     const auto var_num = _field_var.number();
     121             : 
     122      783545 :     _solution.set(_current_elem->dof_number(sys_num, var_num, 0),
     123      783545 :                   std::sqrt(value) / static_cast<Real>(n_flux_faces));
     124             :   }
     125      783545 : }