actuationDiskSource.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2016 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  option,
42  actuationDiskSource,
43  dictionary
44  );
45 }
46 }
47 
48 
49 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
50 
51 void Foam::fv::actuationDiskSource::checkData() const
52 {
53  if (magSqr(diskArea_) <= VSMALL)
54  {
56  << "diskArea is approximately zero"
57  << exit(FatalIOError);
58  }
59  if (Cp_ <= VSMALL || Ct_ <= VSMALL)
60  {
62  << "Cp and Ct must be greater than zero"
63  << exit(FatalIOError);
64  }
65  if (mag(diskDir_) < VSMALL)
66  {
68  << "disk direction vector is approximately zero"
69  << exit(FatalIOError);
70  }
71  if (returnReduce(upstreamCellId_, maxOp<label>()) == -1)
72  {
74  << "upstream location " << upstreamPoint_ << " not found in mesh"
75  << exit(FatalIOError);
76  }
77 }
78 
79 
80 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
81 
82 Foam::fv::actuationDiskSource::actuationDiskSource
83 (
84  const word& name,
85  const word& modelType,
86  const dictionary& dict,
87  const fvMesh& mesh
88 )
89 :
90  cellSetOption(name, modelType, dict, mesh),
91  diskDir_(coeffs_.lookup("diskDir")),
92  Cp_(readScalar(coeffs_.lookup("Cp"))),
93  Ct_(readScalar(coeffs_.lookup("Ct"))),
94  diskArea_(readScalar(coeffs_.lookup("diskArea"))),
95  upstreamPoint_(coeffs_.lookup("upstreamPoint")),
96  upstreamCellId_(-1)
97 {
98  coeffs_.lookup("fields") >> fieldNames_;
99  applied_.setSize(fieldNames_.size(), false);
100 
101  Info<< " - creating actuation disk zone: "
102  << this->name() << endl;
103 
104  upstreamCellId_ = mesh.findCell(upstreamPoint_);
105 
106  checkData();
107 }
108 
109 
110 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
111 
113 (
114  fvMatrix<vector>& eqn,
115  const label fieldi
116 )
117 {
118  const scalarField& cellsV = mesh_.V();
119  vectorField& Usource = eqn.source();
120  const vectorField& U = eqn.psi();
121 
122  if (V() > VSMALL)
123  {
124  addActuationDiskAxialInertialResistance
125  (
126  Usource,
127  cells_,
128  cellsV,
130  U
131  );
132  }
133 }
134 
135 
137 (
138  const volScalarField& rho,
139  fvMatrix<vector>& eqn,
140  const label fieldi
141 )
142 {
143  const scalarField& cellsV = mesh_.V();
144  vectorField& Usource = eqn.source();
145  const vectorField& U = eqn.psi();
146 
147  if (V() > VSMALL)
148  {
149  addActuationDiskAxialInertialResistance
150  (
151  Usource,
152  cells_,
153  cellsV,
154  rho,
155  U
156  );
157  }
158 }
159 
160 
162 {
163  if (cellSetOption::read(dict))
164  {
165  coeffs_.readIfPresent("diskDir", diskDir_);
166  coeffs_.readIfPresent("Cp", Cp_);
167  coeffs_.readIfPresent("Ct", Ct_);
168  coeffs_.readIfPresent("diskArea", diskArea_);
169 
170  checkData();
171 
172  return true;
173  }
174  else
175  {
176  return false;
177  }
178 }
179 
180 
181 // ************************************************************************* //
defineTypeNameAndDebug(option, 0)
U
Definition: pEqn.H:83
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:281
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
addToRunTimeSelectionTable(option, fixedTemperatureConstraint, dictionary)
virtual bool read(const dictionary &dict)
Read source dictionary.
Macros for easy insertion into run-time selection tables.
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)
virtual void addSup(fvMatrix< vector > &eqn, const label fieldi)
Source term to momentum equation.
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:71
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
dimensioned< scalar > magSqr(const dimensioned< Type > &)
Field< Type > & source()
Definition: fvMatrix.H:291
virtual bool read(const dictionary &dict)
Read dictionary.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
label findCell(const point &p, const cellDecomposition=CELL_TETS) const
Find cell enclosing this location and return index.
Definition: polyMesh.C:1408
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
Cell-set options abtract base class. Provides a base set of controls, e.g.:
Definition: cellSetOption.H:72
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
Namespace for OpenFOAM.
IOerror FatalIOError