DiagTensorI.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-2022 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 #include "SphericalTensor.H"
27 #include "SymmTensor.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class Cmpt>
33 {}
34 
35 
36 template<class Cmpt>
38 :
39  VectorSpace<DiagTensor<Cmpt>, Cmpt, 3>(Zero)
40 {}
41 
42 
43 template<class Cmpt>
44 template<class Cmpt2>
46 (
47  const VectorSpace<DiagTensor<Cmpt2>, Cmpt2, 3>& vs
48 )
49 :
50  VectorSpace<DiagTensor<Cmpt>, Cmpt, 3>(vs)
51 {}
52 
53 
54 template<class Cmpt>
56 (
57  const Cmpt& vxx,
58  const Cmpt& vyy,
59  const Cmpt& vzz
60 )
61 {
62  this->v_[XX] = vxx;
63  this->v_[YY] = vyy;
64  this->v_[ZZ] = vzz;
65 }
66 
67 
68 template<class Cmpt>
70 :
71  VectorSpace<DiagTensor<Cmpt>, Cmpt, 3>(v)
72 {}
73 
74 
75 template<class Cmpt>
77 :
78  VectorSpace<DiagTensor<Cmpt>, Cmpt, 3>(is)
79 {}
80 
81 
82 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
83 
84 template<class Cmpt>
85 inline const Cmpt& Foam::DiagTensor<Cmpt>::xx() const
86 {
87  return this->v_[XX];
88 }
89 
90 template<class Cmpt>
91 inline const Cmpt& Foam::DiagTensor<Cmpt>::yy() const
92 {
93  return this->v_[YY];
94 }
95 
96 template<class Cmpt>
97 inline const Cmpt& Foam::DiagTensor<Cmpt>::zz() const
98 {
99  return this->v_[ZZ];
100 }
101 
102 
103 template<class Cmpt>
105 {
106  return this->v_[XX];
107 }
108 
109 template<class Cmpt>
111 {
112  return this->v_[YY];
113 }
114 
115 template<class Cmpt>
117 {
118  return this->v_[ZZ];
119 }
120 
121 template<class Cmpt>
123 {
124  return Vector<Cmpt>(xx(), yy(), zz());
125 }
126 
127 
128 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
129 
130 namespace Foam
131 {
132 
133 // * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
134 
135 template<class Cmpt>
136 inline Tensor<Cmpt>
138 {
139  return Tensor<Cmpt>
140  (
141  dt1.xx() + t2.xx(), t2.xy(), t2.xz(),
142  t2.yx(), dt1.yy() + t2.yy(), t2.yz(),
143  t2.zx(), t2.zy(), dt1.zz() + t2.zz()
144  );
145 }
146 
147 
148 template<class Cmpt>
149 inline Tensor<Cmpt>
151 {
152  return Tensor<Cmpt>
153  (
154  t1.xx() + dt2.xx(), t1.xy(), t1.xz(),
155  t1.yx(), t1.yy() + dt2.yy(), t1.yz(),
156  t1.zx(), t1.zy(), t1.zz() + dt2.zz()
157  );
158 }
159 
160 
161 template<class Cmpt>
162 inline Tensor<Cmpt>
164 {
165  return Tensor<Cmpt>
166  (
167  dt1.xx() - t2.xx(), -t2.xy(), -t2.xz(),
168  -t2.yx(), dt1.yy() - t2.yy(), -t2.yz(),
169  -t2.zx(), -t2.zy(), dt1.zz() - t2.zz()
170  );
171 }
172 
173 
174 template<class Cmpt>
175 inline Tensor<Cmpt>
177 {
178  return Tensor<Cmpt>
179  (
180  t1.xx() - dt2.xx(), t1.xy(), t1.xz(),
181  t1.yx(), t1.yy() - dt2.yy(), t1.yz(),
182  t1.zx(), t1.zy(), t1.zz() - dt2.zz()
183  );
184 }
185 
186 
187 template<class Cmpt>
188 inline SymmTensor<Cmpt>
190 {
191  return SymmTensor<Cmpt>
192  (
193  dt1.xx() + t2.xx(), t2.xy(), t2.xz(),
194  dt1.yy() + t2.yy(), t2.yz(),
195  dt1.zz() + t2.zz()
196  );
197 }
198 
199 
200 template<class Cmpt>
201 inline SymmTensor<Cmpt>
203 {
204  return SymmTensor<Cmpt>
205  (
206  t1.xx() + dt2.xx(), t1.xy(), t1.xz(),
207  t1.yy() + dt2.yy(), t1.yz(),
208  t1.zz() + dt2.zz()
209  );
210 }
211 
212 
213 template<class Cmpt>
214 inline SymmTensor<Cmpt>
216 {
217  return SymmTensor<Cmpt>
218  (
219  dt1.xx() - t2.xx(), -t2.xy(), -t2.xz(),
220  dt1.yy() - t2.yy(), -t2.yz(),
221  dt1.zz() - t2.zz()
222  );
223 }
224 
225 
226 template<class Cmpt>
227 inline SymmTensor<Cmpt>
229 {
230  return SymmTensor<Cmpt>
231  (
232  t1.xx() - dt2.xx(), t1.xy(), t1.xz(),
233  t1.yy() - dt2.yy(), t1.yz(),
234  t1.zz() - dt2.zz()
235  );
236 }
237 
238 
239 //- Inner-product between two diagonal tensors
240 template<class Cmpt>
241 inline DiagTensor<Cmpt>
243 {
244  return DiagTensor<Cmpt>
245  (
246  dt1.xx()*dt2.xx(),
247  dt1.yy()*dt2.yy(),
248  dt1.zz()*dt2.zz()
249  );
250 }
251 
252 
253 //- Inner-product between a diagonal tensor and a tensor
254 template<class Cmpt>
255 inline Tensor<Cmpt>
257 {
258  return Tensor<Cmpt>
259  (
260  dt1.xx()*t2.xx(),
261  dt1.xx()*t2.xy(),
262  dt1.xx()*t2.xz(),
263 
264  dt1.yy()*t2.yx(),
265  dt1.yy()*t2.yy(),
266  dt1.yy()*t2.yz(),
267 
268  dt1.zz()*t2.zx(),
269  dt1.zz()*t2.zy(),
270  dt1.zz()*t2.zz()
271  );
272 }
273 
274 
275 //- Inner-product between a tensor and a diagonal tensor
276 template<class Cmpt>
277 inline Tensor<Cmpt>
279 {
280  return Tensor<Cmpt>
281  (
282  t1.xx()*dt2.xx(),
283  t1.xy()*dt2.yy(),
284  t1.xz()*dt2.zz(),
285 
286  t1.yx()*dt2.xx(),
287  t1.yy()*dt2.yy(),
288  t1.yz()*dt2.zz(),
289 
290  t1.zx()*dt2.xx(),
291  t1.zy()*dt2.yy(),
292  t1.zz()*dt2.zz()
293  );
294 }
295 
296 
297 //- Inner-product between a diagonal tensor and a vector
298 template<class Cmpt>
299 inline Vector<Cmpt>
301 {
302  return Vector<Cmpt>
303  (
304  dt.xx()*v.x(),
305  dt.yy()*v.y(),
306  dt.zz()*v.z()
307  );
308 }
309 
310 
311 //- Inner-product between a vector and a diagonal tensor
312 template<class Cmpt>
313 inline Vector<Cmpt>
315 {
316  return Vector<Cmpt>
317  (
318  v.x()*dt.xx(),
319  v.y()*dt.yy(),
320  v.z()*dt.zz()
321  );
322 }
323 
324 
325 //- Division of a scalar by a diagonalTensor
326 template<class Cmpt>
327 inline DiagTensor<Cmpt>
328 operator/(const scalar s, const DiagTensor<Cmpt>& dt)
329 {
330  return DiagTensor<Cmpt>(s/dt.xx(), s/dt.yy(), s/dt.zz());
331 }
332 
333 
334 //- Division of a vector by a diagonalTensor
335 template<class Cmpt>
336 inline Vector<Cmpt>
338 {
339  return Vector<Cmpt>(v.x()/dt.xx(), v.y()/dt.yy(), v.z()/dt.zz());
340 }
341 
342 
343 //- Return the trace of a diagonal tensor
344 template<class Cmpt>
345 inline Cmpt tr(const DiagTensor<Cmpt>& dt)
346 {
347  return dt.xx() + dt.yy() + dt.zz();
348 }
349 
350 
351 //- Return the spherical part of a diagonal tensor
352 template<class Cmpt>
354 {
355  return 0.5*tr(dt);
356 }
357 
358 
359 //- Return the determinant of a diagonal tensor
360 template<class Cmpt>
361 inline Cmpt det(const DiagTensor<Cmpt>& t)
362 {
363  return t.xx()*t.yy()*t.zz();
364 }
365 
366 
367 //- Return the inverse of a diagonal tensor
368 template<class Cmpt>
370 {
371  return DiagTensor<Cmpt>(1.0/dt.xx(), 1.0/dt.yy(), 1.0/dt.zz());
372 }
373 
374 
375 //- Return the diagonal of a symmetric tensor as a diagonal tensor
376 template<class Cmpt>
378 {
379  return DiagTensor<Cmpt>(t.xx(), t.yy(), t.zz());
380 }
381 
382 
383 //- Return the diagonal of a tensor as a diagonal tensor
384 template<class Cmpt>
386 {
387  return DiagTensor<Cmpt>(t.xx(), t.yy(), t.zz());
388 }
389 
390 
391 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
392 
393 } // End namespace Foam
394 
395 // ************************************************************************* //
Templated 3D DiagTensor derived from VectorSpace.
Definition: DiagTensor.H:56
const Cmpt & xx() const
Definition: DiagTensorI.H:85
DiagTensor()
Construct null.
Definition: DiagTensorI.H:32
const Cmpt & zz() const
Definition: DiagTensorI.H:97
const Cmpt & yy() const
Definition: DiagTensorI.H:91
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
Templated 3D SphericalTensor derived from VectorSpace adding construction from 1 component,...
Templated 3D symmetric tensor derived from VectorSpace adding construction from 6 components,...
Definition: SymmTensor.H:56
const Cmpt & xx() const
Definition: SymmTensorI.H:87
const Cmpt & yz() const
Definition: SymmTensorI.H:111
const Cmpt & xz() const
Definition: SymmTensorI.H:99
const Cmpt & zz() const
Definition: SymmTensorI.H:117
const Cmpt & xy() const
Definition: SymmTensorI.H:93
const Cmpt & yy() const
Definition: SymmTensorI.H:105
Templated 3D tensor derived from MatrixSpace adding construction from 9 components,...
Definition: Tensor.H:67
const Cmpt & xx() const
Definition: TensorI.H:163
const Cmpt & yx() const
Definition: TensorI.H:184
const Cmpt & yz() const
Definition: TensorI.H:198
const Cmpt & xz() const
Definition: TensorI.H:177
const Cmpt & zz() const
Definition: TensorI.H:219
const Cmpt & xy() const
Definition: TensorI.H:170
const Cmpt & zx() const
Definition: TensorI.H:205
const Cmpt & zy() const
Definition: TensorI.H:212
const Cmpt & yy() const
Definition: TensorI.H:191
Templated vector space.
Definition: VectorSpace.H:85
Templated 3D Vector derived from VectorSpace adding construction from 3 components,...
Definition: Vector.H:60
const Cmpt & z() const
Definition: VectorI.H:87
const Cmpt & y() const
Definition: VectorI.H:81
const Cmpt & x() const
Definition: VectorI.H:75
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:50
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.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Namespace for OpenFOAM.
dimensionedScalar det(const dimensionedSphericalTensor &dt)
static const zero Zero
Definition: zero.H:97
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
SphericalTensor< Cmpt > sph(const DiagTensor< Cmpt > &dt)
Return the spherical part of a diagonal tensor.
Definition: DiagTensorI.H:353
void diag(pointPatchField< vector > &, const pointPatchField< tensor > &)
tmp< VolField< Type > > operator&(const fvMatrix< Type > &, const DimensionedField< Type, volMesh > &)
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
tmp< fvMatrix< Type > > operator+(const fvMatrix< Type > &, const fvMatrix< Type > &)
tmp< fvMatrix< Type > > operator-(const fvMatrix< Type > &)
tmp< fvMatrix< Type > > operator/(const fvMatrix< Type > &, const volScalarField::Internal &)