ENH: expose calculation of coordinate rotations as static methods

- improve handling of degenerate cases for the two-axes specification.

src/meshTools/coordinate/rotation/EulerCoordinateRotation.C

@@ -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) 2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2017-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -24,8 +24,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "EulerCoordinateRotation.H"
-
+#include "unitConversion.H"
-#include "mathematicalConstants.H"
 #include "addToRunTimeSelectionTable.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -47,6 +46,43 @@ namespace Foam
     );
 }
 
+// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
+
+Foam::tensor Foam::EulerCoordinateRotation::rotation
+(
+    const vector& angles,
+    bool degrees
+)
+{
+    scalar phi   = angles.component(vector::X); // 1. Rotate about Z
+    scalar theta = angles.component(vector::Y); // 2. Rotate about X
+    scalar psi   = angles.component(vector::Z); // 3. Rotate about Z
+
+    if (degrees)
+    {
+        phi   *= degToRad();
+        theta *= degToRad();
+        psi   *= degToRad();
+    }
+
+    return
+        tensor
+        (
+            cos(phi)*cos(psi) - sin(phi)*sin(psi)*cos(theta),
+            -sin(phi)*cos(psi)*cos(theta) - cos(phi)*sin(psi),
+            sin(phi)*sin(theta),
+
+            cos(phi)*sin(psi)*cos(theta) + sin(phi)*cos(psi),
+            cos(phi)*cos(psi)*cos(theta) - sin(phi)*sin(psi),
+            -cos(phi)*sin(theta),
+
+            sin(psi)*sin(theta),
+            cos(psi)*sin(theta),
+            cos(theta)
+        );
+}
+
+
 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
 Foam::vector Foam::EulerCoordinateRotation::transform(const vector& st) const
@@ -147,49 +183,6 @@ Foam::symmTensor Foam::EulerCoordinateRotation::transformVector
 }
 
 
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-void Foam::EulerCoordinateRotation::calcTransform
-(
-    const scalar phiAngle,
-    const scalar thetaAngle,
-    const scalar psiAngle,
-    const bool inDegrees
-)
-{
-    scalar phi   = phiAngle;
-    scalar theta = thetaAngle;
-    scalar psi   = psiAngle;
-
-    if (inDegrees)
-    {
-        phi   *= constant::mathematical::pi/180.0;
-        theta *= constant::mathematical::pi/180.0;
-        psi   *= constant::mathematical::pi/180.0;
-    }
-
-    R_ =
-    (
-        tensor
-        (
-            cos(phi)*cos(psi) - sin(phi)*sin(psi)*cos(theta),
-            -sin(phi)*cos(psi)*cos(theta) - cos(phi)*sin(psi),
-            sin(phi)*sin(theta),
-
-            cos(phi)*sin(psi)*cos(theta) + sin(phi)*cos(psi),
-            cos(phi)*cos(psi)*cos(theta) - sin(phi)*sin(psi),
-            -cos(phi)*sin(theta),
-
-            sin(psi)*sin(theta),
-            cos(psi)*sin(theta),
-            cos(theta)
-        )
-    );
-
-    Rtr_ = R_.T();
-}
-
-
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::EulerCoordinateRotation::EulerCoordinateRotation()
@@ -212,35 +205,25 @@ Foam::EulerCoordinateRotation::EulerCoordinateRotation
 Foam::EulerCoordinateRotation::EulerCoordinateRotation
 (
     const vector& phiThetaPsi,
-    const bool inDegrees
+    const bool degrees
 )
 :
-    R_(sphericalTensor::I),
+    R_(rotation(phiThetaPsi, degrees)),
-    Rtr_(sphericalTensor::I)
+    Rtr_(R_.T())
-{
+{}
-    calcTransform
-    (
-        phiThetaPsi.component(vector::X),
-        phiThetaPsi.component(vector::Y),
-        phiThetaPsi.component(vector::Z),
-        inDegrees
-    );
-}
 
 
 Foam::EulerCoordinateRotation::EulerCoordinateRotation
 (
-    const scalar phiAngle,
+    const scalar phi,
-    const scalar thetaAngle,
+    const scalar theta,
-    const scalar psiAngle,
+    const scalar psi,
-    const bool inDegrees
+    const bool degrees
 )
 :
-    R_(sphericalTensor::I),
+    R_(rotation(vector(phi, theta, psi), degrees)),
-    Rtr_(sphericalTensor::I)
+    Rtr_(R_.T())
-{
+{}
-    calcTransform(phiAngle, thetaAngle, psiAngle, inDegrees);
-}
 
 
 Foam::EulerCoordinateRotation::EulerCoordinateRotation
@@ -248,19 +231,16 @@ Foam::EulerCoordinateRotation::EulerCoordinateRotation
     const dictionary& dict
 )
 :
-    R_(sphericalTensor::I),
+    R_
-    Rtr_(sphericalTensor::I)
-{
-    const vector rotation(dict.lookup("rotation"));
-
-    calcTransform
     (
-        rotation.component(vector::X),
+        rotation
-        rotation.component(vector::Y),
+        (
-        rotation.component(vector::Z),
+            dict.get<vector>("rotation"),
             dict.lookupOrDefault("degrees", true)
-    );
+        )
-}
+    ),
+    Rtr_(R_.T())
+{}
 
 
 Foam::EulerCoordinateRotation::EulerCoordinateRotation

src/meshTools/coordinate/rotation/EulerCoordinateRotation.H

@@ -73,7 +73,6 @@ class EulerCoordinateRotation
 :
     public coordinateRotation
 {
-
     // Private Member Data
 
         //- Local-to-global transformation tensor
@@ -83,18 +82,6 @@ class EulerCoordinateRotation
         tensor Rtr_;
 
 
-    // Private Member Functions
-
-        //- Calculate transformation tensor
-        void calcTransform
-        (
-            const scalar phiAngle,
-            const scalar thetaAngle,
-            const scalar psiAngle,
-            const bool inDegrees
-        );
-
-
 public:
 
     //- Runtime type information
@@ -113,16 +100,16 @@ public:
         EulerCoordinateRotation
         (
             const vector& phiThetaPsi,
-            const bool inDegrees
+            const bool degrees
         );
 
         //- Construct from components of rotation vector
         EulerCoordinateRotation
         (
-            const scalar phiAngle,
+            const scalar phi,
-            const scalar thetaAngle,
+            const scalar theta,
-            const scalar psiAngle,
+            const scalar psi,
-            const bool inDegrees
+            const bool degrees
         );
 
         //- Construct from dictionary
@@ -144,6 +131,13 @@ public:
         }
 
 
+    // Static Member Functions
+
+        //- The rotation tensor calculated for the specified Euler angles
+        //- interpreted as phi/theta/psi (z-x-z order)
+        static tensor rotation(const vector& angles, bool degrees);
+
+
     // Member Functions
 
         //- Reset rotation to an identity rotation

src/meshTools/coordinate/rotation/STARCDCoordinateRotation.C

@@ -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) 2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2017-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -24,8 +24,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "STARCDCoordinateRotation.H"
-
+#include "unitConversion.H"
-#include "mathematicalConstants.H"
 #include "addToRunTimeSelectionTable.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -148,29 +147,26 @@ Foam::symmTensor Foam::STARCDCoordinateRotation::transformVector
 }
 
 
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
 
-void Foam::STARCDCoordinateRotation::calcTransform
+Foam::tensor Foam::STARCDCoordinateRotation::rotation
 (
-    const scalar rotZ,
+    const vector& angles,
-    const scalar rotX,
+    const bool degrees
-    const scalar rotY,
-    const bool inDegrees
 )
 {
-    scalar x = rotX;
+    scalar z = angles.component(vector::X);    // 1. Rotate about Z
-    scalar y = rotY;
+    scalar x = angles.component(vector::Y);    // 2. Rotate about X
-    scalar z = rotZ;
+    scalar y = angles.component(vector::Z);    // 3. Rotate about Y
 
-    if (inDegrees)
+    if (degrees)
     {
-        x *= constant::mathematical::pi/180.0;
+        x *= degToRad();
-        y *= constant::mathematical::pi/180.0;
+        y *= degToRad();
-        z *= constant::mathematical::pi/180.0;
+        z *= degToRad();
     }
 
-    R_ =
+    return
-    (
         tensor
         (
             cos(y)*cos(z) - sin(x)*sin(y)*sin(z),
@@ -184,10 +180,7 @@ void Foam::STARCDCoordinateRotation::calcTransform
             -cos(x)*sin(y),
             sin(x),
             cos(x)*cos(y)
-        )
         );
-
-    Rtr_ = R_.T();
 }
 
 
@@ -213,20 +206,12 @@ Foam::STARCDCoordinateRotation::STARCDCoordinateRotation
 Foam::STARCDCoordinateRotation::STARCDCoordinateRotation
 (
     const vector& rotZrotXrotY,
-    const bool inDegrees
+    const bool degrees
 )
 :
-    R_(sphericalTensor::I),
+    R_(rotation(rotZrotXrotY, degrees)),
-    Rtr_(sphericalTensor::I)
+    Rtr_(R_.T())
-{
+{}
-    calcTransform
-    (
-        rotZrotXrotY.component(vector::X),
-        rotZrotXrotY.component(vector::Y),
-        rotZrotXrotY.component(vector::Z),
-        inDegrees
-    );
-}
 
 
 Foam::STARCDCoordinateRotation::STARCDCoordinateRotation
@@ -234,14 +219,12 @@ Foam::STARCDCoordinateRotation::STARCDCoordinateRotation
     const scalar rotZ,
     const scalar rotX,
     const scalar rotY,
-    const bool inDegrees
+    const bool degrees
 )
 :
-    R_(sphericalTensor::I),
+    R_(rotation(vector(rotZ, rotX, rotY), degrees)),
-    Rtr_(sphericalTensor::I)
+    Rtr_(R_.T())
-{
+{}
-    calcTransform(rotZ, rotX, rotY, inDegrees);
-}
 
 
 Foam::STARCDCoordinateRotation::STARCDCoordinateRotation
@@ -249,19 +232,16 @@ Foam::STARCDCoordinateRotation::STARCDCoordinateRotation
     const dictionary& dict
 )
 :
-    R_(sphericalTensor::I),
+    R_
-    Rtr_(sphericalTensor::I)
-{
-    const vector rotation(dict.lookup("rotation"));
-
-    calcTransform
     (
-        rotation.component(vector::X),
+        rotation
-        rotation.component(vector::Y),
+        (
-        rotation.component(vector::Z),
+            dict.get<vector>("rotation"),
             dict.lookupOrDefault("degrees", true)
-    );
+        )
-}
+    ),
+    Rtr_(R_.T())
+{}
 
 
 Foam::STARCDCoordinateRotation::STARCDCoordinateRotation

src/meshTools/coordinate/rotation/STARCDCoordinateRotation.H

@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2017-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -65,7 +65,6 @@ class STARCDCoordinateRotation
 :
     public coordinateRotation
 {
-
     // Private Member Data
 
         //- Local-to-Global transformation tensor
@@ -75,18 +74,6 @@ class STARCDCoordinateRotation
         tensor Rtr_;
 
 
-    // Private Member Functions
-
-        //- Calculate transformation tensor
-        void calcTransform
-        (
-            const scalar rotZ,
-            const scalar rotX,
-            const scalar rotY,
-            const bool inDegrees
-        );
-
-
 public:
 
     //- Runtime type information
@@ -105,7 +92,7 @@ public:
         STARCDCoordinateRotation
         (
             const vector& rotZrotXrotY,
-            const bool inDegrees
+            const bool degrees
         );
 
         //- Construct from components of rotation vector
@@ -114,7 +101,7 @@ public:
             const scalar rotZ,
             const scalar rotX,
             const scalar rotY,
-            const bool inDegrees
+            const bool degrees
         );
 
         //- Construct from dictionary
@@ -137,6 +124,13 @@ public:
         }
 
 
+    // Static Member Functions
+
+        //- The rotation tensor calculated for the specified STARCD angles
+        //- interpreted as rotate-Z, rotate-X, rotate-Y
+        static tensor rotation(const vector& angles, bool degrees);
+
+
     // Member Functions
 
         //- Reset rotation to an identity rotation

src/meshTools/coordinate/rotation/axesRotation.C

@@ -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) 2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2017-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -47,6 +47,151 @@ namespace Foam
 }
 
 
+// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
+
+Foam::tensor Foam::axesRotation::rotation
+(
+    const vector& axis1,
+    const vector& axis2,
+    axisOrder order
+)
+{
+    const scalar magAxis1(mag(axis1));
+    scalar magAxis2(mag(axis2));
+
+    if (magAxis1 < ROOTVSMALL)
+    {
+        FatalErrorInFunction
+            << "Dominant coordinate axis cannot have zero length"
+            << nl << endl
+            << abort(FatalError);
+    }
+
+    const vector ax1(axis1 / magAxis1);  // normalise
+    vector ax2(axis2);
+
+    if (magAxis2 < ROOTVSMALL)
+    {
+        // axis2 == Zero : Use best-guess for the second axis.
+        ax2 = findOrthogonal(axis1);
+    }
+
+    // Remove colinear component
+    ax2 -= ((ax1 & ax2) * ax1);
+
+    magAxis2 = mag(ax2);
+
+    if (magAxis2 < SMALL)
+    {
+        WarningInFunction
+            << "axis1, axis2 appear to be co-linear: "
+            << axis1 << ", " << axis2 << "  Revert to guessing axis2"
+            << nl << endl;
+
+        ax2 = findOrthogonal(axis1);
+
+        // Remove colinear component
+        ax2 -= ((ax1 & ax2) * ax1);
+
+        magAxis2 = mag(ax2);
+
+        if (magAxis2 < SMALL)
+        {
+            FatalErrorInFunction
+                << "Could not find an appropriate second axis"
+                << nl << endl
+                << abort(FatalError);
+        }
+    }
+
+    ax2 /= magAxis2;  // normalise
+
+
+    // The local axes are columns of the rotation matrix
+
+    tensor rotTensor;
+
+    switch (order)
+    {
+        case E1_E2:
+        {
+            rotTensor.col<0>(ax1);
+            rotTensor.col<1>(ax2);
+            rotTensor.col<2>(ax1^ax2);
+            break;
+        }
+        case E2_E3:
+        {
+            rotTensor.col<0>(ax1^ax2);
+            rotTensor.col<1>(ax1);
+            rotTensor.col<2>(ax2);
+            break;
+        }
+        case E3_E1:
+        case E3_E1_COMPAT:
+        {
+            rotTensor.col<0>(ax2);
+            rotTensor.col<1>(ax1^ax2);
+            rotTensor.col<2>(ax1);
+            break;
+        }
+    }
+
+    return rotTensor;
+}
+
+
+// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
+
+void Foam::axesRotation::read(const dictionary& dict)
+{
+    vector axis1, axis2;
+    axisOrder order = E3_E1;
+
+    if
+    (
+        dict.readIfPresent("e1", axis1)
+     && dict.readIfPresent("e2", axis2)
+    )
+    {
+        order = E1_E2;
+    }
+    else if
+    (
+        dict.readIfPresent("e2", axis1)
+     && dict.readIfPresent("e3", axis2)
+    )
+    {
+        order = E2_E3;
+    }
+    else if
+    (
+        dict.readIfPresent("e3", axis1)
+     && dict.readIfPresent("e1", axis2)
+    )
+    {
+        order = E3_E1;
+    }
+    else if
+    (
+        dict.readIfPresent("axis", axis1)
+     && dict.readIfPresent("direction", axis2)
+    )
+    {
+        order = E3_E1_COMPAT;
+    }
+    else
+    {
+        FatalErrorInFunction
+            << "No entries of the type (e1, e2) or (e2, e3) or (e3, e1) found"
+            << exit(FatalError);
+    }
+
+    R_ = rotation(axis1, axis2, order);
+    Rtr_ = R_.T();
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::axesRotation::axesRotation()
@@ -77,11 +222,9 @@ Foam::axesRotation::axesRotation
     const axisOrder& order
 )
 :
-    R_(sphericalTensor::I),
+    R_(rotation(axis, dir, order)),
-    Rtr_(sphericalTensor::I)
+    Rtr_(R_.T())
-{
+{}
-    setTransform(axis, dir, order);
-}
 
 
 Foam::axesRotation::axesRotation
@@ -89,37 +232,9 @@ Foam::axesRotation::axesRotation
     const vector& axis
 )
 :
-    R_(sphericalTensor::I),
+    R_(rotation(axis, findOrthogonal(axis), E3_E1)),
-    Rtr_(sphericalTensor::I)
+    Rtr_(R_.T())
-{
+{}
-    direction maxCmpt = 0, dirCmpt = 1;
-
-    scalar maxVal = mag(axis[maxCmpt]);
-    bool negative = (axis[maxCmpt] < 0);
-
-    for (direction cmpt = 1; cmpt < vector::nComponents; ++cmpt)
-    {
-        const scalar val = mag(axis[cmpt]);
-
-        if (maxVal < val)
-        {
-            maxVal  = val;
-            maxCmpt = cmpt;
-            dirCmpt = maxCmpt+1;
-            negative = (axis[cmpt] < 0);
-
-            if (dirCmpt >= vector::nComponents)
-            {
-                dirCmpt = 0;
-            }
-        }
-    }
-
-    vector dir = Zero;
-    dir.component(dirCmpt) = (negative ? -1 : 1);
-
-    setTransform(axis, dir, E3_E1);
-}
 
 
 Foam::axesRotation::axesRotation
@@ -130,7 +245,7 @@ Foam::axesRotation::axesRotation
     R_(sphericalTensor::I),
     Rtr_(sphericalTensor::I)
 {
-    operator=(dict);
+    read(dict);
 }
 
 
@@ -146,62 +261,6 @@ Foam::axesRotation::axesRotation
 
 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
-void Foam::axesRotation::setTransform
-(
-    const vector& axis1,
-    const vector& axis2,
-    const axisOrder& order
-)
-{
-    const vector a = axis1/mag(axis1);
-    vector b = axis2;
-
-    b = b - (b & a)*a;
-
-    if (mag(b) < SMALL)
-    {
-        FatalErrorInFunction
-            << "axis1, axis2 appear to be co-linear: "
-            << axis1 << ", " << axis2 << endl
-            << abort(FatalError);
-    }
-
-    b = b/mag(b);
-    const vector c = a^b;
-
-    // Global->local transformation
-    switch (order)
-    {
-        case E1_E2:
-        {
-            Rtr_ = tensor(a, b, c);
-            break;
-        }
-        case E2_E3:
-        {
-            Rtr_ = tensor(c, a, b);
-            break;
-        }
-        case E3_E1:
-        {
-            Rtr_ = tensor(b, c, a);
-            break;
-        }
-        default:
-        {
-            FatalErrorInFunction
-                << "Unhandled axes specification" << endl
-                << abort(FatalError);
-
-            break;
-        }
-    }
-
-    // Local->global transformation
-    R_ = Rtr_.T();
-}
-
-
 const Foam::tensorField& Foam::axesRotation::Tr() const
 {
     NotImplemented;
@@ -298,36 +357,7 @@ Foam::symmTensor Foam::axesRotation::transformVector
 
 void Foam::axesRotation::operator=(const dictionary& dict)
 {
-    vector axis1, axis2;
+    read(dict);
-
-    if (dict.readIfPresent("e1", axis1) && dict.readIfPresent("e2", axis2))
-    {
-        setTransform(axis1, axis2, E1_E2);
-    }
-    else if (dict.readIfPresent("e2", axis1) && dict.readIfPresent("e3", axis2))
-    {
-        setTransform(axis1, axis2, E2_E3);
-    }
-    else if (dict.readIfPresent("e3", axis1) && dict.readIfPresent("e1", axis2))
-    {
-        setTransform(axis1, axis2, E3_E1);
-    }
-    else if (dict.found("axis") || dict.found("direction"))
-    {
-        // Both "axis" and "direction" are required
-        // If one is missing the appropriate error message will be generated
-        dict.readEntry("axis", axis1);
-        dict.readEntry("direction", axis2);
-
-        setTransform(axis1, axis2, E3_E1);
-    }
-    else
-    {
-        FatalErrorInFunction
-            << "not entry of the type (e1, e2) or (e2, e3) or (e3, e1) "
-            << "found "
-            << exit(FatalError);
-    }
 }
 
 

src/meshTools/coordinate/rotation/axesRotation.H

@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2017-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -73,9 +73,10 @@ public:
     //  Note that these follow the right-hand rule.
     enum axisOrder
     {
-        E1_E2,  //!< The axis is X-dominant, the direction is Y-dominant
+        E1_E2,  //!< The axis1 (dominant) is local X, axis2 is local Y.
-        E2_E3,  //!< The axis is Y-dominant, the direction is Z-dominant
+        E2_E3,  //!< The axis1 (dominant) is local Y, axis2 is local Z.
-        E3_E1   //!< The axis is Z-dominant, the direction is X-dominant
+        E3_E1,  //!< The axis1 (dominant) is local Z, axis2 is local X.
+        E3_E1_COMPAT, //!< E3_E1 specified as axis/direction.
     };
 
 
@@ -90,6 +91,12 @@ private:
         tensor Rtr_;
 
 
+    // Protected Member Functions
+
+        //- Read from dictionary
+        void read(const dictionary& dict);
+
+
 public:
 
     //- Runtime type information
@@ -149,14 +156,19 @@ public:
             Rtr_ = sphericalTensor::I;
         }
 
-        //- Set the transformation tensors from two axes (axis and direction)
+
-        void setTransform
+    // Static Member Functions
+
+        //- The rotation tensor calculated from two axes and their order.
+        //  The input axes will be normalised.
+        static tensor rotation
         (
             const vector& axis1,
             const vector& axis2,
-            const axisOrder& order = E3_E1
+            axisOrder order = E3_E1
         );
 
+
         //- Update the rotation for a list of cells
         virtual void updateCells(const polyMesh&, const labelList&)
         {}

src/meshTools/coordinate/rotation/coordinateRotation.C

@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2017-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -37,41 +37,65 @@ namespace Foam
 }
 
 
-// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
+
+Foam::vector Foam::coordinateRotation::findOrthogonal(const vector& axis)
+{
+    direction maxCmpt = 0;
+    scalar maxVal = mag(axis[maxCmpt]);
+
+    for (direction cmpt=1; cmpt < vector::nComponents; ++cmpt)
+    {
+        const scalar val = mag(axis[cmpt]);
+
+        if (maxVal < val)
+        {
+            maxVal  = val;
+            maxCmpt = cmpt;
+        }
+    }
+
+    direction cmpt = ((maxCmpt == vector::nComponents-1) ? 0 : (maxCmpt+1));
+
+    vector dirn(Zero);
+    dirn.component(cmpt) = ((axis[maxCmpt] < 0) ? -1 : 1);
+
+    return dirn;
+}
+
 
 Foam::symmTensor Foam::coordinateRotation::transformPrincipal
 (
     const tensor& tt,
-    const vector& st
+    const vector& v
-) const
+)
 {
     return symmTensor
     (
-        tt.xx()*st.x()*tt.xx()
+        tt.xx()*v.x()*tt.xx()
-      + tt.xy()*st.y()*tt.xy()
+      + tt.xy()*v.y()*tt.xy()
-      + tt.xz()*st.z()*tt.xz(),
+      + tt.xz()*v.z()*tt.xz(),
 
-        tt.xx()*st.x()*tt.yx()
+        tt.xx()*v.x()*tt.yx()
-      + tt.xy()*st.y()*tt.yy()
+      + tt.xy()*v.y()*tt.yy()
-      + tt.xz()*st.z()*tt.yz(),
+      + tt.xz()*v.z()*tt.yz(),
 
-        tt.xx()*st.x()*tt.zx()
+        tt.xx()*v.x()*tt.zx()
-      + tt.xy()*st.y()*tt.zy()
+      + tt.xy()*v.y()*tt.zy()
-      + tt.xz()*st.z()*tt.zz(),
+      + tt.xz()*v.z()*tt.zz(),
 
-        tt.yx()*st.x()*tt.yx()
+        tt.yx()*v.x()*tt.yx()
-      + tt.yy()*st.y()*tt.yy()
+      + tt.yy()*v.y()*tt.yy()
-      + tt.yz()*st.z()*tt.yz(),
+      + tt.yz()*v.z()*tt.yz(),
 
-        tt.yx()*st.x()*tt.zx()
+        tt.yx()*v.x()*tt.zx()
-      + tt.yy()*st.y()*tt.zy()
+      + tt.yy()*v.y()*tt.zy()
-      + tt.yz()*st.z()*tt.zz(),
+      + tt.yz()*v.z()*tt.zz(),
 
-        tt.zx()*st.x()*tt.zx()
+        tt.zx()*v.x()*tt.zx()
-      + tt.zy()*st.y()*tt.zy()
+      + tt.zy()*v.y()*tt.zy()
-      + tt.zz()*st.z()*tt.zz()
+      + tt.zz()*v.z()*tt.zz()
     );
-
 }
 
 

src/meshTools/coordinate/rotation/coordinateRotation.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.
@@ -75,12 +75,14 @@ protected:
     // Protected member functions
 
         //- Transform principal
-        symmTensor transformPrincipal(const tensor&, const vector&) const;
+        static symmTensor transformPrincipal(const tensor&, const vector&);
+
+        //- Determine best-guess for an orthogonal axis
+        static vector findOrthogonal(const vector& axis);
 
 
 public:
 
-
     //- Runtime type information
     TypeName("coordinateRotation");
 
@@ -135,8 +137,7 @@ public:
 
 
     //- Destructor
-    virtual ~coordinateRotation()
+    virtual ~coordinateRotation() = default;
-    {}
 
 
     // Member Functions