From eef3c7ceb68143831597f02aaadd446266949fac Mon Sep 17 00:00:00 2001 From: andy Date: Wed, 9 Sep 2009 16:30:26 +0100 Subject: [PATCH] added pow025 function - graham showed significant speedup (8-9x) for sqrt(sqrt(x)) compared to pow(x, 0.25) --- src/OpenFOAM/dimensionSet/dimensionSet.C | 33 +++++++++++++++++++ src/OpenFOAM/dimensionSet/dimensionSet.H | 2 ++ .../dimensionedScalar/dimensionedScalar.C | 27 +++++++++++++++ .../dimensionedScalar/dimensionedScalar.H | 1 + .../DimensionedScalarField.C | 1 + .../DimensionedScalarField.H | 1 + .../scalarFieldField/scalarFieldField.C | 1 + .../scalarFieldField/scalarFieldField.H | 3 ++ .../fields/Fields/scalarField/scalarField.C | 1 + .../fields/Fields/scalarField/scalarField.H | 1 + .../GeometricScalarField.C | 1 + .../GeometricScalarField.H | 1 + src/OpenFOAM/primitives/Scalar/Scalar.H | 29 ++++++++++++++++ 13 files changed, 102 insertions(+) diff --git a/src/OpenFOAM/dimensionSet/dimensionSet.C b/src/OpenFOAM/dimensionSet/dimensionSet.C index 9fba472d98..e60a3998d4 100644 --- a/src/OpenFOAM/dimensionSet/dimensionSet.C +++ b/src/OpenFOAM/dimensionSet/dimensionSet.C @@ -110,6 +110,7 @@ Foam::scalar Foam::dimensionSet::operator[](const dimensionType type) const return exponents_[type]; } + Foam::scalar& Foam::dimensionSet::operator[](const dimensionType type) { return exponents_[type]; @@ -130,6 +131,7 @@ bool Foam::dimensionSet::operator==(const dimensionSet& ds) const return equall; } + bool Foam::dimensionSet::operator!=(const dimensionSet& ds) const { return !(operator==(ds)); @@ -163,6 +165,7 @@ bool Foam::dimensionSet::operator+=(const dimensionSet& ds) const return true; } + bool Foam::dimensionSet::operator-=(const dimensionSet& ds) const { if (dimensionSet::debug && *this != ds) @@ -176,6 +179,7 @@ bool Foam::dimensionSet::operator-=(const dimensionSet& ds) const return true; } + bool Foam::dimensionSet::operator*=(const dimensionSet& ds) { reset((*this)*ds); @@ -183,6 +187,7 @@ bool Foam::dimensionSet::operator*=(const dimensionSet& ds) return true; } + bool Foam::dimensionSet::operator/=(const dimensionSet& ds) { reset((*this)/ds); @@ -206,6 +211,7 @@ Foam::dimensionSet Foam::max(const dimensionSet& ds1, const dimensionSet& ds2) return ds1; } + Foam::dimensionSet Foam::min(const dimensionSet& ds1, const dimensionSet& ds2) { if (dimensionSet::debug && ds1 != ds2) @@ -256,6 +262,7 @@ Foam::dimensionSet Foam::pow(const dimensionSet& ds, const scalar p) return dimPow; } + Foam::dimensionSet Foam::pow ( const dimensionSet& ds, @@ -283,6 +290,7 @@ Foam::dimensionSet Foam::pow return dimPow; } + Foam::dimensionSet Foam::pow ( const dimensionedScalar& dS, @@ -309,61 +317,79 @@ Foam::dimensionSet Foam::sqr(const dimensionSet& ds) return pow(ds, 2); } + Foam::dimensionSet Foam::pow3(const dimensionSet& ds) { return pow(ds, 3); } + Foam::dimensionSet Foam::pow4(const dimensionSet& ds) { return pow(ds, 4); } + Foam::dimensionSet Foam::pow5(const dimensionSet& ds) { return pow(ds, 5); } + Foam::dimensionSet Foam::pow6(const dimensionSet& ds) { return pow(ds, 6); } + +Foam::dimensionSet Foam::pow025(const dimensionSet& ds) +{ + return sqrt(sqrt(ds)); +} + + Foam::dimensionSet Foam::sqrt(const dimensionSet& ds) { return pow(ds, 0.5); } + Foam::dimensionSet Foam::magSqr(const dimensionSet& ds) { return pow(ds, 2); } + Foam::dimensionSet Foam::mag(const dimensionSet& ds) { return ds; } + Foam::dimensionSet Foam::sign(const dimensionSet&) { return dimless; } + Foam::dimensionSet Foam::pos(const dimensionSet&) { return dimless; } + Foam::dimensionSet Foam::neg(const dimensionSet&) { return dimless; } + Foam::dimensionSet Foam::inv(const dimensionSet& ds) { return dimless/ds; } + Foam::dimensionSet Foam::trans(const dimensionSet& ds) { if (dimensionSet::debug && !ds.dimensionless()) @@ -376,6 +402,7 @@ Foam::dimensionSet Foam::trans(const dimensionSet& ds) return ds; } + Foam::dimensionSet Foam::transform(const dimensionSet& ds) { return ds; @@ -389,6 +416,7 @@ Foam::dimensionSet Foam::operator-(const dimensionSet& ds) return ds; } + Foam::dimensionSet Foam::operator+ ( const dimensionSet& ds1, @@ -409,6 +437,7 @@ Foam::dimensionSet Foam::operator+ return dimSum; } + Foam::dimensionSet Foam::operator- ( const dimensionSet& ds1, @@ -429,6 +458,7 @@ Foam::dimensionSet Foam::operator- return dimDifference; } + Foam::dimensionSet Foam::operator* ( const dimensionSet& ds1, @@ -445,6 +475,7 @@ Foam::dimensionSet Foam::operator* return dimProduct; } + Foam::dimensionSet Foam::operator/ ( const dimensionSet& ds1, @@ -471,6 +502,7 @@ Foam::dimensionSet Foam::operator& return ds1*ds2; } + Foam::dimensionSet Foam::operator^ ( const dimensionSet& ds1, @@ -480,6 +512,7 @@ Foam::dimensionSet Foam::operator^ return ds1*ds2; } + Foam::dimensionSet Foam::operator&& ( const dimensionSet& ds1, diff --git a/src/OpenFOAM/dimensionSet/dimensionSet.H b/src/OpenFOAM/dimensionSet/dimensionSet.H index e61d3ec4df..13b0a899a9 100644 --- a/src/OpenFOAM/dimensionSet/dimensionSet.H +++ b/src/OpenFOAM/dimensionSet/dimensionSet.H @@ -70,6 +70,7 @@ 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 magSqr(const dimensionSet&); @@ -226,6 +227,7 @@ public: 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&); diff --git a/src/OpenFOAM/dimensionedTypes/dimensionedScalar/dimensionedScalar.C b/src/OpenFOAM/dimensionedTypes/dimensionedScalar/dimensionedScalar.C index 8761396a0f..054d1eb43b 100644 --- a/src/OpenFOAM/dimensionedTypes/dimensionedScalar/dimensionedScalar.C +++ b/src/OpenFOAM/dimensionedTypes/dimensionedScalar/dimensionedScalar.C @@ -38,32 +38,38 @@ dimensionedScalar operator+(const dimensionedScalar& ds1, const scalar s2) return ds1 + dimensionedScalar(s2); } + dimensionedScalar operator+(const scalar s1, const dimensionedScalar& ds2) { return dimensionedScalar(s1) + ds2; } + dimensionedScalar operator-(const dimensionedScalar& ds1, const scalar s2) { return ds1 - dimensionedScalar(s2); } + dimensionedScalar operator-(const scalar s1, const dimensionedScalar& ds2) { return dimensionedScalar(s1) - ds2; } + dimensionedScalar operator*(const dimensionedScalar& ds1, const scalar s2) { return ds1 * dimensionedScalar(s2); } + dimensionedScalar operator/(const scalar s1, const dimensionedScalar& ds2) { return dimensionedScalar(s1)/ds2; } + dimensionedScalar pow ( const dimensionedScalar& ds, @@ -78,6 +84,7 @@ dimensionedScalar pow ); } + dimensionedScalar pow3(const dimensionedScalar& ds) { return dimensionedScalar @@ -88,6 +95,7 @@ dimensionedScalar pow3(const dimensionedScalar& ds) ); } + dimensionedScalar pow4(const dimensionedScalar& ds) { return dimensionedScalar @@ -98,6 +106,7 @@ dimensionedScalar pow4(const dimensionedScalar& ds) ); } + dimensionedScalar pow5(const dimensionedScalar& ds) { return dimensionedScalar @@ -108,6 +117,7 @@ dimensionedScalar pow5(const dimensionedScalar& ds) ); } + dimensionedScalar pow6(const dimensionedScalar& ds) { return dimensionedScalar @@ -118,6 +128,18 @@ dimensionedScalar pow6(const dimensionedScalar& ds) ); } + +dimensionedScalar pow025(const dimensionedScalar& ds) +{ + return dimensionedScalar + ( + "pow025(" + ds.name() + ')', + pow025(ds.dimensions()), + pow025(ds.value()) + ); +} + + dimensionedScalar sqrt(const dimensionedScalar& ds) { return dimensionedScalar @@ -128,6 +150,7 @@ dimensionedScalar sqrt(const dimensionedScalar& ds) ); } + dimensionedScalar cbrt(const dimensionedScalar& ds) { return dimensionedScalar @@ -138,6 +161,7 @@ dimensionedScalar cbrt(const dimensionedScalar& ds) ); } + dimensionedScalar hypot ( const dimensionedScalar& x, @@ -152,6 +176,7 @@ dimensionedScalar hypot ); } + dimensionedScalar sign(const dimensionedScalar& ds) { return dimensionedScalar @@ -162,6 +187,7 @@ dimensionedScalar sign(const dimensionedScalar& ds) ); } + dimensionedScalar pos(const dimensionedScalar& ds) { return dimensionedScalar @@ -172,6 +198,7 @@ dimensionedScalar pos(const dimensionedScalar& ds) ); } + dimensionedScalar neg(const dimensionedScalar& ds) { return dimensionedScalar diff --git a/src/OpenFOAM/dimensionedTypes/dimensionedScalar/dimensionedScalar.H b/src/OpenFOAM/dimensionedTypes/dimensionedScalar/dimensionedScalar.H index 492d5064ef..c6e46168c9 100644 --- a/src/OpenFOAM/dimensionedTypes/dimensionedScalar/dimensionedScalar.H +++ b/src/OpenFOAM/dimensionedTypes/dimensionedScalar/dimensionedScalar.H @@ -61,6 +61,7 @@ dimensionedScalar pow3(const dimensionedScalar&); dimensionedScalar pow4(const dimensionedScalar&); dimensionedScalar pow5(const dimensionedScalar&); dimensionedScalar pow6(const dimensionedScalar&); +dimensionedScalar pow025(const dimensionedScalar&); dimensionedScalar sqrt(const dimensionedScalar&); dimensionedScalar cbrt(const dimensionedScalar&); diff --git a/src/OpenFOAM/fields/DimensionedFields/DimensionedScalarField/DimensionedScalarField.C b/src/OpenFOAM/fields/DimensionedFields/DimensionedScalarField/DimensionedScalarField.C index ade367e8ff..4aad204f79 100644 --- a/src/OpenFOAM/fields/DimensionedFields/DimensionedScalarField/DimensionedScalarField.C +++ b/src/OpenFOAM/fields/DimensionedFields/DimensionedScalarField/DimensionedScalarField.C @@ -376,6 +376,7 @@ UNARY_FUNCTION(scalar, scalar, pow3, pow3) UNARY_FUNCTION(scalar, scalar, pow4, pow4) UNARY_FUNCTION(scalar, scalar, pow5, pow5) UNARY_FUNCTION(scalar, scalar, pow6, pow6) +UNARY_FUNCTION(scalar, scalar, pow025, pow025) UNARY_FUNCTION(scalar, scalar, sqrt, sqrt) UNARY_FUNCTION(scalar, scalar, sign, sign) UNARY_FUNCTION(scalar, scalar, pos, pos) diff --git a/src/OpenFOAM/fields/DimensionedFields/DimensionedScalarField/DimensionedScalarField.H b/src/OpenFOAM/fields/DimensionedFields/DimensionedScalarField/DimensionedScalarField.H index 1b5e5d22bf..e2d565283e 100644 --- a/src/OpenFOAM/fields/DimensionedFields/DimensionedScalarField/DimensionedScalarField.H +++ b/src/OpenFOAM/fields/DimensionedFields/DimensionedScalarField/DimensionedScalarField.H @@ -84,6 +84,7 @@ UNARY_FUNCTION(scalar, scalar, pow3, pow3) UNARY_FUNCTION(scalar, scalar, pow4, pow4) UNARY_FUNCTION(scalar, scalar, pow5, pow5) UNARY_FUNCTION(scalar, scalar, pow6, pow6) +UNARY_FUNCTION(scalar, scalar, pow025, pow025) UNARY_FUNCTION(scalar, scalar, sqrt, sqrt) UNARY_FUNCTION(scalar, scalar, sign, sign) UNARY_FUNCTION(scalar, scalar, pos, pos) diff --git a/src/OpenFOAM/fields/FieldFields/scalarFieldField/scalarFieldField.C b/src/OpenFOAM/fields/FieldFields/scalarFieldField/scalarFieldField.C index fe9067a8f5..8a308adbb5 100644 --- a/src/OpenFOAM/fields/FieldFields/scalarFieldField/scalarFieldField.C +++ b/src/OpenFOAM/fields/FieldFields/scalarFieldField/scalarFieldField.C @@ -104,6 +104,7 @@ UNARY_FUNCTION(scalar, scalar, pow3) UNARY_FUNCTION(scalar, scalar, pow4) UNARY_FUNCTION(scalar, scalar, pow5) UNARY_FUNCTION(scalar, scalar, pow6) +UNARY_FUNCTION(scalar, scalar, pow025) UNARY_FUNCTION(scalar, scalar, sqrt) UNARY_FUNCTION(scalar, scalar, sign) UNARY_FUNCTION(scalar, scalar, pos) diff --git a/src/OpenFOAM/fields/FieldFields/scalarFieldField/scalarFieldField.H b/src/OpenFOAM/fields/FieldFields/scalarFieldField/scalarFieldField.H index 45f1e86ea0..019b731979 100644 --- a/src/OpenFOAM/fields/FieldFields/scalarFieldField/scalarFieldField.H +++ b/src/OpenFOAM/fields/FieldFields/scalarFieldField/scalarFieldField.H @@ -57,6 +57,7 @@ void stabilise const scalar s ); + template class Field> tmp > stabilise ( @@ -64,6 +65,7 @@ tmp > stabilise const scalar s ); + template class Field> tmp > stabilise ( @@ -95,6 +97,7 @@ UNARY_FUNCTION(scalar, scalar, pow3) UNARY_FUNCTION(scalar, scalar, pow4) UNARY_FUNCTION(scalar, scalar, pow5) UNARY_FUNCTION(scalar, scalar, pow6) +UNARY_FUNCTION(scalar, scalar, pow025) UNARY_FUNCTION(scalar, scalar, sqrt) UNARY_FUNCTION(scalar, scalar, sign) UNARY_FUNCTION(scalar, scalar, pos) diff --git a/src/OpenFOAM/fields/Fields/scalarField/scalarField.C b/src/OpenFOAM/fields/Fields/scalarField/scalarField.C index 4f2eccfa9c..ed8e2da5e3 100644 --- a/src/OpenFOAM/fields/Fields/scalarField/scalarField.C +++ b/src/OpenFOAM/fields/Fields/scalarField/scalarField.C @@ -109,6 +109,7 @@ UNARY_FUNCTION(scalar, scalar, pow3) UNARY_FUNCTION(scalar, scalar, pow4) UNARY_FUNCTION(scalar, scalar, pow5) UNARY_FUNCTION(scalar, scalar, pow6) +UNARY_FUNCTION(scalar, scalar, pow025) UNARY_FUNCTION(scalar, scalar, sqrt) UNARY_FUNCTION(scalar, scalar, sign) UNARY_FUNCTION(scalar, scalar, pos) diff --git a/src/OpenFOAM/fields/Fields/scalarField/scalarField.H b/src/OpenFOAM/fields/Fields/scalarField/scalarField.H index a1d0865d20..bd25d64a42 100644 --- a/src/OpenFOAM/fields/Fields/scalarField/scalarField.H +++ b/src/OpenFOAM/fields/Fields/scalarField/scalarField.H @@ -96,6 +96,7 @@ UNARY_FUNCTION(scalar, scalar, pow3) UNARY_FUNCTION(scalar, scalar, pow4) UNARY_FUNCTION(scalar, scalar, pow5) UNARY_FUNCTION(scalar, scalar, pow6) +UNARY_FUNCTION(scalar, scalar, pow025) UNARY_FUNCTION(scalar, scalar, sqrt) UNARY_FUNCTION(scalar, scalar, sign) UNARY_FUNCTION(scalar, scalar, pos) diff --git a/src/OpenFOAM/fields/GeometricFields/GeometricScalarField/GeometricScalarField.C b/src/OpenFOAM/fields/GeometricFields/GeometricScalarField/GeometricScalarField.C index fa6a2a58d0..3aaf5ddf3b 100644 --- a/src/OpenFOAM/fields/GeometricFields/GeometricScalarField/GeometricScalarField.C +++ b/src/OpenFOAM/fields/GeometricFields/GeometricScalarField/GeometricScalarField.C @@ -447,6 +447,7 @@ UNARY_FUNCTION(scalar, scalar, pow3, pow3) UNARY_FUNCTION(scalar, scalar, pow4, pow4) UNARY_FUNCTION(scalar, scalar, pow5, pow5) UNARY_FUNCTION(scalar, scalar, pow6, pow6) +UNARY_FUNCTION(scalar, scalar, pow025, pow025) UNARY_FUNCTION(scalar, scalar, sqrt, sqrt) UNARY_FUNCTION(scalar, scalar, sign, sign) UNARY_FUNCTION(scalar, scalar, pos, pos) diff --git a/src/OpenFOAM/fields/GeometricFields/GeometricScalarField/GeometricScalarField.H b/src/OpenFOAM/fields/GeometricFields/GeometricScalarField/GeometricScalarField.H index 4a2d9166b4..3b50684200 100644 --- a/src/OpenFOAM/fields/GeometricFields/GeometricScalarField/GeometricScalarField.H +++ b/src/OpenFOAM/fields/GeometricFields/GeometricScalarField/GeometricScalarField.H @@ -92,6 +92,7 @@ UNARY_FUNCTION(scalar, scalar, pow3, pow3) UNARY_FUNCTION(scalar, scalar, pow4, pow4) UNARY_FUNCTION(scalar, scalar, pow5, pow5) UNARY_FUNCTION(scalar, scalar, pow6, pow6) +UNARY_FUNCTION(scalar, scalar, pow025, pow025) UNARY_FUNCTION(scalar, scalar, sqrt, sqrt) UNARY_FUNCTION(scalar, scalar, sign, sign) UNARY_FUNCTION(scalar, scalar, pos, pos) diff --git a/src/OpenFOAM/primitives/Scalar/Scalar.H b/src/OpenFOAM/primitives/Scalar/Scalar.H index cf5d9cb31d..69961ab2cc 100644 --- a/src/OpenFOAM/primitives/Scalar/Scalar.H +++ b/src/OpenFOAM/primitives/Scalar/Scalar.H @@ -90,56 +90,67 @@ inline Scalar& setComponent(Scalar& s, const direction) return s; } + inline Scalar component(const Scalar s, const direction) { return s; } + inline Scalar sign(const Scalar s) { return (s >= 0)? 1: -1; } + inline Scalar pos(const Scalar s) { return (s >= 0)? 1: 0; } + inline Scalar neg(const Scalar s) { return (s < 0)? 1: 0; } + inline bool equal(const Scalar& s1, const Scalar& s2) { return mag(s1 - s2) <= ScalarVSMALL; } + inline bool notEqual(const Scalar s1, const Scalar s2) { return mag(s1 - s2) > ScalarVSMALL; } + inline Scalar limit(const Scalar s1, const Scalar s2) { return (mag(s1) < mag(s2)) ? s1: 0.0; } + inline Scalar minMod(const Scalar s1, const Scalar s2) { return (mag(s1) < mag(s2)) ? s1: s2; } + inline Scalar magSqr(const Scalar s) { return s*s; } + inline Scalar sqr(const Scalar s) { return s*s; } + inline Scalar sqrtSumSqr(const Scalar a, const Scalar b) { Scalar maga = mag(a); @@ -155,61 +166,79 @@ inline Scalar sqrtSumSqr(const Scalar a, const Scalar b) } } + inline Scalar pow3(const Scalar s) { return s*sqr(s); } + inline Scalar pow4(const Scalar s) { return sqr(sqr(s)); } + inline Scalar pow5(const Scalar s) { return s*pow4(s); } + inline Scalar pow6(const Scalar s) { return pow3(sqr(s)); } + +inline Scalar pow025(const Scalar s) +{ + return sqrt(sqrt(s)); +} + + inline Scalar inv(const Scalar s) { return 1.0/s; } + inline Scalar dot(const Scalar s1, const Scalar s2) { return s1*s2; } + inline Scalar cmptMultiply(const Scalar s1, const Scalar s2) { return s1*s2; } + inline Scalar cmptDivide(const Scalar s1, const Scalar s2) { return s1/s2; } + inline Scalar cmptMax(const Scalar s) { return s; } + inline Scalar cmptMin(const Scalar s) { return s; } + inline Scalar cmptAv(const Scalar s) { return s; } + inline Scalar cmptMag(const Scalar s) { return mag(s);