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-2020 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 true;
85 }
86 
87 
89 {
90  word rhoNMeanName = "rhoNMean";
91  word rhoMMeanName = "rhoMMean";
92  word momentumMeanName = "momentumMean";
93  word linearKEMeanName = "linearKEMean";
94  word internalEMeanName = "internalEMean";
95  word iDofMeanName = "iDofMean";
96  word fDMeanName = "fDMean";
97 
98  const volScalarField& rhoNMean = obr_.lookupObject<volScalarField>
99  (
100  rhoNMeanName
101  );
102 
103  const volScalarField& rhoMMean = obr_.lookupObject<volScalarField>
104  (
105  rhoMMeanName
106  );
107 
108  const volVectorField& momentumMean = obr_.lookupObject<volVectorField>
109  (
110  momentumMeanName
111  );
112 
113  const volScalarField& linearKEMean = obr_.lookupObject<volScalarField>
114  (
115  linearKEMeanName
116  );
117 
118  const volScalarField& internalEMean = obr_.lookupObject<volScalarField>
119  (
120  internalEMeanName
121  );
122 
123  const volScalarField& iDofMean = obr_.lookupObject<volScalarField>
124  (
125  iDofMeanName
126  );
127 
128  const volVectorField& fDMean = obr_.lookupObject<volVectorField>
129  (
130  fDMeanName
131  );
132 
133  if (min(mag(rhoNMean)).value() > vSmall)
134  {
135  Info<< "Calculating dsmcFields." << endl;
136 
137  Info<< " Calculating UMean field." << endl;
139  (
140  IOobject
141  (
142  "UMean",
143  obr_.time().timeName(),
144  obr_,
146  ),
147  momentumMean/rhoMMean
148  );
149 
150  Info<< " Calculating translationalT field." << endl;
151  volScalarField translationalT
152  (
153  IOobject
154  (
155  "translationalT",
156  obr_.time().timeName(),
157  obr_,
159  ),
160 
161  2.0/(3.0*physicoChemical::k.value()*rhoNMean)
162  *(linearKEMean - 0.5*rhoMMean*(UMean & UMean))
163  );
164 
165  Info<< " Calculating internalT field." << endl;
166  volScalarField internalT
167  (
168  IOobject
169  (
170  "internalT",
171  obr_.time().timeName(),
172  obr_,
174  ),
175  (2.0/physicoChemical::k.value())*(internalEMean/iDofMean)
176  );
177 
178  Info<< " Calculating overallT field." << endl;
179  volScalarField overallT
180  (
181  IOobject
182  (
183  "overallT",
184  obr_.time().timeName(),
185  obr_,
187  ),
188  2.0/(physicoChemical::k.value()*(3.0*rhoNMean + iDofMean))
189  *(linearKEMean - 0.5*rhoMMean*(UMean & UMean) + internalEMean)
190  );
191 
192  Info<< " Calculating pressure field." << endl;
194  (
195  IOobject
196  (
197  "p",
198  obr_.time().timeName(),
199  obr_,
201  ),
202  physicoChemical::k.value()*rhoNMean*translationalT
203  );
204 
205  volScalarField::Boundary& pBf = p.boundaryFieldRef();
206 
207  forAll(mesh_.boundaryMesh(), i)
208  {
209  const polyPatch& patch = mesh_.boundaryMesh()[i];
210 
211  if (isA<wallPolyPatch>(patch))
212  {
213  pBf[i] =
214  fDMean.boundaryField()[i]
215  & (patch.faceAreas()/patch.magFaceAreas());
216  }
217  }
218 
219  Info<< " mag(UMean) max/min : "
220  << max(mag(UMean)).value() << " "
221  << min(mag(UMean)).value() << endl;
222 
223  Info<< " translationalT max/min : "
224  << max(translationalT).value() << " "
225  << min(translationalT).value() << endl;
226 
227  Info<< " internalT max/min : "
228  << max(internalT).value() << " "
229  << min(internalT).value() << endl;
230 
231  Info<< " overallT max/min : "
232  << max(overallT).value() << " "
233  << min(overallT).value() << endl;
234 
235  Info<< " p max/min : "
236  << max(p).value() << " "
237  << min(p).value() << endl;
238 
239  UMean.write();
240 
241  translationalT.write();
242 
243  internalT.write();
244 
245  overallT.write();
246 
247  p.write();
248 
249  Info<< "dsmcFields written." << nl << endl;
250 
251  return true;
252  }
253  else
254  {
255  Info<< "Small value (" << min(mag(rhoNMean))
256  << ") found in rhoNMean field. "
257  << "Not calculating dsmcFields to avoid division by zero."
258  << endl;
259 
260  return false;
261  }
262 }
263 
264 
265 // ************************************************************************* //
Collection of constants.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:291
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
volVectorField UMean(UMeanHeader, mesh)
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:88
virtual ~dsmcFields()
Destructor.
Definition: dsmcFields.C:70
Macros for easy insertion into run-time selection tables.
virtual bool execute()
Do nothing.
Definition: dsmcFields.C:82
bool read(const char *, int32_t &)
Definition: int32IO.C:85
A class for handling words, derived from string.
Definition: word.H:59
const vectorField::subField faceAreas() const
Return face areas.
Definition: polyPatch.C:303
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
static const char nl
Definition: Ostream.H:260
defineTypeNameAndDebug(combustionModel, 0)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const scalarField::subField magFaceAreas() const
Return face area magnitudes.
Definition: polyPatch.C:309
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:92
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.