complexI.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2014 OpenFOAM Foundation
9 Copyright (C) 2019-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30
31inline constexpr Foam::complex::complex() noexcept
32:
33 re(0),
34 im(0)
35{}
36
37
39:
40 re(0),
41 im(0)
42{}
43
44
45inline constexpr Foam::complex::complex(const scalar r) noexcept
46:
47 re(r),
48 im(0)
49{}
50
51
52inline constexpr Foam::complex::complex(const scalar r, const scalar i) noexcept
53:
54 re(r),
55 im(i)
56{}
57
58
59inline Foam::complex::complex(const std::complex<float>& c)
60:
61 re(c.real()),
62 im(c.imag())
63{}
64
65
66inline Foam::complex::complex(const std::complex<double>& c)
67:
68 re(c.real()),
69 im(c.imag())
70{}
71
72
73// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
74
75inline void Foam::complex::real(scalar val)
76{
77 re = val;
78}
79
80
81inline void Foam::complex::imag(scalar val)
82{
83 im = val;
84}
85
86
87inline Foam::scalar Foam::complex::Re() const
88{
89 return re;
90}
91
92
93inline Foam::scalar Foam::complex::Im() const
94{
95 return im;
96}
97
98
99inline Foam::scalar& Foam::complex::Re()
100{
101 return re;
102}
103
104
105inline Foam::scalar& Foam::complex::Im()
106{
107 return im;
108}
109
110
112{
113 return complex(re, -im);
114}
115
116
117// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
118
120{
121 re = 0;
122 im = 0;
123}
124
125
126inline void Foam::complex::operator=(const scalar s)
127{
128 re = s;
129 im = 0;
130}
131
132
134{
135 re += c.re;
136 im += c.im;
137}
138
139
140inline void Foam::complex::operator+=(const scalar s)
141{
142 re += s;
143}
144
145
147{
148 re -= c.re;
149 im -= c.im;
150}
151
152
153inline void Foam::complex::operator-=(const scalar s)
154{
155 re -= s;
156}
157
158
160{
161 *this = (*this)*c;
162}
163
164
165inline void Foam::complex::operator*=(const scalar s)
166{
167 re *= s;
168 im *= s;
169}
170
171
173{
174 *this = *this/c;
175}
176
177
178inline void Foam::complex::operator/=(const scalar s)
179{
180 re /= s;
181 im /= s;
182}
183
184
185inline bool Foam::complex::operator==(const complex& c) const
186{
187 return (equal(re, c.re) && equal(im, c.im));
188}
189
190
191inline bool Foam::complex::operator!=(const complex& c) const
192{
193 return !operator==(c);
194}
195
196
197// * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
198
200{
201 return c.conjugate();
202}
203
204
205// * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * //
206
207namespace Foam
208{
209
210inline scalar magSqr(const complex& c)
211{
212 return (c.re*c.re + c.im*c.im);
213}
214
215
216inline scalar mag(const complex& c)
217{
218 return std::hypot(c.re, c.im);
219}
220
221
222inline complex sqr(const complex& c)
223{
224 return c*c;
225}
226
227
228inline complex sign(const complex& c)
229{
230 const scalar s(mag(c));
231 return s < ROOTVSMALL ? Zero : c/s;
232}
233
234
235inline scalar csign(const complex& c)
236{
237 return equal(c.Re(), 0) ? sign(c.Im()) : sign(c.Re());
238}
239
240
241inline const complex& min(const complex& c1, const complex& c2)
242{
243 if (magSqr(c1) < magSqr(c2))
244 {
245 return c1;
246 }
247
248 return c2;
249}
250
251
252inline const complex& max(const complex& c1, const complex& c2)
253{
254 if (magSqr(c1) < magSqr(c2))
255 {
256 return c2;
257 }
258
259 return c1;
260}
261
262
263inline complex limit(const complex& c1, const complex& c2)
264{
265 return complex(limit(c1.re, c2.re), limit(c1.im, c2.im));
266}
267
268
269inline const complex& sum(const complex& c)
270{
271 return c;
272}
273
274
275template<class Cmpt> class Tensor;
276
278{
279 return c;
280}
281
282
283// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
284
285inline complex operator-(const complex& c)
286{
287 return complex(-c.re, -c.im);
288}
289
290
291inline complex operator+(const complex& c1, const complex& c2)
292{
293 return complex
294 (
295 c1.re + c2.re,
296 c1.im + c2.im
297 );
298}
299
300
301inline complex operator+(const complex& c, const scalar s)
302{
303 return complex(c.re + s, c.im);
304}
305
306
307inline complex operator+(const scalar s, const complex& c)
308{
309 return complex(c.re + s, c.im);
310}
311
312
313inline complex operator-(const complex& c1, const complex& c2)
314{
315 return complex
316 (
317 c1.re - c2.re,
318 c1.im - c2.im
319 );
320}
321
322
323inline complex operator-(const complex& c, const scalar s)
324{
325 return complex(c.re - s, c.im);
326}
327
328
329inline complex operator-(const scalar s, const complex& c)
330{
331 return complex(s - c.re, -c.im);
332}
333
334
335inline complex operator*(const complex& c1, const complex& c2)
336{
337 return complex
338 (
339 c1.re*c2.re - c1.im*c2.im,
340 c1.im*c2.re + c1.re*c2.im
341 );
342}
343
344
345inline complex operator*(const complex& c, const scalar s)
346{
347 return complex(s*c.re, s*c.im);
348}
349
350
351inline complex operator*(const scalar s, const complex& c)
352{
353 return complex(s*c.re, s*c.im);
354}
355
356
357inline complex operator/(const complex& c1, const complex& c2)
358{
359 const scalar sqrC2 = magSqr(c2);
360
361 return complex
362 (
363 (c1.re*c2.re + c1.im*c2.im)/sqrC2,
364 (c1.im*c2.re - c1.re*c2.im)/sqrC2
365 );
366}
367
368
369inline complex operator/(const complex& c, const scalar s)
370{
371 return complex(c.re/s, c.im/s);
372}
373
374
375inline complex operator/(const scalar s, const complex& c)
376{
377 return complex(s)/c;
378}
379
380
381// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
382
383} // End namespace Foam
384
385
386// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
387
388// Complex transcendental functions
389namespace Foam
390{
391 #define transFunc(func) \
392 inline complex func(const complex& c) \
393 { \
394 return std:: func (static_cast<std::complex<scalar>>(c)); \
395 }
396
413
414// Special treatment for pow()
415inline complex pow(const complex& x, const complex& y)
416{
417 return std::pow
418 (
419 static_cast<std::complex<scalar>>(x),
420 static_cast<std::complex<scalar>>(y)
421 );
422}
423
424
425// Combinations of complex and integral/float
426#define powFuncs(type2) \
427 inline complex pow(const complex& x, const type2 y) \
428 { \
429 return std::pow(static_cast<std::complex<scalar>>(x), y); \
430 } \
431 \
432 inline complex pow(const type2 x, const complex& y) \
433 { \
434 return std::pow \
435 ( \
436 static_cast<std::complex<scalar>>(x), \
437 static_cast<std::complex<scalar>>(y) \
438 ); \
439 }
440
443#if defined(__APPLE__) && WM_LABEL_SIZE == 64
444powFuncs(int64_t)
445#endif
447powFuncs(double)
448
449
450inline complex pow3(const complex& c)
451{
452 return c*sqr(c);
453}
454
455
456inline complex pow4(const complex& c)
457{
458 return sqr(sqr(c));
459}
460
461
462inline complex pow5(const complex& c)
463{
464 return c*pow4(c);
465}
466
467
468inline complex pow6(const complex& c)
469{
470 return pow3(sqr(c));
471}
472
473
474inline complex pow025(const complex& c)
475{
476 return sqrt(sqrt(c));
477}
478
479} // End namespace Foam
480
481#undef transFunc
482#undef powFuncs
483
484// ************************************************************************* //
scalar y
A templated (3 x 3) tensor of objects of <T> derived from MatrixSpace.
Definition: Tensor.H:64
A complex number, similar to the C++ complex type.
Definition: complex.H:83
complex conjugate() const
Complex conjugate.
Definition: complexI.H:111
complex & operator=(const complex &)=default
Copy assignment.
void operator+=(const complex &c)
Definition: complexI.H:133
constexpr scalar imag() const
Imaginary part of complex number - STL naming.
Definition: complex.H:142
constexpr complex() noexcept
Default construct, as zero-initialized.
Definition: complexI.H:31
scalar Im() const
Imaginary part of complex number.
Definition: complexI.H:93
void operator-=(const complex &c)
Definition: complexI.H:146
scalar Re() const
Real part of complex number.
Definition: complexI.H:87
void operator*=(const complex &c)
Definition: complexI.H:159
void operator/=(const complex &c)
Definition: complexI.H:172
constexpr scalar real() const
Real part of complex number - STL naming.
Definition: complex.H:136
friend bool operator!=(const refineCell &rc1, const refineCell &rc2)
Definition: refineCell.H:106
friend bool operator==(const refineCell &rc1, const refineCell &rc2)
Definition: refineCell.H:97
A class representing the concept of 0 (zero) that can be used to avoid manipulating objects known to ...
Definition: zero.H:63
#define powFuncs(type2)
Definition: complexI.H:426
#define transFunc(func)
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Namespace for OpenFOAM.
tmp< faMatrix< Type > > operator-(const faMatrix< Type > &)
Unary negation.
bitSet operator~(const bitSet &bitset)
Bitset complement, returns a copy of the bitset with all its bits flipped.
Definition: bitSetI.H:762
dimensionedScalar pow6(const dimensionedScalar &ds)
dimensionedScalar pow5(const dimensionedScalar &ds)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dimensionedScalar asin(const dimensionedScalar &ds)
dimensionedScalar exp(const dimensionedScalar &ds)
dimensionedScalar tan(const dimensionedScalar &ds)
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
Definition: dimensionSet.C:536
dimensionedScalar sign(const dimensionedScalar &ds)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< faMatrix< Type > > operator+(const faMatrix< Type > &, const faMatrix< Type > &)
dimensionedScalar pow3(const dimensionedScalar &ds)
dimensionedScalar cosh(const dimensionedScalar &ds)
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
dimensionedScalar sin(const dimensionedScalar &ds)
dimensionedScalar tanh(const dimensionedScalar &ds)
dimensionedScalar sinh(const dimensionedScalar &ds)
dimensionedScalar log10(const dimensionedScalar &ds)
complex limit(const complex &, const complex &)
Definition: complexI.H:263
tmp< faMatrix< Type > > operator*(const areaScalarField::Internal &, const faMatrix< Type > &)
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar operator/(const scalar s1, const dimensionedScalar &ds2)
bool equal(const T &s1, const T &s2)
Compare two values for equality.
Definition: doubleFloat.H:46
dimensionedScalar acosh(const dimensionedScalar &ds)
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
dimensionedScalar sqrt(const dimensionedScalar &ds)
scalar csign(const complex &c)
Definition: complexI.H:235
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
dimensionedScalar pow4(const dimensionedScalar &ds)
dimensionedScalar atanh(const dimensionedScalar &ds)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
const direction noexcept
Definition: Scalar.H:223
dimensionedScalar atan(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar acos(const dimensionedScalar &ds)
dimensionedScalar pow025(const dimensionedScalar &ds)
dimensionedScalar asinh(const dimensionedScalar &ds)