complexI.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
27 
28 namespace Foam
29 {
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
34 {}
35 
36 
37 inline complex::complex(const scalar Re, const scalar Im)
38 :
39  re(Re),
40  im(Im)
41 {}
42 
43 
44 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
45 
46 inline scalar complex::Re() const
47 {
48  return re;
49 }
50 
51 
52 inline scalar complex::Im() const
53 {
54  return im;
55 }
56 
57 
58 inline scalar& complex::Re()
59 {
60  return re;
61 }
62 
63 
64 inline scalar& complex::Im()
65 {
66  return im;
67 }
68 
69 
71 {
72  return complex(re, -im);
73 }
74 
75 
76 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
77 
78 inline void complex::operator+=(const complex& c)
79 {
80  re += c.re;
81  im += c.im;
82 }
83 
84 
85 inline void complex::operator-=(const complex& c)
86 {
87  re -= c.re;
88  im -= c.im;
89 }
90 
91 
92 inline void complex::operator*=(const complex& c)
93 {
94  *this = (*this)*c;
95 }
96 
97 
98 inline void complex::operator/=(const complex& c)
99 {
100  *this = *this/c;
101 }
102 
103 
104 inline void complex::operator=(const scalar s)
105 {
106  re = s;
107  im = 0.0;
108 }
109 
110 
111 inline void complex::operator+=(const scalar s)
112 {
113  re += s;
114 }
115 
116 
117 inline void complex::operator-=(const scalar s)
118 {
119  re -= s;
120 }
121 
122 
123 inline void complex::operator*=(const scalar s)
124 {
125  re *= s;
126  im *= s;
127 }
128 
129 
130 inline void complex::operator/=(const scalar s)
131 {
132  re /= s;
133  im /= s;
134 }
135 
136 
138 {
139  return conjugate();
140 }
141 
142 
143 inline bool complex::operator==(const complex& c) const
144 {
145  return (equal(re, c.re) && equal(im, c.im));
146 }
147 
148 
149 inline bool complex::operator!=(const complex& c) const
150 {
151  return !operator==(c);
152 }
153 
154 
155 // * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * //
156 
157 
158 inline scalar magSqr(const complex& c)
159 {
160  return (c.re*c.re + c.im*c.im);
161 }
162 
163 
164 inline complex sqr(const complex& c)
165 {
166  return c * c;
167 }
168 
169 
170 inline scalar mag(const complex& c)
171 {
172  return sqrt(magSqr(c));
173 }
174 
175 
176 inline const complex& max(const complex& c1, const complex& c2)
177 {
178  if (mag(c1) > mag(c2))
179  {
180  return c1;
181  }
182  else
183  {
184  return c2;
185  }
186 }
187 
188 
189 inline const complex& min(const complex& c1, const complex& c2)
190 {
191  if (mag(c1) < mag(c2))
192  {
193  return c1;
194  }
195  else
196  {
197  return c2;
198  }
199 }
200 
201 
202 inline complex limit(const complex& c1, const complex& c2)
203 {
204  return complex(limit(c1.re, c2.re), limit(c1.im, c2.im));
205 }
206 
207 
208 inline const complex& sum(const complex& c)
209 {
210  return c;
211 }
212 
213 
214 template<class Cmpt>
215 class Tensor;
216 
217 inline complex transform(const Tensor<scalar>&, const complex c)
218 {
219  return c;
220 }
221 
222 
223 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
224 
225 inline complex operator+(const complex& c1, const complex& c2)
226 {
227  return complex
228  (
229  c1.re + c2.re,
230  c1.im + c2.im
231  );
232 }
233 
234 
235 inline complex operator-(const complex& c)
236 {
237  return complex
238  (
239  -c.re,
240  -c.im
241  );
242 }
243 
244 
245 inline complex operator-(const complex& c1, const complex& c2)
246 {
247  return complex
248  (
249  c1.re - c2.re,
250  c1.im - c2.im
251  );
252 }
253 
254 
255 inline complex operator*(const complex& c1, const complex& c2)
256 {
257  return complex
258  (
259  c1.re*c2.re - c1.im*c2.im,
260  c1.im*c2.re + c1.re*c2.im
261  );
262 }
263 
264 
265 inline complex operator/(const complex& c1, const complex& c2)
266 {
267  scalar sqrC2 = magSqr(c2);
268 
269  return complex
270  (
271  (c1.re*c2.re + c1.im*c2.im)/sqrC2,
272  (c1.im*c2.re - c1.re*c2.im)/sqrC2
273  );
274 }
275 
276 
277 inline complex operator*(const scalar s, const complex& c)
278 {
279  return complex(s*c.re, s*c.im);
280 }
281 
282 
283 inline complex operator*(const complex& c, const scalar s)
284 {
285  return complex(s*c.re, s*c.im);
286 }
287 
288 
289 inline complex operator/(const complex& c, const scalar s)
290 {
291  return complex(c.re/s, c.im/s);
292 }
293 
294 
295 inline complex operator/(const scalar s, const complex& c)
296 {
297  return complex(s/c.re, s/c.im);
298 }
299 
300 
301 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
302 
303 } // End namespace Foam
304 
305 // ************************************************************************* //
complex()
Construct null.
Definition: complexI.H:33
bool operator==(const complex &) const
Definition: complexI.H:143
friend complex operator+(const complex &, const complex &)
Definition: complexI.H:225
friend const complex & sum(const complex &)
Definition: complexI.H:208
const dimensionedScalar & c1
First radiation constant: default SI units: [W/m^2].
scalar Im() const
Definition: complexI.H:52
void operator-=(const complex &)
Definition: complexI.H:85
dimensionedScalar sqrt(const dimensionedScalar &ds)
friend complex sqr(const complex &c)
Definition: complexI.H:164
friend scalar mag(const complex &c)
Definition: complexI.H:170
friend complex operator-(const complex &)
Definition: complexI.H:235
const dimensionedScalar & c
Speed of light in a vacuum.
friend const complex & min(const complex &, const complex &)
Definition: complexI.H:189
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().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))
void operator=(const scalar)
Definition: complexI.H:104
complex conjugate() const
Definition: complexI.H:70
complex operator!() const
Definition: complexI.H:137
bool operator!=(const complex &) const
Definition: complexI.H:149
void operator/=(const complex &)
Definition: complexI.H:98
friend const complex & max(const complex &, const complex &)
Definition: complexI.H:176
friend scalar magSqr(const complex &c)
Definition: complexI.H:158
friend complex limit(const complex &, const complex &)
Definition: complexI.H:202
friend complex operator/(const complex &, const complex &)
Definition: complexI.H:265
bool equal(const T &s1, const T &s2)
Definition: doubleFloat.H:62
friend complex operator*(const complex &, const complex &)
Definition: complexI.H:255
const dimensionedScalar & c2
Second radiation constant: default SI units: [m K].
Extension to the c++ complex library type.
Definition: complex.H:76
Templated 3D tensor derived from MatrixSpace adding construction from 9 components, element access using xx(), xy() etc. member functions and the inner-product (dot-product) and outer-product of two Vectors (tensor-product) operators.
Definition: complexI.H:215
scalar Re() const
Definition: complexI.H:46
void operator+=(const complex &)
Definition: complexI.H:78
Namespace for OpenFOAM.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:477
void operator*=(const complex &)
Definition: complexI.H:92