SchaefferFrictionalStress.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-2018 OpenFOAM Foundation
9  Copyright (C) 2019 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 
31 #include "unitConversion.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace kineticTheoryModels
38 {
39 namespace frictionalStressModels
40 {
41  defineTypeNameAndDebug(Schaeffer, 0);
42 
44  (
45  frictionalStressModel,
46  Schaeffer,
47  dictionary
48  );
49 }
50 }
51 }
52 
53 
54 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
55 
57 (
58  const dictionary& dict
59 )
60 :
62  coeffDict_(dict.optionalSubDict(typeName + "Coeffs")),
63  phi_("phi", dimless, coeffDict_)
64 {
65  phi_ *= degToRad();
66 }
67 
68 
69 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
70 
72 {}
73 
74 
75 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
76 
80 (
81  const phaseModel& phase,
82  const dimensionedScalar& alphaMinFriction,
84 ) const
85 {
86  const volScalarField& alpha = phase;
87 
88  return
89  dimensionedScalar("", dimensionSet(1, -1, -2, 0, 0), 1e24)
90  *pow(Foam::max(alpha - alphaMinFriction, scalar(0)), 10.0);
91 }
92 
93 
97 (
98  const phaseModel& phase,
99  const dimensionedScalar& alphaMinFriction,
101 ) const
102 {
103  const volScalarField& alpha = phase;
104 
105  return
106  dimensionedScalar("", dimensionSet(1, -1, -2, 0, 0), 1e25)
107  *pow(Foam::max(alpha - alphaMinFriction, scalar(0)), 9.0);
108 }
109 
110 
113 (
114  const phaseModel& phase,
115  const dimensionedScalar& alphaMinFriction,
117  const volScalarField& pf,
118  const volSymmTensorField& D
119 ) const
120 {
121  const volScalarField& alpha = phase;
122 
124  (
125  new volScalarField
126  (
127  IOobject
128  (
129  "Schaeffer:nu",
130  phase.mesh().time().timeName(),
131  phase.mesh(),
134  false
135  ),
136  phase.mesh(),
137  dimensionedScalar(dimensionSet(0, 2, -1, 0, 0))
138  )
139  );
140 
141  volScalarField& nuf = tnu.ref();
142 
143  forAll(D, celli)
144  {
145  if (alpha[celli] > alphaMinFriction.value())
146  {
147  nuf[celli] =
148  0.5*pf[celli]*sin(phi_.value())
149  /(
150  sqrt((1.0/3.0)*sqr(tr(D[celli])) - invariantII(D[celli]))
151  + SMALL
152  );
153  }
154  }
155 
156  const fvPatchList& patches = phase.mesh().boundary();
157  const volVectorField& U = phase.U();
158 
159  volScalarField::Boundary& nufBf = nuf.boundaryFieldRef();
160 
161  forAll(patches, patchi)
162  {
163  if (!patches[patchi].coupled())
164  {
165  nufBf[patchi] =
166  (
167  pf.boundaryField()[patchi]*sin(phi_.value())
168  /(
169  mag(U.boundaryField()[patchi].snGrad())
170  + SMALL
171  )
172  );
173  }
174  }
175 
176  // Correct coupled BCs
178 
179  return tnu;
180 }
181 
182 
184 {
185  coeffDict_ <<= dict_.optionalSubDict(typeName + "Coeffs");
186 
187  phi_.read(coeffDict_);
188  phi_ *= degToRad();
189 
190  return true;
191 }
192 
193 
194 // ************************************************************************* //
Foam::kineticTheoryModels::frictionalStressModel
Definition: frictionalStressModel.H:54
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
Foam::phaseModel
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:54
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::kineticTheoryModels::frictionalStressModels::Schaeffer::frictionalPressure
virtual tmp< volScalarField > frictionalPressure(const phaseModel &phase, const dimensionedScalar &alphaMinFriction, const dimensionedScalar &alphaMax) const
Definition: SchaefferFrictionalStress.C:80
Foam::phase
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phase.H:54
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::constant::atomic::alpha
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Definition: readThermalProperties.H:212
Foam::sin
dimensionedScalar sin(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:264
Foam::kineticTheoryModels::frictionalStressModels::Schaeffer::~Schaeffer
virtual ~Schaeffer()
Destructor.
Definition: SchaefferFrictionalStress.C:71
unitConversion.H
Unit conversion functions.
Foam::dimensioned::value
const Type & value() const
Return const reference to value.
Definition: dimensionedType.C:434
Foam::dimensionSet
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
Definition: dimensionSet.H:108
SchaefferFrictionalStress.H
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::tmp::ref
T & ref() const
Definition: tmpI.H:227
Foam::kineticTheoryModels::frictionalStressModels::Schaeffer::read
virtual bool read()
Definition: SchaefferFrictionalStress.C:183
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
Foam::PtrList
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: List.H:59
Foam::pow
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Definition: dimensionedScalar.C:75
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
alphaMax
dimensionedScalar alphaMax("alphaMax", dimless/dimTime, laminarTransport)
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::kineticTheoryModels::frictionalStressModels::defineTypeNameAndDebug
defineTypeNameAndDebug(JohnsonJackson, 0)
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
Foam::kineticTheoryModels::frictionalStressModels::addToRunTimeSelectionTable
addToRunTimeSelectionTable(frictionalStressModel, JohnsonJackson, dictionary)
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::dimensioned< scalar >
Foam::invariantII
Cmpt invariantII(const SymmTensor< Cmpt > &st)
Return the 2nd invariant of a SymmTensor.
Definition: SymmTensorI.H:350
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::degToRad
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
Definition: unitConversion.H:48
Foam::kineticTheoryModels::frictionalStressModels::Schaeffer::frictionalPressurePrime
virtual tmp< volScalarField > frictionalPressurePrime(const phaseModel &phase, const dimensionedScalar &alphaMinFriction, const dimensionedScalar &alphaMax) const
Definition: SchaefferFrictionalStress.C:97
Foam::GeometricField::correctBoundaryConditions
void correctBoundaryConditions()
Correct boundary field.
Definition: GeometricField.C:940
Foam::kineticTheoryModels::frictionalStressModels::Schaeffer::Schaeffer
Schaeffer(const dictionary &dict)
Construct from components.
Definition: SchaefferFrictionalStress.C:57
U
U
Definition: pEqn.H:72
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
Foam::GeometricField::boundaryFieldRef
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
Definition: GeometricField.C:783
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
D
const dimensionedScalar & D
Definition: solveBulkSurfactant.H:4
patches
const polyBoundaryMesh & patches
Definition: convertProcessorPatches.H:65
Foam::tr
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
Definition: dimensionedSphericalTensor.C:51
coupled
bool coupled(solutionDict.getOrDefault("coupledEnergyField", false))
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::IOobject::NO_READ
Definition: IOobject.H:188
Foam::kineticTheoryModels::frictionalStressModels::Schaeffer::nu
virtual tmp< volScalarField > nu(const phaseModel &phase, const dimensionedScalar &alphaMinFriction, const dimensionedScalar &alphaMax, const volScalarField &pf, const volSymmTensorField &D) const
Definition: SchaefferFrictionalStress.C:113
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
Foam::dimless
const dimensionSet dimless
Dimensionless.
Definition: dimensionSets.C:189