MixedDiffuseSpecular.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 "MixedDiffuseSpecular.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class CloudType>
32 (
33  const dictionary& dict,
34  CloudType& cloud
35 )
36 :
37  WallInteractionModel<CloudType>(dict, cloud, typeName),
38  diffuseFraction_(readScalar(this->coeffDict().lookup("diffuseFraction")))
39 {}
40 
41 
42 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
43 
44 template<class CloudType>
46 {}
47 
48 
49 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
50 
51 template<class CloudType>
53 (
54  typename CloudType::parcelType& p,
55  const wallPolyPatch& wpp
56 )
57 {
58  vector& U = p.U();
59 
60  scalar& Ei = p.Ei();
61 
62  label typeId = p.typeId();
63 
64  label wppIndex = wpp.index();
65 
66  label wppLocalFace = wpp.whichFace(p.face());
67 
68  vector nw = p.normal();
69  nw /= mag(nw);
70 
71  // Normal velocity magnitude
72  scalar U_dot_nw = U & nw;
73 
74  CloudType& cloud(this->owner());
75 
76  Random& rndGen(cloud.rndGen());
77 
78  if (diffuseFraction_ > rndGen.scalar01())
79  {
80  // Diffuse reflection
81 
82  // Wall tangential velocity (flow direction)
83  vector Ut = U - U_dot_nw*nw;
84 
85  while (mag(Ut) < SMALL)
86  {
87  // If the incident velocity is parallel to the face normal, no
88  // tangential direction can be chosen. Add a perturbation to the
89  // incoming velocity and recalculate.
90 
91  U = vector
92  (
93  U.x()*(0.8 + 0.2*rndGen.scalar01()),
94  U.y()*(0.8 + 0.2*rndGen.scalar01()),
95  U.z()*(0.8 + 0.2*rndGen.scalar01())
96  );
97 
98  U_dot_nw = U & nw;
99 
100  Ut = U - U_dot_nw*nw;
101  }
102 
103  // Wall tangential unit vector
104  vector tw1 = Ut/mag(Ut);
105 
106  // Other tangential unit vector
107  vector tw2 = nw^tw1;
108 
109  scalar T = cloud.boundaryT().boundaryField()[wppIndex][wppLocalFace];
110 
111  scalar mass = cloud.constProps(typeId).mass();
112 
113  direction iDof = cloud.constProps(typeId).internalDegreesOfFreedom();
114 
115  U =
116  sqrt(physicoChemical::k.value()*T/mass)
117  *(
118  rndGen.GaussNormal()*tw1
119  + rndGen.GaussNormal()*tw2
120  - sqrt(-2.0*log(max(1 - rndGen.scalar01(), VSMALL)))*nw
121  );
122 
123  U += cloud.boundaryU().boundaryField()[wppIndex][wppLocalFace];
124 
125  Ei = cloud.equipartitionInternalEnergy(T, iDof);
126  }
127  else
128  {
129  // Specular reflection
130 
131  if (U_dot_nw > 0.0)
132  {
133  U -= 2.0*U_dot_nw*nw;
134  }
135  }
136 
137 }
138 
139 
140 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
cachedRandom rndGen(label(0),-1)
const Cmpt & z() const
Definition: VectorI.H:87
dictionary dict
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 > &)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
label k
Boltzmann constant.
const volScalarField & boundaryT() const
Return macroscopic temperature.
Definition: DSMCCloudI.H:193
virtual ~MixedDiffuseSpecular()
Destructor.
stressControl lookup("compactNormalStress") >> compactNormalStress
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
virtual void correct(typename CloudType::parcelType &p, const wallPolyPatch &wpp)
Apply wall correction.
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
const List< typename ParcelType::constantProperties > & constProps() const
Return all of the constant properties.
Definition: DSMCCloudI.H:95
const volScalarField & T
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
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
MixedDiffuseSpecular(const dictionary &dict, CloudType &cloud)
Construct from dictionary.