ENH: add sign(), csign() methods for complex

- use std::hypot for complex mag() instead of long-hand version

- Detail::conj() function for complex or non-complex

applications/test/complex/Test-complex.C

@@ -24,11 +24,12 @@ License
 Application
 
 Description
-    Some tests for complex numbers
+    Tests for complex numbers
 
 \*---------------------------------------------------------------------------*/
 
 #include "argList.H"
+#include "complex.H"
 #include "complexFields.H"
 #include "ops.H"
 #include "ListOps.H"
@@ -41,8 +42,7 @@ void print1(const complex& z)
 }
 
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-//  Main program:
+// * * * * * * * * * * * * * * * Main Program  * * * * * * * * * * * * * * * //
 
 int main(int argc, char *argv[])
 {
@@ -113,11 +113,81 @@ int main(int argc, char *argv[])
 
     // Info<< "operator * : " << (fld1 * fld2) << nl;
     // Info<< "operator / : " << (fld1 / fld2) << nl;
-    Info<< "operator / : " << (fld1 / 2) << nl;
+    // Info<< "operator / : " << (fld1 / 2) << nl;
     // Info<< "operator / : " << (fld1 / fld2) << nl;
-    Info<< "sqrt   : " << sqrt(fld1) << nl;
+    // Info<< "sqrt   : " << sqrt(fld1) << nl;
     // Info<< "pow(2) : " << pow(fld1, 2) << nl;
 
+    Info<< nl << "Elementary complex arithmetic operations:" << nl;
+    {
+        complex a (6, 1);
+        complex b = a;
+
+        Info << "Compound assignment operations:" << nl;
+
+        // Multiplication
+        b *= a;
+        Info<< "b *= a:" << tab << "b=" << b << nl;
+
+        // Addition
+        b += a;
+        Info<< "b += a:" << tab << "b=" << b << nl;
+
+        // Subtraction
+        b -= a;
+        Info<< "b -= a:" << tab << "b=" << b << nl;
+
+        // Division
+        b /= a;
+        Info<< "b /= a:" << tab << "b=" << b << nl;
+
+        Info << "Operations with scalars:" << nl;
+
+        Info<< "b=" << b << nl;
+
+        // Scalar multiplication
+        b *= 2.0;
+        Info<< "b*2 (elementwise multiplication):" << tab << b << nl;
+
+        // Scalar addition
+        b += 1.0;
+        Info<< "b + 1 (only real part):" << tab << b << nl;
+
+        // Scalar subtraction
+        b -= 1.0;
+        Info<< "b - 1 (only real part):" << tab << b << nl;
+
+        // Scalar division
+        b = 1.0/b;
+        Info<< "1/b (elementwise division):" << tab << b << nl;
+
+   }
+
+    Info<< nl << "Other mathematical expressions:" << nl;
+    {
+        complex a (4.3, -3.14);
+        complex b (0, -4.3);
+
+        Info<< "a=" << a << tab << "b=" << b << nl;
+
+        // Square-root
+        //Info<< "sqrt(a)=" << sqrt(a) << tab << "sqrt(b)=" << sqrt(b) << nl;
+
+        // Square
+        Info<< "sqr(a)=" << sqr(a) << tab << "sqr(b)=" << sqr(b) << nl;
+
+        // n^th power
+        //Info<< "pow(a,-1)=" << pow(a,-1) << tab
+        //      << "pow(b,-1)=" << pow(b,-1) << nl;
+
+        // Exponential
+        //Info<< "exp(a)=" << exp(a) << tab << "exp(b)=" << exp(b) << nl;
+
+        // Natural logarithm
+        //Info<< "log(a)=" << log(a) << tab << "log(b)=" << log(b) << nl;
+    }
+
+    Info<< nl << "End" << nl;
 
     // Make some changes
     {
@@ -150,6 +220,7 @@ int main(int argc, char *argv[])
     // Info<< "min/max = " << MinMax<complex>(fld1) << nl;
 
     Info<< "\nEnd\n" << endl;
+
     return 0;
 }
 

src/OpenFOAM/primitives/complex/complex.C

@@ -30,7 +30,6 @@ License
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
-const char* const Foam::complex::typeName = "complex";
 const Foam::complex Foam::complex::zero(0, 0);
 const Foam::complex Foam::complex::one(1, 0);
 
@@ -69,7 +68,6 @@ Foam::pTraits<Foam::complex>::pTraits(Istream& is)
 }
 
 
-
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::complex::complex(Istream& is)
@@ -78,7 +76,7 @@ Foam::complex::complex(Istream& is)
 }
 
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+// * * * * * * * * * * * * * * * Global Functions  * * * * * * * * * * * * * //
 
 Foam::word Foam::name(const complex& c)
 {

src/OpenFOAM/primitives/complex/complex.H

@@ -39,6 +39,7 @@ SourceFiles
 #define complex_H
 
 #include <complex>
+#include <type_traits>
 #include "scalar.H"
 #include "word.H"
 #include "zero.H"
@@ -91,7 +92,7 @@ public:
     // Static Data Members
 
         //- The type name is "complex"
-        static const char* const typeName;
+        static constexpr const char* const typeName = "complex";
 
         //- A complex zero (0,0)
         static const complex zero;
@@ -209,12 +210,17 @@ public:
         friend scalar magSqr(const complex& c);
         friend scalar mag(const complex& c);
         friend complex sqr(const complex& c);
-        friend const complex& min(const complex&, const complex&);
-        friend const complex& max(const complex&, const complex&);
 
-        friend complex limit(const complex&, const complex&);
+        //- sgn() https://en.wikipedia.org/wiki/Sign_function#Complex_signum
+        friend complex sign(const complex& c);
 
-        friend const complex& sum(const complex&);
+        //- csgn() https://en.wikipedia.org/wiki/Sign_function#Complex_signum
+        friend scalar csign(const complex& c);
+
+        friend const complex& min(const complex& c1, const complex& c2);
+        friend const complex& max(const complex& c1, const complex& c2);
+        friend complex limit(const complex& c1, const complex& c2);
+        friend const complex& sum(const complex& c);
 
 
     // Friend Operators
@@ -232,6 +238,10 @@ public:
 };
 
 
+/*---------------------------------------------------------------------------*\
+                      Class pTraits<complex> Declaration
+\*---------------------------------------------------------------------------*/
+
 // Template specialisation for pTraits<complex>
 template<>
 class pTraits<complex>
@@ -247,7 +257,7 @@ public:
     typedef label labelType;
 
 
-    // Member constants
+    // Member Constants
 
         //- Dimensionality of space
         static constexpr direction dim = 3;
@@ -296,6 +306,41 @@ public:
 };
 
 
+/*---------------------------------------------------------------------------*\
+                               Namespace Detail
+\*---------------------------------------------------------------------------*/
+
+namespace Detail
+{
+    // Helper functions for complex, in Detail namespace to avoid possible
+    // name collisions (could change in the future)
+
+    //- The 'conjugate' of non-complex returns itself (pass-through)
+    //- it does not return a complex!
+    template<class T>
+    typename std::enable_if
+    <
+        !std::is_same<complex, T>::value,
+        const T&
+    >::type conj(const T& val)
+    {
+        return val;
+    }
+
+    //- The conjugate of a complex number
+    template<class T>
+    typename std::enable_if
+    <
+        std::is_same<complex, T>::value,
+        complex
+    >::type conj(const T& val)
+    {
+        return val.conjugate();
+    }
+
+} // End namespace Detail
+
+
 // * * * * * * * * * * * * * * * Global Operators  * * * * * * * * * * * * * //
 
 Istream& operator>>(Istream& is, complex& c);

src/OpenFOAM/primitives/complex/complexI.H

@@ -221,7 +221,7 @@ inline scalar magSqr(const complex& c)
 
 inline scalar mag(const complex& c)
 {
-    return sqrt(magSqr(c)); // OR std::hypot(c.re, c.im);
+    return std::hypot(c.re, c.im);
 }
 
 
@@ -231,6 +231,19 @@ inline complex sqr(const complex& c)
 }
 
 
+inline complex sign(const complex& c)
+{
+    const scalar s(mag(c));
+    return s < ROOTVSMALL ? Zero : c/s;
+}
+
+
+inline scalar csign(const complex& c)
+{
+    return equal(c.Re(), 0) ? sign(c.Im()) : sign(c.Re());
+}
+
+
 inline const complex& min(const complex& c1, const complex& c2)
 {
     if (magSqr(c1) < magSqr(c2))