FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

Pe.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------* \
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 Application
00025     Pe
00026 
00027 Description
00028     Calculates and writes the Pe number as a surfaceScalarField obtained from
00029     field phi.
00030 
00031     The -noWrite option just outputs the max/min values without writing
00032     the field.
00033 
00034 Usage
00035 
00036     - Pe [OPTIONS]
00037 
00038     @param -noWrite \n
00039     Suppress output to files.
00040 
00041     @param -dict <dictionary name>\n
00042     Use named dictionary instead of system/controlDict.
00043 
00044     @param -noZero \n
00045     Ignore timestep 0.
00046 
00047     @param -constant \n
00048     Include the constant directory.
00049 
00050     @param -time <time>\n
00051     Apply only to specific time.
00052 
00053     @param -latestTime \n
00054     Only apply to latest time step.
00055 
00056     @param -case <dir>\n
00057     Case directory.
00058 
00059     @param -parallel \n
00060     Run in parallel.
00061 
00062     @param -help \n
00063     Display help message.
00064 
00065     @param -doc \n
00066     Display Doxygen API documentation page for this application.
00067 
00068     @param -srcDoc \n
00069     Display Doxygen source documentation page for this application.
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         // can also check how many cells exceed a particular Pe limit
00372         /*
00373         {
00374             label count = 0;
00375             label PeLimit = 200;
00376             forAll(PePtr(), i)
00377             {
00378                 if (PePtr()[i] > PeLimit)
00379                 {
00380                     count++;
00381                 }
00382 
00383             }
00384 
00385             Info<< "Fraction > " << PeLimit << " = "
00386                 << scalar(count)/Pe.size() << endl;
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 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines