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

Q.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     Q
00026 
00027 Description
00028     Calculates and writes the second invariant of the velocity gradient tensor.
00029 
00030     The -noWrite option just outputs the max/min values without writing
00031     the field.
00032 
00033 Usage
00034 
00035     - Q [OPTIONS]
00036 
00037     @param -noWrite \n
00038     Suppress output to files.
00039 
00040     @param -dict <dictionary name>\n
00041     Use named dictionary instead of system/controlDict.
00042 
00043     @param -noZero \n
00044     Ignore timestep 0.
00045 
00046     @param -constant \n
00047     Include the constant directory.
00048 
00049     @param -time <time>\n
00050     Apply only to specific time.
00051 
00052     @param -latestTime \n
00053     Only apply to latest time step.
00054 
00055     @param -case <dir>\n
00056     Case directory.
00057 
00058     @param -parallel \n
00059     Run in parallel.
00060 
00061     @param -help \n
00062     Display help message.
00063 
00064     @param -doc \n
00065     Display Doxygen API documentation page for this application.
00066 
00067     @param -srcDoc \n
00068     Display Doxygen source documentation page for this application.
00069 
00070 \*---------------------------------------------------------------------------*/
00071 
00072 #include <postCalc/calc.H>
00073 #include <finiteVolume/fvc.H>
00074 
00075 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00076 
00077 void Foam::calc(const argList& args, const Time& runTime, const fvMesh& mesh)
00078 {
00079     bool writeResults = !args.optionFound("noWrite");
00080 
00081     IOobject Uheader
00082     (
00083         "U",
00084         runTime.timeName(),
00085         mesh,
00086         IOobject::MUST_READ
00087     );
00088 
00089     if (Uheader.headerOk())
00090     {
00091         Info<< "    Reading U" << endl;
00092         volVectorField U(Uheader, mesh);
00093         volTensorField gradU = fvc::grad(U);
00094 
00095         volScalarField Q
00096         (
00097             IOobject
00098             (
00099                 "Q",
00100                 runTime.timeName(),
00101                 mesh,
00102                 IOobject::NO_READ,
00103                 IOobject::NO_WRITE
00104             ),
00105             0.5*(sqr(tr(gradU)) - tr(((gradU)&(gradU))))
00106         );
00107 
00108         /*
00109         // This is a second way of calculating Q, that delivers results
00110         // very close, but not identical to the first approach.
00111 
00112         volSymmTensorField S = symm(gradU);  // symmetric part of tensor
00113         volTensorField W = skew(gradU);  // anti-symmetric part
00114 
00115         volScalarField SS =  S&&S;
00116         volScalarField WW =  W&&W;
00117 
00118         volScalarField Q
00119         (
00120             IOobject
00121             (
00122                 "Q",
00123                 runTime.timeName(),
00124                 mesh,
00125                 IOobject::NO_READ,
00126                 IOobject::NO_WRITE
00127             ),
00128             0.5*(WW - SS)
00129         );
00130         */
00131 
00132         Info<< "mag(Q) max/min : "
00133             << max(Q).value() << " "
00134             << min(Q).value() << endl;
00135 
00136         if (writeResults)
00137         {
00138             Q.write();
00139         }
00140     }
00141     else
00142     {
00143         Info<< "    No U" << endl;
00144     }
00145 
00146     Info<< "\nEnd\n" << endl;
00147 }
00148 
00149 
00150 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines