JohnsonJacksonFrictionalStress.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 {
42 
44  (
48  );
49 }
50 }
51 }
52 
53 
54 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
55 
58 (
59  const dictionary& dict
60 )
61 :
63  coeffDict_(dict.optionalSubDict(typeName + "Coeffs")),
64  Fr_("Fr", dimensionSet(1, -1, -2, 0, 0), coeffDict_),
65  eta_("eta", dimless, coeffDict_),
66  p_("p", dimless, coeffDict_),
67  phi_("phi", dimless, coeffDict_),
68  alphaDeltaMin_("alphaDeltaMin", dimless, coeffDict_)
69 {
70  phi_ *= degToRad();
71 }
72 
73 
74 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
75 
78 {}
79 
80 
81 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
82 
86 (
87  const phaseModel& phase,
88  const dimensionedScalar& alphaMinFriction,
90 ) const
91 {
92  const volScalarField& alpha = phase;
93 
94  return
95  Fr_*pow(max(alpha - alphaMinFriction, scalar(0)), eta_)
96  /pow(max(alphaMax - alpha, alphaDeltaMin_), p_);
97 }
98 
99 
103 (
104  const phaseModel& phase,
105  const dimensionedScalar& alphaMinFriction,
107 ) const
108 {
109  const volScalarField& alpha = phase;
110 
111  return Fr_*
112  (
113  eta_*pow(max(alpha - alphaMinFriction, scalar(0)), eta_ - 1)
114  *(alphaMax - alpha)
115  + p_*pow(max(alpha - alphaMinFriction, scalar(0)), eta_)
116  )/pow(max(alphaMax - alpha, alphaDeltaMin_), p_ + 1);
117 }
118 
119 
122 (
123  const phaseModel& phase,
124  const dimensionedScalar& alphaMinFriction,
126  const volScalarField& pf,
127  const volSymmTensorField& D
128 ) const
129 {
130  return dimensionedScalar("1/2",dimTime, 0.5)*pf*sin(phi_);
131 }
132 
133 
135 {
136  coeffDict_ <<= dict_.optionalSubDict(typeName + "Coeffs");
137 
138  Fr_.read(coeffDict_);
139  eta_.read(coeffDict_);
140  p_.read(coeffDict_);
141 
142  phi_.read(coeffDict_);
143  phi_ *= degToRad();
144 
145  alphaDeltaMin_.read(coeffDict_);
146 
147  return true;
148 }
149 
150 
151 // ************************************************************************* //
Foam::kineticTheoryModels::frictionalStressModel
Definition: frictionalStressModel.H:54
Foam::phaseModel
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:54
Foam::kineticTheoryModels::frictionalStressModels::JohnsonJackson::frictionalPressure
virtual tmp< volScalarField > frictionalPressure(const phaseModel &phase, const dimensionedScalar &alphaMinFriction, const dimensionedScalar &alphaMax) const
Definition: JohnsonJacksonFrictionalStress.C:86
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
unitConversion.H
Unit conversion functions.
Foam::dimensionSet
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
Definition: dimensionSet.H:108
Foam::kineticTheoryModels::frictionalStressModels::JohnsonJackson::read
virtual bool read()
Definition: JohnsonJacksonFrictionalStress.C:134
Foam::dimTime
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
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)
JohnsonJacksonFrictionalStress.H
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::dimensioned< scalar >
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::kineticTheoryModels::frictionalStressModels::JohnsonJackson::~JohnsonJackson
virtual ~JohnsonJackson()
Destructor.
Definition: JohnsonJacksonFrictionalStress.C:77
Foam::degToRad
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
Definition: unitConversion.H:48
Foam::kineticTheoryModels::frictionalStressModels::JohnsonJackson::frictionalPressurePrime
virtual tmp< volScalarField > frictionalPressurePrime(const phaseModel &phase, const dimensionedScalar &alphaMinFriction, const dimensionedScalar &alphaMax) const
Definition: JohnsonJacksonFrictionalStress.C:103
Foam::kineticTheoryModels::frictionalStressModels::JohnsonJackson::JohnsonJackson
JohnsonJackson(const dictionary &dict)
Construct from components.
Definition: JohnsonJacksonFrictionalStress.C:58
D
const dimensionedScalar & D
Definition: solveBulkSurfactant.H:4
Foam::kineticTheoryModels::frictionalStressModels::JohnsonJackson
Definition: JohnsonJacksonFrictionalStress.H:53
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::kineticTheoryModels::frictionalStressModels::JohnsonJackson::nu
virtual tmp< volScalarField > nu(const phaseModel &phase, const dimensionedScalar &alphaMinFriction, const dimensionedScalar &alphaMax, const volScalarField &pf, const volSymmTensorField &D) const
Definition: JohnsonJacksonFrictionalStress.C:122
Foam::dimless
const dimensionSet dimless
Dimensionless.
Definition: dimensionSets.C:189