ATCModel.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) 2007-2021 PCOpt/NTUA
9  Copyright (C) 2013-2021 FOSS GP
10  Copyright (C) 2019 OpenCFD Ltd.
11 -------------------------------------------------------------------------------
12 License
13  This file is part of OpenFOAM.
14 
15  OpenFOAM is free software: you can redistribute it and/or modify it
16  under the terms of the GNU General Public License as published by
17  the Free Software Foundation, either version 3 of the License, or
18  (at your option) any later version.
19 
20  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
21  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
22  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
23  for more details.
24 
25  You should have received a copy of the GNU General Public License
26  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
27 
28 
29 Class
30  Foam::ATCModel
31 
32 Description
33  Base class for selecting the adjoint transpose convection model.
34  Inherits from regIOobject to add lookup functionality
35 
36 
37 SourceFiles
38  ATCModel.C
39 
40 \*---------------------------------------------------------------------------*/
41 
42 #ifndef ATCModel_H
43 #define ATCModel_H
44 
45 #include "regIOobject.H"
46 #include "autoPtr.H"
47 #include "zeroATCcells.H"
48 #include "incompressibleVars.H"
50 #include "runTimeSelectionTables.H"
51 
52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53 
54 namespace Foam
55 {
56 
57 /*---------------------------------------------------------------------------*\
58  Class ATCModel Declaration
59 \*---------------------------------------------------------------------------*/
60 
61 class ATCModel
62 :
63  public regIOobject
64 {
65 private:
66 
67  // Private Member Functions
68 
69  //- No copy construct
70  ATCModel(const ATCModel&) = delete;
71 
72  //- No copy assignment
73  void operator=(const ATCModel&) = delete;
74 
75 
76 protected:
77 
78  // Protected data
79 
80  const fvMesh& mesh_;
83 
84 
85  const dictionary& dict_;
86  const scalar extraConvection_;
87  const scalar extraDiffusion_;
88  const label nSmooth_;
94 
95 
96  // Protected functions
97 
98  //- Compute limiter based on the cells given by zeroATCcells
99  void computeLimiter();
100 
101  //- Limit ATC field using ATClimiter_
102  void smoothATC();
103 
104 
105 public:
106 
107  //- Runtime type information
108  TypeName("ATCModel");
109 
110 
111  // Declare run-time constructor selection table
112 
114  (
115  autoPtr,
116  ATCModel,
117  dictionary,
118  (
119  const fvMesh& mesh,
120  const incompressibleVars& primalVars,
121  const incompressibleAdjointVars& adjointVars,
122  const dictionary& dict
123  ),
124  (mesh, primalVars, adjointVars, dict)
125  );
126 
127 
128  // Constructors
129 
130  //- Construct from components
131  ATCModel
132  (
133  const fvMesh& mesh,
134  const incompressibleVars& primalVars,
135  const incompressibleAdjointVars& adjointVars,
136  const dictionary& dict
137  );
138 
139 
140  // Selectors
141 
142  //- Return a reference to the selected turbulence model
143  static autoPtr<ATCModel> New
144  (
145  const fvMesh& mesh,
146  const incompressibleVars& primalVars,
147  const incompressibleAdjointVars& adjointVars,
148  const dictionary& dict
149  );
150 
151 
152  //- Destructor
153  virtual ~ATCModel() = default;
154 
155 
156  // Member Functions
157 
158  //- Add ATC to the given matrix
159  virtual void addATC(fvVectorMatrix& UaEqn) = 0;
160 
161  //- Get the list of cells on which to zero ATC
162  const labelList& getZeroATCcells();
163 
164  //- Get the extra convection multiplier
166 
167  //- Get the extra diffusion multiplier
169 
170  //- Get the list of cells on which to zero ATC
171  const volScalarField& getLimiter() const;
172 
173  //- Compute limiter based on the prescribed cells and number of
174  //- smoothing iterations
175  static void computeLimiter
176  (
178  const labelList& smoothCells,
179  const label nSmooth
180  );
181 
182  //- Return a limiter based on given cells.
183  //- For use with classes other than ATC which employ the same smoothing
184  //- mechanism
186  (
187  const fvMesh& mesh,
188  const dictionary& dict
189  );
190 
191  //- Smooth an arbitrary field on a given list of cells
192  template<class Type>
194  (
196  const labelList& cells
197  );
198 
199  //- Get the FI sensitivity derivatives term coming from the ATC
200  virtual tmp<volTensorField> getFISensitivityTerm() const = 0;
201 
202  //- Dummy writeData function required from regIOobject
203  virtual bool writeData(Ostream&) const;
204 
205  //- Update quantities related with the primal fields
206  virtual void updatePrimalBasedQuantities();
207 };
208 
209 
210 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
211 
212 } // End namespace Foam
213 
214 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
215 
216 #ifdef NoRepository
217 # include "ATCModelTemplates.C"
218 #endif
219 
220 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
221 
222 #endif
223 
224 // ************************************************************************* //
void computeLimiter()
Compute limiter based on the cells given by zeroATCcells.
Definition: ATCModel.C:38
virtual void updatePrimalBasedQuantities()
Update quantities related with the primal fields.
Definition: ATCModel.C:258
dictionary dict
volScalarField ATClimiter_
Definition: ATCModel.H:91
void smoothATC()
Limit ATC field using ATClimiter_.
Definition: ATCModel.C:44
A list of keyword definitions, which are a keyword followed by a number of values (eg...
Definition: dictionary.H:129
autoPtr< zeroATCcells > zeroATCcells_
Definition: ATCModel.H:90
word adjointSolverName_
Definition: ATCModel.H:89
tmp< areaScalarField > limiter(const areaScalarField &phi)
Definition: faNVDscheme.C:31
scalar getExtraDiffusionMultiplier()
Get the extra diffusion multiplier.
Definition: ATCModel.C:168
const incompressibleVars & primalVars_
Definition: ATCModel.H:80
TypeName("ATCModel")
Runtime type information.
Class including all adjoint fields for incompressible flows.
virtual void addATC(fvVectorMatrix &UaEqn)=0
Add ATC to the given matrix.
const incompressibleAdjointVars & adjointVars_
Definition: ATCModel.H:81
virtual ~ATCModel()=default
Destructor.
void smoothFieldBasedOnCells(GeometricField< Type, fvPatchField, volMesh > &vf, const labelList &cells)
Smooth an arbitrary field on a given list of cells.
Base class for solution control classes.
dynamicFvMesh & mesh
const cellShapeList & cells
const label nSmooth_
Definition: ATCModel.H:87
virtual bool writeData(Ostream &) const
Dummy writeData function required from regIOobject.
Definition: ATCModel.C:251
A class for handling words, derived from Foam::string.
Definition: word.H:63
const fvMesh & mesh_
Definition: ATCModel.H:79
Base class for selecting the adjoint transpose convection model. Inherits from regIOobject to add loo...
Definition: ATCModel.H:56
static tmp< volScalarField > createLimiter(const fvMesh &mesh, const dictionary &dict)
Return a limiter based on given cells. For use with classes other than ATC which employ the same smoo...
Definition: ATCModel.C:216
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:64
volVectorField ATC_
Definition: ATCModel.H:92
const scalar extraConvection_
Definition: ATCModel.H:85
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:56
scalar getExtraConvectionMultiplier()
Get the extra convection multiplier.
Definition: ATCModel.C:162
static autoPtr< ATCModel > New(const fvMesh &mesh, const incompressibleVars &primalVars, const incompressibleAdjointVars &adjointVars, const dictionary &dict)
Return a reference to the selected turbulence model.
Definition: ATCModel.C:123
const labelList & getZeroATCcells()
Get the list of cells on which to zero ATC.
Definition: ATCModel.C:156
const scalar extraDiffusion_
Definition: ATCModel.H:86
const bool reconstructGradients_
Definition: ATCModel.H:88
const volScalarField & getLimiter() const
Get the list of cells on which to zero ATC.
Definition: ATCModel.C:174
virtual tmp< volTensorField > getFISensitivityTerm() const =0
Get the FI sensitivity derivatives term coming from the ATC.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Definition: regIOobject.H:66
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
Definition: HashPtrTable.H:48
Macros to ease declaration of run-time selection tables.
A class for managing temporary objects.
Definition: HashPtrTable.H:50
declareRunTimeSelectionTable(autoPtr, ATCModel, dictionary,(const fvMesh &mesh, const incompressibleVars &primalVars, const incompressibleAdjointVars &adjointVars, const dictionary &dict),(mesh, primalVars, adjointVars, dict))
Namespace for OpenFOAM.
const dictionary & dict_
Definition: ATCModel.H:84