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 "ADCubicTransition.h"
      11             : #include "MooseError.h"
      12             : #include "ADReal.h"
      13             : 
      14             : #include "libmesh/dense_vector.h"
      15             : #include "DenseMatrix.h"
      16             : 
      17          72 : ADCubicTransition::ADCubicTransition(const ADReal & x_center, const ADReal & transition_width)
      18             :   : ADSmoothTransition(x_center, transition_width),
      19             : 
      20         144 :     _A(0.0),
      21          72 :     _B(0.0),
      22          72 :     _C(0.0),
      23          72 :     _D(0.0),
      24             : 
      25          72 :     _initialized(false)
      26             : {
      27          72 : }
      28             : 
      29             : void
      30          72 : ADCubicTransition::initialize(const ADReal & f1_end_value,
      31             :                               const ADReal & f2_end_value,
      32             :                               const ADReal & df1dx_end_value,
      33             :                               const ADReal & df2dx_end_value)
      34             : {
      35             :   // compute cubic polynomial coefficients
      36             : 
      37          72 :   DenseMatrix<ADReal> mat(4, 4);
      38             : 
      39          72 :   mat(0, 0) = std::pow(_x1, 3);
      40          72 :   mat(0, 1) = std::pow(_x1, 2);
      41          72 :   mat(0, 2) = _x1;
      42          72 :   mat(0, 3) = 1.0;
      43             : 
      44          72 :   mat(1, 0) = std::pow(_x2, 3);
      45          72 :   mat(1, 1) = std::pow(_x2, 2);
      46          72 :   mat(1, 2) = _x2;
      47          72 :   mat(1, 3) = 1.0;
      48             : 
      49         144 :   mat(2, 0) = 3.0 * std::pow(_x1, 2);
      50         144 :   mat(2, 1) = 2.0 * _x1;
      51          72 :   mat(2, 2) = 1.0;
      52          72 :   mat(2, 3) = 0.0;
      53             : 
      54         144 :   mat(3, 0) = 3.0 * std::pow(_x2, 2);
      55         144 :   mat(3, 1) = 2.0 * _x2;
      56          72 :   mat(3, 2) = 1.0;
      57          72 :   mat(3, 3) = 0.0;
      58             : 
      59          72 :   DenseVector<ADReal> rhs(4);
      60          72 :   rhs(0) = f1_end_value;
      61          72 :   rhs(1) = f2_end_value;
      62          72 :   rhs(2) = df1dx_end_value;
      63          72 :   rhs(3) = df2dx_end_value;
      64             : 
      65          72 :   DenseVector<ADReal> coefs(4);
      66          72 :   mat.lu_solve(rhs, coefs);
      67             : 
      68          72 :   _A = coefs(0);
      69          72 :   _B = coefs(1);
      70          72 :   _C = coefs(2);
      71          72 :   _D = coefs(3);
      72             : 
      73          72 :   _initialized = true;
      74          72 : }
      75             : 
      76             : ADReal
      77        1200 : ADCubicTransition::value(const ADReal & x, const ADReal & f1, const ADReal & f2) const
      78             : {
      79             :   mooseAssert(_initialized, "initialize() must be called.");
      80             : 
      81        1200 :   if (x <= _x1)
      82         450 :     return f1;
      83         750 :   else if (x >= _x2)
      84         450 :     return f2;
      85             :   else
      86        1200 :     return _A * std::pow(x, 3) + _B * std::pow(x, 2) + _C * x + _D;
      87             : }