FreeFOAM The Cross-Platform CFD Toolkit
HerschelBulkley.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 | Copyright (C) 1991-2011 OpenCFD Ltd.
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 "HerschelBulkley.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace viscosityModels
35 {
36  defineTypeNameAndDebug(HerschelBulkley, 0);
37 
39  (
40  viscosityModel,
41  HerschelBulkley,
42  dictionary
43  );
44 }
45 }
46 
47 
48 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
49 
51 Foam::viscosityModels::HerschelBulkley::calcNu() const
52 {
53  dimensionedScalar tone("tone", dimTime, 1.0);
54  dimensionedScalar rtone("rtone", dimless/dimTime, 1.0);
55 
56  tmp<volScalarField> sr(strainRate());
57 
58  // return
59  // (
60  // min
61  // (
62  // nu0_,
63  // (tau0_ + k_*rtone*(pow(tone*sr(), n_) - pow(tone*tau0_/nu0_, n_)))
64  // /max(sr(), dimensionedScalar("VSMALL", dimless/dimTime, VSMALL))
65  // )
66  // );
67 
68  return
69  (
70  min
71  (
72  nu0_,
73  (tau0_ + k_*rtone*pow(tone*sr(), n_))
74  /(max(sr(), dimensionedScalar ("VSMALL", dimless/dimTime, VSMALL)))
75  )
76  );
77 }
78 
79 
80 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
81 
83 (
84  const word& name,
85  const dictionary& viscosityProperties,
86  const volVectorField& U,
87  const surfaceScalarField& phi
88 )
89 :
90  viscosityModel(name, viscosityProperties, U, phi),
91  HerschelBulkleyCoeffs_(viscosityProperties.subDict(typeName + "Coeffs")),
92  k_(HerschelBulkleyCoeffs_.lookup("k")),
93  n_(HerschelBulkleyCoeffs_.lookup("n")),
94  tau0_(HerschelBulkleyCoeffs_.lookup("tau0")),
95  nu0_(HerschelBulkleyCoeffs_.lookup("nu0")),
96  nu_
97  (
98  IOobject
99  (
100  name,
101  U_.time().timeName(),
102  U_.db(),
105  ),
106  calcNu()
107  )
108 {}
109 
110 
111 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
112 
114 (
115  const dictionary& viscosityProperties
116 )
117 {
118  viscosityModel::read(viscosityProperties);
119 
120  HerschelBulkleyCoeffs_ = viscosityProperties.subDict(typeName + "Coeffs");
121 
122  HerschelBulkleyCoeffs_.lookup("k") >> k_;
123  HerschelBulkleyCoeffs_.lookup("n") >> n_;
124  HerschelBulkleyCoeffs_.lookup("tau0") >> tau0_;
125  HerschelBulkleyCoeffs_.lookup("nu0") >> nu0_;
126 
127  return true;
128 }
129 
130 
131 // ************************ vim: set sw=4 sts=4 et: ************************ //