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