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-2018 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  filmViscosityModel,
51  thixotropicViscosity,
52  dictionary
53 );
54 
55 
56 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
57 
58 thixotropicViscosity::thixotropicViscosity
59 (
61  const dictionary& dict,
62  volScalarField& mu
63 )
64 :
65  filmViscosityModel(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  K_(1 - sqrt(muInf_/mu0_)),
73  lambda_
74  (
75  IOobject
76  (
77  typeName + ":lambda",
78  film.regionMesh().time().timeName(),
79  film.regionMesh(),
82  ),
83  film.regionMesh()
84  )
85 {
86  lambda_.min(1);
87  lambda_.max(0);
88 
89  // Initialise viscosity to inf value because it cannot be evaluated yet
90  mu_ = muInf_;
91  mu_.correctBoundaryConditions();
92 }
93 
94 
95 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
96 
98 {}
99 
100 
101 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
102 
104 (
105  const volScalarField& p,
106  const volScalarField& T
107 )
108 {
109  const kinematicSingleLayer& film = filmType<kinematicSingleLayer>();
110 
111  const volVectorField& U = film.U();
112  const volVectorField& Uw = film.Uw();
113  const volScalarField& delta = film.delta();
114  const volScalarField& deltaRho = film.deltaRho();
115  const surfaceScalarField& phi = film.phi();
116  const volScalarField& alpha = film.alpha();
117  const Time& runTime = this->film().regionMesh().time();
118 
119  // Shear rate
120  const volScalarField gDot
121  (
122  "gDot",
123  alpha*mag(U - Uw)/(delta + film.deltaSmall())
124  );
125 
126  if (debug && runTime.writeTime())
127  {
128  gDot.write();
129  }
130 
131  const dimensionedScalar deltaRho0
132  (
133  "deltaRho0",
134  deltaRho.dimensions(),
135  rootVSmall
136  );
137 
138  const surfaceScalarField phiU(phi/fvc::interpolate(deltaRho + deltaRho0));
139 
140  const dimensionedScalar c0("c0", dimless/dimTime, rootVSmall);
141  const volScalarField coeff("coeff", -c_*pow(gDot, d_) + c0);
142 
143  fvScalarMatrix lambdaEqn
144  (
146  + fvm::div(phiU, lambda_)
147  - fvm::Sp(fvc::div(phiU), lambda_)
148  ==
149  a_*pow((1 - lambda_), b_)
150  + fvm::SuSp(coeff, lambda_)
151 
152  // Include the effect of the impinging droplets added with lambda = 0
153  - fvm::Sp
154  (
155  max
156  (
157  -film.rhoSp(),
158  dimensionedScalar("zero", film.rhoSp().dimensions(), 0)
159  )/(deltaRho + deltaRho0),
160  lambda_
161  )
162  );
163 
164  lambdaEqn.relax();
165  lambdaEqn.solve();
166 
167  lambda_.min(1);
168  lambda_.max(0);
169 
170  mu_ = muInf_/(sqr(1 - K_*lambda_) + rootVSmall);
172 }
173 
174 
175 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
176 
177 } // End namespace surfaceFilmModels
178 } // End namespace regionModels
179 } // End namespace Foam
180 
181 // ************************************************************************* //
scalar delta
virtual const volVectorField & U() const
Return the film velocity [m/s].
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
virtual const volScalarField & deltaRho() const
Return the film thickness*density (helper field) [kg/m3].
surfaceScalarField & phi
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:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
engineTime & runTime
dimensionedScalar sqrt(const dimensionedScalar &ds)
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh)
Macros for easy insertion into run-time selection tables.
const surfaceFilmRegionModel & film() const
Return const access to the film surface film model.
bool writeTime() const
Return true if this is a write time.
Definition: TimeStateI.H:65
const dimensionSet & dimensions() const
Return dimensions.
Base class for surface film viscosity models.
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.
virtual const surfaceScalarField & phi() const
Return the film flux [kg.m/s].
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:61
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
void min(const dimensioned< Type > &)
virtual void correct(const volScalarField &p, const volScalarField &T)
Correct.
const dimensionSet dimPressure
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:523
volScalarField & mu_
Reference to the viscosity field.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const volScalarField & delta() const
Return const access to the film thickness [m].
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Calculate the matrix for the divergence of the given field and flux.
U
Definition: pEqn.H:72
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
virtual const volVectorField & Uw() const
Return the film wall velocity [m/s].
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
void max(const dimensioned< Type > &)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
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 > &)
volScalarField & p
virtual bool write(const bool valid=true) const
Write using setting from DB.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
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.
const volScalarField & alpha() const
Return the film coverage, 1 = covered, 0 = uncovered [].
Namespace for OpenFOAM.