All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
actuationDiskSource.C
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-2021 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 "actuationDiskSource.H"
27 #include "fvMesh.H"
28 #include "fvMatrix.H"
29 #include "geometricOneField.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace fv
37 {
38  defineTypeNameAndDebug(actuationDiskSource, 0);
40  (
41  fvModel,
42  actuationDiskSource,
43  dictionary
44  );
45 }
46 }
47 
48 
49 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
50 
51 void Foam::fv::actuationDiskSource::readCoeffs()
52 {
53  UName_ = coeffs().lookupOrDefault<word>("U", "U");
54 
55  diskDir_ = coeffs().lookup<vector>("diskDir");
56  if (mag(diskDir_) < vSmall)
57  {
59  << "disk direction vector is approximately zero"
60  << exit(FatalIOError);
61  }
62 
63  Cp_ = coeffs().lookup<scalar>("Cp");
64  Ct_ = coeffs().lookup<scalar>("Ct");
65  if (Cp_ <= vSmall || Ct_ <= vSmall)
66  {
68  << "Cp and Ct must be greater than zero"
69  << exit(FatalIOError);
70  }
71 
72  diskArea_ = coeffs().lookup<scalar>("diskArea");
73  if (magSqr(diskArea_) <= vSmall)
74  {
76  << "diskArea is approximately zero"
77  << exit(FatalIOError);
78  }
79 
80  upstreamPoint_ = coeffs().lookup<point>("upstreamPoint");
81  upstreamCellId_ = mesh().findCell(upstreamPoint_);
82  if (returnReduce(upstreamCellId_, maxOp<label>()) == -1)
83  {
85  << "upstream location " << upstreamPoint_ << " not found in mesh"
86  << exit(FatalIOError);
87  }
88 }
89 
90 
91 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
92 
94 (
95  const word& name,
96  const word& modelType,
97  const dictionary& dict,
98  const fvMesh& mesh
99 )
100 :
101  fvModel(name, modelType, dict, mesh),
102  set_(coeffs(), mesh),
103  UName_(word::null),
104  diskDir_(vector::uniform(NaN)),
105  Cp_(NaN),
106  Ct_(NaN),
107  diskArea_(NaN),
108  upstreamPoint_(vector::uniform(NaN)),
109  upstreamCellId_(-1)
110 {
111  readCoeffs();
112 }
113 
114 
115 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
116 
118 {
119  return wordList(1, UName_);
120 }
121 
122 
124 (
125  fvMatrix<vector>& eqn,
126  const word& fieldName
127 ) const
128 {
129  const scalarField& cellsV = mesh().V();
130  vectorField& Usource = eqn.source();
131  const vectorField& U = eqn.psi();
132 
133  if (set_.V() > vSmall)
134  {
135  addActuationDiskAxialInertialResistance
136  (
137  Usource,
138  set_.cells(),
139  cellsV,
141  U
142  );
143  }
144 }
145 
146 
148 (
149  const volScalarField& rho,
150  fvMatrix<vector>& eqn,
151  const word& fieldName
152 ) const
153 {
154  const scalarField& cellsV = mesh().V();
155  vectorField& Usource = eqn.source();
156  const vectorField& U = eqn.psi();
157 
158  if (set_.V() > vSmall)
159  {
160  addActuationDiskAxialInertialResistance
161  (
162  Usource,
163  set_.cells(),
164  cellsV,
165  rho,
166  U
167  );
168  }
169 }
170 
171 
173 {
174  set_.updateMesh(mpm);
175 }
176 
177 
179 {
180  if (fvModel::read(dict))
181  {
182  set_.read(coeffs());
183  readCoeffs();
184  return true;
185  }
186  else
187  {
188  return false;
189  }
190 }
191 
192 
193 // ************************************************************************* //
defineTypeNameAndDebug(fixedTemperatureConstraint, 0)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:282
virtual bool read(const dictionary &dict)
Read source dictionary.
Definition: fvModel.C:158
actuationDiskSource(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from components.
Finite volume model abstract base class.
Definition: fvModel.H:55
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
Macros for easy insertion into run-time selection tables.
virtual void updateMesh(const mapPolyMesh &)
Update for mesh changes.
virtual wordList addSupFields() const
Return the list of fields for which the fvModel adds source term.
virtual void addSup(fvMatrix< vector > &eqn, const word &fieldName) const
Source term to momentum equation.
dynamicFvMesh & mesh
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
static const word null
An empty word.
Definition: word.H:77
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
dimensioned< scalar > magSqr(const dimensioned< Type > &)
Field< Type > & source()
Definition: fvMatrix.H:292
virtual bool read(const dictionary &dict)
Read dictionary.
addToRunTimeSelectionTable(fvConstraint, fixedTemperatureConstraint, dictionary)
List< word > wordList
A List of words.
Definition: fileName.H:54
U
Definition: pEqn.H:72
vector point
Point is a vector.
Definition: point.H:41
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
static Vector< scalar > uniform(const scalar &s)
Return a VectorSpace with all elements = s.
Definition: VectorSpaceI.H:168
dimensioned< scalar > mag(const dimensioned< Type > &)
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
Namespace for OpenFOAM.
IOerror FatalIOError