diff --git a/applications/test/dimensionSet/Make/files b/applications/test/dimensionSet/Make/files
new file mode 100644
index 0000000000..43d84d1cbc
--- /dev/null
+++ b/applications/test/dimensionSet/Make/files
@@ -0,0 +1,3 @@
+Test-dimensionSet.C
+
+EXE = $(FOAM_USER_APPBIN)/Test-dimensionSet
diff --git a/applications/test/dimensionSet/Make/options b/applications/test/dimensionSet/Make/options
new file mode 100644
index 0000000000..18e6fe47af
--- /dev/null
+++ b/applications/test/dimensionSet/Make/options
@@ -0,0 +1,2 @@
+/* EXE_INC = */
+/* EXE_LIBS = */
diff --git a/applications/test/dimensionSet/Test-dimensionSet.C b/applications/test/dimensionSet/Test-dimensionSet.C
new file mode 100644
index 0000000000..1c48defacc
--- /dev/null
+++ b/applications/test/dimensionSet/Test-dimensionSet.C
@@ -0,0 +1,101 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 OpenCFD Ltd.
+     \\/     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 <http://www.gnu.org/licenses/>.
+
+Description
+    Print values of predefined dimensionSets, and some other tests
+
+\*---------------------------------------------------------------------------*/
+
+#include "dimensionSet.H"
+#include "IOmanip.H"
+
+using namespace Foam;
+
+// Macros to stringify macro contents.
+#define STRINGIFY(content)      #content
+#define STRING_QUOTE(input)     STRINGIFY(input)
+
+#define PRINT_DIMSET(arg)       \
+    Info<< STRING_QUOTE(arg) << "  " << arg << nl
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+// Main program:
+
+int main(int argc, char *argv[])
+{
+    Info<< "Predefined dimensionSets" << nl << nl;
+
+    PRINT_DIMSET(dimless);
+    PRINT_DIMSET(dimMass);
+    PRINT_DIMSET(dimLength);
+    PRINT_DIMSET(dimTime);
+    PRINT_DIMSET(dimTemperature);
+    PRINT_DIMSET(dimMoles);
+    PRINT_DIMSET(dimCurrent);
+    PRINT_DIMSET(dimLuminousIntensity);
+
+    PRINT_DIMSET(dimArea);
+    PRINT_DIMSET(dimVolume);
+    PRINT_DIMSET(dimVol);
+
+    PRINT_DIMSET(dimVelocity);
+    PRINT_DIMSET(dimAcceleration);
+
+    PRINT_DIMSET(dimDensity);
+    PRINT_DIMSET(dimForce);
+    PRINT_DIMSET(dimEnergy);
+    PRINT_DIMSET(dimPower);
+
+    PRINT_DIMSET(dimPressure);
+    PRINT_DIMSET(dimCompressibility);
+    PRINT_DIMSET(dimGasConstant);
+    PRINT_DIMSET(dimSpecificHeatCapacity);
+    PRINT_DIMSET(dimViscosity);
+    PRINT_DIMSET(dimDynamicViscosity);
+
+
+    Info<< nl << "Operators" << nl << nl;
+
+    // Operators
+    {
+        Info<< "dimLength/dimTime  " << (dimLength/dimTime) << nl;
+        Info<< "1/dimTime  " << (dimless/dimTime) << nl;
+        Info<< "-dimTime  " << (-dimTime) << nl;
+        Info<< "~dimTime  " << (~dimTime) << nl;
+
+        Info<< "sqr(dimVelocity) " << sqr(dimVelocity) << nl;
+        Info<< "pow2(dimVelocity) " << pow2(dimVelocity) << nl;
+        Info<< "sqrt(dimVelocity) " << sqrt(dimVelocity) << nl;
+        Info<< "pow025(dimVelocity) " << pow025(dimVelocity) << nl;
+        Info<< "sqrt(sqrt(dimVelocity)) " << sqrt(sqrt(dimVelocity)) << nl;
+    }
+
+
+    Info<< nl << "End\n" << endl;
+
+    return 0;
+}
+
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/dimensionSet/dimensionSet.C b/src/OpenFOAM/dimensionSet/dimensionSet.C
index d2ec48e3f8..3386fd370b 100644
--- a/src/OpenFOAM/dimensionSet/dimensionSet.C
+++ b/src/OpenFOAM/dimensionSet/dimensionSet.C
@@ -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.
@@ -83,14 +83,10 @@ Foam::dimensionSet::dimensionSet(const dimensionSet& ds)
 
 bool Foam::dimensionSet::dimensionless() const
 {
-    for (int d=0; d<nDimensions; ++d)
+    for (const scalar& val : exponents_)
     {
-        // ie, mag(exponents_[d]) > smallExponent
-        if
-        (
-            exponents_[d] > smallExponent
-         || exponents_[d] < -smallExponent
-        )
+        // ie, mag(val) > smallExponent
+        if ((val > smallExponent) || (val < -smallExponent))
         {
             return false;
         }
@@ -173,7 +169,7 @@ bool Foam::dimensionSet::operator=(const dimensionSet& ds) const
     if (dimensionSet::debug && *this != ds)
     {
         FatalErrorInFunction
-            << "Different dimensions for =" << endl
+            << "Different dimensions for =" << nl
             << "     dimensions : " << *this << " = " << ds << endl
             << abort(FatalError);
     }
@@ -187,7 +183,7 @@ bool Foam::dimensionSet::operator+=(const dimensionSet& ds) const
     if (dimensionSet::debug && *this != ds)
     {
         FatalErrorInFunction
-            << "Different dimensions for +=" << endl
+            << "Different dimensions for +=" << nl
             << "     dimensions : " << *this << " = " << ds << endl
             << abort(FatalError);
     }
@@ -201,7 +197,7 @@ bool Foam::dimensionSet::operator-=(const dimensionSet& ds) const
     if (dimensionSet::debug && *this != ds)
     {
         FatalErrorInFunction
-            << "Different dimensions for -=" << endl
+            << "Different dimensions for -=" << nl
             << "     dimensions : " << *this << " = " << ds << endl
             << abort(FatalError);
     }
@@ -226,14 +222,14 @@ bool Foam::dimensionSet::operator/=(const dimensionSet& ds)
 }
 
 
-// * * * * * * * * * * * * * * * Friend functions * * * * * * * * * * * * * * //
+// * * * * * * * * * * * * * * Global Functions  * * * * * * * * * * * * * * //
 
-Foam::dimensionSet Foam::max(const dimensionSet& ds1, const dimensionSet& ds2)
+Foam::dimensionSet Foam::min(const dimensionSet& ds1, const dimensionSet& ds2)
 {
     if (dimensionSet::debug && ds1 != ds2)
     {
         FatalErrorInFunction
-            << "Arguments of max have different dimensions" << endl
+            << "Arguments of min have different dimensions" << nl
             << "     dimensions : " << ds1 << " and " << ds2 << endl
             << abort(FatalError);
     }
@@ -242,12 +238,12 @@ Foam::dimensionSet Foam::max(const dimensionSet& ds1, const dimensionSet& ds2)
 }
 
 
-Foam::dimensionSet Foam::min(const dimensionSet& ds1, const dimensionSet& ds2)
+Foam::dimensionSet Foam::max(const dimensionSet& ds1, const dimensionSet& ds2)
 {
     if (dimensionSet::debug && ds1 != ds2)
     {
         FatalErrorInFunction
-            << "Arguments of min have different dimensions" << endl
+            << "Arguments of max have different dimensions" << nl
             << "     dimensions : " << ds1 << " and " << ds2 << endl
             << abort(FatalError);
     }
@@ -278,7 +274,7 @@ Foam::dimensionSet Foam::cmptDivide
 
 Foam::dimensionSet Foam::pow(const dimensionSet& ds, const scalar p)
 {
-    dimensionSet dimPow
+    return dimensionSet
     (
         ds[dimensionSet::MASS]*p,
         ds[dimensionSet::LENGTH]*p,
@@ -288,8 +284,6 @@ Foam::dimensionSet Foam::pow(const dimensionSet& ds, const scalar p)
         ds[dimensionSet::CURRENT]*p,
         ds[dimensionSet::LUMINOUS_INTENSITY]*p
     );
-
-    return dimPow;
 }
 
 
@@ -302,22 +296,11 @@ Foam::dimensionSet Foam::pow
     if (dimensionSet::debug && !dS.dimensions().dimensionless())
     {
         FatalErrorInFunction
-            << "Exponent of pow is not dimensionless"
+            << "Exponent of pow is not dimensionless" << endl
             << abort(FatalError);
     }
 
-    dimensionSet dimPow
-    (
-        ds[dimensionSet::MASS]*dS.value(),
-        ds[dimensionSet::LENGTH]*dS.value(),
-        ds[dimensionSet::TIME]*dS.value(),
-        ds[dimensionSet::TEMPERATURE]*dS.value(),
-        ds[dimensionSet::MOLES]*dS.value(),
-        ds[dimensionSet::CURRENT]*dS.value(),
-        ds[dimensionSet::LUMINOUS_INTENSITY]*dS.value()
-    );
-
-    return dimPow;
+    return pow(ds, dS.value());
 }
 
 
@@ -349,6 +332,12 @@ Foam::dimensionSet Foam::sqr(const dimensionSet& ds)
 }
 
 
+Foam::dimensionSet Foam::pow2(const dimensionSet& ds)
+{
+    return pow(ds, 2);
+}
+
+
 Foam::dimensionSet Foam::pow3(const dimensionSet& ds)
 {
     return pow(ds, 3);
@@ -375,7 +364,7 @@ Foam::dimensionSet Foam::pow6(const dimensionSet& ds)
 
 Foam::dimensionSet Foam::pow025(const dimensionSet& ds)
 {
-    return sqrt(sqrt(ds));
+    return pow(ds, 0.25);
 }
 
 
@@ -447,7 +436,16 @@ Foam::dimensionSet Foam::negPart(const dimensionSet& ds)
 
 Foam::dimensionSet Foam::inv(const dimensionSet& ds)
 {
-    return dimless/ds;
+    return dimensionSet
+    (
+        -ds[dimensionSet::MASS],
+        -ds[dimensionSet::LENGTH],
+        -ds[dimensionSet::TIME],
+        -ds[dimensionSet::TEMPERATURE],
+        -ds[dimensionSet::MOLES],
+        -ds[dimensionSet::CURRENT],
+        -ds[dimensionSet::LUMINOUS_INTENSITY]
+    );
 }
 
 
@@ -456,7 +454,7 @@ Foam::dimensionSet Foam::trans(const dimensionSet& ds)
     if (dimensionSet::debug && !ds.dimensionless())
     {
         FatalErrorInFunction
-            << "Argument of trancendental function not dimensionless"
+            << "Argument of trancendental function not dimensionless" << nl
             << abort(FatalError);
     }
 
@@ -469,7 +467,7 @@ Foam::dimensionSet Foam::atan2(const dimensionSet& ds1, const dimensionSet& ds2)
     if (dimensionSet::debug && ds1 != ds2)
     {
         FatalErrorInFunction
-            << "Arguments of atan2 have different dimensions" << endl
+            << "Arguments of atan2 have different dimensions" << nl
             << "     dimensions : " << ds1 << " and " << ds2 << endl
             << abort(FatalError);
     }
@@ -484,7 +482,19 @@ Foam::dimensionSet Foam::transform(const dimensionSet& ds)
 }
 
 
-// * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
+Foam::dimensionSet Foam::invTransform(const dimensionSet& ds)
+{
+    return ds;
+}
+
+
+// * * * * * * * * * * * * * * Global Operators  * * * * * * * * * * * * * * //
+
+Foam::dimensionSet Foam::operator~(const dimensionSet& ds)
+{
+    return inv(ds);
+}
+
 
 Foam::dimensionSet Foam::operator-(const dimensionSet& ds)
 {
@@ -498,17 +508,15 @@ Foam::dimensionSet Foam::operator+
     const dimensionSet& ds2
 )
 {
-    dimensionSet dimSum(ds1);
-
     if (dimensionSet::debug && ds1 != ds2)
     {
         FatalErrorInFunction
-            << "LHS and RHS of + have different dimensions" << endl
+            << "LHS and RHS of '+' have different dimensions" << nl
             << "     dimensions : " << ds1 << " + " << ds2 << endl
             << abort(FatalError);
     }
 
-    return dimSum;
+    return ds1;
 }
 
 
@@ -518,17 +526,15 @@ Foam::dimensionSet Foam::operator-
     const dimensionSet& ds2
 )
 {
-    dimensionSet dimDifference(ds1);
-
     if (dimensionSet::debug && ds1 != ds2)
     {
         FatalErrorInFunction
-            << "LHS and RHS of - have different dimensions" << endl
+            << "LHS and RHS of '-' have different dimensions" << nl
             << "     dimensions : " << ds1 << " - " << ds2 << endl
             << abort(FatalError);
     }
 
-    return dimDifference;
+    return ds1;
 }
 
 
@@ -538,14 +544,17 @@ Foam::dimensionSet Foam::operator*
     const dimensionSet& ds2
 )
 {
-    dimensionSet dimProduct(ds1);
+    dimensionSet result(ds1);
 
-    for (int d=0; d<dimensionSet::nDimensions; ++d)
+    auto rhs = ds2.values().begin();
+
+    for (scalar& val : result.values())
     {
-        dimProduct.exponents_[d] += ds2.exponents_[d];
+        val += *rhs;
+        ++rhs;
     }
 
-    return dimProduct;
+    return result;
 }
 
 
@@ -555,14 +564,17 @@ Foam::dimensionSet Foam::operator/
     const dimensionSet& ds2
 )
 {
-    dimensionSet dimQuotient(ds1);
+    dimensionSet result(ds1);
 
-    for (int d=0; d<dimensionSet::nDimensions; ++d)
+    auto rhs = ds2.values().begin();
+
+    for (scalar& val : result.values())
     {
-        dimQuotient.exponents_[d] -= ds2.exponents_[d];
+        val -= *rhs;
+        ++rhs;
     }
 
-    return dimQuotient;
+    return result;
 }
 
 
diff --git a/src/OpenFOAM/dimensionSet/dimensionSet.H b/src/OpenFOAM/dimensionSet/dimensionSet.H
index 7c102b852c..4c2dc86ebb 100644
--- a/src/OpenFOAM/dimensionSet/dimensionSet.H
+++ b/src/OpenFOAM/dimensionSet/dimensionSet.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.
@@ -33,7 +33,6 @@ Description
 SourceFiles
     dimensionSet.C
     dimensionSetIO.C
-    dimensionSets.C
 
 \*---------------------------------------------------------------------------*/
 
@@ -57,63 +56,6 @@ namespace Foam
 class dimensionSet;
 class dimensionSets;
 
-// Friend functions
-
-dimensionSet max(const dimensionSet&, const dimensionSet&);
-dimensionSet min(const dimensionSet&, const dimensionSet&);
-dimensionSet cmptMultiply(const dimensionSet&, const dimensionSet&);
-dimensionSet cmptDivide(const dimensionSet&, const dimensionSet&);
-
-dimensionSet pow(const dimensionSet&, const scalar);
-dimensionSet pow(const dimensionSet&, const dimensionedScalar&);
-dimensionSet pow(const dimensionedScalar&, const dimensionSet&);
-
-dimensionSet sqr(const dimensionSet&);
-dimensionSet pow3(const dimensionSet&);
-dimensionSet pow4(const dimensionSet&);
-dimensionSet pow5(const dimensionSet&);
-dimensionSet pow6(const dimensionSet&);
-dimensionSet pow025(const dimensionSet&);
-
-dimensionSet sqrt(const dimensionSet&);
-dimensionSet cbrt(const dimensionSet&);
-dimensionSet magSqr(const dimensionSet&);
-dimensionSet mag(const dimensionSet&);
-dimensionSet sign(const dimensionSet&);
-dimensionSet pos(const dimensionSet&);
-dimensionSet pos0(const dimensionSet&);
-dimensionSet neg(const dimensionSet&);
-dimensionSet neg0(const dimensionSet&);
-dimensionSet posPart(const dimensionSet&);
-dimensionSet negPart(const dimensionSet&);
-dimensionSet inv(const dimensionSet&);
-
-// Function to check the argument is dimensionless
-//  for transcendental functions
-dimensionSet trans(const dimensionSet&);
-
-dimensionSet atan2(const dimensionSet&, const dimensionSet&);
-
-// Return the argument; transformations do not change the dimensions
-dimensionSet transform(const dimensionSet&);
-
-// Friend operators
-
-dimensionSet operator-(const dimensionSet&);
-dimensionSet operator+(const dimensionSet&, const dimensionSet&);
-dimensionSet operator-(const dimensionSet&, const dimensionSet&);
-dimensionSet operator*(const dimensionSet&, const dimensionSet&);
-dimensionSet operator/(const dimensionSet&, const dimensionSet&);
-dimensionSet operator&(const dimensionSet&, const dimensionSet&);
-dimensionSet operator^(const dimensionSet&, const dimensionSet&);
-dimensionSet operator&&(const dimensionSet&, const dimensionSet&);
-
-// IOstream operators
-
-Istream& operator>>(Istream&, dimensionSet&);
-Ostream& operator<<(Ostream&, const dimensionSet&);
-
-
 /*---------------------------------------------------------------------------*\
                         Class dimensionSet Declaration
 \*---------------------------------------------------------------------------*/
@@ -260,7 +202,7 @@ public:
         explicit dimensionSet(Istream& is);
 
 
-    // Member functions
+    // Member Functions
 
         //- Return true if it is dimensionless
         bool dimensionless() const;
@@ -274,7 +216,7 @@ public:
         void reset(const dimensionSet& ds);
 
 
-    // I/O
+    // IO
 
         //- Read using provided units. Used only in initial parsing
         Istream& read
@@ -315,7 +257,7 @@ public:
         ) const;
 
 
-    // Member operators
+    // Member Operators
 
         scalar operator[](const dimensionType) const;
         scalar& operator[](const dimensionType);
@@ -332,105 +274,82 @@ public:
         bool operator-=(const dimensionSet&) const;
         bool operator*=(const dimensionSet&);
         bool operator/=(const dimensionSet&);
-
-
-    // Friend functions
-
-        friend dimensionSet max(const dimensionSet&, const dimensionSet&);
-        friend dimensionSet min(const dimensionSet&, const dimensionSet&);
-        friend dimensionSet cmptMultiply
-        (
-            const dimensionSet&,
-            const dimensionSet&
-        );
-        friend dimensionSet cmptDivide
-        (
-            const dimensionSet&,
-            const dimensionSet&
-        );
-
-        friend dimensionSet pow(const dimensionSet&, const scalar);
-        friend dimensionSet pow(const dimensionSet&, const dimensionedScalar&);
-        friend dimensionSet pow(const dimensionedScalar&, const dimensionSet&);
-
-        friend dimensionSet sqr(const dimensionSet&);
-        friend dimensionSet pow3(const dimensionSet&);
-        friend dimensionSet pow4(const dimensionSet&);
-        friend dimensionSet pow5(const dimensionSet&);
-        friend dimensionSet pow6(const dimensionSet&);
-        friend dimensionSet pow025(const dimensionSet&);
-
-        friend dimensionSet sqrt(const dimensionSet&);
-        friend dimensionSet magSqr(const dimensionSet&);
-        friend dimensionSet mag(const dimensionSet&);
-        friend dimensionSet sign(const dimensionSet&);
-        friend dimensionSet pos0(const dimensionSet&);
-        friend dimensionSet neg(const dimensionSet&);
-        friend dimensionSet inv(const dimensionSet&);
-
-        //- Function to check the argument is dimensionless
-        //  for transcendental functions
-        friend dimensionSet trans(const dimensionSet&);
-
-        friend dimensionSet atan2(const dimensionSet&, const dimensionSet&);
-
-        //- Return the argument; transformations do not change the dimensions
-        friend dimensionSet transform(const dimensionSet&);
-
-
-    // Friend operators
-
-        friend dimensionSet operator-(const dimensionSet& ds);
-
-        friend dimensionSet operator+
-        (
-            const dimensionSet& ds1,
-            const dimensionSet& ds2
-        );
-
-        friend dimensionSet operator-
-        (
-            const dimensionSet& ds1,
-            const dimensionSet& ds2
-        );
-
-        friend dimensionSet operator*
-        (
-            const dimensionSet& ds1,
-            const dimensionSet& ds2
-        );
-
-        friend dimensionSet operator/
-        (
-            const dimensionSet& ds1,
-            const dimensionSet& ds2
-        );
-
-        friend dimensionSet operator&
-        (
-            const dimensionSet& ds1,
-            const dimensionSet& ds2
-        );
-
-        friend dimensionSet operator^
-        (
-            const dimensionSet& ds1,
-            const dimensionSet& ds2
-        );
-
-        friend dimensionSet operator&&
-        (
-            const dimensionSet& ds1,
-            const dimensionSet& ds2
-        );
-
-
-    // IOstream operators
-
-        friend Istream& operator>>(Istream&, dimensionSet&);
-        friend Ostream& operator<<(Ostream&, const dimensionSet&);
 };
 
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// IOstream Operators
+
+Istream& operator>>(Istream& is, dimensionSet& ds);
+Ostream& operator<<(Ostream& os, const dimensionSet& ds);
+
+
+// Global Functions
+
+dimensionSet min(const dimensionSet& ds1, const dimensionSet& ds2);
+dimensionSet max(const dimensionSet& ds1, const dimensionSet& ds2);
+
+dimensionSet cmptMultiply(const dimensionSet& ds1, const dimensionSet& ds2);
+dimensionSet cmptDivide(const dimensionSet& ds1, const dimensionSet& ds2);
+
+
+dimensionSet pow(const dimensionSet& ds, const scalar p);
+dimensionSet pow(const dimensionSet& ds, const dimensionedScalar& dS);
+dimensionSet pow(const dimensionedScalar& dS, const dimensionSet& ds);
+
+dimensionSet sqr(const dimensionSet& ds);
+dimensionSet pow2(const dimensionSet& ds);
+dimensionSet pow3(const dimensionSet& ds);
+dimensionSet pow4(const dimensionSet& ds);
+dimensionSet pow5(const dimensionSet& ds);
+dimensionSet pow6(const dimensionSet& ds);
+dimensionSet pow025(const dimensionSet& ds);
+
+dimensionSet sqrt(const dimensionSet& ds);
+dimensionSet cbrt(const dimensionSet& ds);
+
+dimensionSet magSqr(const dimensionSet& ds);
+dimensionSet mag(const dimensionSet& ds);
+dimensionSet sign(const dimensionSet&);
+dimensionSet pos(const dimensionSet&);
+dimensionSet pos0(const dimensionSet&);
+dimensionSet neg(const dimensionSet&);
+dimensionSet neg0(const dimensionSet&);
+dimensionSet posPart(const dimensionSet&);
+dimensionSet negPart(const dimensionSet&);
+
+//- The dimensionSet inverted
+dimensionSet inv(const dimensionSet& ds);
+
+//- Check the argument is dimensionless (for transcendental functions)
+dimensionSet trans(const dimensionSet& ds);
+
+//- Check that the arguments have the same dimensions
+dimensionSet atan2(const dimensionSet& ds1, const dimensionSet& ds2);
+
+//- Return the argument; transformations do not change the dimensions
+dimensionSet transform(const dimensionSet& ds);
+
+//- Return the argument; transformations do not change the dimensions
+dimensionSet invTransform(const dimensionSet& ds);
+
+
+// Global Operators
+
+//- The dimensionSet inverted
+dimensionSet operator~(const dimensionSet& ds);
+
+//- Unary negation.
+dimensionSet operator-(const dimensionSet& ds);
+
+dimensionSet operator+(const dimensionSet& ds1, const dimensionSet& ds2);
+dimensionSet operator-(const dimensionSet& ds1, const dimensionSet& ds2);
+dimensionSet operator*(const dimensionSet& ds1, const dimensionSet& ds2);
+dimensionSet operator/(const dimensionSet& ds1, const dimensionSet& ds2);
+dimensionSet operator&(const dimensionSet& ds1, const dimensionSet& ds2);
+dimensionSet operator^(const dimensionSet& ds1, const dimensionSet& ds2);
+dimensionSet operator&&(const dimensionSet& ds1, const dimensionSet& ds2);
+
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
diff --git a/src/OpenFOAM/dimensionSet/dimensionSetIO.C b/src/OpenFOAM/dimensionSet/dimensionSetIO.C
index fc7117ca8c..b1e0d296c5 100644
--- a/src/OpenFOAM/dimensionSet/dimensionSetIO.C
+++ b/src/OpenFOAM/dimensionSet/dimensionSetIO.C
@@ -218,22 +218,22 @@ void Foam::dimensionSet::tokeniser::putBack(const token& t)
 
 void Foam::dimensionSet::round(const scalar tol)
 {
-    for (int i=0; i < dimensionSet::nDimensions; ++i)
+    scalar integralPart;
+    for (scalar& val : exponents_)
     {
-        scalar integralPart;
-        scalar fractionalPart = std::modf(exponents_[i], &integralPart);
+        const scalar fractionalPart = std::modf(val, &integralPart);
 
         if (mag(fractionalPart-1.0) <= tol)
         {
-            exponents_[i] = 1.0+integralPart;
+            val = 1.0+integralPart;
         }
         else if (mag(fractionalPart+1.0) <= tol)
         {
-            exponents_[i] = -1.0+integralPart;
+            val = -1.0+integralPart;
         }
         else if (mag(fractionalPart) <= tol)
         {
-            exponents_[i] = integralPart;
+            val = integralPart;
         }
     }
 }
@@ -427,10 +427,7 @@ Foam::Istream& Foam::dimensionSet::read
         dimensionedScalar ds(parse(0, tis, readSet));
 
         multiplier = ds.value();
-        for (int i=0; i < dimensionSet::nDimensions; ++i)
-        {
-            exponents_[i] = ds.dimensions()[i];
-        }
+        exponents_ = ds.dimensions().values();
     }
     else
     {
@@ -524,17 +521,17 @@ Foam::Istream& Foam::dimensionSet::read
             // Parse unit
             dimensionSet symbolSet(dimless);
 
-            size_t index = symbolPow.find('^');
-            if (index != string::npos)
+            const auto index = symbolPow.find('^');
+            if (index != std::string::npos)
             {
                 const word symbol = symbolPow.substr(0, index);
-                const word exp = symbolPow.substr(index+1);
-                scalar exponent = readScalar(exp);
+                const scalar exponent = readScalar(symbolPow.substr(index+1));
 
                 dimensionedScalar s;
                 s.read(readSet.lookup(symbol), readSet);
 
                 symbolSet.reset(pow(s.dimensions(), exponent));
+
                 // Round to nearest integer if close to it
                 symbolSet.round(10*smallExponent);
                 multiplier *= Foam::pow(s.value(), exponent);
@@ -598,8 +595,8 @@ Foam::Istream& Foam::dimensionSet::read
         if (nextToken != token::END_SQR)
         {
             FatalIOErrorInFunction(is)
-                << "Expected a " << token::END_SQR << " in dimensionSet "
-                << endl << "in stream " << is.info()
+                << "Expected a " << token::END_SQR << " in dimensionSet " << nl
+                << "in stream " << is.info() << endl
                 << exit(FatalIOError);
         }
     }
@@ -693,12 +690,12 @@ Foam::Ostream& Foam::dimensionSet::write
 
 // * * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * //
 
-Foam::Istream& Foam::operator>>(Istream& is, dimensionSet& dset)
+Foam::Istream& Foam::operator>>(Istream& is, dimensionSet& ds)
 {
-    scalar multiplier;
-    dset.read(is, multiplier);
+    scalar mult(1.0);
+    ds.read(is, mult);
 
-    if (mag(multiplier-1.0) > dimensionSet::smallExponent)
+    if (mag(mult-1.0) > dimensionSet::smallExponent)
     {
         FatalIOErrorInFunction(is)
             << "Cannot use scaled units in dimensionSet"
@@ -710,10 +707,10 @@ Foam::Istream& Foam::operator>>(Istream& is, dimensionSet& dset)
 }
 
 
-Foam::Ostream& Foam::operator<<(Ostream& os, const dimensionSet& dset)
+Foam::Ostream& Foam::operator<<(Ostream& os, const dimensionSet& ds)
 {
-    scalar multiplier;
-    dset.write(os, multiplier);
+    scalar mult(1.0);
+    ds.write(os, mult);
 
     os.check(FUNCTION_NAME);
     return os;
diff --git a/src/OpenFOAM/dimensionSet/dimensionSets.C b/src/OpenFOAM/dimensionSet/dimensionSets.C
index 176b746a55..43d89065d5 100644
--- a/src/OpenFOAM/dimensionSet/dimensionSets.C
+++ b/src/OpenFOAM/dimensionSet/dimensionSets.C
@@ -109,20 +109,29 @@ const HashTable<dimensionedScalar>& unitSet()
 
         const word unitSetCoeffs(dict.get<word>("unitSet") + "Coeffs");
 
-        if (!dict.found(unitSetCoeffs))
+        const dictionary* unitDictPtr = dict.findDict(unitSetCoeffs);
+
+        if (!unitDictPtr)
         {
             FatalIOErrorInFunction(dict)
                 << "Cannot find " << unitSetCoeffs << " in dictionary "
-                << dict.name() << exit(FatalIOError);
+                << dict.name() << nl
+                << exit(FatalIOError);
         }
 
-        const dictionary& unitDict = dict.subDict(unitSetCoeffs);
+        const dictionary& unitDict = *unitDictPtr;
 
-        unitSetPtr_ = new HashTable<dimensionedScalar>(unitDict.size());
+        unitSetPtr_ = new HashTable<dimensionedScalar>(2*unitDict.size());
+
+        wordList writeUnitNames;
 
         for (const entry& dEntry : unitDict)
         {
-            if (dEntry.keyword() != "writeUnits")
+            if ("writeUnits" == dEntry.keyword())
+            {
+                dEntry.readEntry(writeUnitNames);
+            }
+            else
             {
                 dimensionedScalar dt;
                 dt.read(dEntry.stream(), unitDict);
@@ -138,17 +147,6 @@ const HashTable<dimensionedScalar>& unitSet()
             }
         }
 
-        wordList writeUnitNames
-        (
-            unitDict.lookupOrDefault<wordList>
-            (
-                "writeUnits",
-                wordList(0)
-            )
-        );
-
-        writeUnitSetPtr_ = new dimensionSets(*unitSetPtr_, writeUnitNames);
-
         if (writeUnitNames.size() != 0 && writeUnitNames.size() != 7)
         {
             FatalIOErrorInFunction(dict)
@@ -156,6 +154,9 @@ const HashTable<dimensionedScalar>& unitSet()
                 << " or it is not a wordList of size 7"
                 << exit(FatalIOError);
         }
+
+        writeUnitSetPtr_ = new dimensionSets(*unitSetPtr_, writeUnitNames);
+
     }
     return *unitSetPtr_;
 }
diff --git a/src/OpenFOAM/dimensionSet/dimensionSets.H b/src/OpenFOAM/dimensionSet/dimensionSets.H
index ecaa3cf8d5..643c8c4bb8 100644
--- a/src/OpenFOAM/dimensionSet/dimensionSets.H
+++ b/src/OpenFOAM/dimensionSet/dimensionSets.H
@@ -96,15 +96,15 @@ public:
 
     // Constructors
 
-        //- Construct from all units and set of units to use for inversion
-        //  (writing)
+        //- Construct from all units and set of units to use
+        //- for inversion (writing)
         dimensionSets
         (
             const HashTable<dimensionedScalar>&,
             const wordList& unitNames
         );
 
-    // Member functions
+    // Member Functions
 
         //- Return the units
         const PtrList<dimensionedScalar>& units() const