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          48 : ADCubicTransition::ADCubicTransition(const ADReal & x_center, const ADReal & transition_width)
      18             :   : ADSmoothTransition(x_center, transition_width),
      19             : 
      20          96 :     _A(0.0),
      21          48 :     _B(0.0),
      22          48 :     _C(0.0),
      23          48 :     _D(0.0),
      24             : 
      25          48 :     _initialized(false)
      26             : {
      27          48 : }
      28             : 
      29             : void
      30          48 : 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             :   using std::pow;
      37             : 
      38          48 :   DenseMatrix<ADReal> mat(4, 4);
      39             : 
      40          48 :   mat(0, 0) = pow(_x1, 3);
      41          48 :   mat(0, 1) = pow(_x1, 2);
      42          48 :   mat(0, 2) = _x1;
      43          48 :   mat(0, 3) = 1.0;
      44             : 
      45          48 :   mat(1, 0) = pow(_x2, 3);
      46          48 :   mat(1, 1) = pow(_x2, 2);
      47          48 :   mat(1, 2) = _x2;
      48          48 :   mat(1, 3) = 1.0;
      49             : 
      50          96 :   mat(2, 0) = 3.0 * pow(_x1, 2);
      51          96 :   mat(2, 1) = 2.0 * _x1;
      52          48 :   mat(2, 2) = 1.0;
      53          48 :   mat(2, 3) = 0.0;
      54             : 
      55          96 :   mat(3, 0) = 3.0 * pow(_x2, 2);
      56          96 :   mat(3, 1) = 2.0 * _x2;
      57          48 :   mat(3, 2) = 1.0;
      58          48 :   mat(3, 3) = 0.0;
      59             : 
      60          48 :   DenseVector<ADReal> rhs(4);
      61          48 :   rhs(0) = f1_end_value;
      62          48 :   rhs(1) = f2_end_value;
      63          48 :   rhs(2) = df1dx_end_value;
      64          48 :   rhs(3) = df2dx_end_value;
      65             : 
      66          48 :   DenseVector<ADReal> coefs(4);
      67          48 :   mat.lu_solve(rhs, coefs);
      68             : 
      69          48 :   _A = coefs(0);
      70          48 :   _B = coefs(1);
      71          48 :   _C = coefs(2);
      72          48 :   _D = coefs(3);
      73             : 
      74          48 :   _initialized = true;
      75          48 : }
      76             : 
      77             : ADReal
      78        1000 : ADCubicTransition::value(const ADReal & x, const ADReal & f1, const ADReal & f2) const
      79             : {
      80             :   mooseAssert(_initialized, "initialize() must be called.");
      81             :   using std::pow;
      82             : 
      83        1000 :   if (x <= _x1)
      84         375 :     return f1;
      85         625 :   else if (x >= _x2)
      86         375 :     return f2;
      87             :   else
      88        1000 :     return _A * pow(x, 3) + _B * pow(x, 2) + _C * x + _D;
      89             : }