NonInertialFrameForce.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  Copyright (C) 2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "NonInertialFrameForce.H"
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
34 template<class CloudType>
36 (
37  CloudType& owner,
38  const fvMesh& mesh,
39  const dictionary& dict
40 )
41 :
42  ParticleForce<CloudType>(owner, mesh, dict, typeName, true),
43  WName_
44  (
45  this->coeffs().template getOrDefault<word>
46  (
47  "linearAcceleration",
48  "linearAcceleration"
49  )
50  ),
51  W_(Zero),
52  omegaName_
53  (
54  this->coeffs().template getOrDefault<word>
55  (
56  "angularVelocity",
57  "angularVelocity"
58  )
59  ),
60  omega_(Zero),
61  omegaDotName_
62  (
63  this->coeffs().template getOrDefault<word>
64  (
65  "angularAcceleration",
66  "angularAcceleration"
67  )
68  ),
69  omegaDot_(Zero),
70  centreOfRotationName_
71  (
72  this->coeffs().template getOrDefault<word>
73  (
74  "centreOfRotation",
75  "centreOfRotation"
76  )
77  ),
78  centreOfRotation_(Zero)
79 {}
80 
81 
82 template<class CloudType>
84 (
85  const NonInertialFrameForce& niff
86 )
87 :
89  WName_(niff.WName_),
90  W_(niff.W_),
91  omegaName_(niff.omegaName_),
92  omega_(niff.omega_),
93  omegaDotName_(niff.omegaDotName_),
94  omegaDot_(niff.omegaDot_),
95  centreOfRotationName_(niff.centreOfRotationName_),
96  centreOfRotation_(niff.centreOfRotation_)
97 {}
98 
99 
100 // * * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * //
101 
102 template<class CloudType>
104 {}
105 
106 
107 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
108 
109 template<class CloudType>
111 {
112  W_ = Zero;
113  omega_ = Zero;
114  omegaDot_ = Zero;
115  centreOfRotation_ = Zero;
116 
117  if (store)
118  {
119  if
120  (
121  this->mesh().template foundObject<uniformDimensionedVectorField>
122  (
123  WName_
124  )
125  )
126  {
127  const uniformDimensionedVectorField& W = this->mesh().template
128  lookupObject<uniformDimensionedVectorField>(WName_);
129 
130  W_ = W.value();
131  }
132 
133  if
134  (
135  this->mesh().template foundObject<uniformDimensionedVectorField>
136  (
137  omegaName_
138  )
139  )
140  {
141  const uniformDimensionedVectorField& omega = this->mesh().template
142  lookupObject<uniformDimensionedVectorField>(omegaName_);
143 
144  omega_ = omega.value();
145  }
146 
147  if
148  (
149  this->mesh().template foundObject<uniformDimensionedVectorField>
150  (
151  omegaDotName_
152  )
153  )
154  {
155  const uniformDimensionedVectorField& omegaDot =
156  this->mesh().template
157  lookupObject<uniformDimensionedVectorField>(omegaDotName_);
158 
159  omegaDot_ = omegaDot.value();
160  }
161 
162  if
163  (
164  this->mesh().template foundObject<uniformDimensionedVectorField>
165  (
166  centreOfRotationName_
167  )
168  )
169  {
170  const uniformDimensionedVectorField& centreOfRotation =
171  this->mesh().template
172  lookupObject<uniformDimensionedVectorField>
173  (
174  centreOfRotationName_
175  );
176 
177  centreOfRotation_ = centreOfRotation.value();
178  }
179  }
180 }
181 
182 
183 template<class CloudType>
185 (
186  const typename CloudType::parcelType& p,
187  const typename CloudType::parcelType::trackingData& td,
188  const scalar dt,
189  const scalar mass,
190  const scalar Re,
191  const scalar muc
192 ) const
193 {
194  forceSuSp value(Zero);
195 
196  const vector r = p.position() - centreOfRotation_;
197 
198  value.Su() =
199  mass
200  *(
201  -W_
202  + (r ^ omegaDot_)
203  + 2.0*(p.U() ^ omega_)
204  + (omega_ ^ (r ^ omega_))
205  );
206 
207  return value;
208 }
209 
210 
211 // ************************************************************************* //
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::NonInertialFrameForce
Calculates particle non-inertial reference frame force. Variable names as from Landau and Lifshitz,...
Definition: NonInertialFrameForce.H:57
Foam::NonInertialFrameForce::cacheFields
virtual void cacheFields(const bool store)
Cache fields.
Definition: NonInertialFrameForce.C:110
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
NonInertialFrameForce.H
Foam::NonInertialFrameForce::~NonInertialFrameForce
virtual ~NonInertialFrameForce()
Destructor.
Definition: NonInertialFrameForce.C:103
Foam::dimensioned::value
const Type & value() const
Return const reference to value.
Definition: dimensionedType.C:434
Foam::forceSuSp::Su
const vector & Su() const
Return const access to the explicit contribution [kg.m/s2].
Definition: forceSuSpI.H:61
Foam::UniformDimensionedField< vector >
Foam::forceSuSp
Helper container for force Su and Sp terms.
Definition: forceSuSp.H:64
Foam::ParticleForce
Abstract base class for particle forces.
Definition: ParticleForce.H:59
Foam::DSMCCloud
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:71
Foam::NonInertialFrameForce::NonInertialFrameForce
NonInertialFrameForce(CloudType &owner, const fvMesh &mesh, const dictionary &dict)
Construct from mesh.
Definition: NonInertialFrameForce.C:36
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::NonInertialFrameForce::calcNonCoupled
virtual forceSuSp calcNonCoupled(const typename CloudType::parcelType &p, const typename CloudType::parcelType::trackingData &td, const scalar dt, const scalar mass, const scalar Re, const scalar muc) const
Calculate the non-coupled force.
Definition: NonInertialFrameForce.C:185
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:85
uniformDimensionedFields.H
Foam::Vector< scalar >
Foam::Re
scalarField Re(const UList< complex > &cf)
Extract real component.
Definition: complexField.C:159
Foam::DSMCCloud::parcelType
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:220