codedFvModel.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) 2012-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 Class
25  Foam::fv::codedFvModel
26 
27 Description
28  Constructs on-the-fly fvModel source from user-supplied code
29 
30 Usage
31  Example usage in constant/fvModels:
32  \verbatim
33  energySource
34  {
35  type coded;
36 
37  selectionMode all;
38 
39  field h;
40 
41  codeInclude
42  #{
43  #};
44 
45  codeAddSup
46  #{
47  Pout<< "**codeAddSup**" << endl;
48  const Time& time = mesh().time();
49  const scalarField& V = mesh().V();
50  scalarField& heSource = eqn.source();
51  heSource -= 0.1*sqr(time.value())*V;
52  #};
53 
54  codeAddRhoSup
55  #{
56  Pout<< "**codeAddRhoSup**" << endl;
57  #};
58 
59  codeAddAlphaRhoSup
60  #{
61  Pout<< "**codeAddAlphaRhoSup**" << endl;
62  #};
63  }
64  \endverbatim
65 
66 
67 SourceFiles
68  codedFvModel.C
69 
70 \*---------------------------------------------------------------------------*/
71 
72 #ifndef codedFvModel_H
73 #define codedFvModel_H
74 
75 #include "fvModel.H"
76 #include "codedBase.H"
77 
78 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
79 
80 namespace Foam
81 {
82 namespace fv
83 {
84 
85 /*---------------------------------------------------------------------------*\
86  Class codedFvModel Declaration
87 \*---------------------------------------------------------------------------*/
88 
89 class codedFvModel
90 :
91  public fvModel,
92  public codedBase
93 {
94  // Private data
95 
96  //- The name of the field that this fvModel applies to
97  word fieldName_;
98 
99  //- Underlying functionObject
100  mutable autoPtr<fvModel> redirectFvModelPtr_;
101 
102 
103  // Private Member Functions
104 
105  //- Non-virtual read
106  void readCoeffs();
107 
108  //- Return the name of the field's primitive type
109  word fieldPrimitiveTypeName() const;
110 
111  //- Adapt the context for the current object
112  virtual void prepare(dynamicCode&, const dynamicCodeContext&) const;
113 
114  //- Name of the dynamically generated CodedType
115  virtual const word& codeName() const;
116 
117  //- Return a description (type + name) for the output
118  virtual string description() const;
119 
120  //- Clear any redirected objects
121  virtual void clearRedirect() const;
122 
123  //- Get the dictionary to initialise the codeContext
124  virtual const dictionary& codeDict() const;
125 
126  //- Get the keywords associated with source code
127  virtual wordList codeKeys() const;
128 
129  //- Dynamically compiled fvModel
130  fvModel& redirectFvModel() const;
131 
132 
133  // Sources
134 
135  //- Add a source term to an equation
136  template<class Type>
137  void addSupType
138  (
139  fvMatrix<Type>& eqn,
140  const word& fieldName
141  ) const;
142 
143  //- Add a source term to a compressible equation
144  template<class Type>
145  void addSupType
146  (
147  const volScalarField& rho,
148  fvMatrix<Type>& eqn,
149  const word& fieldName
150  ) const;
151 
152  //- Add a source term to a phase equation
153  template<class Type>
154  void addSupType
155  (
156  const volScalarField& alpha,
157  const volScalarField& rho,
158  fvMatrix<Type>& eqn,
159  const word& fieldName
160  ) const;
161 
162 
163 public:
164 
165  //- Runtime type information
166  TypeName("coded");
167 
168 
169  // Constructors
170 
171  //- Construct from components
173  (
174  const word& name,
175  const word& modelType,
176  const dictionary& dict,
177  const fvMesh& mesh
178  );
179 
180 
181  // Member Functions
182 
183  // Checks
184 
185  //- Return the list of fields for which the fvModel adds source term
186  // to the transport equation
187  virtual wordList addSupFields() const;
188 
189 
190  // Sources
191 
192  //- Add a source term to an equation
194 
195  //- Add a source term to a compressible equation
197 
198  //- Add a source term to a phase equation
200 
201 
202  // Mesh changes
203 
204  //- Update for mesh motion
205  virtual bool movePoints();
206 
207  //- Update topology using the given map
208  virtual void topoChange(const polyTopoChangeMap&);
209 
210  //- Update from another mesh using the given map
211  virtual void mapMesh(const polyMeshMap&);
212 
213  //- Redistribute or update using the given distribution map
214  virtual void distribute(const polyDistributionMap&);
215 
216 
217  // IO
218 
219  //- Read source dictionary
220  virtual bool read(const dictionary& dict);
221 };
222 
223 
224 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
225 
226 } // End namespace fv
227 } // End namespace Foam
228 
229 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
230 
231 #endif
232 
233 // ************************************************************************* //
dictionary dict
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
const word & name() const
Return const access to the source name.
Definition: fvModelI.H:28
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define DEFINE_FV_MODEL_ADD_SUP(Type, nullArg)
Definition: fvModelM.H:26
virtual bool read(const dictionary &dict)
Read source dictionary.
Definition: codedFvModel.C:283
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
virtual bool movePoints()
Update for mesh motion.
Definition: codedFvModel.C:259
codedFvModel(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from components.
Definition: codedFvModel.C:228
#define DEFINE_FV_MODEL_ADD_ALPHA_RHO_SUP(Type, nullArg)
Definition: fvModelM.H:62
Finite volume model abstract base class.
Definition: fvModel.H:57
FOR_ALL_FIELD_TYPES(DEFINE_FV_MODEL_ADD_SUP)
Add a source term to an equation.
const fvMesh & mesh() const
Return const access to the mesh database.
Definition: fvModelI.H:34
A class for handling words, derived from string.
Definition: word.H:59
#define DEFINE_FV_MODEL_ADD_RHO_SUP(Type, nullArg)
Definition: fvModelM.H:43
labelList fv(nPoints)
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
Definition: codedFvModel.C:265
Base class for function objects and boundary conditions using dynamic code.
Definition: codedBase.H:53
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
virtual void distribute(const polyDistributionMap &)
Redistribute or update using the given distribution map.
Definition: codedFvModel.C:277
Constructs on-the-fly fvModel source from user-supplied code.
Definition: codedFvModel.H:88
Tools for handling dynamic code compilation.
Definition: dynamicCode.H:56
Encapsulation of dynamic code dictionaries.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
virtual wordList addSupFields() const
Return the list of fields for which the fvModel adds source term.
Definition: codedFvModel.C:244
TypeName("coded")
Runtime type information.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
Definition: codedFvModel.C:271
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:50
Namespace for OpenFOAM.