From 28b492bd54b163bb2904dbaf21d92a78b02d499e Mon Sep 17 00:00:00 2001
From: Mark Olesen <Mark.Olesen@esi-group.com>
Date: Fri, 24 Feb 2023 11:20:39 +0100
Subject: [PATCH] ENH: unroll loops for complexVector functions
STYLE: prefer STL naming real(), imag() instead of Re(), Im()
---
applications/test/complex/Test-complex.C | 20 ++
.../test/matrices/Matrix/Test-Matrix.C | 12 +-
.../fields/Fields/complex/complexField.C | 173 +++++++++-----
.../fields/Fields/complex/complexField.H | 77 +++++--
.../fields/Fields/complex/complexFields.H | 4 +-
.../Fields/complex/complexVectorField.C | 212 +++++++++++++-----
.../Fields/complex/complexVectorField.H | 89 ++++++--
.../matrices/EigenMatrix/EigenMatrix.C | 16 +-
.../primitives/Vector/complex/complexVector.C | 3 -
.../primitives/Vector/complex/complexVector.H | 4 +-
.../Vector/complex/complexVectorI.H | 15 ++
src/OpenFOAM/primitives/complex/complex.H | 36 +--
src/OpenFOAM/primitives/complex/complexI.H | 6 +
.../energySpectrum/energySpectrum.C | 2 +-
src/randomProcesses/fft/calcEk.C | 2 +-
src/randomProcesses/fft/fft.C | 18 +-
src/randomProcesses/turbulence/turbGen.C | 7 +-
17 files changed, 495 insertions(+), 201 deletions(-)
diff --git a/applications/test/complex/Test-complex.C b/applications/test/complex/Test-complex.C
index 7e76b09324f..e8940e07a78 100644
--- a/applications/test/complex/Test-complex.C
+++ b/applications/test/complex/Test-complex.C
@@ -145,6 +145,8 @@ int main(int argc, char *argv[])
complexField cmplx(4);
+ zip(cmplx, reals, zero{});
+ zip(cmplx, 1, imags);
zip(cmplx, reals, imags);
Info<< nl
<< "zip " << reals << nl
@@ -360,6 +362,24 @@ int main(int argc, char *argv[])
// MinMax fails since there is no less comparison operator
// Info<< "min/max = " << MinMax<complex>(fld1) << nl;
+
+ // Cross-product
+ {
+ const vector vec(1, 2, 3);
+ const vector realValue(4, 5, 6);
+ const vector imagValue(7, 8, 9);
+
+ complexVector cmplxVec(zip(realValue, imagValue));
+
+ Info<< "complexVector: " << cmplxVec << nl;
+ Info<< "cross: " << (vec ^ cmplxVec) << nl;
+
+ Info<< "cross real: " << (vec ^ realValue) << nl
+ << "cross imag: " << (vec ^ imagValue) << nl
+ << "cross : "
+ << zip((vec ^ realValue), (vec ^ imagValue)) << nl;
+ }
+
Info<< "\nEnd\n" << endl;
return 0;
diff --git a/applications/test/matrices/Matrix/Test-Matrix.C b/applications/test/matrices/Matrix/Test-Matrix.C
index 73d23b5fced..4c87f078970 100644
--- a/applications/test/matrices/Matrix/Test-Matrix.C
+++ b/applications/test/matrices/Matrix/Test-Matrix.C
@@ -468,13 +468,13 @@ int main(int argc, char *argv[])
for (auto& val : A)
{
- val.Re() = rndGen.GaussNormal<scalar>();
- val.Im() = rndGen.GaussNormal<scalar>();
+ val.real(rndGen.GaussNormal<scalar>());
+ val.imag(rndGen.GaussNormal<scalar>());
}
for (auto& val : B)
{
- val.Re() = rndGen.GaussNormal<scalar>();
- val.Im() = rndGen.GaussNormal<scalar>();
+ val.real(rndGen.GaussNormal<scalar>());
+ val.imag(rndGen.GaussNormal<scalar>());
}
Info<< nl << "# (A.T()).T() = A:" << nl;
@@ -893,8 +893,8 @@ int main(int argc, char *argv[])
for (auto& val : A)
{
- val.Re() = rndGen.GaussNormal<scalar>();
- val.Im() = rndGen.GaussNormal<scalar>();
+ val.real(rndGen.GaussNormal<scalar>());
+ val.imag(rndGen.GaussNormal<scalar>());
}
Info<< "# SquareMatrix<complex>.symmetric():" << nl
diff --git a/src/OpenFOAM/fields/Fields/complex/complexField.C b/src/OpenFOAM/fields/Fields/complex/complexField.C
index f24fafbf658..e2c08e889d4 100644
--- a/src/OpenFOAM/fields/Fields/complex/complexField.C
+++ b/src/OpenFOAM/fields/Fields/complex/complexField.C
@@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011 OpenFOAM Foundation
- Copyright (C) 2019 OpenCFD Ltd.
+ Copyright (C) 2019-2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -46,27 +46,80 @@ namespace Foam
void Foam::zip
(
complexField& result,
- const UList<scalar>& re,
- const UList<scalar>& im
+ const UList<scalar>& realValues,
+ const UList<scalar>& imagValues
)
{
const label len = result.size();
#ifdef FULLDEBUG
- if (len != re.size() || len != im.size())
+ if (len != realValues.size() || len != imagValues.size())
{
FatalErrorInFunction
<< "Components sizes do not match: " << len << " ("
- << re.size() << ' ' << im.size() << ')'
- << nl
+ << realValues.size() << ' ' << imagValues.size() << ')' << nl
<< abort(FatalError);
}
#endif
for (label i=0; i < len; ++i)
{
- result[i].Re() = re[i];
- result[i].Im() = im[i];
+ result[i].real(realValues[i]);
+ result[i].imag(imagValues[i]);
+ }
+}
+
+
+void Foam::zip
+(
+ complexField& result,
+ const UList<scalar>& realValues,
+ const scalar imagValue
+)
+{
+ const label len = result.size();
+
+ #ifdef FULLDEBUG
+ if (len != realValues.size())
+ {
+ FatalErrorInFunction
+ << "Components sizes do not match: " << len
+ << " != " << realValues.size() << nl
+ << abort(FatalError);
+ }
+ #endif
+
+ for (label i=0; i < len; ++i)
+ {
+ result[i].real(realValues[i]);
+ result[i].imag(imagValue);
+ }
+}
+
+
+void Foam::zip
+(
+ complexField& result,
+ const scalar realValue,
+ const UList<scalar>& imagValues
+)
+{
+ const label len = result.size();
+
+ #ifdef FULLDEBUG
+ if (len != imagValues.size())
+ {
+ FatalErrorInFunction
+ << "Components sizes do not match: " << len
+ << " != " << imagValues.size() << nl
+ << abort(FatalError);
+ }
+ #endif
+
+ for (label i=0; i < len; ++i)
+ {
+ result[i].real(realValue);
+ result[i].imag(imagValues[i]);
}
}
@@ -74,27 +127,26 @@ void Foam::zip
void Foam::unzip
(
const UList<complex>& input,
- scalarField& re,
- scalarField& im
+ scalarField& realValues,
+ scalarField& imagValues
)
{
const label len = input.size();
#ifdef FULLDEBUG
- if (len != re.size() || len != im.size())
+ if (len != realValues.size() || len != imagValues.size())
{
FatalErrorInFunction
<< "Components sizes do not match: " << len << " ("
- << re.size() << ' ' << im.size() << ')'
- << nl
+ << realValues.size() << ' ' << imagValues.size() << ')' << nl
<< abort(FatalError);
}
#endif
for (label i=0; i < len; ++i)
{
- re[i] = input[i].Re();
- im[i] = input[i].Im();
+ realValues[i] = input[i].real();
+ imagValues[i] = input[i].imag();
}
}
@@ -103,82 +155,91 @@ void Foam::unzip
Foam::complexField Foam::ComplexField
(
- const UList<scalar>& re,
- const UList<scalar>& im
+ const UList<scalar>& realValues,
+ const UList<scalar>& imagValues
)
{
- complexField result(re.size());
+ complexField result(realValues.size());
- Foam::zip(result, re, im);
+ Foam::zip(result, realValues, imagValues);
return result;
}
-Foam::complexField Foam::ReComplexField(const UList<scalar>& re)
+Foam::complexField Foam::ComplexField
+(
+ const UList<scalar>& realValues,
+ const scalar imagValue
+)
{
- complexField cf(re.size());
+ complexField result(realValues.size());
- forAll(cf, i)
- {
- cf[i].Re() = re[i];
- cf[i].Im() = Zero;
- }
+ Foam::zip(result, realValues, imagValue);
- return cf;
+ return result;
}
-Foam::complexField Foam::ImComplexField(const UList<scalar>& im)
+Foam::complexField Foam::ComplexField
+(
+ const scalar realValue,
+ const UList<scalar>& imagValues
+)
{
- complexField cf(im.size());
+ complexField result(imagValues.size());
- forAll(cf, i)
- {
- cf[i].Re() = Zero;
- cf[i].Im() = im[i];
- }
+ Foam::zip(result, realValue, imagValues);
- return cf;
+ return result;
}
-Foam::scalarField Foam::ReImSum(const UList<complex>& cf)
+Foam::scalarField Foam::ReImSum(const UList<complex>& cmplx)
{
- scalarField sf(cf.size());
+ scalarField result(cmplx.size());
- forAll(sf, i)
- {
- sf[i] = cf[i].Re() + cf[i].Im();
- }
+ std::transform
+ (
+ cmplx.cbegin(),
+ cmplx.cend(),
+ result.begin(),
+ [](const complex& c) { return c.cmptSum(); }
+ );
- return sf;
+ return result;
}
-Foam::scalarField Foam::Re(const UList<complex>& cf)
+Foam::scalarField Foam::Re(const UList<complex>& cmplx)
{
- scalarField sf(cf.size());
+ scalarField result(cmplx.size());
- forAll(sf, i)
- {
- sf[i] = cf[i].Re();
- }
+ std::transform
+ (
+ cmplx.cbegin(),
+ cmplx.cend(),
+ result.begin(),
+ [](const complex& c) { return c.real(); }
+ );
- return sf;
+ return result;
}
-Foam::scalarField Foam::Im(const UList<complex>& cf)
+Foam::scalarField Foam::Im(const UList<complex>& cmplx)
{
- scalarField sf(cf.size());
+ scalarField result(cmplx.size());
- forAll(sf, i)
- {
- sf[i] = cf[i].Im();
- }
+ std::transform
+ (
+ cmplx.cbegin(),
+ cmplx.cend(),
+ result.begin(),
+ [](const complex& c) { return c.imag(); }
+ );
- return sf;
+ return result;
}
diff --git a/src/OpenFOAM/fields/Fields/complex/complexField.H b/src/OpenFOAM/fields/Fields/complex/complexField.H
index 70a2d44e905..8de655acae9 100644
--- a/src/OpenFOAM/fields/Fields/complex/complexField.H
+++ b/src/OpenFOAM/fields/Fields/complex/complexField.H
@@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011 OpenFOAM Foundation
- Copyright (C) 2019 OpenCFD Ltd.
+ Copyright (C) 2019-2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -35,8 +35,8 @@ SourceFiles
\*---------------------------------------------------------------------------*/
-#ifndef complexField_H
-#define complexField_H
+#ifndef Foam_complexField_H
+#define Foam_complexField_H
#include "complex.H"
#include "scalarField.H"
@@ -53,46 +53,87 @@ namespace Foam
typedef Field<complex> complexField;
+
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
-//- Zip together complex field from components
+//- Zip together complex field from real/imag components
+void zip
+(
+ complexField& result,
+ const UList<scalar>& realValues,
+ const UList<scalar>& imagValues
+);
+
+//- Zip together complex field from real components and constant imag component
+void zip
+(
+ complexField& result,
+ const UList<scalar>& realValues,
+ const scalar imagValue
+);
+
+//- Zip together complex field from constant real component and imag components
void zip
(
complexField& result,
- const UList<scalar>& re,
- const UList<scalar>& im
+ const scalar realValue,
+ const UList<scalar>& imagValues
);
//- Unzip complex field into components
void unzip
(
const UList<complex>& input,
- scalarField& re,
- scalarField& im
+ scalarField& realValues,
+ scalarField& imagValues
);
-//- Zip up two lists of values into a list of complex
+//- Create complex field by zipping two lists of real/imag values
complexField ComplexField
(
- const UList<scalar>& re,
- const UList<scalar>& im
+ const UList<scalar>& realValues,
+ const UList<scalar>& imagValues
);
-//- Create complex field from a list of real (using imag == 0)
-complexField ReComplexField(const UList<scalar>& re);
+//- Create complex field by zipping a list of real components
+//- and a constant imag component
+complexField ComplexField
+(
+ const UList<scalar>& realValues,
+ const scalar imagValue
+);
+
+//- Create complex field by zipping a constant real component
+//- and a list of imag components
+complexField ComplexField
+(
+ const scalar realValue,
+ const UList<scalar>& imagValues
+);
-//- Create complex field from a list of imag (using real == 0)
-complexField ImComplexField(const UList<scalar>& im);
//- Extract real component
-scalarField Re(const UList<complex>& cf);
+scalarField Re(const UList<complex>& cmplx);
//- Extract imag component
-scalarField Im(const UList<complex>& cf);
+scalarField Im(const UList<complex>& cmplx);
//- Sum real and imag components
-scalarField ReImSum(const UList<complex>& cf);
+scalarField ReImSum(const UList<complex>& cmplx);
+
+
+//- Create complex field from a list of real (using imag == 0)
+inline complexField ReComplexField(const UList<scalar>& realValues)
+{
+ return ComplexField(realValues, scalar(0));
+}
+
+//- Create complex field from a list of imag (using real == 0)
+inline complexField ImComplexField(const UList<scalar>& imagValues)
+{
+ return ComplexField(scalar(0), imagValues);
+}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/OpenFOAM/fields/Fields/complex/complexFields.H b/src/OpenFOAM/fields/Fields/complex/complexFields.H
index 4a9ac3f25d4..0eabcaf2bce 100644
--- a/src/OpenFOAM/fields/Fields/complex/complexFields.H
+++ b/src/OpenFOAM/fields/Fields/complex/complexFields.H
@@ -31,8 +31,8 @@ Description
\*---------------------------------------------------------------------------*/
-#ifndef complexFields_H
-#define complexFields_H
+#ifndef Foam_complexFields_H
+#define Foam_complexFields_H
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/OpenFOAM/fields/Fields/complex/complexVectorField.C b/src/OpenFOAM/fields/Fields/complex/complexVectorField.C
index 99ce5ae5444..dfbc2a08746 100644
--- a/src/OpenFOAM/fields/Fields/complex/complexVectorField.C
+++ b/src/OpenFOAM/fields/Fields/complex/complexVectorField.C
@@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011 OpenFOAM Foundation
- Copyright (C) 2019 OpenCFD Ltd.
+ Copyright (C) 2019-2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -40,117 +40,209 @@ namespace Foam
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
-Foam::complexVectorField Foam::ComplexField
+void Foam::zip
(
- const UList<vector>& re,
- const UList<vector>& im
+ complexVectorField& result,
+ const UList<vector>& realValues,
+ const UList<vector>& imagValues
)
{
- complexVectorField cvf(re.size());
+ const label len = result.size();
- for (direction cmpt=0; cmpt<vector::nComponents; ++cmpt)
+ #ifdef FULLDEBUG
+ if (len != realValues.size() || len != imagValues.size())
{
- forAll(cvf, i)
- {
- cvf[i].component(cmpt).Re() = re[i].component(cmpt);
- cvf[i].component(cmpt).Im() = im[i].component(cmpt);
- }
+ FatalErrorInFunction
+ << "Components sizes do not match: " << len << " ("
+ << realValues.size() << ' ' << imagValues.size() << ')' << nl
+ << abort(FatalError);
}
+ #endif
- return cvf;
+ for (label i=0; i < len; ++i)
+ {
+ result[i] = Foam::zip(realValues[i], imagValues[i]);
+ }
}
-Foam::complexVectorField Foam::ReComplexField(const UList<vector>& re)
+void Foam::zip
+(
+ complexVectorField& result,
+ const UList<vector>& realValues,
+ const vector& imagValue
+)
{
- complexVectorField cvf(re.size());
+ const label len = result.size();
- for (direction cmpt=0; cmpt<vector::nComponents; ++cmpt)
+ #ifdef FULLDEBUG
+ if (len != realValues.size())
{
- forAll(cvf, i)
- {
- cvf[i].component(cmpt).Re() = re[i].component(cmpt);
- cvf[i].component(cmpt).Im() = 0.0;
- }
+ FatalErrorInFunction
+ << "Components sizes do not match: " << len << " != "
+ << realValues.size() << nl
+ << abort(FatalError);
}
+ #endif
- return cvf;
+ for (label i=0; i < len; ++i)
+ {
+ result[i] = Foam::zip(realValues[i], imagValue);
+ }
}
-Foam::complexVectorField Foam::ImComplexField(const UList<vector>& im)
+void Foam::zip
+(
+ complexVectorField& result,
+ const vector& realValue,
+ const UList<vector>& imagValues
+)
{
- complexVectorField cvf(im.size());
+ const label len = result.size();
- for (direction cmpt=0; cmpt<vector::nComponents; ++cmpt)
+ #ifdef FULLDEBUG
+ if (len != imagValues.size())
{
- forAll(cvf, i)
- {
- cvf[i].component(cmpt).Re() = 0.0;
- cvf[i].component(cmpt).Im() = im[i].component(cmpt);
- }
+ FatalErrorInFunction
+ << "Components sizes do not match: " << len << " != "
+ << imagValues.size() << nl
+ << abort(FatalError);
}
+ #endif
- return cvf;
+ for (label i=0; i < len; ++i)
+ {
+ result[i] = Foam::zip(realValue, imagValues[i]);
+ }
}
-Foam::vectorField Foam::ReImSum(const UList<complexVector>& cvf)
+Foam::complexVectorField Foam::ComplexField
+(
+ const UList<vector>& realValues,
+ const UList<vector>& imagValues
+)
{
- vectorField vf(cvf.size());
+ complexVectorField result(realValues.size());
- for (direction cmpt=0; cmpt<vector::nComponents; ++cmpt)
- {
- forAll(cvf, i)
- {
- vf[i].component(cmpt) =
- cvf[i].component(cmpt).Re() + cvf[i].component(cmpt).Im();
- }
- }
+ Foam::zip(result, realValues, imagValues);
- return vf;
+ return result;
}
-Foam::vectorField Foam::Re(const UList<complexVector>& cvf)
+Foam::complexVectorField Foam::ComplexField
+(
+ const UList<vector>& realValues,
+ const vector& imagValue
+)
{
- vectorField vf(cvf.size());
+ complexVectorField result(realValues.size());
- for (direction cmpt=0; cmpt<vector::nComponents; ++cmpt)
- {
- forAll(cvf, i)
+ Foam::zip(result, realValues, imagValue);
+
+ return result;
+}
+
+
+Foam::complexVectorField Foam::ComplexField
+(
+ const vector& realValue,
+ const UList<vector>& imagValues
+)
+{
+ complexVectorField result(imagValues.size());
+
+ Foam::zip(result, realValue, imagValues);
+
+ return result;
+}
+
+
+Foam::vectorField Foam::ReImSum(const UList<complexVector>& cmplx)
+{
+ vectorField result(cmplx.size());
+
+ std::transform
+ (
+ cmplx.cbegin(),
+ cmplx.cend(),
+ result.begin(),
+ [](const complexVector& c)
{
- vf[i].component(cmpt) = cvf[i].component(cmpt).Re();
+ return vector(c.x().cmptSum(), c.y().cmptSum(), c.z().cmptSum());
}
- }
+ );
- return vf;
+ return result;
}
-Foam::vectorField Foam::Im(const UList<complexVector>& cvf)
+Foam::vectorField Foam::Re(const UList<complexVector>& cmplx)
{
- vectorField vf(cvf.size());
+ vectorField result(cmplx.size());
+
+ std::transform
+ (
+ cmplx.cbegin(),
+ cmplx.cend(),
+ result.begin(),
+ [](const complexVector& c)
+ {
+ return vector(c.x().real(), c.y().real(), c.z().real());
+ }
+ );
- for (direction cmpt=0; cmpt<vector::nComponents; ++cmpt)
- {
- forAll(cvf, i)
+ return result;
+}
+
+
+Foam::vectorField Foam::Im(const UList<complexVector>& cmplx)
+{
+ vectorField result(cmplx.size());
+
+ std::transform
+ (
+ cmplx.cbegin(),
+ cmplx.cend(),
+ result.begin(),
+ [](const complexVector& c)
{
- vf[i].component(cmpt) = cvf[i].component(cmpt).Im();
+ return vector(c.x().imag(), c.y().imag(), c.z().imag());
}
- }
+ );
- return vf;
+ return result;
}
Foam::complexVectorField Foam::operator^
(
- const UList<vector>& vf,
- const UList<complexVector>& cvf
+ const UList<vector>& vec,
+ const UList<complexVector>& cmplx
)
{
- return ComplexField(vf^Re(cvf), vf^Im(cvf));
+ const label len = cmplx.size();
+
+ #ifdef FULLDEBUG
+ if (len != vec.size())
+ {
+ FatalErrorInFunction
+ << "Parameter sizes do not match: " << vec.size()
+ << " != " << len << nl
+ << abort(FatalError);
+ }
+ #endif
+
+ complexVectorField result(len);
+
+ for (label i=0; i < len; ++i)
+ {
+ result[i] = (vec[i] ^ cmplx[i]);
+ }
+
+ return result;
}
diff --git a/src/OpenFOAM/fields/Fields/complex/complexVectorField.H b/src/OpenFOAM/fields/Fields/complex/complexVectorField.H
index 15b8978581f..ab8e33215ef 100644
--- a/src/OpenFOAM/fields/Fields/complex/complexVectorField.H
+++ b/src/OpenFOAM/fields/Fields/complex/complexVectorField.H
@@ -6,7 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011 OpenFOAM Foundation
- Copyright (C) 2019 OpenCFD Ltd.
+ Copyright (C) 2019-2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -35,8 +35,8 @@ SourceFiles
\*---------------------------------------------------------------------------*/
-#ifndef complexVectorField_H
-#define complexVectorField_H
+#ifndef Foam_complexVectorField_H
+#define Foam_complexVectorField_H
#include "complex.H"
#include "complexVector.H"
@@ -51,33 +51,92 @@ namespace Foam
typedef Field<complexVector> complexVectorField;
-//- Zip up two lists of values into a list of complex
+
+// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
+
+//- Zip together complex vector field from real/imag components
+void zip
+(
+ complexVectorField& result,
+ const UList<vector>& realValues,
+ const UList<vector>& imagValues
+);
+
+
+//- Zip together complex vector field from real components
+//- and constant imag component
+void zip
+(
+ complexVectorField& result,
+ const UList<vector>& realValues,
+ const vector& imagValue
+);
+
+
+//- Zip together complex vector field from constant real component
+//- and imag components
+void zip
+(
+ complexVectorField& result,
+ const vector& realValue,
+ const UList<vector>& imagValues
+);
+
+
+//- Create complexVector field by zipping two lists of real/imag values
complexVectorField ComplexField
(
- const UList<vector>& re,
- const UList<vector>& im
+ const UList<vector>& realValues,
+ const UList<vector>& imagValues
);
+//- Create complexVector field by zipping a list of real components
+//- and a constant imag component
+complexVectorField ComplexField
+(
+ const UList<vector>& realValues,
+ const vector& imagValue
+);
-//- Create complex field from a list of real (using imag == 0)
-complexVectorField ReComplexField(const UList<vector>& re);
+//- Create complexVector field by zipping a constant real component
+//- and a list of imag components
+complexVectorField ComplexField
+(
+ const vector& realValue,
+ const UList<vector>& imagValues
+);
-//- Create complex field from a list of imag (using real == 0)
-complexVectorField ImComplexField(const UList<vector>& im);
//- Extract real component
-vectorField Re(const UList<complexVector>& cvf);
+vectorField Re(const UList<complexVector>& cmplx);
//- Extract imag component
-vectorField Im(const UList<complexVector>& cvf);
+vectorField Im(const UList<complexVector>& cmplx);
//- Sum real and imag components
-vectorField ReImSum(const UList<complexVector>& cvf);
+vectorField ReImSum(const UList<complexVector>& cmplx);
+
+
+//- Create complexVector field from a list of real (using imag == 0)
+inline complexVectorField ReComplexField(const UList<vector>& realValues)
+{
+ return ComplexField(realValues, vector::zero);
+}
+
+//- Create complexVector field from a list of imag (using real == 0)
+inline complexVectorField ImComplexField(const UList<vector>& imagValues)
+{
+ return ComplexField(vector::zero, imagValues);
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+//- Cross-product
complexVectorField operator^
(
- const UList<vector>& vf,
- const UList<complexVector>& cvf
+ const UList<vector>& vec,
+ const UList<complexVector>& cmplx
);
diff --git a/src/OpenFOAM/matrices/EigenMatrix/EigenMatrix.C b/src/OpenFOAM/matrices/EigenMatrix/EigenMatrix.C
index 30ddb10c36d..30f28d72106 100644
--- a/src/OpenFOAM/matrices/EigenMatrix/EigenMatrix.C
+++ b/src/OpenFOAM/matrices/EigenMatrix/EigenMatrix.C
@@ -873,8 +873,8 @@ void Foam::EigenMatrix<cmptType>::realSchur()
complex cDiv = complex(0.0, -H_[n - 1][n])
/complex(H_[n - 1][n-1]-p, q);
- H_[n - 1][n - 1] = cDiv.Re();
- H_[n - 1][n] = cDiv.Im();
+ H_[n - 1][n - 1] = cDiv.real();
+ H_[n - 1][n] = cDiv.imag();
}
H_[n][n - 1] = Zero;
@@ -907,8 +907,8 @@ void Foam::EigenMatrix<cmptType>::realSchur()
if (EValsIm_[i] == 0)
{
complex cDiv = complex(-ra, -sa)/complex(w, q);
- H_[i][n - 1] = cDiv.Re();
- H_[i][n] = cDiv.Im();
+ H_[i][n - 1] = cDiv.real();
+ H_[i][n] = cDiv.imag();
}
else
{
@@ -930,8 +930,8 @@ void Foam::EigenMatrix<cmptType>::realSchur()
complex(x*r - z*ra + q*sa, x*s - z*sa - q*ra)
/complex(vr, vi);
- H_[i][n - 1] = cDiv.Re();
- H_[i][n] = cDiv.Im();
+ H_[i][n - 1] = cDiv.real();
+ H_[i][n] = cDiv.imag();
if (mag(x) > (mag(z) + mag(q)))
{
@@ -946,8 +946,8 @@ void Foam::EigenMatrix<cmptType>::realSchur()
complex cDiv = complex(-r - y*H_[i][n - 1], -s
- y*H_[i][n])/complex(z, q);
- H_[i + 1][n - 1] = cDiv.Re();
- H_[i + 1][n] = cDiv.Im();
+ H_[i + 1][n - 1] = cDiv.real();
+ H_[i + 1][n] = cDiv.imag();
}
}
diff --git a/src/OpenFOAM/primitives/Vector/complex/complexVector.C b/src/OpenFOAM/primitives/Vector/complex/complexVector.C
index b1838bd1b5c..819b54b40c6 100644
--- a/src/OpenFOAM/primitives/Vector/complex/complexVector.C
+++ b/src/OpenFOAM/primitives/Vector/complex/complexVector.C
@@ -24,9 +24,6 @@ License
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
-Description
- Vector of complex numbers.
-
\*---------------------------------------------------------------------------*/
#include "complexVector.H"
diff --git a/src/OpenFOAM/primitives/Vector/complex/complexVector.H b/src/OpenFOAM/primitives/Vector/complex/complexVector.H
index fbba707ef45..7016b19109d 100644
--- a/src/OpenFOAM/primitives/Vector/complex/complexVector.H
+++ b/src/OpenFOAM/primitives/Vector/complex/complexVector.H
@@ -36,8 +36,8 @@ SourceFiles
\*---------------------------------------------------------------------------*/
-#ifndef complexVector_H
-#define complexVector_H
+#ifndef Foam_complexVector_H
+#define Foam_complexVector_H
#include "complex.H"
#include "vector.H"
diff --git a/src/OpenFOAM/primitives/Vector/complex/complexVectorI.H b/src/OpenFOAM/primitives/Vector/complex/complexVectorI.H
index fbd0a769262..bf01a86b674 100644
--- a/src/OpenFOAM/primitives/Vector/complex/complexVectorI.H
+++ b/src/OpenFOAM/primitives/Vector/complex/complexVectorI.H
@@ -6,6 +6,7 @@
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2011 OpenFOAM Foundation
+ Copyright (C) 2023 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
@@ -30,6 +31,20 @@ License
namespace Foam
{
+// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
+
+//- Zip together complexVector from real/imag vectors
+inline complexVector zip(const vector& realValue, const vector& imagValue)
+{
+ return complexVector
+ (
+ complex(realValue.x(), imagValue.x()),
+ complex(realValue.y(), imagValue.y()),
+ complex(realValue.z(), imagValue.z())
+ );
+}
+
+
// * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
inline complexVector operator*(const complex& v1, const complexVector& v2)
diff --git a/src/OpenFOAM/primitives/complex/complex.H b/src/OpenFOAM/primitives/complex/complex.H
index 6fa4bc5bf04..cc54994d23e 100644
--- a/src/OpenFOAM/primitives/complex/complex.H
+++ b/src/OpenFOAM/primitives/complex/complex.H
@@ -138,24 +138,6 @@ public:
void imag(scalar val) noexcept { im = val; }
- // Access
-
- //- Real part of complex number
- scalar Re() const noexcept { return re; }
-
- //- Imaginary part of complex number
- scalar Im() const noexcept { return im; }
-
-
- // Edit
-
- //- Real part of complex number
- scalar& Re() noexcept { return re; }
-
- //- Imaginary part of complex number
- scalar& Im() noexcept { return im; }
-
-
// Methods
//- Complex conjugate
@@ -168,6 +150,9 @@ public:
//- The L2-norm squared of complex
inline scalar magSqr() const;
+ //- The sum of real/imag components
+ inline scalar cmptSum() const noexcept;
+
// Member Operators
@@ -219,6 +204,21 @@ public:
friend complex operator/(const complex& c1, const complex& c2);
friend complex operator/(const complex& c, const scalar s);
friend complex operator/(const scalar s, const complex& c);
+
+
+ // Housekeeping
+
+ //- Get real part of complex number. Same as real()
+ scalar Re() const noexcept { return re; }
+
+ //- Get imaginary part of complex number. Same as imag()
+ scalar Im() const noexcept { return im; }
+
+ //- Non-const access to real part. Prefer real() setter method
+ scalar& Re() noexcept { return re; }
+
+ //- Non-const access to imaginary part. Prefer imag() setter method
+ scalar& Im() noexcept { return im; }
};
diff --git a/src/OpenFOAM/primitives/complex/complexI.H b/src/OpenFOAM/primitives/complex/complexI.H
index cfb57ae2ffc..ceda93cb89e 100644
--- a/src/OpenFOAM/primitives/complex/complexI.H
+++ b/src/OpenFOAM/primitives/complex/complexI.H
@@ -90,6 +90,12 @@ inline Foam::scalar Foam::complex::magSqr() const
}
+inline Foam::scalar Foam::complex::cmptSum() const noexcept
+{
+ return (re + im);
+}
+
+
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
inline void Foam::complex::operator=(const Foam::zero)
diff --git a/src/functionObjects/randomProcesses/energySpectrum/energySpectrum.C b/src/functionObjects/randomProcesses/energySpectrum/energySpectrum.C
index 8e1d24af700..0a8a159b287 100644
--- a/src/functionObjects/randomProcesses/energySpectrum/energySpectrum.C
+++ b/src/functionObjects/randomProcesses/energySpectrum/energySpectrum.C
@@ -73,7 +73,7 @@ void Foam::functionObjects::energySpectrum::calcAndWriteSpectrum
(
fft::forwardTransform
(
- ReComplexField(U),
+ ComplexField(U, vector::zero),
List<int>({N.x(), N.y(), N.z()})
)
/scalar(cmptProduct(N))
diff --git a/src/randomProcesses/fft/calcEk.C b/src/randomProcesses/fft/calcEk.C
index ebc234cfb26..5baf01cb475 100644
--- a/src/randomProcesses/fft/calcEk.C
+++ b/src/randomProcesses/fft/calcEk.C
@@ -41,7 +41,7 @@ graph calcEk
(
fft::forwardTransform
(
- ReComplexField(U.primitiveField()),
+ ComplexField(U.primitiveField(), vector::zero),
K.nn()
)*recRootN,
K
diff --git a/src/randomProcesses/fft/fft.C b/src/randomProcesses/fft/fft.C
index 8d483b60310..27c4dc5299b 100644
--- a/src/randomProcesses/fft/fft.C
+++ b/src/randomProcesses/fft/fft.C
@@ -47,7 +47,7 @@ void Foam::fft::fftRenumberRecurse
// We've worked out the renumbering scheme. Now copy
// the components across
- data[l1] = complex(renumData[l2].Re(), renumData[l2].Im());
+ data[l1] = renumData[l2];
}
else
{
@@ -150,12 +150,12 @@ Foam::fft::realTransform1D(const scalarField& field)
auto tresult = tmp<complexField>::New(nBy2 + 1);
auto& result = tresult.ref();
- result[0].Re() = out[0];
- result[nBy2].Re() = out[nBy2];
+ result[0].real(out[0]);
+ result[nBy2].real(out[nBy2]);
for (label i = 1; i < nBy2; ++i)
{
- result[i].Re() = out[i];
- result[i].Im() = out[n - i];
+ result[i].real(out[i]);
+ result[i].imag(out[n - i]);
}
fftw_destroy_plan(plan);
@@ -198,8 +198,8 @@ void Foam::fft::transform
forAll(field, i)
{
- inPtr[i][0] = field[i].Re();
- inPtr[i][1] = field[i].Im();
+ inPtr[i][0] = field[i].real();
+ inPtr[i][1] = field[i].imag();
}
// Create the plan
@@ -220,8 +220,8 @@ void Foam::fft::transform
forAll(field, i)
{
- field[i].Re() = outPtr[i][0];
- field[i].Im() = outPtr[i][1];
+ field[i].real(outPtr[i][0]);
+ field[i].imag(outPtr[i][1]);
}
fftw_destroy_plan(plan);
diff --git a/src/randomProcesses/turbulence/turbGen.C b/src/randomProcesses/turbulence/turbGen.C
index 23b119f87a1..dd50dea39a7 100644
--- a/src/randomProcesses/turbulence/turbGen.C
+++ b/src/randomProcesses/turbulence/turbGen.C
@@ -74,8 +74,11 @@ Foam::vectorField Foam::turbGen::U()
(
fft::reverseTransform
(
- ComplexField(cos(constant::mathematical::twoPi*rndPhases)*s,
- sin(constant::mathematical::twoPi*rndPhases)*s),
+ ComplexField
+ (
+ cos(constant::mathematical::twoPi*rndPhases)*s,
+ sin(constant::mathematical::twoPi*rndPhases)*s
+ ),
K.nn()
)*recRootN
);
--
GitLab