MaxwellianThermal.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-2015 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 "MaxwellianThermal.H"
27 #include "constants.H"
28 
29 using namespace Foam::constant;
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
33 template<class CloudType>
35 (
36  const dictionary& dict,
37  CloudType& cloud
38 )
39 :
41 {}
42 
43 
44 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
45 
46 template<class CloudType>
48 {}
49 
50 
51 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
52 
53 template<class CloudType>
55 (
56  typename CloudType::parcelType& p,
57  const wallPolyPatch& wpp
58 )
59 {
60  vector& U = p.U();
61 
62  scalar& Ei = p.Ei();
63 
64  label typeId = p.typeId();
65 
66  label wppIndex = wpp.index();
67 
68  label wppLocalFace = wpp.whichFace(p.face());
69 
70  vector nw = p.normal();
71  nw /= mag(nw);
72 
73  // Normal velocity magnitude
74  scalar U_dot_nw = U & nw;
75 
76  // Wall tangential velocity (flow direction)
77  vector Ut = U - U_dot_nw*nw;
78 
79  CloudType& cloud(this->owner());
80 
81  Random& rndGen(cloud.rndGen());
82 
83  while (mag(Ut) < SMALL)
84  {
85  // If the incident velocity is parallel to the face normal, no
86  // tangential direction can be chosen. Add a perturbation to the
87  // incoming velocity and recalculate.
88 
89  U = vector
90  (
91  U.x()*(0.8 + 0.2*rndGen.scalar01()),
92  U.y()*(0.8 + 0.2*rndGen.scalar01()),
93  U.z()*(0.8 + 0.2*rndGen.scalar01())
94  );
95 
96  U_dot_nw = U & nw;
97 
98  Ut = U - U_dot_nw*nw;
99  }
100 
101  // Wall tangential unit vector
102  vector tw1 = Ut/mag(Ut);
103 
104  // Other tangential unit vector
105  vector tw2 = nw^tw1;
106 
107  scalar T = cloud.boundaryT().boundaryField()[wppIndex][wppLocalFace];
108 
109  scalar mass = cloud.constProps(typeId).mass();
110 
111  direction iDof = cloud.constProps(typeId).internalDegreesOfFreedom();
112 
113  U =
114  sqrt(physicoChemical::k.value()*T/mass)
115  *(
116  rndGen.GaussNormal()*tw1
117  + rndGen.GaussNormal()*tw2
118  - sqrt(-2.0*log(max(1 - rndGen.scalar01(), VSMALL)))*nw
119  );
120 
121  U += cloud.boundaryU().boundaryField()[wppIndex][wppLocalFace];
122 
123  Ei = cloud.equipartitionInternalEnergy(T, iDof);
124 }
125 
126 
127 // ************************************************************************* //
Collection of constants.
cachedRandom rndGen(label(0),-1)
const Cmpt & z() const
Definition: VectorI.H:87
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
uint8_t direction
Definition: direction.H:46
label whichFace(const label l) const
Return label of face in patch from global face label.
Definition: polyPatch.H:377
dimensionedScalar log(const dimensionedScalar &ds)
const Cmpt & x() const
Definition: VectorI.H:75
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 > &)
virtual void correct(typename CloudType::parcelType &p, const wallPolyPatch &wpp)
Apply wall correction.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
const volScalarField & boundaryT() const
Return macroscopic temperature.
Definition: DSMCCloudI.H:193
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Random & rndGen()
Return refernce to the random object.
Definition: DSMCCloudI.H:121
Foam::wallPolyPatch.
Definition: wallPolyPatch.H:48
const Cmpt & y() const
Definition: VectorI.H:81
Simple random number generator.
Definition: Random.H:49
scalar equipartitionInternalEnergy(scalar temperature, direction internalDegreesOfFreedom)
Generate a random internal energy, sampled from the.
Definition: DSMCCloud.C:1053
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:217
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const List< typename ParcelType::constantProperties > & constProps() const
Return all of the constant properties.
Definition: DSMCCloudI.H:95
MaxwellianThermal(const dictionary &dict, CloudType &cloud)
Construct from dictionary.
const dimensionedScalar k
Boltzmann constant.
const volVectorField & boundaryU() const
Return macroscopic velocity.
Definition: DSMCCloudI.H:201
dimensioned< scalar > mag(const dimensioned< Type > &)
Templated wall interaction model class.
Definition: DSMCCloud.H:58
virtual ~MaxwellianThermal()
Destructor.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:68
label index() const
Return the index of this patch in the boundaryMesh.
label nw
Definition: createFields.H:25