Merge branch 'master' of github.com:OpenFOAM/OpenFOAM-2.2.x
[OpenFOAM-2.2.x.git] / src / postProcessing / functionObjects / utilities / Peclet / Peclet.C
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
4    \\    /   O peration     |
5     \\  /    A nd           | Copyright (C) 2013 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 "Peclet.H"
27 #include "volFields.H"
28 #include "dictionary.H"
29 #include "surfaceFields.H"
30 #include "incompressible/turbulenceModel/turbulenceModel.H"
31 #include "compressible/turbulenceModel/turbulenceModel.H"
32 #include "surfaceInterpolate.H"
33
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35
36 namespace Foam
37 {
38 defineTypeNameAndDebug(Peclet, 0);
39 }
40
41
42 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
43
44 Foam::Peclet::Peclet
45 (
46     const word& name,
47     const objectRegistry& obr,
48     const dictionary& dict,
49     const bool loadFromFiles
50 )
51 :
52     name_(name),
53     obr_(obr),
54     active_(true),
55     phiName_("phi"),
56     rhoName_("rho")
57 {
58     // Check if the available mesh is an fvMesh, otherwise deactivate
59     if (!isA<fvMesh>(obr_))
60     {
61         active_ = false;
62         WarningIn
63         (
64             "Peclet::Peclet"
65             "("
66                 "const word&, "
67                 "const objectRegistry&, "
68                 "const dictionary&, "
69                 "const bool"
70             ")"
71         )   << "No fvMesh available, deactivating " << name_ << nl
72             << endl;
73     }
74
75     read(dict);
76
77     if (active_)
78     {
79         const fvMesh& mesh = refCast<const fvMesh>(obr_);
80
81         surfaceScalarField* PecletPtr
82         (
83             new surfaceScalarField
84             (
85                 IOobject
86                 (
87                     type(),
88                     mesh.time().timeName(),
89                     mesh,
90                     IOobject::NO_READ,
91                     IOobject::NO_WRITE
92                 ),
93                 mesh,
94                 dimensionedScalar("0", dimless, 0.0)
95             )
96         );
97
98         mesh.objectRegistry::store(PecletPtr);
99     }
100 }
101
102
103 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
104
105 Foam::Peclet::~Peclet()
106 {}
107
108
109 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
110
111 void Foam::Peclet::read(const dictionary& dict)
112 {
113     if (active_)
114     {
115         phiName_ = dict.lookupOrDefault<word>("phiName", "phi");
116         rhoName_ = dict.lookupOrDefault<word>("rhoName", "rho");
117     }
118 }
119
120
121 void Foam::Peclet::execute()
122 {
123     // Do nothing - only valid on write
124 }
125
126
127 void Foam::Peclet::end()
128 {
129     // Do nothing - only valid on write
130 }
131
132 void Foam::Peclet::timeSet()
133 {
134     // Do nothing - only valid on write
135 }
136
137
138 void Foam::Peclet::write()
139 {
140     typedef compressible::turbulenceModel cmpTurbModel;
141     typedef incompressible::turbulenceModel icoTurbModel;
142
143     if (active_)
144     {
145         const fvMesh& mesh = refCast<const fvMesh>(obr_);
146
147         tmp<volScalarField> nuEff;
148         if (mesh.foundObject<cmpTurbModel>("turbulenceModel"))
149         {
150             const cmpTurbModel& model =
151                 mesh.lookupObject<cmpTurbModel>("turbulenceModel");
152
153             const volScalarField& rho =
154                 mesh.lookupObject<volScalarField>(rhoName_);
155
156             nuEff = model.muEff()/rho;
157         }
158         else if (mesh.foundObject<icoTurbModel>("turbulenceModel"))
159         {
160             const icoTurbModel& model =
161                 mesh.lookupObject<icoTurbModel>("turbulenceModel");
162
163             nuEff = model.nuEff();
164         }
165         else if (mesh.foundObject<dictionary>("transportProperties"))
166         {
167             const dictionary& model =
168                 mesh.lookupObject<dictionary>("transportProperties");
169
170             nuEff =
171                 tmp<volScalarField>
172                 (
173                     new volScalarField
174                     (
175                         IOobject
176                         (
177                             "nuEff",
178                             mesh.time().timeName(),
179                             mesh,
180                             IOobject::NO_READ,
181                             IOobject::NO_WRITE
182                         ),
183                         mesh,
184                         dimensionedScalar(model.lookup("nu"))
185                     )
186                 );
187         }
188         else
189         {
190             FatalErrorIn("void Foam::Peclet::write()")
191                 << "Unable to determine the viscosity"
192                 << exit(FatalError);
193         }
194
195         const surfaceScalarField& phi =
196             mesh.lookupObject<surfaceScalarField>(phiName_);
197
198         surfaceScalarField& Peclet =
199             const_cast<surfaceScalarField&>
200             (
201                 mesh.lookupObject<surfaceScalarField>(type())
202             );
203
204         Peclet =
205             mag(phi)
206            /(
207                 mesh.magSf()
208                *mesh.surfaceInterpolation::deltaCoeffs()
209                *fvc::interpolate(nuEff)
210             );
211
212         Info<< type() << " " << name_ << " output:" << nl
213             << "    writing field " << Peclet.name() << nl
214             << endl;
215
216         Peclet.write();
217     }
218 }
219
220
221 // ************************************************************************* //