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 "DGDiffusion.h"
      11             : #include "MooseVariableFE.h"
      12             : 
      13             : #include "libmesh/utility.h"
      14             : 
      15             : registerMooseObject("MooseApp", DGDiffusion);
      16             : 
      17             : InputParameters
      18       15783 : DGDiffusion::validParams()
      19             : {
      20       15783 :   InputParameters params = DGKernel::validParams();
      21       15783 :   params.addClassDescription("Computes residual contribution for the diffusion operator using "
      22             :                              "discontinous Galerkin method.");
      23             :   // See header file for sigma and epsilon
      24       15783 :   params.addRequiredParam<Real>("sigma", "sigma");
      25       15783 :   params.addRequiredParam<Real>("epsilon", "epsilon");
      26       47349 :   params.addParam<MaterialPropertyName>(
      27       31566 :       "diff", 1., "The diffusion (or thermal conductivity or viscosity) coefficient.");
      28       15783 :   return params;
      29           0 : }
      30             : 
      31         787 : DGDiffusion::DGDiffusion(const InputParameters & parameters)
      32             :   : DGKernel(parameters),
      33         787 :     _epsilon(getParam<Real>("epsilon")),
      34         787 :     _sigma(getParam<Real>("sigma")),
      35         787 :     _diff(getMaterialProperty<Real>("diff")),
      36        1574 :     _diff_neighbor(getNeighborMaterialProperty<Real>("diff"))
      37             : {
      38         787 : }
      39             : 
      40             : Real
      41   152575856 : DGDiffusion::computeQpResidual(Moose::DGResidualType type)
      42             : {
      43   152575856 :   Real r = 0.0;
      44             : 
      45   152575856 :   const int elem_b_order = std::max(libMesh::Order(1), _var.order());
      46             :   const Real h_elem =
      47   152575856 :       _current_elem_volume / _current_side_volume * 1.0 / Utility::pow<2>(elem_b_order);
      48             : 
      49   152575856 :   switch (type)
      50             :   {
      51    76329760 :     case Moose::Element:
      52   152659520 :       r -= 0.5 *
      53    76329760 :            (_diff[_qp] * _grad_u[_qp] * _normals[_qp] +
      54    76329760 :             _diff_neighbor[_qp] * _grad_u_neighbor[_qp] * _normals[_qp]) *
      55    76329760 :            _test[_i][_qp];
      56    76329760 :       r += _epsilon * 0.5 * (_u[_qp] - _u_neighbor[_qp]) * _diff[_qp] * _grad_test[_i][_qp] *
      57    76329760 :            _normals[_qp];
      58    76329760 :       r += _sigma / h_elem * (_u[_qp] - _u_neighbor[_qp]) * _test[_i][_qp];
      59    76329760 :       break;
      60             : 
      61    76246096 :     case Moose::Neighbor:
      62   152492192 :       r += 0.5 *
      63    76246096 :            (_diff[_qp] * _grad_u[_qp] * _normals[_qp] +
      64    76246096 :             _diff_neighbor[_qp] * _grad_u_neighbor[_qp] * _normals[_qp]) *
      65    76246096 :            _test_neighbor[_i][_qp];
      66    76246096 :       r += _epsilon * 0.5 * (_u[_qp] - _u_neighbor[_qp]) * _diff_neighbor[_qp] *
      67   152492192 :            _grad_test_neighbor[_i][_qp] * _normals[_qp];
      68    76246096 :       r -= _sigma / h_elem * (_u[_qp] - _u_neighbor[_qp]) * _test_neighbor[_i][_qp];
      69    76246096 :       break;
      70             :   }
      71             : 
      72   152575856 :   return r;
      73             : }
      74             : 
      75             : Real
      76   223135712 : DGDiffusion::computeQpJacobian(Moose::DGJacobianType type)
      77             : {
      78   223135712 :   Real r = 0.0;
      79             : 
      80   223135712 :   const int elem_b_order = std::max(libMesh::Order(1), _var.order());
      81             :   const Real h_elem =
      82   223135712 :       _current_elem_volume / _current_side_volume * 1.0 / Utility::pow<2>(elem_b_order);
      83             : 
      84   223135712 :   switch (type)
      85             :   {
      86    55820216 :     case Moose::ElementElement:
      87    55820216 :       r -= 0.5 * _diff[_qp] * _grad_phi[_j][_qp] * _normals[_qp] * _test[_i][_qp];
      88    55820216 :       r += _epsilon * 0.5 * _phi[_j][_qp] * _diff[_qp] * _grad_test[_i][_qp] * _normals[_qp];
      89    55820216 :       r += _sigma / h_elem * _phi[_j][_qp] * _test[_i][_qp];
      90    55820216 :       break;
      91             : 
      92    55778744 :     case Moose::ElementNeighbor:
      93    55778744 :       r -= 0.5 * _diff_neighbor[_qp] * _grad_phi_neighbor[_j][_qp] * _normals[_qp] * _test[_i][_qp];
      94    55778744 :       r += _epsilon * 0.5 * -_phi_neighbor[_j][_qp] * _diff[_qp] * _grad_test[_i][_qp] *
      95    55778744 :            _normals[_qp];
      96    55778744 :       r += _sigma / h_elem * -_phi_neighbor[_j][_qp] * _test[_i][_qp];
      97    55778744 :       break;
      98             : 
      99    55778744 :     case Moose::NeighborElement:
     100    55778744 :       r += 0.5 * _diff[_qp] * _grad_phi[_j][_qp] * _normals[_qp] * _test_neighbor[_i][_qp];
     101    55778744 :       r += _epsilon * 0.5 * _phi[_j][_qp] * _diff_neighbor[_qp] * _grad_test_neighbor[_i][_qp] *
     102    55778744 :            _normals[_qp];
     103    55778744 :       r -= _sigma / h_elem * _phi[_j][_qp] * _test_neighbor[_i][_qp];
     104    55778744 :       break;
     105             : 
     106    55758008 :     case Moose::NeighborNeighbor:
     107    55758008 :       r += 0.5 * _diff_neighbor[_qp] * _grad_phi_neighbor[_j][_qp] * _normals[_qp] *
     108    55758008 :            _test_neighbor[_i][_qp];
     109    55758008 :       r += _epsilon * 0.5 * -_phi_neighbor[_j][_qp] * _diff_neighbor[_qp] *
     110   111516016 :            _grad_test_neighbor[_i][_qp] * _normals[_qp];
     111    55758008 :       r -= _sigma / h_elem * -_phi_neighbor[_j][_qp] * _test_neighbor[_i][_qp];
     112    55758008 :       break;
     113             :   }
     114             : 
     115   223135712 :   return r;
     116             : }