dsmcFields.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 "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 
56 (
57  const word& name,
58  const Time& runTime,
59  const dictionary& dict
60 )
61 :
62  fvMeshFunctionObject(name, runTime, dict)
63 {
64  read(dict);
65 }
66 
67 
68 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
69 
71 {}
72 
73 
74 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
75 
77 {
78  return true;
79 }
80 
81 
83 {
84  return wordList
85  {
86  "rhoNMean",
87  "rhoMMean",
88  "momentumMean",
89  "linearKEMean",
90  "internalEMean",
91  "iDofMean",
92  "fDMean"
93  };
94 }
95 
96 
98 {
99  return true;
100 }
101 
102 
104 {
105  word rhoNMeanName = "rhoNMean";
106  word rhoMMeanName = "rhoMMean";
107  word momentumMeanName = "momentumMean";
108  word linearKEMeanName = "linearKEMean";
109  word internalEMeanName = "internalEMean";
110  word iDofMeanName = "iDofMean";
111  word fDMeanName = "fDMean";
112 
113  const volScalarField& rhoNMean = obr_.lookupObject<volScalarField>
114  (
115  rhoNMeanName
116  );
117 
118  const volScalarField& rhoMMean = obr_.lookupObject<volScalarField>
119  (
120  rhoMMeanName
121  );
122 
123  const volVectorField& momentumMean = obr_.lookupObject<volVectorField>
124  (
125  momentumMeanName
126  );
127 
128  const volScalarField& linearKEMean = obr_.lookupObject<volScalarField>
129  (
130  linearKEMeanName
131  );
132 
133  const volScalarField& internalEMean = obr_.lookupObject<volScalarField>
134  (
135  internalEMeanName
136  );
137 
138  const volScalarField& iDofMean = obr_.lookupObject<volScalarField>
139  (
140  iDofMeanName
141  );
142 
143  const volVectorField& fDMean = obr_.lookupObject<volVectorField>
144  (
145  fDMeanName
146  );
147 
148  if (min(mag(rhoNMean)).value() > vSmall)
149  {
150  Info<< "Calculating dsmcFields." << endl;
151 
152  Info<< " Calculating UMean field." << endl;
153  volVectorField UMean
154  (
155  IOobject
156  (
157  "UMean",
158  obr_.time().timeName(),
159  obr_,
161  ),
162  momentumMean/rhoMMean
163  );
164 
165  Info<< " Calculating translationalT field." << endl;
166  volScalarField translationalT
167  (
168  IOobject
169  (
170  "translationalT",
171  obr_.time().timeName(),
172  obr_,
174  ),
175 
176  2.0/(3.0*physicoChemical::k.value()*rhoNMean)
177  *(linearKEMean - 0.5*rhoMMean*(UMean & UMean))
178  );
179 
180  Info<< " Calculating internalT field." << endl;
181  volScalarField internalT
182  (
183  IOobject
184  (
185  "internalT",
186  obr_.time().timeName(),
187  obr_,
189  ),
190  (2.0/physicoChemical::k.value())*(internalEMean/iDofMean)
191  );
192 
193  Info<< " Calculating overallT field." << endl;
194  volScalarField overallT
195  (
196  IOobject
197  (
198  "overallT",
199  obr_.time().timeName(),
200  obr_,
202  ),
203  2.0/(physicoChemical::k.value()*(3.0*rhoNMean + iDofMean))
204  *(linearKEMean - 0.5*rhoMMean*(UMean & UMean) + internalEMean)
205  );
206 
207  Info<< " Calculating pressure field." << endl;
209  (
210  IOobject
211  (
212  "p",
213  obr_.time().timeName(),
214  obr_,
216  ),
217  physicoChemical::k.value()*rhoNMean*translationalT
218  );
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()/patch.magFaceAreas());
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.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
virtual Ostream & write(const char)=0
Write character.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:297
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:69
virtual bool write()
Calculate and write the DSMC fields.
Definition: dsmcFields.C:103
virtual ~dsmcFields()
Destructor.
Definition: dsmcFields.C:70
Macros for easy insertion into run-time selection tables.
virtual bool execute()
Do nothing.
Definition: dsmcFields.C:97
bool read(const char *, int32_t &)
Definition: int32IO.C:85
virtual wordList fields() const
Return the list of fields required.
Definition: dsmcFields.C:82
A class for handling words, derived from string.
Definition: word.H:59
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
const vectorField::subField faceAreas() const
Return face areas.
Definition: polyPatch.C:309
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
static const char nl
Definition: Ostream.H:260
defineTypeNameAndDebug(combustionModel, 0)
const scalarField::subField magFaceAreas() const
Return face area magnitudes.
Definition: polyPatch.C:315
const dimensionedScalar k
Boltzmann constant.
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
Specialisation of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
virtual bool write(const bool write=true) const
Write using setting from DB.
volScalarField & p
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:98
dsmcFields(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Definition: dsmcFields.C:56
virtual bool read(const dictionary &)
Read the dsmcFields data.
Definition: dsmcFields.C:76
Namespace for OpenFOAM.