Commit 4a61042f authored by Mark Olesen's avatar Mark Olesen
Browse files

ENH: add column access and other methods for Tensor

- Can now retrieve or set a column/row of a tensor.
  Either compile-time or run-time checks.

  Get
     t.col<1>();   t.col(1);
     t.row<1>();   t.row(1);

  Set
     t.col<1>(vec);   t.col(1,vec);
     t.row<1>(vec);   t.row(1,vec);

  The templated versions are compile-time checked

     t.col<3>();
     t.col<3>(vec);

  The parameter versions are run-time checked

     t.col(3);
     t.col(3,vec);

ENH: provide named access to tensor/tensor inner product as inner()
parent 2d5f77f2
#include "tensor.H"
#include "triad.H"
#include "symmTensor.H"
#include "transform.H"
#include "stringList.H"
......@@ -13,7 +14,40 @@ int main()
tensor t3 = t1 + t2;
Info<< t3 << endl;
Info<< "tensor " << t1 << " + " << t2 << nl
<< " = " << t3 << nl << nl;
{
Info<< "rows:" << nl;
for (direction i=0; i < 3; ++i)
{
Info<< " (" << i << ") = " << t3.row(i) << nl;
}
}
{
Info<< "cols:" << nl;
for (direction i=0; i < 3; ++i)
{
Info<< " (" << i << ") = " << t3.col(i) << nl;
}
Info<< "col<0> = " << t3.col<0>() << nl;
Info<< "col<1> = " << t3.col<1>() << nl;
Info<< "col<2> = " << t3.col<2>() << nl;
// Compilation error: Info << "col<3> = " << t3.col<3>() << nl;
t3.col<1>({0, 2, 4});
Info<< "replaced col<1> = " << t3.col<1>() << nl;
Info<< "tensor " << t3 << nl;
t3.row<2>(Zero);
Info<< "replaced row<1> = " << t3.row<2>() << nl;
Info<< "tensor " << t3 << nl;
triad tr3(t3);
Info<< "triad " << tr3 << " :: T() " << tr3.T() << nl;
}
Info<< nl;
tensor t4(3,-2,1,-2,2,0,1, 0, 4);
......
......@@ -3,7 +3,7 @@
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
\\/ M anipulation | Copyright (C) 2018 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
......@@ -95,22 +95,22 @@ public:
//- Construct given MatrixSpace of the same rank
template<class Cmpt2>
inline Tensor(const MatrixSpace<Tensor<Cmpt2>, Cmpt2, 3, 3>&);
inline Tensor(const MatrixSpace<Tensor<Cmpt2>, Cmpt2, 3, 3>& vs);
//- Construct given VectorSpace of the same rank
template<class Cmpt2>
inline Tensor(const VectorSpace<Tensor<Cmpt2>, Cmpt2, 9>&);
inline Tensor(const VectorSpace<Tensor<Cmpt2>, Cmpt2, 9>& vs);
//- Construct given SphericalTensor
inline Tensor(const SphericalTensor<Cmpt>&);
inline Tensor(const SphericalTensor<Cmpt>& st);
//- Construct given SymmTensor
inline Tensor(const SymmTensor<Cmpt>&);
inline Tensor(const SymmTensor<Cmpt>& st);
//- Construct given triad
inline Tensor(const Vector<Vector<Cmpt>>&);
//- Construct given triad of row vectors
inline Tensor(const Vector<Vector<Cmpt>>& tr);
//- Construct given the three vector components
//- Construct given the three row vectors
inline Tensor
(
const Vector<Cmpt>& x,
......@@ -139,7 +139,7 @@ public:
);
//- Construct from Istream
inline Tensor(Istream&);
inline Tensor(Istream& is);
// Member Functions
......@@ -166,12 +166,72 @@ public:
inline Cmpt& zy();
inline Cmpt& zz();
// Column-vector access.
//- Extract vector for column 0
inline Vector<Cmpt> cx() const;
//- Extract vector for column 1
inline Vector<Cmpt> cy() const;
//- Extract vector for column 2
inline Vector<Cmpt> cz() const;
//- Extract vector for given column.
// Compile-time check of column index.
template<direction Col>
inline Vector<Cmpt> col() const;
//- Extract vector for given column (0,1,2).
// Runtime check of column index.
inline Vector<Cmpt> col(const direction c) const;
//- Set values of given column
// Compile-time check of column index.
template<direction Col>
inline void col(const Vector<Cmpt>& v);
//- Set values of given column (0,1,2)
// Runtime check of column index.
inline void col(const direction c, const Vector<Cmpt>& v);
// Row-vector access.
//- Extract vector for row 0
inline Vector<Cmpt> x() const;
//- Extract vector for row 1
inline Vector<Cmpt> y() const;
//- Extract vector for row 2
inline Vector<Cmpt> z() const;
inline Vector<Cmpt> vectorComponent(const direction) const;
//- Extract vector for given row.
// Compile-time check of row index.
template<direction Row>
inline Vector<Cmpt> row() const;
//- Extract vector for given row (0,1,2)
// Runtime check of row index.
inline Vector<Cmpt> row(const direction r) const;
//- Set values of given row
// Compile-time check of row index.
template<direction Row>
inline void row(const Vector<Cmpt>& v);
//- Set values of given row (0,1,2)
// Runtime check of row index.
inline void row(const direction r, const Vector<Cmpt>& v);
//- Return vector for given row (0,1,2)
// Runtime check of row index.
inline Vector<Cmpt> vectorComponent(const direction cmpt) const;
// Tensor Operations
//- Return transpose
inline Tensor<Cmpt> T() const;
......@@ -179,11 +239,22 @@ public:
//- Return inverse
inline Tensor<Cmpt> inv() const;
//- Inner-product of this with another Tensor.
inline Tensor<Cmpt> inner(const Tensor<Cmpt>& t2) const;
//- Inner-product of this with transpose of another Tensor.
// Primarily useful for coordinate transformations
// (where transpose is the same as the inverse).
inline Tensor<Cmpt> innerT(const Tensor<Cmpt>& t2) const;
//- Schur-product of this with another Tensor.
inline Tensor<Cmpt> schur(const Tensor<Cmpt>& t2) const;
// Member Operators
//- Inner-product with a Tensor
inline void operator&=(const Tensor<Cmpt>&);
//- Assign inner-product of this with another Tensor.
inline void operator&=(const Tensor<Cmpt>& t);
//- Inherit MatrixSpace assignment operators
using Tensor::msType::operator=;
......@@ -198,7 +269,7 @@ public:
//- Assign to a SymmTensor
inline void operator=(const SymmTensor<Cmpt>&);
//- Assign to a triad
//- Assign to a triad of row vectors
inline void operator=(const Vector<Vector<Cmpt>>&);
};
......
......@@ -3,7 +3,7 @@
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation | Copyright (C) 2016 OpenCFD Ltd.
\\/ M anipulation | Copyright (C) 2016-2018 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
......@@ -297,31 +297,178 @@ inline Foam::Vector<Cmpt> Foam::Tensor<Cmpt>::z() const
template<class Cmpt>
inline Foam::Vector<Cmpt> Foam::Tensor<Cmpt>::vectorComponent
(
const direction cmpt
) const
inline Foam::Vector<Cmpt> Foam::Tensor<Cmpt>::cx() const
{
return Vector<Cmpt>(this->v_[XX], this->v_[YX], this->v_[ZX]);
}
template<class Cmpt>
inline Foam::Vector<Cmpt> Foam::Tensor<Cmpt>::cy() const
{
return Vector<Cmpt>(this->v_[XY], this->v_[YY], this->v_[ZY]);
}
template<class Cmpt>
inline Foam::Vector<Cmpt> Foam::Tensor<Cmpt>::cz() const
{
return Vector<Cmpt>(this->v_[XZ], this->v_[YZ], this->v_[ZZ]);
}
template<class Cmpt>
template<Foam::direction Col>
inline Foam::Vector<Cmpt> Foam::Tensor<Cmpt>::col() const
{
switch (cmpt)
if (Col == 0) return cx();
else if (Col == 1) return cy();
else if (Col == 2) return cz();
static_assert(Col < 3, "Invalid column access");
return Zero;
}
template<class Cmpt>
inline Foam::Vector<Cmpt> Foam::Tensor<Cmpt>::col(const direction c) const
{
switch (c)
{
case 0:
return x();
break;
case 1:
return y();
case 0: return cx(); break;
case 1: return cy(); break;
case 2: return cz(); break;
}
FatalErrorInFunction
<< "Invalid column access " << c << abort(FatalError);
return Zero;
}
template<class Cmpt>
template<Foam::direction Row>
inline Foam::Vector<Cmpt> Foam::Tensor<Cmpt>::row() const
{
if (Row == 0) return x();
else if (Row == 1) return y();
else if (Row == 2) return z();
static_assert(Row < 3, "Invalid row access");
return Zero;
}
template<class Cmpt>
inline Foam::Vector<Cmpt> Foam::Tensor<Cmpt>::row(const direction r) const
{
switch (r)
{
case 0: return x(); break;
case 1: return y(); break;
case 2: return z(); break;
}
FatalErrorInFunction
<< "Invalid row access " << r << abort(FatalError);
return Zero;
}
template<class Cmpt>
template<Foam::direction Col>
inline void Foam::Tensor<Cmpt>::col(const Vector<Cmpt>& v)
{
if (Col == 0)
{
this->v_[XX] = v.x();
this->v_[YX] = v.y();
this->v_[ZX] = v.z();
}
else if (Col == 1)
{
this->v_[XY] = v.x();
this->v_[YY] = v.y();
this->v_[ZY] = v.z();
}
else if (Col == 2)
{
this->v_[XZ] = v.x();
this->v_[YZ] = v.y();
this->v_[ZZ] = v.z();
}
static_assert(Col < 3, "Invalid column access");
}
template<class Cmpt>
template<Foam::direction Row>
inline void Foam::Tensor<Cmpt>::row(const Vector<Cmpt>& v)
{
if (Row == 0)
{
this->v_[XX] = v.x(); this->v_[XY] = v.y(); this->v_[XZ] = v.z();
}
else if (Row == 1)
{
this->v_[YX] = v.x(); this->v_[YY] = v.y(); this->v_[YZ] = v.z();
}
else if (Row == 2)
{
this->v_[ZX] = v.x(); this->v_[ZY] = v.y(); this->v_[ZZ] = v.z();
}
static_assert(Row < 3, "Invalid row access");
}
template<class Cmpt>
inline void Foam::Tensor<Cmpt>::col(const direction c, const Vector<Cmpt>& v)
{
switch (c)
{
case 0: col<0>(v); break;
case 1: col<1>(v); break;
case 2: col<2>(v); break;
default:
FatalErrorInFunction
<< "Invalid column access " << c << abort(FatalError);
break;
case 2:
return z();
}
}
template<class Cmpt>
inline void Foam::Tensor<Cmpt>::row(const direction r, const Vector<Cmpt>& v)
{
switch (r)
{
case 0: row<0>(v); break;
case 1: row<1>(v); break;
case 2: row<2>(v); break;
default:
FatalErrorInFunction
<< "Invalid row access " << r << abort(FatalError);
break;
}
}
FatalErrorInFunction
<< "Unhandled component " << cmpt << abort(FatalError);
return x();
template<class Cmpt>
inline Foam::Vector<Cmpt> Foam::Tensor<Cmpt>::vectorComponent
(
const direction cmpt
) const
{
return row(cmpt);
}
// * * * * * * * * * * * * * * * Member Operations * * * * * * * * * * * * * //
template<class Cmpt>
inline Foam::Tensor<Cmpt> Foam::Tensor<Cmpt>::T() const
{
......@@ -334,31 +481,76 @@ inline Foam::Tensor<Cmpt> Foam::Tensor<Cmpt>::T() const
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
template<class Cmpt>
inline Foam::Tensor<Cmpt>
Foam::Tensor<Cmpt>::inner(const Tensor<Cmpt>& t2) const
{
const Tensor<Cmpt>& t1 = *this;
return Tensor<Cmpt>
(
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<class Cmpt>
inline void Foam::Tensor<Cmpt>::operator&=(const Tensor<Cmpt>& t)
inline Foam::Tensor<Cmpt>
Foam::Tensor<Cmpt>::innerT(const Tensor<Cmpt>& t2) const
{
*this =
const Tensor<Cmpt>& t1 = *this;
return Tensor<Cmpt>
(
Tensor<Cmpt>
(
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()
)
t1.xx()*t2.xx() + t1.xy()*t2.xy() + t1.xz()*t2.xz(),
t1.xx()*t2.yx() + t1.xy()*t2.yy() + t1.xz()*t2.yz(),
t1.xx()*t2.zx() + t1.xy()*t2.zy() + t1.xz()*t2.zz(),
t1.yx()*t2.xx() + t1.yy()*t2.xy() + t1.yz()*t2.xz(),
t1.yx()*t2.yx() + t1.yy()*t2.yy() + t1.yz()*t2.yz(),
t1.yx()*t2.zx() + t1.yy()*t2.zy() + t1.yz()*t2.zz(),
t1.zx()*t2.xx() + t1.zy()*t2.xy() + t1.zz()*t2.xz(),
t1.zx()*t2.yx() + t1.zy()*t2.yy() + t1.zz()*t2.yz(),
t1.zx()*t2.zx() + t1.zy()*t2.zy() + t1.zz()*t2.zz()
);
}
template<class Cmpt>
inline Foam::Tensor<Cmpt>
Foam::Tensor<Cmpt>::schur(const Tensor<Cmpt>& t2) const
{
const Tensor<Cmpt>& t1 = *this;
return Tensor<Cmpt>
(
t1.xx()*t2.xx(), t1.xy()*t2.xy(), t1.xz()*t2.xz(),
t1.yx()*t2.yx(), t1.yy()*t2.yy(), t1.yz()*t2.yz(),
t1.zx()*t2.zx(), t1.zy()*t2.zy(), t1.zz()*t2.zz()
);
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
template<class Cmpt>
inline void Foam::Tensor<Cmpt>::operator&=(const Tensor<Cmpt>& t)
{
*this = this->inner(t);
}
template<class Cmpt>
template<class Cmpt2>
inline void Foam::Tensor<Cmpt>::operator=
......@@ -435,20 +627,7 @@ template<class Cmpt>
inline typename innerProduct<Tensor<Cmpt>, Tensor<Cmpt>>::type
operator&(const Tensor<Cmpt>& t1, const Tensor<Cmpt>& t2)
{
return Tensor<Cmpt>
(
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()
);
return t1.inner(t2);
}
......
......@@ -3,7 +3,7 @@
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
\\/ M anipulation | Copyright (C) 2018 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
......@@ -97,7 +97,7 @@ public:
//- Construct given SphericalTensor2D
inline Tensor2D(const SphericalTensor2D<Cmpt>&);
//- Construct given the two vectors
//- Construct given the two row vectors
inline Tensor2D
(
const Vector2D<Cmpt>& x,
......@@ -117,7 +117,7 @@ public:
// Member Functions
// Access
// Component access
inline const Cmpt& xx() const;
inline const Cmpt& xy() const;
......@@ -129,14 +129,41 @@ public:
inline Cmpt& yx();
inline Cmpt& yy();
// Access vector components.
// Column-vector access.
//- Extract vector for column 0
inline Vector2D<Cmpt> cx() const;
//- Extract vector for column 1
inline Vector2D<Cmpt> cy() const;
// Row-vector access.
//- Extract vector for row 0
inline Vector2D<Cmpt> x() const;
//- Extract vector for row 1
inline Vector2D<Cmpt> y() const;
// Tensor Operations
//- Transpose
inline Tensor2D<Cmpt> T() const;
//- Inner-product of this with another Tensor2D.
inline Tensor2D<Cmpt> inner(const Tensor2D<Cmpt>& t2) const;
//- Inner-product of this with transpose of another Tensor2D.
// Primarily useful for coordinate transformations
// (where transpose is the same as the inverse).
inline Tensor2D<Cmpt> innerT(const Tensor2D<Cmpt>& t2) const;
//- Schur-product of this with another Tensor2D.
inline Tensor2D<Cmpt> schur(const Tensor2D<Cmpt>& t2) const;
// Member Operators
......
......@@ -3,7 +3,7 @@
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
\\/ M anipulation | Copyright (C) 2018 OpenCFD Ltd.