00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073 #include <postCalc/calc.H>
00074 #include <finiteVolume/fvc.H>
00075
00076 #include <incompressibleTransportModels/singlePhaseTransportModel.H>
00077 #include <incompressibleRASModels/RASModel.H>
00078 #include <incompressibleLESModels/LESModel.H>
00079 #include <basicThermophysicalModels/basicPsiThermo.H>
00080 #include <compressibleRASModels/RASModel.H>
00081 #include <compressibleLESModels/LESModel.H>
00082
00083
00084
00085 void Foam::calc(const argList& args, const Time& runTime, const fvMesh& mesh)
00086 {
00087 bool writeResults = !args.optionFound("noWrite");
00088
00089 IOobject phiHeader
00090 (
00091 "phi",
00092 runTime.timeName(),
00093 mesh,
00094 IOobject::MUST_READ
00095 );
00096
00097 if (phiHeader.headerOk())
00098 {
00099 autoPtr<surfaceScalarField> PePtr;
00100
00101 Info<< " Reading phi" << endl;
00102 surfaceScalarField phi(phiHeader, mesh);
00103
00104 volVectorField U
00105 (
00106 IOobject
00107 (
00108 "U",
00109 runTime.timeName(),
00110 mesh,
00111 IOobject::MUST_READ
00112 ),
00113 mesh
00114 );
00115
00116 IOobject RASPropertiesHeader
00117 (
00118 "RASProperties",
00119 runTime.constant(),
00120 mesh,
00121 IOobject::MUST_READ,
00122 IOobject::NO_WRITE
00123 );
00124
00125 IOobject LESPropertiesHeader
00126 (
00127 "LESProperties",
00128 runTime.constant(),
00129 mesh,
00130 IOobject::MUST_READ,
00131 IOobject::NO_WRITE
00132 );
00133
00134 Info<< " Calculating Pe" << endl;
00135
00136 if (phi.dimensions() == dimensionSet(0, 3, -1, 0, 0))
00137 {
00138 if (RASPropertiesHeader.headerOk())
00139 {
00140 IOdictionary RASProperties(RASPropertiesHeader);
00141
00142 singlePhaseTransportModel laminarTransport(U, phi);
00143
00144 autoPtr<incompressible::RASModel> RASModel
00145 (
00146 incompressible::RASModel::New
00147 (
00148 U,
00149 phi,
00150 laminarTransport
00151 )
00152 );
00153
00154 PePtr.set
00155 (
00156 new surfaceScalarField
00157 (
00158 IOobject
00159 (
00160 "Pe",
00161 runTime.timeName(),
00162 mesh,
00163 IOobject::NO_READ
00164 ),
00165 mag(phi)
00166 /(
00167 mesh.magSf()
00168 * mesh.surfaceInterpolation::deltaCoeffs()
00169 * fvc::interpolate(RASModel->nuEff())
00170 )
00171 )
00172 );
00173 }
00174 else if (LESPropertiesHeader.headerOk())
00175 {
00176 IOdictionary LESProperties(LESPropertiesHeader);
00177
00178 singlePhaseTransportModel laminarTransport(U, phi);
00179
00180 autoPtr<incompressible::LESModel> sgsModel
00181 (
00182 incompressible::LESModel::New(U, phi, laminarTransport)
00183 );
00184
00185 PePtr.set
00186 (
00187 new surfaceScalarField
00188 (
00189 IOobject
00190 (
00191 "Pe",
00192 runTime.timeName(),
00193 mesh,
00194 IOobject::NO_READ
00195 ),
00196 mag(phi)
00197 /(
00198 mesh.magSf()
00199 * mesh.surfaceInterpolation::deltaCoeffs()
00200 * fvc::interpolate(sgsModel->nuEff())
00201 )
00202 )
00203 );
00204 }
00205 else
00206 {
00207 IOdictionary transportProperties
00208 (
00209 IOobject
00210 (
00211 "transportProperties",
00212 runTime.constant(),
00213 mesh,
00214 IOobject::MUST_READ,
00215 IOobject::NO_WRITE
00216 )
00217 );
00218
00219 dimensionedScalar nu(transportProperties.lookup("nu"));
00220
00221 PePtr.set
00222 (
00223 new surfaceScalarField
00224 (
00225 IOobject
00226 (
00227 "Pe",
00228 runTime.timeName(),
00229 mesh,
00230 IOobject::NO_READ
00231 ),
00232 mesh.surfaceInterpolation::deltaCoeffs()
00233 * (mag(phi)/mesh.magSf())*(runTime.deltaT()/nu)
00234 )
00235 );
00236 }
00237 }
00238 else if (phi.dimensions() == dimensionSet(1, 0, -1, 0, 0))
00239 {
00240 if (RASPropertiesHeader.headerOk())
00241 {
00242 IOdictionary RASProperties(RASPropertiesHeader);
00243
00244 autoPtr<basicPsiThermo> thermo(basicPsiThermo::New(mesh));
00245
00246 volScalarField rho
00247 (
00248 IOobject
00249 (
00250 "rho",
00251 runTime.timeName(),
00252 mesh
00253 ),
00254 thermo->rho()
00255 );
00256
00257 autoPtr<compressible::RASModel> RASModel
00258 (
00259 compressible::RASModel::New
00260 (
00261 rho,
00262 U,
00263 phi,
00264 thermo()
00265 )
00266 );
00267
00268 PePtr.set
00269 (
00270 new surfaceScalarField
00271 (
00272 IOobject
00273 (
00274 "Pe",
00275 runTime.timeName(),
00276 mesh,
00277 IOobject::NO_READ
00278 ),
00279 mag(phi)
00280 /(
00281 mesh.magSf()
00282 * mesh.surfaceInterpolation::deltaCoeffs()
00283 * fvc::interpolate(RASModel->muEff())
00284 )
00285 )
00286 );
00287 }
00288 else if (LESPropertiesHeader.headerOk())
00289 {
00290 IOdictionary LESProperties(LESPropertiesHeader);
00291
00292 autoPtr<basicPsiThermo> thermo(basicPsiThermo::New(mesh));
00293
00294 volScalarField rho
00295 (
00296 IOobject
00297 (
00298 "rho",
00299 runTime.timeName(),
00300 mesh
00301 ),
00302 thermo->rho()
00303 );
00304
00305 autoPtr<compressible::LESModel> sgsModel
00306 (
00307 compressible::LESModel::New(rho, U, phi, thermo())
00308 );
00309
00310 PePtr.set
00311 (
00312 new surfaceScalarField
00313 (
00314 IOobject
00315 (
00316 "Pe",
00317 runTime.timeName(),
00318 mesh,
00319 IOobject::NO_READ
00320 ),
00321 mag(phi)
00322 /(
00323 mesh.magSf()
00324 * mesh.surfaceInterpolation::deltaCoeffs()
00325 * fvc::interpolate(sgsModel->muEff())
00326 )
00327 )
00328 );
00329 }
00330 else
00331 {
00332 IOdictionary transportProperties
00333 (
00334 IOobject
00335 (
00336 "transportProperties",
00337 runTime.constant(),
00338 mesh,
00339 IOobject::MUST_READ,
00340 IOobject::NO_WRITE
00341 )
00342 );
00343
00344 dimensionedScalar mu(transportProperties.lookup("mu"));
00345
00346 PePtr.set
00347 (
00348 new surfaceScalarField
00349 (
00350 IOobject
00351 (
00352 "Pe",
00353 runTime.timeName(),
00354 mesh,
00355 IOobject::NO_READ
00356 ),
00357 mesh.surfaceInterpolation::deltaCoeffs()
00358 * (mag(phi)/(mesh.magSf()))*(runTime.deltaT()/mu)
00359 )
00360 );
00361 }
00362 }
00363 else
00364 {
00365 FatalErrorIn(args.executable())
00366 << "Incorrect dimensions of phi: " << phi.dimensions()
00367 << abort(FatalError);
00368 }
00369
00370
00371
00372
00373
00374
00375
00376
00377
00378
00379
00380
00381
00382
00383
00384
00385
00386
00387
00388
00389
00390 Info << "Pe max : " << max(PePtr()).value() << endl;
00391
00392 if (writeResults)
00393 {
00394 PePtr().write();
00395 }
00396 }
00397 else
00398 {
00399 Info<< " No phi" << endl;
00400 }
00401
00402 Info<< "\nEnd\n" << endl;
00403 }
00404
00405
00406