LagrangianEqn.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) 2025 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 Class
25  Foam::LagrangianEqn
26 
27 Description
28  This class stores the coefficients of a Lagrangian equation, and
29  facilitates solving that equation and updating the associated field. It is
30  designed to behave and be used similarly to fvMatrix.
31 
32 SourceFiles
33  LagrangianEqn.C
34 
35 \*---------------------------------------------------------------------------*/
36 
37 #ifndef LagrangianEqn_H
38 #define LagrangianEqn_H
39 
40 #include "LagrangianSp.H"
41 #include "LagrangianEqnBase.H"
42 
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 
45 namespace Foam
46 {
47 
48 /*---------------------------------------------------------------------------*\
49  Class LagrangianEqn Declaration
50 \*---------------------------------------------------------------------------*/
51 
52 template<class Type>
53 class LagrangianEqn
54 :
55  public tmp<LagrangianEqn<Type>>::refCount,
56  public LagrangianEqnBase
57 {
58  // Private Data
59 
60  //- The time-step field
62 
63  //- Reference to the sub-sub-field. One of this and psiSub will be
64  // valid and the other will be null.
65  const LagrangianSubSubField<Type>& psiSubSub_;
66 
67  //- Reference to the sub-field. One of this and psiSubSub will be
68  // valid and the other will be null.
69  const LagrangianSubField<Type>& psiSub_;
70 
71  //- Non-const pointer to the field. Might be null.
73 
74  //- Previous implicit time coefficient field into which to cache the
75  // implicit time coefficient on destruction. Might be null.
76  LagrangianDynamicField<scalar>* deltaTSp0Ptr_;
77 
78  //- Previous source field into which to cache the source on
79  // destruction. Might be null.
81 
82 
83  // Private Classes
84 
85  //- Helper class to create the psiSub/psiSubSub references
86  template<template<class> class PrimitiveField>
87  struct PsiRef;
88 
89 
90 public:
91 
92  // Public Data
93 
94  //- Explicit time-coefficient
96 
97  //- Implicit time-coefficient
99 
100  //- Explicit coefficient
102 
103  //- Implicit coefficient
105 
106 
107  // Constructors
108 
109  //- Construct for a const field and a tmp time-step with a name
110  template<template<class> class PrimitiveField>
112  (
113  const word& name,
114  const tmp<LagrangianSubScalarField>& tDeltaT,
116  LagrangianDynamicField<scalar>& deltaTSp0 =
120  );
121 
122  //- Construct for a const field and a time-step with a name
123  template<template<class> class PrimitiveField>
125  (
126  const word& name,
127  const LagrangianSubScalarField& deltaT,
129  LagrangianDynamicField<scalar>& deltaTSp0 =
133  );
134 
135  //- Construct for a const field with a name
136  template<template<class> class PrimitiveField>
138  (
139  const word& name,
141  LagrangianDynamicField<scalar>& deltaTSp0 =
145  );
146 
147  //- Construct for a const field and a tmp time-step. Name will be null.
148  template<template<class> class PrimitiveField>
150  (
151  const tmp<LagrangianSubScalarField>& tDeltaT,
153  LagrangianDynamicField<scalar>& deltaTSp0 =
157  );
158 
159  //- Construct for a const field and a time-step. Name will be null.
160  template<template<class> class PrimitiveField>
162  (
163  const LagrangianSubScalarField& deltaT,
165  LagrangianDynamicField<scalar>& deltaTSp0 =
169  );
170 
171  //- Construct for a const field. Name will be null.
172  template<template<class> class PrimitiveField>
174  (
176  LagrangianDynamicField<scalar>& deltaTSp0 =
180  );
181 
182  //- Construct for a non-const field and a tmp time-step. Will be named
183  // automatically as the field name plus the suffix "Eqn".
185  (
186  const tmp<LagrangianSubScalarField>& tDeltaT,
188  LagrangianDynamicField<scalar>& deltaTSp0 =
192  );
193 
194  //- Construct for a non-const field and a time-step. Will be named
195  // automatically as the field name plus the suffix "Eqn".
197  (
198  const LagrangianSubScalarField& deltaT,
200  LagrangianDynamicField<scalar>& deltaTSp0 =
204  );
205 
206  //- Construct for a non-const field. Will be named
207  // automatically as the field name plus the suffix "Eqn".
209  (
211  LagrangianDynamicField<scalar>& deltaTSp0 =
215  );
216 
217  //- Copy construct
219 
220  //- Move construct
222 
223  //- Construct from tmp
224  LagrangianEqn(const tmp<LagrangianEqn<Type>>& tEqn);
225 
226 
227  //- Destructor
228  ~LagrangianEqn();
229 
230 
231  // Member Functions
232 
233  //- Construct from another equation, with empty coefficients
235 
236  //- Return the field
238 
239  //- Return the field IO
240  const regIOobject& psiIo() const;
241 
242  //- Return the field name
243  const word& psiName() const;
244 
245  //- Return whether the given field is that of the equation
246  template<template<class> class PrimitiveField>
248 
249  //- Check the operation with another equation
250  void op(const LagrangianEqn<Type>& other);
251 
252  //- Set the previous field to that stored by another tmp equation
253  void setPrevious(const tmp<LagrangianEqn<Type>>& tOther);
254 
255  //- Return the combined time and non-time explicit coefficient
257 
258  //- Return the combined time and non-time implicit coefficient
259  tmp<LagrangianSp<Type>> allSp() const;
260 
261  //- Return the combined time and non-time explicit diagonal
262  // coefficient; i.e., with the off-diagonal parts of any tensor Sp
263  // coefficients converted into explicit values
265 
266  //- Return the combined time and non-time implicit diagonal
267  // coefficient; i.e., with the off-diagonal parts of any tensor Sp
268  // coefficients removed
270 
271  //- Solve
272  void solve(const bool final);
273 
274 
275  // Member Operators
276 
277  //- Addition assignment
278  void operator+=(const LagrangianEqn<Type>& other);
279 
280  //- Addition assignment
281  void operator+=(const tmp<LagrangianEqn<Type>>& tOther);
282 
283  //- Subtraction assignment
284  void operator-=(const LagrangianEqn<Type>& other);
285 
286  //- Subtraction assignment
287  void operator-=(const tmp<LagrangianEqn<Type>>& tOther);
288 };
289 
290 
291 //- Negation
292 template<class Type>
294 
295 //- Negation
296 template<class Type>
298 
299 #define LAGRANGIAN_EQN_EQN_OPERATOR(Op) \
300  \
301  template<class Type> \
302  tmp<LagrangianEqn<Type>> operator Op \
303  ( \
304  const LagrangianEqn<Type>& a, \
305  const LagrangianEqn<Type>& b \
306  ); \
307  \
308  template<class Type> \
309  tmp<LagrangianEqn<Type>> operator Op \
310  ( \
311  const tmp<LagrangianEqn<Type>>& tA, \
312  const LagrangianEqn<Type>& b \
313  ); \
314  \
315  template<class Type> \
316  tmp<LagrangianEqn<Type>> operator Op \
317  ( \
318  const tmp<LagrangianEqn<Type>>& tA, \
319  const tmp<LagrangianEqn<Type>>& tB \
320  ); \
321  \
322  template<class Type> \
323  tmp<LagrangianEqn<Type>> operator Op \
324  ( \
325  const LagrangianEqn<Type>& a, \
326  const tmp<LagrangianEqn<Type>>& tB \
327  );
328 
329 //- Addition
331 
332 //- Subtraction
334 
335 //- Set-equal-to
337 
338 #undef LAGRANGIAN_EQN_EQN_OPERATOR
339 
340 #define LAGRANGIAN_EQN_FIELD_OPERATOR(Op, LagrangianSubField) \
341  \
342  template<class Type> \
343  tmp<LagrangianEqn<Type>> operator Op \
344  ( \
345  const LagrangianEqn<Type>& eqn, \
346  const LagrangianSubField<Type>& field \
347  ); \
348  \
349  template<class Type> \
350  tmp<LagrangianEqn<Type>> operator Op \
351  ( \
352  const tmp<LagrangianEqn<Type>>& tEqn, \
353  const LagrangianSubField<Type>& field \
354  ); \
355  \
356  template<class Type> \
357  tmp<LagrangianEqn<Type>> operator Op \
358  ( \
359  const tmp<LagrangianEqn<Type>>& tEqn, \
360  const tmp<LagrangianSubField<Type>>& tField \
361  ); \
362  \
363  template<class Type> \
364  tmp<LagrangianEqn<Type>> operator Op \
365  ( \
366  const LagrangianEqn<Type>& eqn, \
367  const tmp<LagrangianSubField<Type>>& tField \
368  );
369 
370 #define LAGRANGIAN_FIELD_EQN_OPERATOR(Op, LagrangianSubField) \
371  \
372  template<class Type> \
373  tmp<LagrangianEqn<Type>> operator Op \
374  ( \
375  const LagrangianSubField<Type>& field, \
376  const LagrangianEqn<Type>& eqn \
377  ); \
378  \
379  template<class Type> \
380  tmp<LagrangianEqn<Type>> operator Op \
381  ( \
382  const LagrangianSubField<Type>& field, \
383  const tmp<LagrangianEqn<Type>>& tEqn \
384  ); \
385  \
386  template<class Type> \
387  tmp<LagrangianEqn<Type>> operator Op \
388  ( \
389  const tmp<LagrangianSubField<Type>>& tField, \
390  const tmp<LagrangianEqn<Type>>& tEqn \
391  ); \
392  \
393  template<class Type> \
394  tmp<LagrangianEqn<Type>> operator Op \
395  ( \
396  const tmp<LagrangianSubField<Type>>& tField, \
397  const LagrangianEqn<Type>& eqn \
398  );
399 
400 //- Addition
405 
406 //- Subtraction
411 
412 //- Set-equal-to
417 
418 #undef LAGRANGIAN_EQN_FIELD_OPERATOR
419 #undef LAGRANGIAN_FIELD_EQN_OPERATOR
420 
421 #define LAGRANGIAN_EQN_SCALAR_FIELD_OPERATOR(Op, LagrangianSubField) \
422  \
423  template<class Type> \
424  tmp<LagrangianEqn<Type>> operator Op \
425  ( \
426  const LagrangianEqn<Type>& eqn, \
427  const LagrangianSubField<scalar>& field \
428  ); \
429  \
430  template<class Type> \
431  tmp<LagrangianEqn<Type>> operator Op \
432  ( \
433  const tmp<LagrangianEqn<Type>>& tEqn, \
434  const LagrangianSubField<scalar>& field \
435  ); \
436  \
437  template<class Type> \
438  tmp<LagrangianEqn<Type>> operator Op \
439  ( \
440  const tmp<LagrangianEqn<Type>>& tEqn, \
441  const tmp<LagrangianSubField<scalar>>& tField \
442  ); \
443  \
444  template<class Type> \
445  tmp<LagrangianEqn<Type>> operator Op \
446  ( \
447  const LagrangianEqn<Type>& eqn, \
448  const tmp<LagrangianSubField<scalar>>& tField \
449  );
450 
451 #define LAGRANGIAN_SCALAR_FIELD_EQN_OPERATOR(Op, LagrangianSubField) \
452  \
453  template<class Type> \
454  tmp<LagrangianEqn<Type>> operator Op \
455  ( \
456  const LagrangianSubField<scalar>& field, \
457  const LagrangianEqn<Type>& eqn \
458  ); \
459  \
460  template<class Type> \
461  tmp<LagrangianEqn<Type>> operator Op \
462  ( \
463  const LagrangianSubField<scalar>& field, \
464  const tmp<LagrangianEqn<Type>>& tEqn \
465  ); \
466  \
467  template<class Type> \
468  tmp<LagrangianEqn<Type>> operator Op \
469  ( \
470  const tmp<LagrangianSubField<scalar>>& tField, \
471  const tmp<LagrangianEqn<Type>>& tEqn \
472  ); \
473  \
474  template<class Type> \
475  tmp<LagrangianEqn<Type>> operator Op \
476  ( \
477  const tmp<LagrangianSubField<scalar>>& tField, \
478  const LagrangianEqn<Type>& eqn \
479  );
480 
481 //- Multiplication
486 
487 //- Division
490 
491 #undef LAGRANGIAN_EQN_SCALAR_FIELD_OPERATOR
492 #undef LAGRANGIAN_SCALAR_FIELD_EQN_OPERATOR
493 
494 /*---------------------------------------------------------------------------*\
495  Class LagrangianEqn::PsiRef Declaration
496 \*---------------------------------------------------------------------------*/
497 
498 template<class Type>
499 template<template<class> class PrimitiveField>
500 struct Foam::LagrangianEqn<Type>::PsiRef
501 {
502  //- Primitive field types match. Return a valid reference.
503  inline static const LagrangianSubField<Type, PrimitiveField>& ref
504  (
506  )
507  {
508  return psi;
509  };
510 
511  //- Primitive field types do not match. Return a null reference.
512  template<template<class> class OtherPrimitiveField>
513  inline static const LagrangianSubField<Type, PrimitiveField>& ref
514  (
516  )
517  {
518  return NullObjectRef<LagrangianSubField<Type, PrimitiveField>>();
519  };
520 };
521 
522 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
523 
524 } // End namespace Foam
525 
526 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
527 
528 #ifdef NoRepository
529  #include "LagrangianEqn.C"
530 #endif
531 
532 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
533 
534 #endif
535 
536 // ************************************************************************* //
#define LAGRANGIAN_SCALAR_FIELD_EQN_OPERATOR(Op, LagrangianSubField)
#define LAGRANGIAN_EQN_SCALAR_FIELD_OPERATOR(Op, LagrangianSubField)
#define LAGRANGIAN_FIELD_EQN_OPERATOR(Op, LagrangianSubField)
#define LAGRANGIAN_EQN_EQN_OPERATOR(Op)
#define LAGRANGIAN_EQN_FIELD_OPERATOR(Op, LagrangianSubField)
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricField class.
const word & name() const
Return the equation name.
This class stores the coefficients of a Lagrangian equation, and facilitates solving that equation an...
Definition: LagrangianEqn.H:56
LagrangianCoeff< scalar, true > deltaTSp
Implicit time-coefficient.
Definition: LagrangianEqn.H:97
LagrangianSp< Type > Sp
Implicit coefficient.
void op(const LagrangianEqn< Type > &other)
Check the operation with another equation.
LagrangianCoeff< Type, false > Su
Explicit coefficient.
void solve(const bool final)
Solve.
tmp< LagrangianSp< Type > > allSp() const
Return the combined time and non-time implicit coefficient.
tmp< LagrangianCoeff< scalar, true > > allDiagonalSp() const
Return the combined time and non-time implicit diagonal.
const regIOobject & psiIo() const
Return the field IO.
~LagrangianEqn()
Destructor.
void operator+=(const LagrangianEqn< Type > &other)
Addition assignment.
tmp< LagrangianCoeff< Type, false > > allDiagonalSu() const
Return the combined time and non-time explicit diagonal.
tmp< LagrangianSubSubField< Type > > psi() const
Return the field.
static tmp< LagrangianEqn< Type > > NewEmpty(const LagrangianEqn< Type > &)
Construct from another equation, with empty coefficients.
LagrangianEqn(const word &name, const tmp< LagrangianSubScalarField > &tDeltaT, const LagrangianSubField< Type, PrimitiveField > &psi, LagrangianDynamicField< scalar > &deltaTSp0=NullObjectNonConstRef< LagrangianDynamicField< scalar >>(), LagrangianDynamicField< Type > &S0=NullObjectNonConstRef< LagrangianDynamicField< Type >>())
Construct for a const field and a tmp time-step with a name.
Definition: LagrangianEqn.C:35
void setPrevious(const tmp< LagrangianEqn< Type >> &tOther)
Set the previous field to that stored by another tmp equation.
LagrangianCoeff< Type, false > deltaTSu
Explicit time-coefficient.
Definition: LagrangianEqn.H:86
const word & psiName() const
Return the field name.
void operator-=(const LagrangianEqn< Type > &other)
Subtraction assignment.
bool isPsi(const LagrangianSubField< Type, PrimitiveField > &psi) const
Return whether the given field is that of the equation.
tmp< LagrangianCoeff< Type, false > > allSu() const
Return the combined time and non-time explicit coefficient.
Wrapper around LagrangianCoeff to specialise for the implicit coefficient. Trivial at present....
Definition: LagrangianSp.H:71
Reference counter for various OpenFOAM components.
Definition: refCount.H:50
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Definition: regIOobject.H:55
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
trAU ref().rename("rAU")
Namespace for OpenFOAM.
T & NullObjectNonConstRef()
Return non-const reference to the nullObject of type T.
Definition: nullObjectI.H:33
tmp< fvMatrix< Type > > operator-(const fvMatrix< Type > &)