diff --git a/applications/test/tensor/Test-tensor.C b/applications/test/tensor/Test-tensor.C
index d597985da59972c97da37a5bbf8953894dc405f2..cd72200975132e1e41f89fd780dd9c9c136409e0 100644
--- a/applications/test/tensor/Test-tensor.C
+++ b/applications/test/tensor/Test-tensor.C
@@ -1,4 +1,5 @@
 #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);
 
diff --git a/src/OpenFOAM/primitives/Tensor/Tensor.H b/src/OpenFOAM/primitives/Tensor/Tensor.H
index ee670f07b9e06748f51a358987ef1c5b215f9d8f..d39f8b10e136d607aa71230d1cde320b7124b873 100644
--- a/src/OpenFOAM/primitives/Tensor/Tensor.H
+++ b/src/OpenFOAM/primitives/Tensor/Tensor.H
@@ -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>>&);
 };
 
diff --git a/src/OpenFOAM/primitives/Tensor/TensorI.H b/src/OpenFOAM/primitives/Tensor/TensorI.H
index c8ecbc7ce664cc365059239e175bbd88d6300b82..933b39a5c2e9c6efbcb95ae6e4c671beec6ee50f 100644
--- a/src/OpenFOAM/primitives/Tensor/TensorI.H
+++ b/src/OpenFOAM/primitives/Tensor/TensorI.H
@@ -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);
 }
 
 
diff --git a/src/OpenFOAM/primitives/Tensor2D/Tensor2D.H b/src/OpenFOAM/primitives/Tensor2D/Tensor2D.H
index 5ee3db5827c5a4323d40b8f2dc49c3ebfafd2aaf..38de0da43392abee0ae5cd9762813ab6cb5ff1b7 100644
--- a/src/OpenFOAM/primitives/Tensor2D/Tensor2D.H
+++ b/src/OpenFOAM/primitives/Tensor2D/Tensor2D.H
@@ -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
 
diff --git a/src/OpenFOAM/primitives/Tensor2D/Tensor2DI.H b/src/OpenFOAM/primitives/Tensor2D/Tensor2DI.H
index 7f6be0d324de5f1a6302454e4f70ad6bac620ecb..147155035f122d449f93ca2fac3883ae5669e02c 100644
--- a/src/OpenFOAM/primitives/Tensor2D/Tensor2DI.H
+++ b/src/OpenFOAM/primitives/Tensor2D/Tensor2DI.H
@@ -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.
@@ -96,19 +96,6 @@ inline Foam::Tensor2D<Cmpt>::Tensor2D(Istream& is)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-template<class Cmpt>
-inline Foam::Vector2D<Cmpt> Foam::Tensor2D<Cmpt>::x() const
-{
-    return Vector2D<Cmpt>(this->v_[XX], this->v_[XY]);
-}
-
-template<class Cmpt>
-inline Foam::Vector2D<Cmpt> Foam::Tensor2D<Cmpt>::y() const
-{
-    return Vector2D<Cmpt>(this->v_[YX], this->v_[YY]);
-}
-
-
 template<class Cmpt>
 inline const Cmpt& Foam::Tensor2D<Cmpt>::xx() const
 {
@@ -159,6 +146,34 @@ inline Cmpt& Foam::Tensor2D<Cmpt>::yy()
 }
 
 
+template<class Cmpt>
+inline Foam::Vector2D<Cmpt> Foam::Tensor2D<Cmpt>::x() const
+{
+    return Vector2D<Cmpt>(this->v_[XX], this->v_[XY]);
+}
+
+template<class Cmpt>
+inline Foam::Vector2D<Cmpt> Foam::Tensor2D<Cmpt>::y() const
+{
+    return Vector2D<Cmpt>(this->v_[YX], this->v_[YY]);
+}
+
+
+template<class Cmpt>
+inline Foam::Vector2D<Cmpt> Foam::Tensor2D<Cmpt>::cx() const
+{
+    return Vector2D<Cmpt>(this->v_[XX], this->v_[YX]);
+}
+
+template<class Cmpt>
+inline Foam::Vector2D<Cmpt> Foam::Tensor2D<Cmpt>::cy() const
+{
+    return Vector2D<Cmpt>(this->v_[XY], this->v_[YY]);
+}
+
+
+// * * * * * * * * * * * * * * * Member Operations * * * * * * * * * * * * * //
+
 template<class Cmpt>
 inline Foam::Tensor2D<Cmpt> Foam::Tensor2D<Cmpt>::T() const
 {
@@ -170,6 +185,54 @@ inline Foam::Tensor2D<Cmpt> Foam::Tensor2D<Cmpt>::T() const
 }
 
 
+template<class Cmpt>
+inline Foam::Tensor2D<Cmpt>
+Foam::Tensor2D<Cmpt>::inner(const Tensor2D<Cmpt>& t2) const
+{
+    const Tensor2D<Cmpt>& t1 = *this;
+
+    return Tensor2D<Cmpt>
+    (
+        t1.xx()*t2.xx() + t1.xy()*t2.yx(),
+        t1.xx()*t2.xy() + t1.xy()*t2.yy(),
+
+        t1.yx()*t2.xx() + t1.yy()*t2.yx(),
+        t1.yx()*t2.xy() + t1.yy()*t2.yy()
+    );
+}
+
+
+template<class Cmpt>
+inline Foam::Tensor2D<Cmpt>
+Foam::Tensor2D<Cmpt>::innerT(const Tensor2D<Cmpt>& t2) const
+{
+    const Tensor2D<Cmpt>& t1 = *this;
+
+    return Tensor2D<Cmpt>
+    (
+        t1.xx()*t2.xx() + t1.xy()*t2.xy(),
+        t1.xx()*t2.yx() + t1.xy()*t2.yy(),
+
+        t1.yx()*t2.xx() + t1.yy()*t2.xy(),
+        t1.yx()*t2.yx() + t1.yy()*t2.yy()
+    );
+}
+
+
+template<class Cmpt>
+inline Foam::Tensor2D<Cmpt>
+Foam::Tensor2D<Cmpt>::schur(const Tensor2D<Cmpt>& t2) const
+{
+    const Tensor2D<Cmpt>& t1 = *this;
+
+    return Tensor2D<Cmpt>
+    (
+        t1.xx()*t2.xx(), t1.xy()*t2.xy(),
+        t1.yx()*t2.yx(), t1.yy()*t2.yy()
+    );
+}
+
+
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
 template<class Cmpt>
@@ -200,16 +263,10 @@ template<class Cmpt>
 inline typename innerProduct<Tensor2D<Cmpt>, Tensor2D<Cmpt>>::type
 operator&(const Tensor2D<Cmpt>& t1, const Tensor2D<Cmpt>& t2)
 {
-    return Tensor2D<Cmpt>
-    (
-        t1.xx()*t2.xx() + t1.xy()*t2.yx(),
-        t1.xx()*t2.xy() + t1.xy()*t2.yy(),
-
-        t1.yx()*t2.xx() + t1.yy()*t2.yx(),
-        t1.yx()*t2.xy() + t1.yy()*t2.yy()
-    );
+    return t1.inner(t2);
 }
 
+
 //- Inner-product between a tensor and a vector
 template<class Cmpt>
 inline typename innerProduct<Tensor2D<Cmpt>, Vector2D<Cmpt>>::type
diff --git a/src/OpenFOAM/primitives/triad/triad.C b/src/OpenFOAM/primitives/triad/triad.C
index 0a54dc4aca26d49db99ddda0682398695a84a191..ec800e9f3dce43c6105d1bdc8db6aff0261ab0db 100644
--- a/src/OpenFOAM/primitives/triad/triad.C
+++ b/src/OpenFOAM/primitives/triad/triad.C
@@ -95,14 +95,6 @@ Foam::triad::triad(const quaternion& q)
 }
 
 
-Foam::triad::triad(const tensor& t)
-{
-    x() = t.x();
-    y() = t.y();
-    z() = t.z();
-}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 void Foam::triad::orthogonalize()
@@ -379,16 +371,6 @@ Foam::triad::operator Foam::quaternion() const
 }
 
 
-// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
-
-void Foam::triad::operator=(const tensor& t)
-{
-    x() = t.x();
-    y() = t.y();
-    z() = t.z();
-}
-
-
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 Foam::scalar Foam::diff(const triad& A, const triad& B)
diff --git a/src/OpenFOAM/primitives/triad/triad.H b/src/OpenFOAM/primitives/triad/triad.H
index 30eaba81c12d5b28ee30d96e6edcc4266849ee64..9a1acd1ede6b93ff5e96b0e3272fbf42310e58d7 100644
--- a/src/OpenFOAM/primitives/triad/triad.H
+++ b/src/OpenFOAM/primitives/triad/triad.H
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2012-2016 OpenFOAM Foundation
-     \\/     M anipulation  |
+     \\/     M anipulation  | Copyright (C) 2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -25,7 +25,8 @@ Class
     Foam::triad
 
 Description
-    Representation of a 3D Cartesian coordinate system as a Vector of vectors.
+    Representation of a 3D Cartesian coordinate system as a Vector of
+    row vectors.
 
 See also
     Foam::quaternion
@@ -48,12 +49,10 @@ SourceFiles
 namespace Foam
 {
 
-// Forward declaration of classes
+// Forward declarations
+class triad;
 class Istream;
 class Ostream;
-
-// Forward declaration of friend functions and operators
-class triad;
 Istream& operator>>(Istream&, triad&);
 Ostream& operator<<(Ostream&, const triad&);
 
@@ -78,18 +77,18 @@ public:
         //- Construct from components
         inline triad(const Vector<vector>& vv);
 
-        //- Construct from coordinate axes
+        //- Construct from coordinate axes (row vectors)
         inline triad(const vector& x, const vector& y, const vector& z);
 
+        //- Construct from a tensor
+        inline triad(const tensor& t);
+
         //- Construct from a primary axis with the other two unset
         inline triad(const vector& pa);
 
         //- Construct from a quaternion
         triad(const quaternion& q);
 
-        //- Construct from a tensor
-        triad(const tensor& t);
-
         //- Construct from Istream
         inline triad(Istream&);
 
@@ -130,12 +129,27 @@ public:
         //- Convert to a quaternion
         operator quaternion() const;
 
+        //- Return transpose
+        inline triad T() const;
+
+
+    // Column-vector access.
+
+        //- Extract vector for column 0
+        inline vector cx() const;
+
+        //- Extract vector for column 1
+        inline vector cy() const;
+
+        //- Extract vector for column 2
+        inline vector cz() const;
+
 
     // Member Operators
 
-        inline void operator=(const Vector<vector>&);
+        inline void operator=(const Vector<vector>& vv);
 
-        void operator=(const tensor& t);
+        inline void operator=(const tensor& t);
 
         //- Add the triad t2 to this triad
         //  without normalizing or orthogonalizing
diff --git a/src/OpenFOAM/primitives/triad/triadI.H b/src/OpenFOAM/primitives/triad/triadI.H
index d2592581780a6e2608f9e34c33bed270db65e281..249f5f9157e2dfce3d89f3da41bdfb99ab58411b 100644
--- a/src/OpenFOAM/primitives/triad/triadI.H
+++ b/src/OpenFOAM/primitives/triad/triadI.H
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2012-2016 OpenFOAM Foundation
-     \\/     M anipulation  |
+     \\/     M anipulation  | Copyright (C) 2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -43,6 +43,12 @@ inline Foam::triad::triad(const vector& x, const vector& y, const vector& z)
 {}
 
 
+inline Foam::triad::triad(const tensor& t)
+:
+    Vector<vector>(t.x(), t.y(), t.z())
+{}
+
+
 inline Foam::triad::triad(const vector& pa)
 {
     operator=(triad::unset);
@@ -116,6 +122,30 @@ inline void Foam::triad::normalize()
 }
 
 
+inline Foam::vector Foam::triad::cx() const
+{
+    return vector(x().x(), y().x(), z().x());
+}
+
+
+inline Foam::vector Foam::triad::cy() const
+{
+    return vector(x().y(), y().y(), z().y());
+}
+
+
+inline Foam::vector Foam::triad::cz() const
+{
+    return vector(x().z(), y().z(), z().z());
+}
+
+
+inline Foam::triad Foam::triad::T() const
+{
+    return triad(cx(), cy(), cz());
+}
+
+
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
 inline void Foam::triad::operator=(const Vector<vector>& vv)
@@ -124,6 +154,14 @@ inline void Foam::triad::operator=(const Vector<vector>& vv)
 }
 
 
+inline void Foam::triad::operator=(const tensor& t)
+{
+    x() = t.x();
+    y() = t.y();
+    z() = t.z();
+}
+
+
 // * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * * //
 
 inline Foam::Istream& Foam::operator>>(Istream& is, triad& t)