54 this->v_[XX] = st.
ii(); this->v_[XY] =
Zero; this->v_[XZ] =
Zero;
55 this->v_[YY] = st.
ii(); this->v_[YZ] =
Zero;
56 this->v_[ZZ] = st.
ii();
76 const Cmpt txx,
const Cmpt txy,
const Cmpt txz,
77 const Cmpt tyy,
const Cmpt tyz,
81 this->v_[XX] = txx; this->v_[XY] = txy; this->v_[XZ] = txz;
82 this->v_[YY] = tyy; this->v_[YZ] = tyz;
99 return Vector<Cmpt>(this->v_[XX], this->v_[XY], this->v_[XZ]);
106 return Vector<Cmpt>(this->v_[XY], this->v_[YY], this->v_[YZ]);
113 return Vector<Cmpt>(this->v_[XZ], this->v_[YZ], this->v_[ZZ]);
118 template<Foam::direction Idx>
121 if (Idx == 0)
return x();
122 else if (Idx == 1)
return y();
123 else if (Idx == 2)
return z();
125 static_assert(Idx < 3,
"Invalid row access");
135 case 0:
return x();
break;
136 case 1:
return y();
break;
137 case 2:
return z();
break;
148 template<Foam::direction Idx>
153 this->v_[XX] = v.
x(); this->v_[XY] = v.
y(); this->v_[XZ] = v.
z();
157 this->v_[XY] = v.
x(); this->v_[YY] = v.
y(); this->v_[YZ] = v.
z();
161 this->v_[XZ] = v.
x(); this->v_[YZ] = v.
y(); this->v_[ZZ] = v.
z();
164 static_assert(Idx < 3,
"Invalid row access");
176 this->v_[XX] =
x.x(); this->v_[XY] =
x.y(); this->v_[XZ] =
x.z();
177 this->v_[YY] =
y.y(); this->v_[YZ] =
y.z();
178 this->v_[ZZ] = z.
z();
191 case 0: row<0>(v);
break;
192 case 1: row<1>(v);
break;
193 case 2: row<2>(v);
break;
204 return Vector<Cmpt>(this->v_[XX], this->v_[YY], this->v_[ZZ]);
211 this->v_[XX] = v.
x(); this->v_[YY] = v.
y(); this->v_[ZZ] = v.
z();
218 this->v_[XX] += v.
x(); this->v_[YY] += v.
y(); this->v_[ZZ] += v.
z();
225 this->v_[XX] -= v.
x(); this->v_[YY] -= v.
y(); this->v_[ZZ] -= v.
z();
246 xx()*yy()*zz() + xy()*yz()*xz()
247 + xz()*xy()*yz() - xx()*yz()*yz()
248 - xy()*xy()*zz() - xz()*yy()*xz()
260 return (yy()*zz() - yz()*yz());
265 return (xx()*zz() - xz()*xz());
270 return (xx()*yy() - xy()*xy());
280 yy()*zz() - yz()*yz(), xz()*yz() - xy()*zz(), xy()*yz() - xz()*yy(),
281 xx()*zz() - xz()*xz(), xy()*xz() - xx()*yz(),
282 xx()*yy() - xy()*xy()
291 return this->adjunct();
315 return SymmTensor<Cmpt>
325 return SymmTensor<Cmpt>
340 const Cmpt detval = this->det2D(excludeCmpt);
342 return this->adjunct2D(excludeCmpt)/detval;
350 const Cmpt detval = this->
det();
353 if (
mag(detval) < VSMALL)
356 <<
"Tensor not properly invertible, determinant:" 357 << detval <<
" tensor:" << *
this <<
nl 362 return this->adjunct()/detval;
379 const scalar threshold = SMALL * (magSqr_xx + magSqr_yy + magSqr_zz);
381 const bool small_xx = (magSqr_xx < threshold);
382 const bool small_yy = (magSqr_yy < threshold);
383 const bool small_zz = (magSqr_zz < threshold);
385 if (small_xx || small_yy || small_zz)
393 const Cmpt detval = work.
det();
395 if (
mag(detval) < ROOTVSMALL)
411 const Cmpt detval = this->
det();
413 if (
mag(detval) < ROOTVSMALL)
416 return SymmTensor<Cmpt>(
Zero);
419 return this->adjunct()/detval;
428 this->v_[XX] = st.
ii(); this->v_[XY] =
Zero; this->v_[XZ] =
Zero;
429 this->v_[YY] = st.
ii(); this->v_[YZ] =
Zero;
430 this->v_[ZZ] = st.
ii();
445 return st.
xx() + st.
yy() + st.
zz();
470 inline SymmTensor<Cmpt>
twoSymm(
const SymmTensor<Cmpt>& st)
478 inline SymmTensor<Cmpt>
devSymm(
const SymmTensor<Cmpt>& st)
486 inline SymmTensor<Cmpt>
devTwoSymm(
const SymmTensor<Cmpt>& st)
504 return st - 2*
sph(st);
510 inline Cmpt
det(
const SymmTensor<Cmpt>& st)
518 inline SymmTensor<Cmpt>
cof(
const SymmTensor<Cmpt>& st)
526 inline SymmTensor<Cmpt>
inv(
const SymmTensor<Cmpt>& st,
const Cmpt detval)
529 if (
mag(detval) < VSMALL)
532 <<
"SymmTensor not properly invertible, determinant:" 533 << detval <<
" tensor:" << st <<
nl 538 return st.adjunct()/detval;
552 inline Cmpt
invariantI(
const SymmTensor<Cmpt>& st)
560 inline Cmpt
invariantII(
const SymmTensor<Cmpt>& st)
564 st.xx()*st.yy() + st.yy()*st.zz() + st.xx()*st.zz()
565 - st.xy()*st.xy() - st.yz()*st.yz() - st.xz()*st.xz()
580 inline SymmTensor<Cmpt>
612 inline SymmTensor<Cmpt>
sqr(
const Vector<Cmpt>& v)
614 return SymmTensor<Cmpt>
616 v.x()*v.x(), v.x()*v.y(), v.x()*v.z(),
617 v.y()*v.y(), v.y()*v.z(),
625 inline SymmTensor<Cmpt>
lerp 632 const scalar onet = (1-t);
636 onet*a.
xx() + t*
b.xx(),
637 onet*a.
xy() + t*
b.xy(),
638 onet*a.
xz() + t*
b.xz(),
639 onet*a.
yy() + t*
b.yy(),
640 onet*a.
yz() + t*
b.yz(),
650 inline SymmTensor<Cmpt>
651 operator+(
const SphericalTensor<Cmpt>& spt1,
const SymmTensor<Cmpt>& st2)
653 return SymmTensor<Cmpt>
655 spt1.ii() + st2.xx(), st2.xy(), st2.xz(),
656 spt1.ii() + st2.yy(), st2.yz(),
664 inline SymmTensor<Cmpt>
665 operator+(
const SymmTensor<Cmpt>& st1,
const SphericalTensor<Cmpt>& spt2)
667 return SymmTensor<Cmpt>
669 st1.xx() + spt2.ii(), st1.xy(), st1.xz(),
670 st1.yy() + spt2.ii(), st1.yz(),
678 inline SymmTensor<Cmpt>
679 operator-(
const SphericalTensor<Cmpt>& spt1,
const SymmTensor<Cmpt>& st2)
681 return SymmTensor<Cmpt>
683 spt1.ii() - st2.xx(), -st2.xy(), -st2.xz(),
684 spt1.ii() - st2.yy(), -st2.yz(),
692 inline SymmTensor<Cmpt>
693 operator-(
const SymmTensor<Cmpt>& st1,
const SphericalTensor<Cmpt>& spt2)
695 return SymmTensor<Cmpt>
697 st1.xx() - spt2.ii(), st1.xy(), st1.xz(),
698 st1.yy() - spt2.ii(), st1.yz(),
706 inline Vector<Cmpt>
operator*(
const SymmTensor<Cmpt>& st)
708 return Vector<Cmpt>(st.yz(), -st.xz(), st.xy());
714 inline SymmTensor<Cmpt>
715 operator/(
const SymmTensor<Cmpt>& st,
const Cmpt
s)
719 st.xx()/
s, st.xy()/
s, st.xz()/
s,
720 st.yy()/
s, st.yz()/
s,
729 operator&(
const SymmTensor<Cmpt>& st1,
const SymmTensor<Cmpt>& st2)
733 st1.xx()*st2.xx() + st1.xy()*st2.xy() + st1.xz()*st2.xz(),
734 st1.xx()*st2.xy() + st1.xy()*st2.yy() + st1.xz()*st2.yz(),
735 st1.xx()*st2.xz() + st1.xy()*st2.yz() + st1.xz()*st2.zz(),
737 st1.xy()*st2.xx() + st1.yy()*st2.xy() + st1.yz()*st2.xz(),
738 st1.xy()*st2.xy() + st1.yy()*st2.yy() + st1.yz()*st2.yz(),
739 st1.xy()*st2.xz() + st1.yy()*st2.yz() + st1.yz()*st2.zz(),
741 st1.xz()*st2.xx() + st1.yz()*st2.xy() + st1.zz()*st2.xz(),
742 st1.xz()*st2.xy() + st1.yz()*st2.yy() + st1.zz()*st2.yz(),
743 st1.xz()*st2.xz() + st1.yz()*st2.yz() + st1.zz()*st2.zz()
755 spt1.
ii()*st2.
xx(), spt1.
ii()*st2.
xy(), spt1.
ii()*st2.
xz(),
756 spt1.
ii()*st2.
yy(), spt1.
ii()*st2.
yz(),
764 inline SymmTensor<Cmpt>
765 operator&(
const SymmTensor<Cmpt>& st1,
const SphericalTensor<Cmpt>& spt2)
767 return SymmTensor<Cmpt>
769 st1.xx()*spt2.ii(), st1.xy()*spt2.ii(), st1.xz()*spt2.ii(),
770 st1.yy()*spt2.ii(), st1.yz()*spt2.ii(),
783 st.
xx()*v.
x() + st.
xy()*v.
y() + st.
xz()*v.
z(),
784 st.
xy()*v.
x() + st.
yy()*v.
y() + st.
yz()*v.
z(),
785 st.
xz()*v.
x() + st.
yz()*v.
y() + st.
zz()*v.
z()
793 operator&(
const Vector<Cmpt>& v,
const SymmTensor<Cmpt>& st)
797 v.x()*st.xx() + v.y()*st.xy() + v.z()*st.xz(),
798 v.x()*st.xy() + v.y()*st.yy() + v.z()*st.yz(),
799 v.x()*st.xz() + v.y()*st.yz() + v.z()*st.zz()
807 operator&&(
const SymmTensor<Cmpt>& st1,
const SymmTensor<Cmpt>& st2)
811 st1.xx()*st2.xx() + 2*st1.xy()*st2.xy() + 2*st1.xz()*st2.xz()
812 + st1.yy()*st2.yy() + 2*st1.yz()*st2.yz()
821 operator&&(
const SphericalTensor<Cmpt>& spt1,
const SymmTensor<Cmpt>& st2)
823 return (spt1.ii()*st2.xx() + spt1.ii()*st2.yy() + spt1.ii()*st2.zz());
832 return (st1.
xx()*spt2.
ii() + st1.
yy()*spt2.
ii() + st1.
zz()*spt2.
ii());
839 class outerProduct<SymmTensor<Cmpt>, Cmpt>
843 typedef SymmTensor<Cmpt>
type;
855 class innerProduct<SymmTensor<Cmpt>, SymmTensor<Cmpt>>
859 typedef Tensor<Cmpt>
type;
871 class innerProduct<Vector<Cmpt>, SymmTensor<Cmpt>>
875 typedef Vector<Cmpt>
type;
888 class typeOfSum<SymmTensor<Cmpt>, SphericalTensor<Cmpt>>
896 class innerProduct<SphericalTensor<Cmpt>, SymmTensor<Cmpt>>
Vector< Cmpt > x() const
Extract vector for row 0.
Cmpt invariantII(const SymmTensor< Cmpt > &st)
Return the 2nd invariant of a SymmTensor.
Cmpt det(const SymmTensor< Cmpt > &st)
Return the determinant of a SymmTensor.
const Cmpt & xz() const noexcept
A templated (3 x 3) symmetric tensor of objects of <T>, effectively containing 6 elements, derived from VectorSpace.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
constexpr char nl
The newline '\n' character (0x0a)
A traits class, which is primarily used for primitives and vector-space.
dimensionSet operator &&(const dimensionSet &ds1, const dimensionSet &ds2)
SymmTensor< Cmpt > devSymm(const SymmTensor< Cmpt > &st)
Return the deviatoric part of the symmetric part of a SymmTensor.
dimensionedScalar operator/(const scalar s1, const dimensionedScalar &ds2)
const Cmpt & y() const noexcept
Access to the vector y component.
typeOfRank< typename pTraits< arg1 >::cmptType, direction(pTraits< arg1 >::rank)+direction(pTraits< arg2 >::rank) >::type type
dimensionedScalar det(const dimensionedSphericalTensor &dt)
SphericalTensor< Cmpt > sph(const DiagTensor< Cmpt > &dt)
Return the spherical part of a DiagTensor as a SphericalTensor.
tmp< faMatrix< Type > > operator+(const faMatrix< Type > &, const faMatrix< Type > &)
SymmTensor< Cmpt > safeInv() const
Return inverse, with (ad hoc) failsafe handling of 2D tensors.
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
Cmpt det() const
The determinate.
const Cmpt & zz() const noexcept
SymmTensor< Cmpt > adjunct() const
Return adjunct matrix (transpose of cofactor matrix)
void subtractDiag(const Vector< Cmpt > &v)
Subtract from the diagonal.
Vector< Cmpt > diag() const
Extract the diagonal as a vector.
A templated (3 x 3) diagonal tensor of objects of <T>, effectively containing 1 element, derived from VectorSpace.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
tmp< faMatrix< Type > > operator*(const areaScalarField::Internal &, const faMatrix< Type > &)
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
SymmTensor< Cmpt > inv2D(const direction excludeCmpt) const
Return inverse of 2D tensor (by excluding given direction)
scalar diagSqr() const
The L2-norm squared of the diagonal.
Cmpt invariantIII(const SymmTensor< Cmpt > &st)
Return the 3rd invariant of a SymmTensor.
Vector< Cmpt > y() const
Extract vector for row 1.
SymmTensor< Cmpt > adjunct2D(const direction excludeCmpt) const
Return 2D adjunct matrix by excluding given direction.
const Cmpt & ii() const noexcept
SymmTensor< Cmpt > cof() const
Return cofactor matrix (transpose of adjunct matrix)
const Cmpt & xx() const noexcept
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
errorManip< error > abort(error &err)
tmp< faMatrix< Type > > operator-(const faMatrix< Type > &)
Unary negation.
SymmTensor & operator=(const SymmTensor &)=default
Copy assignment.
Templated 3D Vector derived from VectorSpace adding construction from 3 components, element access using x(), y() and z() member functions and the inner-product (dot-product) and cross-product operators.
Cmpt det2D(const direction excludeCmpt) const
The 2D determinant by excluding given direction.
dimensionedSymmTensor innerSqr(const dimensionedSymmTensor &dt)
const Cmpt & x() const noexcept
Access to the vector x component.
dimensioned< Type > lerp(const dimensioned< Type > &a, const dimensioned< Type > &b, const scalar t)
const Cmpt & z() const noexcept
Access to the vector z component.
dimensionedSymmTensor dev2(const dimensionedSymmTensor &dt)
void addDiag(const Vector< Cmpt > &v)
Add to the diagonal.
const Cmpt & yy() const noexcept
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
tmp< GeometricField< Type, faPatchField, areaMesh > > operator &(const faMatrix< Type > &, const DimensionedField< Type, areaMesh > &)
Cmpt invariantI(const SymmTensor< Cmpt > &st)
Return the 1st invariant of a SymmTensor.
A class representing the concept of 0 (zero) that can be used to avoid manipulating objects known to ...
SymmTensor< Cmpt > inv() const
Return inverse.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
typeOfRank< typename pTraits< arg1 >::cmptType, direction(pTraits< arg1 >::rank)+direction(pTraits< arg2 >::rank) - 2 >::type type
dimensionedSymmTensor cof(const dimensionedSymmTensor &dt)
Vector< Cmpt > row() const
Extract vector for given row: compile-time check of index.
const Cmpt & xy() const noexcept
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
void rows(const Vector< Cmpt > &x, const Vector< Cmpt > &y, const Vector< Cmpt > &z)
Set row values.
Vector< Cmpt > z() const
Extract vector for row 2.
SymmTensor< Cmpt > devTwoSymm(const SymmTensor< Cmpt > &st)
Return the deviatoric part of twice the symmetric part of a SymmTensor.
static constexpr const zero Zero
Global zero (0)
const Cmpt & yz() const noexcept
SymmTensor()=default
Default construct.