MathUtils Namespace Reference

Enumerations
enum	ComputeType { ComputeType::value, ComputeType::derivative }

Functions
Real	poly1Log (Real x, Real tol, unsigned int derivative_order)
Real	poly2Log (Real x, Real tol, unsigned int derivative_order)
Real	poly3Log (Real x, Real tol, unsigned int derivative_order)
Real	poly4Log (Real x, Real tol, unsigned int derivative_order)
Real	taylorLog (Real x)
Point	barycentricToCartesian2D (const Point &p0, const Point &p1, const Point &p2, const Real b0, const Real b1, const Real b2)
	Evaluate Cartesian coordinates of any center point of a triangle given Barycentric coordinates of center point and Cartesian coordinates of triangle's vertices. More...
Point	barycentricToCartesian3D (const Point &p0, const Point &p1, const Point &p2, const Point &p3, const Real b0, const Real b1, const Real b2, const Real b3)
	Evaluate Cartesian coordinates of any center point of a tetrahedron given Barycentric coordinates of center point and Cartesian coordinates of tetrahedon's vertices. More...
Point	circumcenter2D (const Point &p0, const Point &p1, const Point &p2)
	Evaluate circumcenter of a triangle given three arbitrary points. More...
Point	circumcenter3D (const Point &p0, const Point &p1, const Point &p2, const Point &p3)
	Evaluate circumcenter of a tetrahedrom given four arbitrary points. More...
template<typename T >
T	round (T x)
template<typename T >
T	sign (T x)
template<typename T >
T	pow (T x, int e)
template<typename T >
T	heavyside (T x)
template<typename T >
T	regularizedHeavyside (T x, Real smoothing_length)
template<typename T >
T	regularizedHeavysideDerivative (T x, Real smoothing_length)
template<typename T >
T	positivePart (T x)
template<typename T >
T	negativePart (T x)
template<typename T , typename T2 , typename T3 , typename std::enable_if< libMesh::ScalarTraits< T >::value &&libMesh::ScalarTraits< T2 >::value &&libMesh::ScalarTraits< T3 >::value, int >::type = 0>
void	addScaled (const T &a, const T2 &b, T3 &result)
template<typename T , typename T2 , typename T3 , typename std::enable_if< libMesh::ScalarTraits< T >::value, int >::type = 0>
void	addScaled (const T &scalar, const libMesh::NumericVector< T2 > &numeric_vector, libMesh::NumericVector< T3 > &result)
template<typename T , typename T2 , template< typename > class W, template< typename > class W2, typename std::enable_if< std::is_same< typename W< T >::index_type, unsigned int >::value &&std::is_same< typename W2< T2 >::index_type, unsigned int >::value, int >::type = 0>
libMesh::CompareTypes< T, T2 >::supertype	dotProduct (const W< T > &a, const W2< T2 > &b)
template<typename C , typename T , typename R = typename libMesh::CompareTypes<typename C::value_type, T>::supertype>
R	poly (const C &c, const T x, const bool derivative=false)
	Evaluate a polynomial with the coefficients c at x. More...
template<typename C , typename T , typename R = typename libMesh::CompareTypes<typename C::value_type, T>::supertype>
R	polynomial (const C &c, const T x)
	Evaluate a polynomial with the coefficients c at x. More...
template<typename C , typename T , typename R = typename libMesh::CompareTypes<typename C::value_type, T>::supertype>
R	polynomialDerivative (const C &c, const T x)
	Returns the derivative of polynomial(c, x) with respect to x. More...
template<typename T , typename T2 >
T	clamp (const T &x, T2 lowerlimit, T2 upperlimit)
template<typename T , typename T2 >
T	smootherStep (T x, T2 start, T2 end, bool derivative=false)
template<ComputeType compute_type, typename X , typename S , typename E >
auto	smootherStep (const X &x, const S &start, const E &end)
template<typename T >
void	mooseSetToZero (T &v)
	Helper function templates to set a variable to zero. More...
template<typename T >
void	mooseSetToZero (T *&)
template<>
void	mooseSetToZero (std::vector< Real > &vec)
std::vector< std::vector< unsigned int > >	multiIndex (unsigned int dim, unsigned int order)
	generate a complete multi index table for given dimension and order i.e. More...
template<ComputeType compute_type, typename X , typename X1 , typename X2 , typename Y1 , typename Y2 >
auto	linearInterpolation (const X &x, const X1 &x1, const X2 &x2, const Y1 &y1, const Y2 &y2)
template<typename T1 , typename T2 >
std::size_t	euclideanMod (T1 dividend, T2 divisor)
	perform modulo operator for Euclidean division that ensures a non-negative result More...
template<typename T >
T	gradName (const T &base_prop_name)
	automatic prefixing for naming material properties based on gradients of coupled variables/functors More...
template<typename T >
T	timeDerivName (const T &base_prop_name)
	automatic prefixing for naming material properties based on time derivatives of coupled variables/functors More...
void	kron (RealEigenMatrix &product, const RealEigenMatrix &mat_A, const RealEigenMatrix &mat_B)
	Computes the Kronecker product of two matrices. More...
template<>
void	mooseSetToZero< RankFourTensor > (RankFourTensor &v)
	Helper function template specialization to set an object to zero. More...
template<>
void	mooseSetToZero< ADRankFourTensor > (ADRankFourTensor &v)
template<>
void	mooseSetToZero< RankFourTensorTempl< Real > > (RankFourTensorTempl< Real > &v)
template<>
void	mooseSetToZero< RankFourTensorTempl< ADReal > > (RankFourTensorTempl< ADReal > &v)
template<>
void	mooseSetToZero< RankThreeTensorTempl< Real > > (RankThreeTensorTempl< Real > &v)
	Helper function template specialization to set an object to zero. More...
template<>
void	mooseSetToZero< RankThreeTensorTempl< ADReal > > (RankThreeTensorTempl< ADReal > &v)
template<>
void	mooseSetToZero< RankThreeTensor > (RankThreeTensor &v)
template<>
void	mooseSetToZero< ADRankThreeTensor > (ADRankThreeTensor &v)
template<>
void	mooseSetToZero< RankTwoTensor > (RankTwoTensor &v)
	Helper function template specialization to set an object to zero. More...
template<>
void	mooseSetToZero< ADRankTwoTensor > (ADRankTwoTensor &v)
	Helper function template specialization to set an object to zero. More...
template<>
void	mooseSetToZero< SymmetricRankFourTensor > (SymmetricRankFourTensor &v)
	Helper function template specialization to set an object to zero. More...
template<>
void	mooseSetToZero< ADSymmetricRankFourTensor > (ADSymmetricRankFourTensor &v)
template<>
void	mooseSetToZero< SymmetricRankFourTensorTempl< Real > > (SymmetricRankFourTensorTempl< Real > &v)
template<>
void	mooseSetToZero< SymmetricRankFourTensorTempl< ADReal > > (SymmetricRankFourTensorTempl< ADReal > &v)
template<>
void	mooseSetToZero< SymmetricRankTwoTensor > (SymmetricRankTwoTensor &v)
	Helper function template specialization to set an object to zero. More...
template<>
void	mooseSetToZero< ADSymmetricRankTwoTensor > (ADSymmetricRankTwoTensor &v)
	Helper function template specialization to set an object to zero. More...
Real	plainLog (Real x, unsigned int derivative_order)
	FactorizedRankTwoTensorOperatorMapUnary (log, std::log(eigval))
	FactorizedRankTwoTensorOperatorMapDerivativeUnary (dlog, std::log(eigval), 1/eigval)
	FactorizedRankTwoTensorOperatorMapUnary (exp, std::exp(eigval))
	FactorizedRankTwoTensorOperatorMapDerivativeUnary (dexp, std::exp(eigval), std::exp(eigval))
	FactorizedRankTwoTensorOperatorMapUnary (sqrt, std::sqrt(eigval))
	FactorizedRankTwoTensorOperatorMapDerivativeUnary (dsqrt, std::sqrt(eigval), std::pow(eigval, -1./2.)/2.)
	FactorizedRankTwoTensorOperatorMapUnary (cbrt, std::cbrt(eigval))
	FactorizedRankTwoTensorOperatorMapDerivativeUnary (dcbrt, std::cbrt(eigval), std::pow(eigval, -2./3.)/3.)
	FactorizedRankTwoTensorOperatorMapBinary (pow, std::pow(eigval, arg))
	FactorizedRankTwoTensorOperatorMapDerivativeBinary (dpow, std::pow(eigval, arg), arg *std::pow(eigval, arg - 1))

Variables
static constexpr Real	sqrt2 = 1.4142135623730951
	std::sqrt is not constexpr, so we add sqrt(2) as a constant (used in Mandel notation) More...

Enumeration Type Documentation

ComputeType

enum MathUtils::ComputeType

strong

Enumerator
value
derivative

Definition at line 338 of file MathUtils.h.

{
  value,
  derivative
};

Function Documentation

addScaled() [1/2]

template<typename T , typename T2 , typename T3 , typename std::enable_if< libMesh::ScalarTraits< T >::value &&libMesh::ScalarTraits< T2 >::value &&libMesh::ScalarTraits< T3 >::value, int >::type = 0>

void MathUtils::addScaled	(	const T &	a,
		const T2 &	b,
		T3 &	result
	)

Definition at line 170 of file MathUtils.h.

Referenced by NewmarkBeta::computeTimeDerivativeHelper(), and BDF2::computeTimeDerivativeHelper().

{
  result += a * b;
}

addScaled() [2/2]

template<typename T , typename T2 , typename T3 , typename std::enable_if< libMesh::ScalarTraits< T >::value, int >::type = 0>

void MathUtils::addScaled	(	const T &	scalar,
		const libMesh::NumericVector< T2 > &	numeric_vector,
		libMesh::NumericVector< T3 > &	result
	)

Definition at line 180 of file MathUtils.h.

{
  result.add(scalar, numeric_vector);
}

barycentricToCartesian2D()

Point MathUtils::barycentricToCartesian2D	(	const Point &	p0,
		const Point &	p1,
		const Point &	p2,
		const Real	b0,
		const Real	b1,
		const Real	b2
	)

Evaluate Cartesian coordinates of any center point of a triangle given Barycentric coordinates of center point and Cartesian coordinates of triangle's vertices.

Parameters

p0,p1,p2	are the three non-collinear vertices in Cartesian coordinates
b0,b1,b2	is the center point in barycentric coordinates with b0+b1+b2=1, e.g. (1/3,1/3,1/3) for a centroid

Returns

the center point of triangle in Cartesian coordinates

Definition at line 217 of file MathUtils.C.

Referenced by circumcenter2D().

{
  mooseAssert(!MooseUtils::isZero(b0 + b1 + b2 - 1.0), "Barycentric coordinates must sum to one!");

  Point center;

  for (unsigned int d = 0; d < 2; ++d)
    center(d) = p0(d) * b0 + p1(d) * b1 + p2(d) * b2;
  // p0, p1, p2 are vertices of triangle
  // b0, b1, b2 are Barycentric coordinates of the triangle center

  return center;
}

barycentricToCartesian3D()

Point MathUtils::barycentricToCartesian3D	(	const Point &	p0,
		const Point &	p1,
		const Point &	p2,
		const Point &	p3,
		const Real	b0,
		const Real	b1,
		const Real	b2,
		const Real	b3
	)

Evaluate Cartesian coordinates of any center point of a tetrahedron given Barycentric coordinates of center point and Cartesian coordinates of tetrahedon's vertices.

Parameters

p0,p1,p2,p3	are the three non-coplanar vertices in Cartesian coordinates
b0,b1,b2,b3	is the center point in barycentric coordinates with b0+b1+b2+b3=1, e.g. (1/4,1/4,1/4,1/4) for a centroid.

Returns

the center point of tetrahedron in Cartesian coordinates

Definition at line 237 of file MathUtils.C.

Referenced by circumcenter3D().

{
  mooseAssert(!MooseUtils::isZero(b0 + b1 + b2 + b3 - 1.0),
              "Barycentric coordinates must sum to one!");

  Point center;

  for (unsigned int d = 0; d < 3; ++d)
    center(d) = p0(d) * b0 + p1(d) * b1 + p2(d) * b2 + p3(d) * b3;
  // p0, p1, p2, p3 are vertices of tetrahedron
  // b0, b1, b2, b3 are Barycentric coordinates of the tetrahedron center

  return center;
}

circumcenter2D()

Point MathUtils::circumcenter2D	(	const Point &	p0,
		const Point &	p1,
		const Point &	p2
	)

Evaluate circumcenter of a triangle given three arbitrary points.

Parameters

p0,p1,p2

are the three non-collinear vertices in Cartesian coordinates

Returns

the circumcenter in Cartesian coordinates

Definition at line 260 of file MathUtils.C.

{
  // Square of triangle edge lengths
  Real edge01 = (p0 - p1).norm_sq();
  Real edge02 = (p0 - p2).norm_sq();
  Real edge12 = (p1 - p2).norm_sq();

  // Barycentric weights for circumcenter
  Real weight0 = edge12 * (edge01 + edge02 - edge12);
  Real weight1 = edge02 * (edge01 + edge12 - edge02);
  Real weight2 = edge01 * (edge02 + edge12 - edge01);

  Real sum_weights = weight0 + weight1 + weight2;

  // Check to make sure vertices are not collinear
  if (MooseUtils::isZero(sum_weights))
    mooseError("Cannot evaluate circumcenter. Points should be non-collinear.");

  Real inv_sum_weights = 1.0 / sum_weights;

  // Barycentric coordinates
  Real b0 = weight0 * inv_sum_weights;
  Real b1 = weight1 * inv_sum_weights;
  Real b2 = weight2 * inv_sum_weights;

  return MathUtils::barycentricToCartesian2D(p0, p1, p2, b0, b1, b2);
}

circumcenter3D()

Point MathUtils::circumcenter3D	(	const Point &	p0,
		const Point &	p1,
		const Point &	p2,
		const Point &	p3
	)

Evaluate circumcenter of a tetrahedrom given four arbitrary points.

Parameters

p0,p1,p2,p3

are the four non-coplanar vertices in Cartesian coordinates

Returns

the circumcenter in Cartesian coordinates

Definition at line 289 of file MathUtils.C.

{
  // Square of tetrahedron edge lengths
  Real edge01 = (p0 - p1).norm_sq();
  Real edge02 = (p0 - p2).norm_sq();
  Real edge03 = (p0 - p3).norm_sq();
  Real edge12 = (p1 - p2).norm_sq();
  Real edge13 = (p1 - p3).norm_sq();
  Real edge23 = (p2 - p3).norm_sq();

  // Barycentric weights for circumcenter
  Real weight0 = -2 * edge12 * edge13 * edge23 + edge01 * edge23 * (edge13 + edge12 - edge23) +
                 edge02 * edge13 * (edge12 + edge23 - edge13) +
                 edge03 * edge12 * (edge13 + edge23 - edge12);
  Real weight1 = -2 * edge02 * edge03 * edge23 + edge01 * edge23 * (edge02 + edge03 - edge23) +
                 edge13 * edge02 * (edge03 + edge23 - edge02) +
                 edge12 * edge03 * (edge02 + edge23 - edge03);
  Real weight2 = -2 * edge01 * edge03 * edge13 + edge02 * edge13 * (edge01 + edge03 - edge13) +
                 edge23 * edge01 * (edge03 + edge13 - edge01) +
                 edge12 * edge03 * (edge01 + edge13 - edge03);
  Real weight3 = -2 * edge01 * edge02 * edge12 + edge03 * edge12 * (edge01 + edge02 - edge12) +
                 edge23 * edge01 * (edge02 + edge12 - edge01) +
                 edge13 * edge02 * (edge01 + edge12 - edge02);

  Real sum_weights = weight0 + weight1 + weight2 + weight3;

  // Check to make sure vertices are not coplanar
  if (MooseUtils::isZero(sum_weights))
    mooseError("Cannot evaluate circumcenter. Points should be non-coplanar.");

  Real inv_sum_weights = 1.0 / sum_weights;

  // Barycentric coordinates
  Real b0 = weight0 * inv_sum_weights;
  Real b1 = weight1 * inv_sum_weights;
  Real b2 = weight2 * inv_sum_weights;
  Real b3 = weight3 * inv_sum_weights;

  return MathUtils::barycentricToCartesian3D(p0, p1, p2, p3, b0, b1, b2, b3);
}

clamp()

template<typename T , typename T2 >

T MathUtils::clamp	(	const T &	x,
		T2	lowerlimit,
		T2	upperlimit
	)

Definition at line 309 of file MathUtils.h.

{
  if (x < lowerlimit)
    return lowerlimit;
  if (x > upperlimit)
    return upperlimit;
  return x;
}

dotProduct()

template<typename T , typename T2 , template< typename > class W, template< typename > class W2, typename std::enable_if< std::is_same< typename W< T >::index_type, unsigned int >::value &&std::is_same< typename W2< T2 >::index_type, unsigned int >::value, int >::type = 0>

libMesh::CompareTypes< T, T2 >::supertype MathUtils::dotProduct	(	const W< T > &	a,
		const W2< T2 > &	b
	)

Definition at line 198 of file MathUtils.h.

Referenced by ADKernelGradTempl< T >::computeResidual(), and ADKernelGradTempl< T >::computeResidualsForJacobian().

{
  return a * b;
}

euclideanMod()

template<typename T1 , typename T2 >

std::size_t MathUtils::euclideanMod	(	T1	dividend,
		T2	divisor
	)

perform modulo operator for Euclidean division that ensures a non-negative result

Parameters

dividend	dividend of the modulo operation
divisor	divisor of the modulo operation

Returns

the non-negative remainder when the dividend is divided by the divisor

Definition at line 434 of file MathUtils.h.

{
  return (dividend % divisor + divisor) % divisor;
}

FactorizedRankTwoTensorOperatorMapBinary()

MathUtils::FactorizedRankTwoTensorOperatorMapBinary	(	pow	,
		std::pow(eigval, arg)
	)

FactorizedRankTwoTensorOperatorMapDerivativeBinary()

MathUtils::FactorizedRankTwoTensorOperatorMapDerivativeBinary	(	dpow	,
		std::pow(eigval, arg)	,
		arg *	std::poweigval, arg - 1
	)

FactorizedRankTwoTensorOperatorMapDerivativeUnary() [1/4]

MathUtils::FactorizedRankTwoTensorOperatorMapDerivativeUnary	(	dlog	,
		std::log(eigval)	,
		1/	eigval
	)

FactorizedRankTwoTensorOperatorMapDerivativeUnary() [2/4]

MathUtils::FactorizedRankTwoTensorOperatorMapDerivativeUnary	(	dexp	,
		std::exp(eigval)	,
		std::exp(eigval)
	)

FactorizedRankTwoTensorOperatorMapDerivativeUnary() [3/4]

MathUtils::FactorizedRankTwoTensorOperatorMapDerivativeUnary	(	dsqrt	,
		std::sqrt(eigval)	,
		std::pow(eigval, -1./2.)/	2.
	)

FactorizedRankTwoTensorOperatorMapDerivativeUnary() [4/4]

MathUtils::FactorizedRankTwoTensorOperatorMapDerivativeUnary	(	dcbrt	,
		std::cbrt(eigval)	,
		std::pow(eigval, -2./3.)/	3.
	)

FactorizedRankTwoTensorOperatorMapUnary() [1/4]

MathUtils::FactorizedRankTwoTensorOperatorMapUnary	(	log	,
		std::log(eigval)
	)

FactorizedRankTwoTensorOperatorMapUnary() [2/4]

MathUtils::FactorizedRankTwoTensorOperatorMapUnary	(	exp	,
		std::exp(eigval)
	)

FactorizedRankTwoTensorOperatorMapUnary() [3/4]

MathUtils::FactorizedRankTwoTensorOperatorMapUnary	(	sqrt	,
		std::sqrt(eigval)
	)

FactorizedRankTwoTensorOperatorMapUnary() [4/4]

MathUtils::FactorizedRankTwoTensorOperatorMapUnary	(	cbrt	,
		std::cbrt(eigval)
	)

gradName()

template<typename T >

T MathUtils::gradName

(

const T &

base_prop_name

)

automatic prefixing for naming material properties based on gradients of coupled variables/functors

Definition at line 445 of file MathUtils.h.

Referenced by GenericFunctorMaterialTempl< is_ad >::GenericFunctorMaterialTempl().

{
  return "grad_" + base_prop_name;
}

heavyside()

template<typename T >

T MathUtils::heavyside

(

T

x

)

Definition at line 120 of file MathUtils.h.

{
  return x < 0.0 ? 0.0 : 1.0;
}

kron()

void MathUtils::kron	(	RealEigenMatrix &	product,
		const RealEigenMatrix &	mat_A,
		const RealEigenMatrix &	mat_B
	)

Computes the Kronecker product of two matrices.

Parameters

product	Reference to the product matrix
mat_A	Reference to the first matrix
mat_B	Reference to the other matrix

Definition at line 17 of file MathUtils.C.

{
  product.resize(mat_A.rows() * mat_B.rows(), mat_A.cols() * mat_B.cols());
  for (unsigned int i = 0; i < mat_A.rows(); i++)
    for (unsigned int j = 0; j < mat_A.cols(); j++)
      for (unsigned int k = 0; k < mat_B.rows(); k++)
        for (unsigned int l = 0; l < mat_B.cols(); l++)
          product(((i * mat_B.rows()) + k), ((j * mat_B.cols()) + l)) = mat_A(i, j) * mat_B(k, l);
}

linearInterpolation()

template<ComputeType compute_type, typename X , typename X1 , typename X2 , typename Y1 , typename Y2 >

auto MathUtils::linearInterpolation	(	const X &	x,
		const X1 &	x1,
		const X2 &	x2,
		const Y1 &	y1,
		const Y2 &	y2
	)

Definition at line 417 of file MathUtils.h.

{
  const auto m = (y2 - y1) / (x2 - x1);
  if constexpr (compute_type == ComputeType::derivative)
    return m;
  if constexpr (compute_type == ComputeType::value)
    return m * (x - x1) + y1;
}

mooseSetToZero() [1/3]

template<typename T >

void MathUtils::mooseSetToZero ( T &  v )

inline

Helper function templates to set a variable to zero.

Specializations may have to be implemented (for examples see RankTwoTensor, RankFourTensor).

The default for non-pointer types is to assign zero. This should either do something sensible, or throw a compiler error. Otherwise the T type is designed badly.

Definition at line 372 of file MathUtils.h.

Referenced by MaterialBase::getGenericZeroMaterialProperty(), MaterialPropertyInterface::getGenericZeroMaterialProperty(), and MaterialBase::getGenericZeroMaterialPropertyByName().

{
  v = 0;
}

mooseSetToZero() [2/3]

template<typename T >

void MathUtils::mooseSetToZero ( T *&  )

inline

Definition at line 383 of file MathUtils.h.

{
  mooseError("mooseSetToZero does not accept pointers");
}

mooseSetToZero() [3/3]

template<>

void MathUtils::mooseSetToZero ( std::vector< Real > &  vec )

inline

Definition at line 390 of file MathUtils.h.

{
  for (auto & v : vec)
    v = 0.;
}

mooseSetToZero< ADRankFourTensor >()

template<>

void MathUtils::mooseSetToZero< ADRankFourTensor >

(

ADRankFourTensor &

v

)

Definition at line 25 of file RankFourTensor.C.

{
  v.zero();
}

mooseSetToZero< ADRankThreeTensor >()

template<>

void MathUtils::mooseSetToZero< ADRankThreeTensor >

(

ADRankThreeTensor &

v

)

Definition at line 26 of file RankThreeTensor.C.

{
  v.zero();
}

mooseSetToZero< ADRankTwoTensor >()

template<>

void MathUtils::mooseSetToZero< ADRankTwoTensor >

(

ADRankTwoTensor &

v

)

Helper function template specialization to set an object to zero.

Needed by DerivativeMaterialInterface

Definition at line 26 of file RankTwoTensor.C.

{
  v.zero();
}

mooseSetToZero< ADSymmetricRankFourTensor >()

template<>

void MathUtils::mooseSetToZero< ADSymmetricRankFourTensor >

(

ADSymmetricRankFourTensor &

v

)

Definition at line 25 of file SymmetricRankFourTensor.C.

{
  v.zero();
}

mooseSetToZero< ADSymmetricRankTwoTensor >()

template<>

void MathUtils::mooseSetToZero< ADSymmetricRankTwoTensor >

(

ADSymmetricRankTwoTensor &

v

)

Helper function template specialization to set an object to zero.

Needed by DerivativeMaterialInterface

Definition at line 26 of file SymmetricRankTwoTensor.C.

{
  v.zero();
}

mooseSetToZero< RankFourTensor >()

template<>

void MathUtils::mooseSetToZero< RankFourTensor >

(

RankFourTensor &

v

)

Helper function template specialization to set an object to zero.

Needed by DerivativeMaterialInterface

Definition at line 19 of file RankFourTensor.C.

{
  v.zero();
}

mooseSetToZero< RankFourTensorTempl< ADReal > >()

template<>

void MathUtils::mooseSetToZero< RankFourTensorTempl< ADReal > >

(

RankFourTensorTempl< ADReal > &

v

)

mooseSetToZero< RankFourTensorTempl< Real > >()

template<>

void MathUtils::mooseSetToZero< RankFourTensorTempl< Real > >

(

RankFourTensorTempl< Real > &

v

)

mooseSetToZero< RankThreeTensor >()

template<>

void MathUtils::mooseSetToZero< RankThreeTensor >

(

RankThreeTensor &

v

)

Definition at line 19 of file RankThreeTensor.C.

{
  v.zero();
}

mooseSetToZero< RankThreeTensorTempl< ADReal > >()

template<>

void MathUtils::mooseSetToZero< RankThreeTensorTempl< ADReal > >

(

RankThreeTensorTempl< ADReal > &

v

)

mooseSetToZero< RankThreeTensorTempl< Real > >()

template<>

void MathUtils::mooseSetToZero< RankThreeTensorTempl< Real > >

(

RankThreeTensorTempl< Real > &

v

)

Helper function template specialization to set an object to zero.

Needed by DerivativeMaterialInterface

mooseSetToZero< RankTwoTensor >()

template<>

void MathUtils::mooseSetToZero< RankTwoTensor >

(

RankTwoTensor &

v

)

Helper function template specialization to set an object to zero.

Needed by DerivativeMaterialInterface

Definition at line 19 of file RankTwoTensor.C.

{
  v.zero();
}

mooseSetToZero< SymmetricRankFourTensor >()

template<>

void MathUtils::mooseSetToZero< SymmetricRankFourTensor >

(

SymmetricRankFourTensor &

v

)

Helper function template specialization to set an object to zero.

Needed by DerivativeMaterialInterface

Definition at line 19 of file SymmetricRankFourTensor.C.

{
  v.zero();
}

mooseSetToZero< SymmetricRankFourTensorTempl< ADReal > >()

template<>

void MathUtils::mooseSetToZero< SymmetricRankFourTensorTempl< ADReal > >

(

SymmetricRankFourTensorTempl< ADReal > &

v

)

mooseSetToZero< SymmetricRankFourTensorTempl< Real > >()

template<>

void MathUtils::mooseSetToZero< SymmetricRankFourTensorTempl< Real > >

(

SymmetricRankFourTensorTempl< Real > &

v

)

mooseSetToZero< SymmetricRankTwoTensor >()

template<>

void MathUtils::mooseSetToZero< SymmetricRankTwoTensor >

(

SymmetricRankTwoTensor &

v

)

Helper function template specialization to set an object to zero.

Needed by DerivativeMaterialInterface

Definition at line 19 of file SymmetricRankTwoTensor.C.

{
  v.zero();
}

multiIndex()

std::vector< std::vector< unsigned int > > MathUtils::multiIndex	(	unsigned int	dim,
		unsigned int	order
	)

generate a complete multi index table for given dimension and order i.e.

given dim = 2, order = 2, generated table will have the following content 0 0 1 0 0 1 2 0 1 1 0 2 The first number in each entry represents the order of the first variable, i.e. x; The second number in each entry represents the order of the second variable, i.e. y. Multiplication is implied between numbers in each entry, i.e. 1 1 represents x^1 * y^1

Parameters

dim	dimension of the multi-index, here dim = mesh dimension
order	generate the multi-index up to certain order

Returns

a data structure holding entries representing the complete multi index

Definition at line 186 of file MathUtils.C.

{
  // first row all zero
  std::vector<std::vector<unsigned int>> multi_index;
  std::vector<std::vector<unsigned int>> n_choose_k;
  std::vector<unsigned int> row(dim, 0);
  multi_index.push_back(row);

  if (dim == 1)
    for (unsigned int q = 1; q <= order; q++)
    {
      row[0] = q;
      multi_index.push_back(row);
    }
  else
    for (unsigned int q = 1; q <= order; q++)
    {
      n_choose_k = multiIndexHelper(dim + q - 1, dim - 1);
      for (unsigned int r = 0; r < n_choose_k.size(); r++)
      {
        row.clear();
        for (unsigned int c = 1; c < n_choose_k[0].size(); c++)
          row.push_back(n_choose_k[r][c] - n_choose_k[r][c - 1] - 1);
        multi_index.push_back(row);
      }
    }

  return multi_index;
}

negativePart()

template<typename T >

T MathUtils::negativePart

(

T

x

)

Definition at line 157 of file MathUtils.h.

{
  return x < 0.0 ? x : 0.0;
}

plainLog()

Real MathUtils::plainLog	(	Real	x,
		unsigned int	derivative_order
	)

Definition at line 28 of file MathUtils.C.

Referenced by poly1Log(), poly2Log(), poly3Log(), and poly4Log().

{
  switch (derivative_order)
  {
    case 0:
      return std::log(x);

    case 1:
      return 1.0 / x;

    case 2:
      return -1.0 / (x * x);

    case 3:
      return 2.0 / (x * x * x);
  }

  mooseError("Unsupported derivative order ", derivative_order);
}

poly()

template<typename C , typename T , typename R = typename libMesh::CompareTypes<typename C::value_type, T>::supertype>

R MathUtils::poly	(	const C &	c,
		const T	x,
		const bool	derivative = false
	)

Evaluate a polynomial with the coefficients c at x.

Note that the Polynomial form is c[0]*x^s + c[1]*x^(s-1) + c[2]*x^(s-2) + ... + c[s-2]*x^2 + c[s-1]*x + c[s] where s = c.size()-1 , which is counter intuitive!

This function will be DEPRECATED soon (10/22/2020)

The coefficient container type can be any container that provides an index operator [] and a .size() method (e.g. std::vector, std::array). The return type is the supertype of the container value type and the argument x. The supertype is the type that can represent both number types.

Definition at line 237 of file MathUtils.h.

Referenced by MortarSegmentHelper::isDisjoint().

{
  const auto size = c.size();
  if (size == 0)
    return 0.0;

  R value = c[0];
  if (derivative)
  {
    value *= size - 1;
    for (std::size_t i = 1; i < size - 1; ++i)
      value = value * x + c[i] * (size - i - 1);
  }
  else
  {
    for (std::size_t i = 1; i < size; ++i)
      value = value * x + c[i];
  }

  return value;
}

poly1Log()

Real MathUtils::poly1Log	(	Real	x,
		Real	tol,
		unsigned int	derivative_order
	)

Definition at line 49 of file MathUtils.C.

{
  if (x >= tol)
    return plainLog(x, derivative_order);

  const auto c1 = [&]() { return 1.0 / tol; };
  const auto c2 = [&]() { return std::log(tol) - 1.0; };

  switch (derivative_order)
  {
    case 0:
      return c1() * x + c2();

    case 1:
      return c1();

    case 2:
      return 0.0;

    case 3:
      return 0.0;
  }

  mooseError("Unsupported derivative order ", derivative_order);
}

poly2Log()

Real MathUtils::poly2Log	(	Real	x,
		Real	tol,
		unsigned int	derivative_order
	)

Definition at line 76 of file MathUtils.C.

{
  if (x >= tol)
    return plainLog(x, derivative_order);

  const auto c1 = [&]() { return -0.5 / (tol * tol); };
  const auto c2 = [&]() { return 2.0 / tol; };
  const auto c3 = [&]() { return std::log(tol) - 3.0 / 2.0; };

  switch (derivative_order)
  {
    case 0:
      return c1() * x * x + c2() * x + c3();

    case 1:
      return 2.0 * c1() * x + c2();

    case 2:
      return 2.0 * c1();

    case 3:
      return 0.0;
  }

  mooseError("Unsupported derivative order ", derivative_order);
}

poly3Log()

Real MathUtils::poly3Log	(	Real	x,
		Real	tol,
		unsigned int	derivative_order
	)

Definition at line 104 of file MathUtils.C.

{
  if (x >= tol)
    return plainLog(x, derivative_order);

  const auto c1 = [&]() { return 1.0 / (3.0 * tol * tol * tol); };
  const auto c2 = [&]() { return -3.0 / (2.0 * tol * tol); };
  const auto c3 = [&]() { return 3.0 / tol; };
  const auto c4 = [&]() { return std::log(tol) - 11.0 / 6.0; };

  switch (derivative_order)
  {
    case 0:
      return c1() * x * x * x + c2() * x * x + c3() * x + c4();

    case 1:
      return 3.0 * c1() * x * x + 2.0 * c2() * x + c3();

    case 2:
      return 6.0 * c1() * x + 2.0 * c2();

    case 3:
      return 6.0 * c1();
  }

  mooseError("Unsupported derivative order ", derivative_order);
}

poly4Log()

Real MathUtils::poly4Log	(	Real	x,
		Real	tol,
		unsigned int	derivative_order
	)

Definition at line 133 of file MathUtils.C.

{
  if (x >= tol)
    return plainLog(x, derivative_order);

  switch (derivative_order)
  {
    case 0:
      return std::log(tol) + (x - tol) / tol -
             Utility::pow<2>(x - tol) / (2.0 * Utility::pow<2>(tol)) +
             Utility::pow<3>(x - tol) / (3.0 * Utility::pow<3>(tol)) -
             Utility::pow<4>(x - tol) / (4.0 * Utility::pow<4>(tol)) +
             Utility::pow<5>(x - tol) / (5.0 * Utility::pow<5>(tol)) -
             Utility::pow<6>(x - tol) / (6.0 * Utility::pow<6>(tol));

    case 1:
      return 1.0 / tol - (x - tol) / Utility::pow<2>(tol) +
             Utility::pow<2>(x - tol) / Utility::pow<3>(tol) -
             Utility::pow<3>(x - tol) / Utility::pow<4>(tol) +
             Utility::pow<4>(x - tol) / Utility::pow<5>(tol) -
             Utility::pow<5>(x - tol) / Utility::pow<6>(tol);

    case 2:
      return -1.0 / Utility::pow<2>(tol) + 2.0 * (x - tol) / Utility::pow<3>(tol) -
             3.0 * Utility::pow<2>(x - tol) / Utility::pow<4>(tol) +
             4.0 * Utility::pow<3>(x - tol) / Utility::pow<5>(tol) -
             5.0 * Utility::pow<4>(x - tol) / Utility::pow<6>(tol);

    case 3:
      return 2.0 / Utility::pow<3>(tol) - 6.0 * (x - tol) / Utility::pow<4>(tol) +
             12.0 * Utility::pow<2>(x - tol) / Utility::pow<5>(tol) -
             20.0 * Utility::pow<3>(x - tol) / Utility::pow<6>(tol);
  }

  mooseError("Unsupported derivative order ", derivative_order);
}

polynomial()

template<typename C , typename T , typename R = typename libMesh::CompareTypes<typename C::value_type, T>::supertype>

R MathUtils::polynomial	(	const C &	c,
		const T	x
	)

Evaluate a polynomial with the coefficients c at x.

Note that the Polynomial form is c[0] + c[1] * x + c[2] * x^2 + ... The coefficient container type can be any container that provides an index operator [] and a .size() method (e.g. std::vector, std::array). The return type is the supertype of the container value type and the argument x. The supertype is the type that can represent both number types.

Definition at line 272 of file MathUtils.h.

Referenced by NodalPatchRecovery::computePVector(), NodalPatchRecoveryBase::evaluateBasisFunctions(), and ProjectedStatefulMaterialNodalPatchRecoveryTempl< T, is_ad >::evaluateBasisFunctions().

{
  auto size = c.size();
  if (size == 0)
    return 0.0;

  size--;
  R value = c[size];
  for (std::size_t i = 1; i <= size; ++i)
    value = value * x + c[size - i];

  return value;
}

polynomialDerivative()

template<typename C , typename T , typename R = typename libMesh::CompareTypes<typename C::value_type, T>::supertype>

R MathUtils::polynomialDerivative	(	const C &	c,
		const T	x
	)

Returns the derivative of polynomial(c, x) with respect to x.

Definition at line 293 of file MathUtils.h.

{
  auto size = c.size();
  if (size <= 1)
    return 0.0;

  size--;
  R value = c[size] * size;
  for (std::size_t i = 1; i < size; ++i)
    value = value * x + c[size - i] * (size - i);

  return value;
}

positivePart()

template<typename T >

T MathUtils::positivePart

(

T

x

)

Definition at line 150 of file MathUtils.h.

{
  return x > 0.0 ? x : 0.0;
}

pow()

template<typename T >

T MathUtils::pow	(	T	x,
		int	e
	)

Definition at line 91 of file MathUtils.h.

Referenced by ProjectedStatefulMaterialNodalPatchRecoveryTempl< T, is_ad >::evaluateBasisFunctions(), PolynomialFit::fillMatrix(), BicubicInterpolation::sample2ndDerivative(), BicubicInterpolation::sampleDerivative(), BicubicInterpolation::sampleInternal(), BicubicInterpolation::sampleValueAndDerivativesInternal(), and PseudoTimestep::timestepEXP().

{
  bool neg = false;
  T result = 1.0;

  if (e < 0)
  {
    neg = true;
    e = -e;
  }

  while (e)
  {
    // if bit 0 is set multiply the current power of two factor of the exponent
    if (e & 1)
      result *= x;

    // x is incrementally set to consecutive powers of powers of two
    x *= x;

    // bit shift the exponent down
    e >>= 1;
  }

  return neg ? 1.0 / result : result;
}

regularizedHeavyside()

template<typename T >

T MathUtils::regularizedHeavyside	(	T	x,
		Real	smoothing_length
	)

Definition at line 127 of file MathUtils.h.

{
  if (x <= -smoothing_length)
    return 0.0;
  else if (x < smoothing_length)
    return 0.5 * (1 + std::sin(libMesh::pi * x / 2 / smoothing_length));
  else
    return 1.0;
}

regularizedHeavysideDerivative()

template<typename T >

T MathUtils::regularizedHeavysideDerivative	(	T	x,
		Real	smoothing_length
	)

Definition at line 139 of file MathUtils.h.

{
  if (x < smoothing_length && x > -smoothing_length)
    return 0.25 * libMesh::pi / smoothing_length *
           (std::cos(libMesh::pi * x / 2 / smoothing_length));
  else
    return 0.0;
}

round()

template<typename T >

T MathUtils::round

(

T

x

)

Definition at line 77 of file MathUtils.h.

Referenced by ParsedSubdomainIDsGenerator::assignElemSubdomainID(), GeneratedMesh::buildMesh(), VariableValueElementSubdomainModifier::computeSubdomainID(), GeneratedMeshGenerator::generate(), ImageSubdomainGenerator::generate(), DistributedRectilinearMeshGenerator::generate(), PointIndexedMap::getRoundedPoint(), FlagElementsThread::onElement(), and ProgressOutput::output().

{
  return ::round(x); // use round from math.h
}

sign()

template<typename T >

T MathUtils::sign

(

T

x

)

Definition at line 84 of file MathUtils.h.

Referenced by AutomaticMortarGeneration::computeNodalGeometry(), HFEMTestJump::computeQpOffDiagJacobian(), HFEMTrialJump::computeQpOffDiagJacobian(), HFEMTrialJump::computeQpResidual(), HFEMTestJump::computeQpResidual(), MooseUnits::parse(), Moose::FV::rF(), VectorFunctionReaction::validParams(), and PiecewiseConstant::value().

{
  return x >= 0.0 ? 1.0 : -1.0;
}

smootherStep() [1/2]

template<typename T , typename T2 >

T MathUtils::smootherStep	(	T	x,
		T2	start,
		T2	end,
		bool	derivative = false
	)

Definition at line 320 of file MathUtils.h.

{
  mooseAssert("start < end", "Start value must be lower than end value for smootherStep");
  if (x <= start)
    return 0.0;
  else if (x >= end)
  {
    if (derivative)
      return 0.0;
    else
      return 1.0;
  }
  x = (x - start) / (end - start);
  if (derivative)
    return 30.0 * libMesh::Utility::pow<2>(x) * (x * (x - 2.0) + 1.0) / (end - start);
  return libMesh::Utility::pow<3>(x) * (x * (x * 6.0 - 15.0) + 10.0);
}

smootherStep() [2/2]

template<ComputeType compute_type, typename X , typename S , typename E >

auto MathUtils::smootherStep	(	const X &	x,
		const S &	start,
		const E &	end
	)

Definition at line 346 of file MathUtils.h.

{
  mooseAssert("start < end", "Start value must be lower than end value for smootherStep");
  if (x <= start)
    return 0.0;
  else if (x >= end)
  {
    if constexpr (compute_type == ComputeType::derivative)
      return 0.0;
    if constexpr (compute_type == ComputeType::value)
      return 1.0;
  }
  const auto u = (x - start) / (end - start);
  if constexpr (compute_type == ComputeType::derivative)
    return 30.0 * libMesh::Utility::pow<2>(u) * (u * (u - 2.0) + 1.0) / (end - start);
  if constexpr (compute_type == ComputeType::value)
    return libMesh::Utility::pow<3>(u) * (u * (u * 6.0 - 15.0) + 10.0);
}

taylorLog()

Real MathUtils::taylorLog

(

Real

x

)

Todo:

This can be done without std::pow!

Definition at line 172 of file MathUtils.C.

{
  Real y = (x - 1.0) / (x + 1.0);
  Real val = 1.0;
  for (unsigned int i = 0; i < 5; ++i)
  {
    Real exponent = i + 2.0;
    val += 1.0 / (2.0 * (i + 1.0) + 1.0) * std::pow(y, exponent);
  }

  return val * 2.0 * y;
}

timeDerivName()

template<typename T >

T MathUtils::timeDerivName

(

const T &

base_prop_name

)

automatic prefixing for naming material properties based on time derivatives of coupled variables/functors

Definition at line 456 of file MathUtils.h.

Referenced by GenericFunctorMaterialTempl< is_ad >::GenericFunctorMaterialTempl(), and GenericVectorFunctorMaterialTempl< is_ad >::GenericVectorFunctorMaterialTempl().

{
  return "d" + base_prop_name + "_dt";
}

Variable Documentation

sqrt2

constexpr Real MathUtils::sqrt2 = 1.4142135623730951

static

std::sqrt is not constexpr, so we add sqrt(2) as a constant (used in Mandel notation)

Definition at line 25 of file MathUtils.h.

Referenced by SymmetricRankTwoTensorTempl< Real >::d2thirdInvariant(), SymmetricRankTwoTensorTempl< Real >::ddet(), SymmetricRankTwoTensorTempl< Real >::det(), SymmetricRankTwoTensorTempl< Real >::dthirdInvariant(), SymmetricRankTwoTensorTempl< Real >::inverse(), SymmetricRankTwoTensorTempl< Real >::mandelFactor(), SymmetricRankFourTensorTempl< T >::rotationMatrix(), SymmetricRankTwoTensorTempl< Real >::thirdInvariant(), and SymmetricRankTwoTensorTempl< Real >::timesTranspose().