thixotropicViscosity.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) 2013-2021 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 "thixotropicViscosity.H"
27 #include "kinematicSingleLayer.H"
29 
30 #include "fvmDdt.H"
31 #include "fvmDiv.H"
32 #include "fvcDiv.H"
33 #include "fvmSup.H"
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 namespace regionModels
40 {
41 namespace surfaceFilmModels
42 {
43 
44 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
45 
46 defineTypeNameAndDebug(thixotropicViscosity, 0);
47 
49 (
50  viscosityModel,
51  thixotropicViscosity,
52  dictionary
53 );
54 
55 
56 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
57 
59 (
61  const dictionary& dict,
62  volScalarField& mu
63 )
64 :
65  viscosityModel(typeName, film, dict, mu),
66  a_("a", dimless/dimTime, coeffDict_),
67  b_("b", dimless, coeffDict_),
68  d_("d", dimless, coeffDict_),
69  c_("c", pow(dimTime, d_.value() - scalar(1)), coeffDict_),
70  mu0_("mu0", dimPressure*dimTime, coeffDict_),
71  muInf_("muInf", mu0_.dimensions(), coeffDict_),
72  BinghamPlastic_(coeffDict_.found("tauy")),
73  tauy_
74  (
75  BinghamPlastic_
76  ? dimensionedScalar("tauy", dimPressure, coeffDict_)
77  : dimensionedScalar("tauy", dimPressure, 0)
78  ),
79  K_(1 - sqrt(muInf_/mu0_)),
80  lambda_
81  (
82  IOobject
83  (
84  IOobject::modelName("lambda", typeName),
85  film.regionMesh().time().timeName(),
86  film.regionMesh(),
87  IOobject::MUST_READ,
88  IOobject::AUTO_WRITE
89  ),
90  film.regionMesh()
91  )
92 {
93  lambda_.min(1);
94  lambda_.max(0);
95 
96  // Initialise viscosity to inf value because it cannot be evaluated yet
97  mu_ = muInf_;
98  mu_.correctBoundaryConditions();
99 }
100 
101 
102 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
103 
105 {}
106 
107 
108 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
109 
111 (
112  const volScalarField& p,
113  const volScalarField& T
114 )
115 {
116  const kinematicSingleLayer& film = filmType<kinematicSingleLayer>();
117 
118  const volScalarField::Internal& delta = film.delta();
119  const volScalarField::Internal alphaRho(film.alpha()()*film.rho()());
120  const volScalarField::Internal& coverage = film.coverage();
121  const surfaceScalarField& phiU = film.phiU();
122 
123  // Shear rate
124  const volScalarField::Internal gDot
125  (
126  "gDot",
127  coverage*mag(film.Us() - film.Uw())/(delta + film.deltaSmall())
128  );
129 
130  const dimensionedScalar alphaRho0
131  (
132  "alphaRho0",
133  alphaRho.dimensions(),
134  rootVSmall
135  );
136 
137  fvScalarMatrix lambdaEqn
138  (
140  + fvm::div(phiU, lambda_)
141  - fvm::Sp(fvc::div(phiU), lambda_)
142  ==
143  a_*pow((1 - lambda_), b_)
144  - fvm::Sp(c_*pow(gDot, d_), lambda_)
145 
146  // Include the effect of the impinging droplets added with lambda = 0
147  - fvm::Sp
148  (
149  max
150  (
151  -film.rhoSp(),
152  dimensionedScalar(film.rhoSp().dimensions(), 0)
153  )/(alphaRho + alphaRho0),
154  lambda_
155  )
156  );
157 
158  lambdaEqn.relax();
159  lambdaEqn.solve();
160 
161  lambda_.min(1);
162  lambda_.max(0);
163 
164  mu_ = muInf_/(sqr(1 - K_*lambda_) + rootVSmall);
165 
166  // Add optional yield stress contribution to the viscosity
167  if (BinghamPlastic_)
168  {
169  dimensionedScalar tauySmall("tauySmall", tauy_.dimensions(), small);
170  dimensionedScalar muMax_("muMax", 100*mu0_);
171 
172  mu_.ref() = min
173  (
174  tauy_/(gDot + 1.0e-4*(tauy_ + tauySmall)/mu0_) + mu_(),
175  muMax_
176  );
177  }
178 
180 }
181 
182 
183 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
184 
185 } // End namespace surfaceFilmModels
186 } // End namespace regionModels
187 } // End namespace Foam
188 
189 // ************************************************************************* //
scalar delta
Base class for surface film viscosity models.
Kinematic form of single-cell layer surface film model.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
volScalarField & mu_
Reference to the viscosity field.
thixotropicViscosity(surfaceFilmRegionModel &film, const dictionary &dict, volScalarField &mu)
Construct from surface film model.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
const dimensionSet dimPressure
dimensionedScalar sqrt(const dimensionedScalar &ds)
const dimensionSet dimless
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
dimensionedScalar mu0_
Limiting viscosity when lambda = 1.
tmp< volVectorField::Internal > Us() const
Return the film surface velocity [m/s].
volScalarField::Internal & rhoSp()
Mass [kg/m^2/s].
Macros for easy insertion into run-time selection tables.
tmp< volVectorField::Internal > Uw() const
Return the film wall velocity [m/s].
const dimensionSet dimTime
const dimensionSet & dimensions() const
Return dimensions.
const surfaceScalarField & phiU() const
Return the film velocity flux [m^3/s].
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
Calculate the matrix for the first temporal derivative.
const dimensionedScalar & deltaSmall() const
Return small delta.
const volScalarField & coverage() const
Return the film coverage, 1 = covered, 0 = uncovered [].
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
word timeName
Definition: getTimeIndex.H:3
void min(const dimensioned< Type > &)
virtual void correct(const volScalarField &p, const volScalarField &T)
Correct.
dimensionedScalar muInf_
Limiting viscosity when lambda = 0.
Calculate the divergence of the given field.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:521
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const volScalarField & delta() const
Return const access to the film thickness [m].
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
Internal & ref()
Return a reference to the dimensioned internal field.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Calculate the matrix for the divergence of the given field and flux.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
void max(const dimensioned< Type > &)
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet & dimensions() const
Return const reference to dimensions.
void correctBoundaryConditions()
Correct boundary field.
dimensionedScalar c_
Model `c&#39; coefficient (read after d since dims depend on d value)
dimensioned< scalar > mag(const dimensioned< Type > &)
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
const volScalarField & rho() const
Return the film density [kg/m^3].
volScalarField & p
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
defineTypeNameAndDebug(kinematicSingleLayer, 0)
Calculate the matrix for implicit and explicit sources.
addToRunTimeSelectionTable(surfaceFilmModel, noFilm, mesh)
const volScalarField & alpha() const
Return const access to the film volume fraction [].
Namespace for OpenFOAM.