DSMCParcel.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "DSMCParcel.H"
29#include "meshTools.H"
30
31// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
32
33template<class ParcelType>
34template<class TrackCloudType>
36(
37 TrackCloudType& cloud,
38 trackingData& td,
39 const scalar trackTime
40)
41{
42 typename TrackCloudType::parcelType& p =
43 static_cast<typename TrackCloudType::parcelType&>(*this);
44
45 td.switchProcessor = false;
46 td.keepParticle = true;
47
48 const polyMesh& mesh = cloud.pMesh();
49
50 // For reduced-D cases, the velocity used to track needs to be
51 // constrained, but the actual U_ of the parcel must not be
52 // altered or used, as it is altered by patch interactions an
53 // needs to retain its 3D value for collision purposes.
54 vector Utracking = U_;
55
56 while (td.keepParticle && !td.switchProcessor && p.stepFraction() < 1)
57 {
58 Utracking = U_;
59
60 // Apply correction to velocity to constrain tracking for
61 // reduced-D cases
62 meshTools::constrainDirection(mesh, mesh.solutionD(), Utracking);
63
64 // Deviation from the mesh centre for reduced-D cases
65 const vector d = p.deviationFromMeshCentre();
66
67 const scalar f = 1 - p.stepFraction();
68 p.trackToAndHitFace(f*trackTime*Utracking - d, f, cloud, td);
69 }
70
71 return td.keepParticle;
72}
73
74
75template<class ParcelType>
76template<class TrackCloudType>
78{
79 return false;
80}
81
82
83template<class ParcelType>
84template<class TrackCloudType>
86(
87 TrackCloudType&,
88 trackingData& td
89)
90{
91 td.switchProcessor = true;
92}
93
94
95template<class ParcelType>
96template<class TrackCloudType>
98(
99 TrackCloudType& cloud,
101)
102{
103 const label wppIndex = this->patch();
104
105 const wallPolyPatch& wpp =
106 static_cast<const wallPolyPatch&>
107 (
108 this->mesh().boundaryMesh()[wppIndex]
109 );
110
111 const label wppLocalFace = wpp.whichFace(this->face());
112
113 const scalar fA = mag(wpp.faceAreas()[wppLocalFace]);
114
115 const scalar deltaT = cloud.pMesh().time().deltaTValue();
116
117 const constantProperties& constProps(cloud.constProps(typeId_));
118
119 scalar m = constProps.mass();
120
121 const vector nw = normalised(wpp.faceAreas()[wppLocalFace]);
122
123 scalar U_dot_nw = U_ & nw;
124
125 vector Ut = U_ - U_dot_nw*nw;
126
127 scalar invMagUnfA = 1/max(mag(U_dot_nw)*fA, VSMALL);
128
129 cloud.rhoNBF()[wppIndex][wppLocalFace] += invMagUnfA;
130
131 cloud.rhoMBF()[wppIndex][wppLocalFace] += m*invMagUnfA;
132
133 cloud.linearKEBF()[wppIndex][wppLocalFace] +=
134 0.5*m*(U_ & U_)*invMagUnfA;
135
136 cloud.internalEBF()[wppIndex][wppLocalFace] += Ei_*invMagUnfA;
137
138 cloud.iDofBF()[wppIndex][wppLocalFace] +=
139 constProps.internalDegreesOfFreedom()*invMagUnfA;
140
141 cloud.momentumBF()[wppIndex][wppLocalFace] += m*Ut*invMagUnfA;
142
143 // pre-interaction energy
144 scalar preIE = 0.5*m*(U_ & U_) + Ei_;
145
146 // pre-interaction momentum
147 vector preIMom = m*U_;
148
149 cloud.wallInteraction().correct(*this);
150
151 U_dot_nw = U_ & nw;
152
153 Ut = U_ - U_dot_nw*nw;
154
155 invMagUnfA = 1/max(mag(U_dot_nw)*fA, VSMALL);
156
157 cloud.rhoNBF()[wppIndex][wppLocalFace] += invMagUnfA;
158
159 cloud.rhoMBF()[wppIndex][wppLocalFace] += m*invMagUnfA;
160
161 cloud.linearKEBF()[wppIndex][wppLocalFace] +=
162 0.5*m*(U_ & U_)*invMagUnfA;
163
164 cloud.internalEBF()[wppIndex][wppLocalFace] += Ei_*invMagUnfA;
165
166 cloud.iDofBF()[wppIndex][wppLocalFace] +=
167 constProps.internalDegreesOfFreedom()*invMagUnfA;
168
169 cloud.momentumBF()[wppIndex][wppLocalFace] += m*Ut*invMagUnfA;
170
171 // post-interaction energy
172 scalar postIE = 0.5*m*(U_ & U_) + Ei_;
173
174 // post-interaction momentum
175 vector postIMom = m*U_;
176
177 scalar deltaQ = cloud.nParticle()*(preIE - postIE)/(deltaT*fA);
178
179 vector deltaFD = cloud.nParticle()*(preIMom - postIMom)/(deltaT*fA);
180
181 cloud.qBF()[wppIndex][wppLocalFace] += deltaQ;
182
183 cloud.fDBF()[wppIndex][wppLocalFace] += deltaFD;
184}
185
186
187template<class ParcelType>
189{
191 U_ = transform(T, U_);
192}
193
194
195template<class ParcelType>
197(
198 const vector& separation
199)
200{
202}
203
204
205// * * * * * * * * * * * * * * * * IOStream operators * * * * * * * * * * * //
206
207#include "DSMCParcelIO.C"
208
209// ************************************************************************* //
Y[inertIndex] max(0.0)
Class to hold DSMC particle constant properties.
Definition: DSMCParcel.H:82
scalar mass() const
Return const access to the particle mass [kg].
Definition: DSMCParcelI.H:95
direction internalDegreesOfFreedom() const
Return the internalDegreesOfFreedom.
Definition: DSMCParcelI.H:118
virtual void transformProperties(const tensor &T)
Transform the physical properties of the particle.
Definition: DSMCParcel.C:188
void hitProcessorPatch(TrackCloudType &, trackingData &)
Overridable function to handle the particle hitting a.
Definition: DSMCParcel.C:86
bool hitPatch(TrackCloudType &, trackingData &)
Overridable function to handle the particle hitting a patch.
Definition: DSMCParcel.C:77
void hitWallPatch(TrackCloudType &, trackingData &)
Overridable function to handle the particle hitting a wallPatch.
Definition: DSMCParcel.C:98
scalar deltaTValue() const noexcept
Return time step value.
Definition: TimeStateI.H:43
A cloud is a registry collection of lagrangian particles.
Definition: cloud.H:60
virtual void move()=0
Class used to pass data into container.
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
const Time & time() const noexcept
Return time registry.
virtual void transformProperties(const tensor &T)
Transform the physical properties of the particle.
Definition: particle.C:1072
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
label whichFace(const label l) const
Return label of face in patch from global face label.
Definition: polyPatch.H:451
const vectorField::subField faceAreas() const
Return face normals.
Definition: polyPatch.C:327
Foam::wallPolyPatch.
Definition: wallPolyPatch.H:53
volScalarField & p
const volScalarField & T
dynamicFvMesh & mesh
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
Definition: dimensionSet.C:536
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
Definition: quaternionI.H:680
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
labelList f(nPoints)