initial commit for version 1.5.x patch release
[OpenFOAM-1.5.x.git] / applications / solvers / multiphase / twoPhaseEulerFoam / kineticTheoryModels / frictionalStressModel / Schaeffer / SchaefferFrictionalStress.C
blob59ced9ee63d772444d29c1d646409bd8a5317356
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
4    \\    /   O peration     |
5     \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
6      \\/     M anipulation  |
7 -------------------------------------------------------------------------------
8 License
9     This file is part of OpenFOAM.
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 the
13     Free Software Foundation; either version 2 of the License, or (at your
14     option) any later version.
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.
21     You should have received a copy of the GNU General Public License
22     along with OpenFOAM; if not, write to the Free Software Foundation,
23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
25 \*---------------------------------------------------------------------------*/
27 #include "SchaefferFrictionalStress.H"
28 #include "addToRunTimeSelectionTable.H"
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 namespace Foam
34     defineTypeNameAndDebug(SchaefferFrictionalStress, 0);
36     addToRunTimeSelectionTable
37     (
38         frictionalStressModel,
39         SchaefferFrictionalStress,
40         dictionary
41     );
45 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
47 Foam::SchaefferFrictionalStress::SchaefferFrictionalStress
49     const dictionary& dict
52     frictionalStressModel(dict)
56 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
58 Foam::SchaefferFrictionalStress::~SchaefferFrictionalStress()
62 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
64 Foam::tmp<Foam::volScalarField> Foam::SchaefferFrictionalStress::
65 frictionalPressure
67     const volScalarField& alpha,
68     const dimensionedScalar& alphaMinFriction,
69     const dimensionedScalar& alphaMax,
70     const dimensionedScalar& Fr,
71     const dimensionedScalar& eta,
72     const dimensionedScalar& p
73 ) const
75     return
76         dimensionedScalar("1e24", dimensionSet(1, -1, -2, 0, 0), 1e24)
77        *pow(Foam::max(alpha - alphaMinFriction, scalar(0)), 10.0);
81 Foam::tmp<Foam::volScalarField> Foam::SchaefferFrictionalStress::
82 frictionalPressurePrime
84     const volScalarField& alpha,
85     const dimensionedScalar& alphaMinFriction,
86     const dimensionedScalar& alphaMax,
87     const dimensionedScalar& Fr,
88     const dimensionedScalar& eta,
89     const dimensionedScalar& p
90 ) const
92     return
93         dimensionedScalar("1e25", dimensionSet(1, -1, -2, 0, 0), 1e25)
94        *pow(Foam::max(alpha - alphaMinFriction, scalar(0)), 9.0);
98 Foam::tmp<Foam::volScalarField> Foam::SchaefferFrictionalStress::muf
100     const volScalarField& alpha,
101     const dimensionedScalar& alphaMax,
102     const volScalarField& pf,
103     const volTensorField& D,
104     const dimensionedScalar& phi
105 ) const
107     const scalar I2Dsmall = 1.0e-15;
109     // Creating muf assuming it should be 0 on the boundary which may not be
110     // true
111     tmp<volScalarField> tmuf
112     (
113         new volScalarField
114         (
115             IOobject
116             (
117                 "muf",
118                 alpha.mesh().time().timeName(),
119                 alpha.mesh()
120             ),
121             alpha.mesh(),
122             dimensionedScalar("muf", dimensionSet(1, -1, -1, 0, 0), 0.0)
123         )
124     );
126     volScalarField& muff = tmuf();
128     forAll (D, celli)
129     {
130         if (alpha[celli] > alphaMax.value()-5e-2)
131         {
132             muff[celli] =
133                 0.5*pf[celli]*sin(phi.value())
134                /(
135                     sqrt(1.0/6.0*(sqr(D[celli].xx() - D[celli].yy())
136                   + sqr(D[celli].yy() - D[celli].zz())
137                   + sqr(D[celli].zz() - D[celli].xx()))
138                   + sqr(D[celli].xy()) + sqr(D[celli].xz())
139                   + sqr(D[celli].yz())) + I2Dsmall
140                 );
141         }
142     }
144     return tmuf;
148 // ************************************************************************* //