/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see .
\*---------------------------------------------------------------------------*/
#include "SymmTensor.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template
inline Tensor::Tensor()
{}
template
template
inline Tensor::Tensor(const VectorSpace, Cmpt2, 9>& vs)
:
VectorSpace, Cmpt, 9>(vs)
{}
template
inline Tensor::Tensor(const SphericalTensor& st)
{
this->v_[XX] = st.ii(); this->v_[XY] = 0; this->v_[XZ] = 0;
this->v_[YX] = 0; this->v_[YY] = st.ii(); this->v_[YZ] = 0;
this->v_[ZX] = 0; this->v_[ZY] = 0; this->v_[ZZ] = st.ii();
}
template
inline Tensor::Tensor(const SymmTensor& st)
{
this->v_[XX] = st.xx(); this->v_[XY] = st.xy(); this->v_[XZ] = st.xz();
this->v_[YX] = st.xy(); this->v_[YY] = st.yy(); this->v_[YZ] = st.yz();
this->v_[ZX] = st.xz(); this->v_[ZY] = st.yz(); this->v_[ZZ] = st.zz();
}
template
inline Tensor::Tensor(const Vector>& tr)
{
this->v_[XX] = tr.x().x();
this->v_[XY] = tr.x().y();
this->v_[XZ] = tr.x().z();
this->v_[YX] = tr.y().x();
this->v_[YY] = tr.y().y();
this->v_[YZ] = tr.y().z();
this->v_[ZX] = tr.z().x();
this->v_[ZY] = tr.z().y();
this->v_[ZZ] = tr.z().z();
}
template
inline Tensor::Tensor
(
const Vector& x,
const Vector& y,
const Vector& z
)
{
this->v_[XX] = x.x(); this->v_[XY] = x.y(); this->v_[XZ] = x.z();
this->v_[YX] = y.x(); this->v_[YY] = y.y(); this->v_[YZ] = y.z();
this->v_[ZX] = z.x(); this->v_[ZY] = z.y(); this->v_[ZZ] = z.z();
}
template
inline Tensor::Tensor
(
const Cmpt txx, const Cmpt txy, const Cmpt txz,
const Cmpt tyx, const Cmpt tyy, const Cmpt tyz,
const Cmpt tzx, const Cmpt tzy, const Cmpt tzz
)
{
this->v_[XX] = txx; this->v_[XY] = txy; this->v_[XZ] = txz;
this->v_[YX] = tyx; this->v_[YY] = tyy; this->v_[YZ] = tyz;
this->v_[ZX] = tzx; this->v_[ZY] = tzy; this->v_[ZZ] = tzz;
}
template
inline Tensor::Tensor(Istream& is)
:
VectorSpace, Cmpt, 9>(is)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template
inline Vector Tensor::x() const
{
return Vector(this->v_[XX], this->v_[XY], this->v_[XZ]);
}
template
inline Vector Tensor::y() const
{
return Vector(this->v_[YX], this->v_[YY], this->v_[YZ]);
}
template
inline Vector Tensor::z() const
{
return Vector(this->v_[ZX], this->v_[ZY], this->v_[ZZ]);
}
template
inline Vector Tensor::vectorComponent(const direction cmpt) const
{
switch (cmpt)
{
case 0:
return x();
break;
case 1:
return y();
break;
case 2:
return z();
break;
}
}
template
inline const Cmpt& Tensor::xx() const
{
return this->v_[XX];
}
template
inline const Cmpt& Tensor::xy() const
{
return this->v_[XY];
}
template
inline const Cmpt& Tensor::xz() const
{
return this->v_[XZ];
}
template
inline const Cmpt& Tensor::yx() const
{
return this->v_[YX];
}
template
inline const Cmpt& Tensor::yy() const
{
return this->v_[YY];
}
template
inline const Cmpt& Tensor::yz() const
{
return this->v_[YZ];
}
template
inline const Cmpt& Tensor::zx() const
{
return this->v_[ZX];
}
template
inline const Cmpt& Tensor::zy() const
{
return this->v_[ZY];
}
template
inline const Cmpt& Tensor::zz() const
{
return this->v_[ZZ];
}
template
inline Cmpt& Tensor::xx()
{
return this->v_[XX];
}
template
inline Cmpt& Tensor::xy()
{
return this->v_[XY];
}
template
inline Cmpt& Tensor::xz()
{
return this->v_[XZ];
}
template
inline Cmpt& Tensor::yx()
{
return this->v_[YX];
}
template
inline Cmpt& Tensor::yy()
{
return this->v_[YY];
}
template
inline Cmpt& Tensor::yz()
{
return this->v_[YZ];
}
template
inline Cmpt& Tensor::zx()
{
return this->v_[ZX];
}
template
inline Cmpt& Tensor::zy()
{
return this->v_[ZY];
}
template
inline Cmpt& Tensor::zz()
{
return this->v_[ZZ];
}
template
inline Tensor Tensor::T() const
{
return Tensor
(
xx(), yx(), zx(),
xy(), yy(), zy(),
xz(), yz(), zz()
);
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template
inline void Tensor::operator&=(const Tensor& t)
{
*this =
(
Tensor
(
this->xx()*t.xx() + this->xy()*t.yx() + this->xz()*t.zx(),
this->xx()*t.xy() + this->xy()*t.yy() + this->xz()*t.zy(),
this->xx()*t.xz() + this->xy()*t.yz() + this->xz()*t.zz(),
this->yx()*t.xx() + this->yy()*t.yx() + this->yz()*t.zx(),
this->yx()*t.xy() + this->yy()*t.yy() + this->yz()*t.zy(),
this->yx()*t.xz() + this->yy()*t.yz() + this->yz()*t.zz(),
this->zx()*t.xx() + this->zy()*t.yx() + this->zz()*t.zx(),
this->zx()*t.xy() + this->zy()*t.yy() + this->zz()*t.zy(),
this->zx()*t.xz() + this->zy()*t.yz() + this->zz()*t.zz()
)
);
}
template
inline void Tensor::operator=(const SphericalTensor& st)
{
this->v_[XX] = st.ii(); this->v_[XY] = 0; this->v_[XZ] = 0;
this->v_[YX] = 0; this->v_[YY] = st.ii(); this->v_[YZ] = 0;
this->v_[ZX] = 0; this->v_[ZY] = 0; this->v_[ZZ] = st.ii();
}
template
inline void Tensor::operator=(const SymmTensor& st)
{
this->v_[XX] = st.xx(); this->v_[XY] = st.xy(); this->v_[XZ] = st.xz();
this->v_[YX] = st.xy(); this->v_[YY] = st.yy(); this->v_[YZ] = st.yz();
this->v_[ZX] = st.xz(); this->v_[ZY] = st.yz(); this->v_[ZZ] = st.zz();
}
template
inline void Tensor::operator=(const Vector>& tr)
{
this->v_[XX] = tr.x().x();
this->v_[XY] = tr.x().y();
this->v_[XZ] = tr.x().z();
this->v_[YX] = tr.y().x();
this->v_[YY] = tr.y().y();
this->v_[YZ] = tr.y().z();
this->v_[ZX] = tr.z().x();
this->v_[ZY] = tr.z().y();
this->v_[ZZ] = tr.z().z();
}
// * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
template
inline Vector operator*(const Tensor& t)
{
return Vector(t.yz(), -t.xz(), t.xy());
}
template
inline Tensor operator*(const Vector& v)
{
return Tensor
(
0, -v.z(), v.y(),
v.z(), 0, -v.x(),
-v.y(), v.x(), 0
);
}
template
inline typename innerProduct, Tensor>::type
operator&(const Tensor& t1, const Tensor& t2)
{
return Tensor
(
t1.xx()*t2.xx() + t1.xy()*t2.yx() + t1.xz()*t2.zx(),
t1.xx()*t2.xy() + t1.xy()*t2.yy() + t1.xz()*t2.zy(),
t1.xx()*t2.xz() + t1.xy()*t2.yz() + t1.xz()*t2.zz(),
t1.yx()*t2.xx() + t1.yy()*t2.yx() + t1.yz()*t2.zx(),
t1.yx()*t2.xy() + t1.yy()*t2.yy() + t1.yz()*t2.zy(),
t1.yx()*t2.xz() + t1.yy()*t2.yz() + t1.yz()*t2.zz(),
t1.zx()*t2.xx() + t1.zy()*t2.yx() + t1.zz()*t2.zx(),
t1.zx()*t2.xy() + t1.zy()*t2.yy() + t1.zz()*t2.zy(),
t1.zx()*t2.xz() + t1.zy()*t2.yz() + t1.zz()*t2.zz()
);
}
template
inline typename innerProduct, Vector>::type
operator&(const Tensor& t, const Vector& v)
{
return Vector
(
t.xx()*v.x() + t.xy()*v.y() + t.xz()*v.z(),
t.yx()*v.x() + t.yy()*v.y() + t.yz()*v.z(),
t.zx()*v.x() + t.zy()*v.y() + t.zz()*v.z()
);
}
template
inline typename innerProduct, Tensor>::type
operator&(const Vector& v, const Tensor& t)
{
return Vector
(
v.x()*t.xx() + v.y()*t.yx() + v.z()*t.zx(),
v.x()*t.xy() + v.y()*t.yy() + v.z()*t.zy(),
v.x()*t.xz() + v.y()*t.yz() + v.z()*t.zz()
);
}
template
inline typename outerProduct, Vector>::type
operator*(const Vector& v1, const Vector& v2)
{
return Tensor
(
v1.x()*v2.x(), v1.x()*v2.y(), v1.x()*v2.z(),
v1.y()*v2.x(), v1.y()*v2.y(), v1.y()*v2.z(),
v1.z()*v2.x(), v1.z()*v2.y(), v1.z()*v2.z()
);
}
template
inline typename innerProduct, Tensor>::type
operator/(const Vector& v, const Tensor& t)
{
return inv(t) & v;
}
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
//- Return the trace of a tensor
template
inline Cmpt tr(const Tensor& t)
{
return t.xx() + t.yy() + t.zz();
}
//- Return the spherical part of a tensor
template
inline SphericalTensor sph(const Tensor& t)
{
return (1.0/3.0)*tr(t);
}
//- Return the symmetric part of a tensor
template
inline SymmTensor symm(const Tensor& t)
{
return SymmTensor
(
t.xx(), 0.5*(t.xy() + t.yx()), 0.5*(t.xz() + t.zx()),
t.yy(), 0.5*(t.yz() + t.zy()),
t.zz()
);
}
//- Return twice the symmetric part of a tensor
template
inline SymmTensor twoSymm(const Tensor& t)
{
return SymmTensor
(
2*t.xx(), (t.xy() + t.yx()), (t.xz() + t.zx()),
2*t.yy(), (t.yz() + t.zy()),
2*t.zz()
);
}
//- Return the skew-symmetric part of a tensor
template
inline Tensor skew(const Tensor& t)
{
return Tensor
(
0.0, 0.5*(t.xy() - t.yx()), 0.5*(t.xz() - t.zx()),
0.5*(t.yx() - t.xy()), 0.0, 0.5*(t.yz() - t.zy()),
0.5*(t.zx() - t.xz()), 0.5*(t.zy() - t.yz()), 0.0
);
}
//- Return the skew-symmetric part of a symmetric tensor
template
inline const Tensor& skew(const SymmTensor& st)
{
return Tensor::zero;
}
//- Return the deviatoric part of a tensor
template
inline Tensor dev(const Tensor& t)
{
return t - SphericalTensor::oneThirdI*tr(t);
}
//- Return the deviatoric part of a tensor
template
inline Tensor dev2(const Tensor& t)
{
return t - SphericalTensor::twoThirdsI*tr(t);
}
//- Return the determinant of a tensor
template
inline Cmpt det(const Tensor& t)
{
return
(
t.xx()*t.yy()*t.zz() + t.xy()*t.yz()*t.zx()
+ t.xz()*t.yx()*t.zy() - t.xx()*t.yz()*t.zy()
- t.xy()*t.yx()*t.zz() - t.xz()*t.yy()*t.zx()
);
}
//- Return the cofactor tensor of a tensor
template
inline Tensor cof(const Tensor& t)
{
return Tensor
(
t.yy()*t.zz() - t.zy()*t.yz(),
t.zx()*t.yz() - t.yx()*t.zz(),
t.yx()*t.zy() - t.yy()*t.zx(),
t.xz()*t.zy() - t.xy()*t.zz(),
t.xx()*t.zz() - t.xz()*t.zx(),
t.xy()*t.zx() - t.xx()*t.zy(),
t.xy()*t.yz() - t.xz()*t.yy(),
t.yx()*t.xz() - t.xx()*t.yz(),
t.xx()*t.yy() - t.yx()*t.xy()
);
}
//- Return the inverse of a tensor given the determinant
template
inline Tensor inv(const Tensor& t, const Cmpt dett)
{
return Tensor
(
t.yy()*t.zz() - t.zy()*t.yz(),
t.xz()*t.zy() - t.xy()*t.zz(),
t.xy()*t.yz() - t.xz()*t.yy(),
t.zx()*t.yz() - t.yx()*t.zz(),
t.xx()*t.zz() - t.xz()*t.zx(),
t.yx()*t.xz() - t.xx()*t.yz(),
t.yx()*t.zy() - t.yy()*t.zx(),
t.xy()*t.zx() - t.xx()*t.zy(),
t.xx()*t.yy() - t.yx()*t.xy()
)/dett;
}
//- Return the inverse of a tensor
template
inline Tensor inv(const Tensor& t)
{
return inv(t, det(t));
}
//- Return the 1st invariant of a tensor
template
inline Cmpt invariantI(const Tensor& t)
{
return tr(t);
}
//- Return the 2nd invariant of a tensor
template
inline Cmpt invariantII(const Tensor& t)
{
return
(
0.5*sqr(tr(t))
- 0.5*
(
t.xx()*t.xx() + t.xy()*t.xy() + t.xz()*t.xz()
+ t.yx()*t.yx() + t.yy()*t.yy() + t.yz()*t.yz()
+ t.zx()*t.zx() + t.zy()*t.zy() + t.zz()*t.zz()
)
);
}
//- Return the 3rd invariant of a tensor
template
inline Cmpt invariantIII(const Tensor& t)
{
return det(t);
}
// * * * * * * * * * Mixed Tensor SphericalTensor Operators * * * * * * * * //
template
inline Tensor
operator+(const SphericalTensor& st1, const Tensor& t2)
{
return Tensor
(
st1.ii() + t2.xx(), t2.xy(), t2.xz(),
t2.yx(), st1.ii() + t2.yy(), t2.yz(),
t2.zx(), t2.zy(), st1.ii() + t2.zz()
);
}
template
inline Tensor
operator+(const Tensor& t1, const SphericalTensor& st2)
{
return Tensor
(
t1.xx() + st2.ii(), t1.xy(), t1.xz(),
t1.yx(), t1.yy() + st2.ii(), t1.yz(),
t1.zx(), t1.zy(), t1.zz() + st2.ii()
);
}
template
inline Tensor
operator-(const SphericalTensor& st1, const Tensor& t2)
{
return Tensor
(
st1.ii() - t2.xx(), -t2.xy(), -t2.xz(),
-t2.yx(), st1.ii() - t2.yy(), -t2.yz(),
-t2.zx(), -t2.zy(), st1.ii() - t2.zz()
);
}
template
inline Tensor
operator-(const Tensor& t1, const SphericalTensor& st2)
{
return Tensor
(
t1.xx() - st2.ii(), t1.xy(), t1.xz(),
t1.yx(), t1.yy() - st2.ii(), t1.yz(),
t1.zx(), t1.zy(), t1.zz() - st2.ii()
);
}
//- Inner-product between a spherical tensor and a tensor
template
inline Tensor
operator&(const SphericalTensor& st1, const Tensor& t2)
{
return Tensor
(
st1.ii()*t2.xx(), st1.ii()*t2.xy(), st1.ii()*t2.xz(),
st1.ii()*t2.yx(), st1.ii()*t2.yy(), st1.ii()*t2.yz(),
st1.ii()*t2.zx(), st1.ii()*t2.zy(), st1.ii()*t2.zz()
);
}
//- Inner-product between a tensor and a spherical tensor
template
inline Tensor
operator&(const Tensor& t1, const SphericalTensor& st2)
{
return Tensor
(
t1.xx()*st2.ii(), t1.xy()*st2.ii(), t1.xz()*st2.ii(),
t1.yx()*st2.ii(), t1.yy()*st2.ii(), t1.yz()*st2.ii(),
t1.zx()*st2.ii(), t1.zy()*st2.ii(), t1.zz()*st2.ii()
);
}
//- Double-dot-product between a spherical tensor and a tensor
template
inline Cmpt
operator&&(const SphericalTensor& st1, const Tensor& t2)
{
return(st1.ii()*t2.xx() + st1.ii()*t2.yy() + st1.ii()*t2.zz());
}
//- Double-dot-product between a tensor and a spherical tensor
template
inline Cmpt
operator&&(const Tensor& t1, const SphericalTensor& st2)
{
return(t1.xx()*st2.ii() + t1.yy()*st2.ii() + t1.zz()*st2.ii());
}
template
class typeOfSum, Tensor>
{
public:
typedef Tensor type;
};
template
class typeOfSum, SphericalTensor>
{
public:
typedef Tensor type;
};
template
class innerProduct, Tensor>
{
public:
typedef Tensor type;
};
template
class innerProduct, SphericalTensor>
{
public:
typedef Tensor type;
};
// * * * * * * * * * * Mixed Tensor SymmTensor Operators * * * * * * * * * * //
template
inline Tensor
operator+(const SymmTensor& st1, const Tensor& t2)
{
return Tensor
(
st1.xx() + t2.xx(), st1.xy() + t2.xy(), st1.xz() + t2.xz(),
st1.xy() + t2.yx(), st1.yy() + t2.yy(), st1.yz() + t2.yz(),
st1.xz() + t2.zx(), st1.yz() + t2.zy(), st1.zz() + t2.zz()
);
}
template
inline Tensor
operator+(const Tensor& t1, const SymmTensor& st2)
{
return Tensor
(
t1.xx() + st2.xx(), t1.xy() + st2.xy(), t1.xz() + st2.xz(),
t1.yx() + st2.xy(), t1.yy() + st2.yy(), t1.yz() + st2.yz(),
t1.zx() + st2.xz(), t1.zy() + st2.yz(), t1.zz() + st2.zz()
);
}
template
inline Tensor
operator-(const SymmTensor& st1, const Tensor& t2)
{
return Tensor
(
st1.xx() - t2.xx(), st1.xy() - t2.xy(), st1.xz() - t2.xz(),
st1.xy() - t2.yx(), st1.yy() - t2.yy(), st1.yz() - t2.yz(),
st1.xz() - t2.zx(), st1.yz() - t2.zy(), st1.zz() - t2.zz()
);
}
template
inline Tensor
operator-(const Tensor& t1, const SymmTensor& st2)
{
return Tensor
(
t1.xx() - st2.xx(), t1.xy() - st2.xy(), t1.xz() - st2.xz(),
t1.yx() - st2.xy(), t1.yy() - st2.yy(), t1.yz() - st2.yz(),
t1.zx() - st2.xz(), t1.zy() - st2.yz(), t1.zz() - st2.zz()
);
}
//- Inner-product between a symmetric tensor and a tensor
template
inline Tensor
operator&(const SymmTensor