LiftForce.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) 2012-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 "LiftForce.H"
27 #include "fvcCurl.H"
28 
29 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
30 
31 template<class CloudType>
33 (
34  const typename CloudType::parcelType& p,
35  const typename CloudType::parcelType::trackingData& td,
36  const vector& curlUc,
37  const scalar Re,
38  const scalar muc
39 ) const
40 {
41  // dummy
42  return 0.0;
43 }
44 
45 
46 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 
48 template<class CloudType>
50 (
51  CloudType& owner,
52  const fvMesh& mesh,
53  const dictionary& dict,
54  const word& forceType
55 )
56 :
57  ParticleForce<CloudType>(owner, mesh, dict, forceType, true),
58  UName_(this->coeffs().template lookupOrDefault<word>("U", "U")),
59  curlUcInterpPtr_(nullptr)
60 {}
61 
62 
63 template<class CloudType>
65 :
67  UName_(lf.UName_),
68  curlUcInterpPtr_(nullptr)
69 {}
70 
71 
72 // * * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * //
73 
74 template<class CloudType>
76 {}
77 
78 
79 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
80 
81 template<class CloudType>
83 {
84  static word fName("curlUcDt");
85 
86  bool fieldExists = this->mesh().template foundObject<volVectorField>(fName);
87 
88  if (store)
89  {
90  if (!fieldExists)
91  {
92  const volVectorField& Uc = this->mesh().template
93  lookupObject<volVectorField>(UName_);
94 
95  volVectorField* curlUcPtr =
96  new volVectorField(fName, fvc::curl(Uc));
97 
98  curlUcPtr->store();
99  }
100 
101  const volVectorField& curlUc = this->mesh().template
102  lookupObject<volVectorField>(fName);
103 
104  curlUcInterpPtr_.reset
105  (
107  (
108  this->owner().solution().interpolationSchemes(),
109  curlUc
110  ).ptr()
111  );
112  }
113  else
114  {
115  curlUcInterpPtr_.clear();
116 
117  if (fieldExists)
118  {
119  const volVectorField& curlUc = this->mesh().template
120  lookupObject<volVectorField>(fName);
121 
122  const_cast<volVectorField&>(curlUc).checkOut();
123  }
124  }
125 }
126 
127 
128 template<class CloudType>
130 (
131  const typename CloudType::parcelType& p,
132  const typename CloudType::parcelType::trackingData& td,
133  const scalar dt,
134  const scalar mass,
135  const scalar Re,
136  const scalar muc
137 ) const
138 {
139  forceSuSp value(Zero, 0.0);
140 
141  vector curlUc =
142  curlUcInterp().interpolate(p.coordinates(), p.currentTetIndices());
143 
144  scalar Cl = this->Cl(p, td, curlUc, Re, muc);
145 
146  value.Su() = mass/p.rho()*td.rhoc()*Cl*((td.Uc() - p.U())^curlUc);
147 
148  return value;
149 }
150 
151 
152 // ************************************************************************* //
virtual void cacheFields(const bool store)
Cache fields.
Definition: LiftForce.C:82
dictionary dict
LiftForce(CloudType &owner, const fvMesh &mesh, const dictionary &dict, const word &forceType)
Construct from mesh.
Definition: LiftForce.C:50
const vector & Su() const
Return const access to the explicit contribution [kg m/s^2].
Definition: forceSuSpI.H:56
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
const fvMesh & mesh() const
Return the mesh database.
Abstract base class for particle forces.
Definition: ParticleForce.H:53
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
const CloudType & owner() const
Return const access to the cloud owner.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:59
fvMesh & mesh
tmp< GeometricField< Type, fvPatchField, volMesh > > curl(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcCurl.C:45
autoPtr< interpolation< vector > > curlUcInterpPtr_
Curk of carrier phase velocity interpolator.
Definition: LiftForce.H:64
Helper container for force Su and Sp terms.
Definition: forceSuSp.H:61
virtual forceSuSp calcCoupled(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: LiftForce.C:130
const interpolation< vector > & curlUcInterp() const
Return the curl of the carrier phase velocity interpolator.
Definition: LiftForceI.H:30
Calculate the curl of the given volField by constructing the Hodge-dual of the symmetric part of the ...
A class for handling words, derived from string.
Definition: word.H:59
static const zero Zero
Definition: zero.H:97
virtual ~LiftForce()
Destructor.
Definition: LiftForce.C:75
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:221
virtual scalar Cl(const typename CloudType::parcelType &p, const typename CloudType::parcelType::trackingData &td, const vector &curlUc, const scalar Re, const scalar muc) const
Calculate the lift coefficient.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
Abstract base class for interpolation.
Base class for particle lift force models.
Definition: LiftForce.H:52
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:48
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:75
const word UName_
Name of velocity field.
Definition: LiftForce.H:61