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             : #pragma once
      11             : 
      12             : /**
      13             :  * Class outside the Moose hierarchy that contains common
      14             :  * functionality for computing derivatives of the viscous
      15             :  * stress tensor.
      16             :  *
      17             :  * This class is templated so that it can be used by either
      18             :  * a Kernel object or a BC object.
      19             :  */
      20             : template <class T>
      21             : class NSViscStressTensorDerivs
      22             : {
      23             : public:
      24             :   NSViscStressTensorDerivs(T & x);
      25             : 
      26             :   /**
      27             :    * The primary interface for computing viscous stress
      28             :    * tensor derivatives.  Requires access to protected data
      29             :    * from the the underlying data.  Uses index notation from
      30             :    * the notes.
      31             :    */
      32             :   Real dtau(unsigned k, unsigned ell, unsigned m);
      33             : 
      34             : private:
      35             :   T & _data;
      36             : };
      37             : 
      38             : template <class T>
      39           0 : NSViscStressTensorDerivs<T>::NSViscStressTensorDerivs(T & x) : _data(x)
      40             : {
      41             : }
      42             : 
      43             : template <class T>
      44             : Real
      45           0 : NSViscStressTensorDerivs<T>::dtau(unsigned k, unsigned ell, unsigned m)
      46             : {
      47             :   // Try to access underlying data.  Since this class is a friend, we can
      48             :   // directly access _qp and other protected data.  This only works if the
      49             :   // individual variables have the **same names** in all types T which may
      50             :   // be used to construct this class.
      51             : 
      52             :   //
      53             :   // Some error checking on input indices...
      54             :   //
      55             : 
      56             :   // 0 <= k,ell <= 2
      57           0 :   if (k > 2 || ell > 2)
      58           0 :     mooseError("Error, 0 <= k,ell <= 2 violated!");
      59             : 
      60             :   // 0 <= m <= 4
      61           0 :   if (m >= 5)
      62           0 :     mooseError("Error, m <= 4 violated!");
      63             : 
      64             :   //
      65             :   // Convenience variables
      66             :   //
      67             : 
      68           0 :   const Real rho = _data._rho[_data._qp];
      69           0 :   const Real rho2 = rho * rho;
      70           0 :   const Real phij = _data._phi[_data._j][_data._qp];
      71             : 
      72           0 :   const Real mu = _data._dynamic_viscosity[_data._qp];
      73           0 :   const Real nu = mu / rho;
      74             : 
      75             :   const RealVectorValue U(
      76           0 :       _data._rho_u[_data._qp], _data._rho_v[_data._qp], _data._rho_w[_data._qp]);
      77             : 
      78           0 :   const Real divU = _data._grad_rho_u[_data._qp](0) + _data._grad_rho_v[_data._qp](1) +
      79           0 :                     _data._grad_rho_w[_data._qp](2);
      80             : 
      81             :   // This makes a copy...but the resulting code is cleaner
      82           0 :   std::vector<RealVectorValue> gradU(3);
      83           0 :   gradU[0] = _data._grad_rho_u[_data._qp];
      84           0 :   gradU[1] = _data._grad_rho_v[_data._qp];
      85           0 :   gradU[2] = _data._grad_rho_w[_data._qp];
      86             : 
      87             :   // So we can refer to gradients without repeated indexing.
      88           0 :   const RealVectorValue & grad_phij = _data._grad_phi[_data._j][_data._qp];
      89           0 :   const RealVectorValue & grad_rho = _data._grad_rho[_data._qp];
      90             : 
      91           0 :   switch (m)
      92             :   {
      93           0 :     case 0: // density
      94             :     {
      95           0 :       const Real term1 = 2.0 / rho2 * (U(k) * grad_rho(ell) + U(ell) * grad_rho(k)) * phij;
      96           0 :       const Real term2 = -1.0 / rho *
      97           0 :                          ((gradU[k](ell) + gradU[ell](k)) * phij +
      98           0 :                           (U(k) * grad_phij(ell) + U(ell) * grad_phij(k)));
      99             : 
     100             :       // Kronecker delta terms
     101             :       Real term3 = 0.0;
     102             :       Real term4 = 0.0;
     103           0 :       if (k == ell)
     104             :       {
     105           0 :         term3 = -4.0 / 3.0 / rho2 * (U * grad_rho) * phij;
     106           0 :         term4 = 2.0 / 3.0 / rho * (U * grad_phij + divU * phij);
     107             :       }
     108             : 
     109             :       // Sum up result and return
     110           0 :       return nu * (term1 + term2 + term3 + term4);
     111             :     }
     112             : 
     113             :     // momentums
     114           0 :     case 1:
     115             :     case 2:
     116             :     case 3:
     117             :     {
     118             :       // note: when comparing m to k or ell, or indexing into Points,
     119             :       // must map m -> 0, 1, 2 by subtracting 1.
     120           0 :       const unsigned m_local = m - 1;
     121             : 
     122             :       // Kronecker delta terms
     123           0 :       const Real delta_km = (k == m_local ? 1.0 : 0.0);
     124           0 :       const Real delta_ellm = (ell == m_local ? 1.0 : 0.0);
     125           0 :       const Real delta_kell = (k == ell ? 1.0 : 0.0);
     126             : 
     127             :       return nu *
     128             :              (
     129           0 :                  /*     */ delta_km * (grad_phij(ell) - (phij / rho) * grad_rho(ell)) +
     130           0 :                  /*     */ delta_ellm * (grad_phij(k) - (phij / rho) * grad_rho(k)) -
     131           0 :                  (2. / 3.) * delta_kell * (grad_phij(m_local) - (phij / rho) * grad_rho(m_local)));
     132             :     } // end case 1,2,3
     133             : 
     134             :     case 4:
     135             :       // Derivative wrt to energy variable is zero.
     136             :       return 0.;
     137             : 
     138             :     default:
     139             :       mooseError("Invalid variable requested.");
     140             :       break;
     141             :   }
     142             : 
     143             :   return 0.;
     144           0 : }