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 "NSMomentumInviscidBC.h"
      11             : 
      12             : InputParameters
      13          38 : NSMomentumInviscidBC::validParams()
      14             : {
      15          38 :   InputParameters params = NSIntegratedBC::validParams();
      16          38 :   params.addClassDescription("his class corresponds to the inviscid part of the 'natural' boundary "
      17             :                              "condition for the momentum equations.");
      18          76 :   params.addRequiredParam<unsigned>(
      19             :       "component", "(0,1,2) = (x,y,z) for which momentum component this BC is applied to");
      20          38 :   return params;
      21           0 : }
      22             : 
      23          20 : NSMomentumInviscidBC::NSMomentumInviscidBC(const InputParameters & parameters)
      24             :   : NSIntegratedBC(parameters),
      25          20 :     _component(getParam<unsigned>("component")),
      26             :     // Object for computing deriviatives of pressure
      27          20 :     _pressure_derivs(*this)
      28             : {
      29          20 : }
      30             : 
      31             : Real
      32      154112 : NSMomentumInviscidBC::pressureQpResidualHelper(Real pressure)
      33             : {
      34             :   // n . (Ip) . v
      35             : 
      36             :   // The pressure contribution: p * n(component) * phi_i
      37      154112 :   Real press_term = pressure * _normals[_qp](_component) * _test[_i][_qp];
      38             : 
      39             :   // Return value, or print it first if debugging...
      40      154112 :   return press_term;
      41             : }
      42             : 
      43             : Real
      44           0 : NSMomentumInviscidBC::pressureQpJacobianHelper(unsigned var_number)
      45             : {
      46           0 :   return _normals[_qp](_component) * _pressure_derivs.get_grad(var_number) * _phi[_j][_qp] *
      47           0 :          _test[_i][_qp];
      48             : }
      49             : 
      50             : Real
      51      154112 : NSMomentumInviscidBC::convectiveQpResidualHelper(Real rhou_udotn)
      52             : {
      53             :   // n . (rho*uu) . v = rho*(u.n)*(u.v) = (rho*u)(u.n) . v
      54             : 
      55             :   // The "inviscid" contribution: (rho*u)(u.n) . v
      56      154112 :   Real conv_term = rhou_udotn * _test[_i][_qp];
      57             : 
      58             :   // Return value, or print it first if debugging...
      59      154112 :   return conv_term;
      60             : }
      61             : 
      62             : Real
      63      374784 : NSMomentumInviscidBC::convectiveQpJacobianHelper(unsigned var_number)
      64             : {
      65             :   // Velocity vector object
      66      374784 :   RealVectorValue vel(_u_vel[_qp], _v_vel[_qp], _w_vel[_qp]);
      67             : 
      68             :   // Variable to store convective contribution to boundary integral.
      69             :   Real conv_term = 0.0;
      70             : 
      71             :   // Inviscid components
      72      374784 :   switch (var_number)
      73             :   {
      74       93696 :     case 0: // density
      75             :       // Note: the minus sign here is correct, it comes from differentiating wrt U_0
      76             :       // (rho) which is in the denominator.
      77       93696 :       conv_term = -vel(_component) * (vel * _normals[_qp]) * _phi[_j][_qp] * _test[_i][_qp];
      78       93696 :       break;
      79             : 
      80      187392 :     case 1:
      81             :     case 2:
      82             :     case 3: // momentums
      83      187392 :       if (var_number - 1 == _component)
      84             :         // See Eqn. (68) from the notes for the inviscid boundary terms
      85       93696 :         conv_term = ((vel * _normals[_qp]) + vel(_component) * _normals[_qp](_component)) *
      86       93696 :                     _phi[_j][_qp] * _test[_i][_qp];
      87             :       else
      88             :         // off-diagonal
      89       93696 :         conv_term =
      90       93696 :             vel(_component) * _normals[_qp](var_number - 1) * _phi[_j][_qp] * _test[_i][_qp];
      91             :       break;
      92             : 
      93             :     case 4: // energy
      94             :       // No derivative wrt energy
      95             :       conv_term = 0.0;
      96             :       break;
      97             : 
      98           0 :     default:
      99           0 :       mooseError("Shouldn't get here!");
     100             :       break;
     101             :   }
     102             : 
     103             :   // Return the result.  We could return it directly from the switch statement, but this is
     104             :   // convenient for printing...
     105      374784 :   return conv_term;
     106             : }