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