dsmcFields.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 "dsmcFields.H"
27 #include "volFields.H"
28 #include "dictionary.H"
29 #include "dsmcCloud.H"
30 #include "constants.H"
32 
33 using namespace Foam::constant;
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 namespace functionObjects
40 {
41  defineTypeNameAndDebug(dsmcFields, 0);
42 
44  (
45  functionObject,
46  dsmcFields,
47  dictionary
48  );
49 }
50 }
51 
52 
53 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
54 
55 Foam::functionObjects::dsmcFields::dsmcFields
56 (
57  const word& name,
58  const Time& runTime,
59  const dictionary& dict
60 )
61 :
62  functionObject(name),
63  obr_
64  (
66  (
68  )
69  )
70 {
71  if (!isA<fvMesh>(obr_))
72  {
74  << "objectRegistry is not an fvMesh" << exit(FatalError);
75  }
76 
77  read(dict);
78 }
79 
80 
81 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
82 
84 {}
85 
86 
87 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
88 
90 {
91  return true;
92 }
93 
94 
96 {
97  return true;
98 }
99 
100 
102 {
103  word rhoNMeanName = "rhoNMean";
104  word rhoMMeanName = "rhoMMean";
105  word momentumMeanName = "momentumMean";
106  word linearKEMeanName = "linearKEMean";
107  word internalEMeanName = "internalEMean";
108  word iDofMeanName = "iDofMean";
109  word fDMeanName = "fDMean";
110 
111  const volScalarField& rhoNMean = obr_.lookupObject<volScalarField>
112  (
113  rhoNMeanName
114  );
115 
116  const volScalarField& rhoMMean = obr_.lookupObject<volScalarField>
117  (
118  rhoMMeanName
119  );
120 
121  const volVectorField& momentumMean = obr_.lookupObject<volVectorField>
122  (
123  momentumMeanName
124  );
125 
126  const volScalarField& linearKEMean = obr_.lookupObject<volScalarField>
127  (
128  linearKEMeanName
129  );
130 
131  const volScalarField& internalEMean = obr_.lookupObject<volScalarField>
132  (
133  internalEMeanName
134  );
135 
136  const volScalarField& iDofMean = obr_.lookupObject<volScalarField>
137  (
138  iDofMeanName
139  );
140 
141  const volVectorField& fDMean = obr_.lookupObject<volVectorField>
142  (
143  fDMeanName
144  );
145 
146  if (min(mag(rhoNMean)).value() > VSMALL)
147  {
148  Info<< "Calculating dsmcFields." << endl;
149 
150  Info<< " Calculating UMean field." << endl;
152  (
153  IOobject
154  (
155  "UMean",
156  obr_.time().timeName(),
157  obr_,
159  ),
160  momentumMean/rhoMMean
161  );
162 
163  Info<< " Calculating translationalT field." << endl;
164  volScalarField translationalT
165  (
166  IOobject
167  (
168  "translationalT",
169  obr_.time().timeName(),
170  obr_,
172  ),
173 
174  2.0/(3.0*physicoChemical::k.value()*rhoNMean)
175  *(linearKEMean - 0.5*rhoMMean*(UMean & UMean))
176  );
177 
178  Info<< " Calculating internalT field." << endl;
179  volScalarField internalT
180  (
181  IOobject
182  (
183  "internalT",
184  obr_.time().timeName(),
185  obr_,
187  ),
188  (2.0/physicoChemical::k.value())*(internalEMean/iDofMean)
189  );
190 
191  Info<< " Calculating overallT field." << endl;
192  volScalarField overallT
193  (
194  IOobject
195  (
196  "overallT",
197  obr_.time().timeName(),
198  obr_,
200  ),
201  2.0/(physicoChemical::k.value()*(3.0*rhoNMean + iDofMean))
202  *(linearKEMean - 0.5*rhoMMean*(UMean & UMean) + internalEMean)
203  );
204 
205  Info<< " Calculating pressure field." << endl;
207  (
208  IOobject
209  (
210  "p",
211  obr_.time().timeName(),
212  obr_,
214  ),
215  physicoChemical::k.value()*rhoNMean*translationalT
216  );
217 
218  const fvMesh& mesh = fDMean.mesh();
219 
220  volScalarField::Boundary& pBf = p.boundaryFieldRef();
221 
222  forAll(mesh.boundaryMesh(), i)
223  {
224  const polyPatch& patch = mesh.boundaryMesh()[i];
225 
226  if (isA<wallPolyPatch>(patch))
227  {
228  pBf[i] =
229  fDMean.boundaryField()[i]
230  & (patch.faceAreas()/mag(patch.faceAreas()));
231  }
232  }
233 
234  Info<< " mag(UMean) max/min : "
235  << max(mag(UMean)).value() << " "
236  << min(mag(UMean)).value() << endl;
237 
238  Info<< " translationalT max/min : "
239  << max(translationalT).value() << " "
240  << min(translationalT).value() << endl;
241 
242  Info<< " internalT max/min : "
243  << max(internalT).value() << " "
244  << min(internalT).value() << endl;
245 
246  Info<< " overallT max/min : "
247  << max(overallT).value() << " "
248  << min(overallT).value() << endl;
249 
250  Info<< " p max/min : "
251  << max(p).value() << " "
252  << min(p).value() << endl;
253 
254  UMean.write();
255 
256  translationalT.write();
257 
258  internalT.write();
259 
260  overallT.write();
261 
262  p.write();
263 
264  Info<< "dsmcFields written." << nl << endl;
265 
266  return true;
267  }
268  else
269  {
270  Info<< "Small value (" << min(mag(rhoNMean))
271  << ") found in rhoNMean field. "
272  << "Not calculating dsmcFields to avoid division by zero."
273  << endl;
274 
275  return false;
276  }
277 }
278 
279 
280 // ************************************************************************* //
Collection of constants.
const vectorField::subField faceAreas() const
Return face normals.
Definition: polyPatch.C:290
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:306
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
volVectorField UMean(UMeanHeader, mesh)
Abstract base-class for Time/database function objects.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
virtual bool write()
Calculate and write the DSMC fields.
Definition: dsmcFields.C:101
virtual ~dsmcFields()
Destructor.
Definition: dsmcFields.C:83
Macros for easy insertion into run-time selection tables.
virtual bool execute()
Do nothing.
Definition: dsmcFields.C:95
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
bool read(const char *, int32_t &)
Definition: int32IO.C:85
dynamicFvMesh & mesh
A class for handling words, derived from string.
Definition: word.H:59
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
static const char nl
Definition: Ostream.H:262
defineTypeNameAndDebug(combustionModel, 0)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const dimensionedScalar k
Boltzmann constant.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual bool write() const
Write using setting from DB.
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual Ostream & write(const token &)=0
Write next token to stream.
volScalarField & p
Registry of regIOobjects.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:278
virtual bool read(const dictionary &)
Read the dsmcFields data.
Definition: dsmcFields.C:89
Namespace for OpenFOAM.