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 "DiffusionIPHDGAssemblyHelper.h"
      11             : #include "MooseTypes.h"
      12             : #include "MooseVariableDependencyInterface.h"
      13             : #include "MooseVariableScalar.h"
      14             : #include "SystemBase.h"
      15             : #include "MooseMesh.h"
      16             : #include "MooseObject.h"
      17             : #include "MaterialPropertyInterface.h"
      18             : 
      19             : using namespace libMesh;
      20             : 
      21             : InputParameters
      22       44435 : DiffusionIPHDGAssemblyHelper::validParams()
      23             : {
      24       44435 :   auto params = IPHDGAssemblyHelper::validParams();
      25       44435 :   params.addRequiredParam<MaterialPropertyName>("diffusivity", "The diffusivity");
      26      133305 :   params.addParam<Real>("alpha",
      27       88870 :                         1,
      28             :                         "The stabilization coefficient required for discontinuous Galerkin "
      29             :                         "schemes.");
      30       44435 :   return params;
      31           0 : }
      32             : 
      33         826 : DiffusionIPHDGAssemblyHelper::DiffusionIPHDGAssemblyHelper(
      34             :     const MooseObject * const moose_obj,
      35             :     MooseVariableDependencyInterface * const mvdi,
      36             :     const TransientInterface * const ti,
      37             :     SystemBase & sys,
      38             :     const Assembly & assembly,
      39             :     const THREAD_ID tid,
      40             :     const std::set<SubdomainID> & block_ids,
      41         826 :     const std::set<BoundaryID> & boundary_ids)
      42             :   : IPHDGAssemblyHelper(moose_obj, mvdi, ti, sys, assembly, tid, block_ids, boundary_ids),
      43         826 :     _diff(this->getADMaterialProperty<Real>("diffusivity")),
      44         826 :     _face_diff(this->getFaceADMaterialProperty<Real>("diffusivity")),
      45        1652 :     _alpha(moose_obj->getParam<Real>("alpha"))
      46             : {
      47         826 : }
      48             : 
      49             : void
      50      572452 : DiffusionIPHDGAssemblyHelper::scalarVolume()
      51             : {
      52     2893760 :   for (const auto qp : make_range(_ip_qrule->n_points()))
      53    10202512 :     for (const auto i : index_range(_scalar_re))
      54     7881204 :       _scalar_re(i) += _ip_JxW[qp] * (_grad_scalar_phi[i][qp] * _diff[qp] * _grad_u_sol[qp]);
      55      572452 : }
      56             : 
      57             : void
      58     2059160 : DiffusionIPHDGAssemblyHelper::scalarFace()
      59             : {
      60     2059160 :   const auto h_elem = _elem_volume / _side_area;
      61             : 
      62     9132752 :   for (const auto i : index_range(_scalar_re))
      63    21292596 :     for (const auto qp : make_range(_ip_qrule_face->n_points()))
      64             :     {
      65    28438008 :       _scalar_re(i) -= _ip_JxW_face[qp] * _face_diff[qp] * _scalar_phi_face[i][qp] *
      66    42657012 :                        (_grad_u_sol[qp] * _ip_normals[qp]);
      67    42657012 :       _scalar_re(i) -= _ip_JxW_face[qp] * _alpha / h_elem * _face_diff[qp] *
      68    56876016 :                        (_lm_u_sol[qp] - _u_sol[qp]) * _scalar_phi_face[i][qp];
      69    42657012 :       _scalar_re(i) += _ip_JxW_face[qp] * (_lm_u_sol[qp] - _u_sol[qp]) * _face_diff[qp] *
      70    56876016 :                        _grad_scalar_phi_face[i][qp] * _ip_normals[qp];
      71             :     }
      72     2059160 : }
      73             : 
      74             : void
      75     2059160 : DiffusionIPHDGAssemblyHelper::lmFace()
      76             : {
      77     2059160 :   const auto h_elem = _elem_volume / _side_area;
      78             : 
      79    17532504 :   for (const auto i : index_range(_lm_re))
      80    46611552 :     for (const auto qp : make_range(_ip_qrule_face->n_points()))
      81             :     {
      82    62276416 :       _lm_re(i) += _ip_JxW_face[qp] * _face_diff[qp] * _grad_u_sol[qp] * _ip_normals[qp] *
      83    62276416 :                    _lm_phi_face[i][qp];
      84    93414624 :       _lm_re(i) += _ip_JxW_face[qp] * _alpha / h_elem * _face_diff[qp] *
      85   124552832 :                    (_lm_u_sol[qp] - _u_sol[qp]) * _lm_phi_face[i][qp];
      86             :     }
      87     2059160 : }
      88             : 
      89             : void
      90       58784 : DiffusionIPHDGAssemblyHelper::scalarDirichlet(const Moose::Functor<Real> & dirichlet_value)
      91             : {
      92       58784 :   const auto h_elem = _elem_volume / _side_area;
      93             : 
      94      177612 :   for (const auto qp : make_range(_ip_qrule_face->n_points()))
      95             :   {
      96      356484 :     const auto scalar_value = dirichlet_value(
      97      118828 :         Moose::ElemSideQpArg{
      98      118828 :             _ip_current_elem, _ip_current_side, qp, _ip_qrule_face, _ip_q_point_face[qp]},
      99      118828 :         _ti.determineState());
     100             : 
     101      517648 :     for (const auto i : index_range(_u_dof_indices))
     102             :     {
     103      797640 :       _scalar_re(i) -= _ip_JxW_face[qp] * _face_diff[qp] * _scalar_phi_face[i][qp] *
     104     1196460 :                        (_grad_u_sol[qp] * _ip_normals[qp]);
     105     1196460 :       _scalar_re(i) -= _ip_JxW_face[qp] * _alpha / h_elem * _face_diff[qp] *
     106     1595280 :                        (scalar_value - _u_sol[qp]) * _scalar_phi_face[i][qp];
     107     1196460 :       _scalar_re(i) += _ip_JxW_face[qp] * (scalar_value - _u_sol[qp]) * _face_diff[qp] *
     108     1595280 :                        _grad_scalar_phi_face[i][qp] * _ip_normals[qp];
     109             :     }
     110             :   }
     111       58784 : }