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 Lambda2 00026 00027 Description 00028 Calculates and writes the second largest eigenvalue of the sum of the 00029 square of the symmetrical and anti-symmetrical parts of the velocity 00030 gradient tensor. 00031 00032 The -noWrite option has no meaning. 00033 00034 Usage 00035 00036 - Lambda2 [OPTIONS] 00037 00038 @param -dict <dictionary name>\n 00039 Use named dictionary instead of system/controlDict. 00040 00041 @param -noZero \n 00042 Ignore timestep 0. 00043 00044 @param -constant \n 00045 Include the constant directory. 00046 00047 @param -time <time>\n 00048 Apply only to specific time. 00049 00050 @param -latestTime \n 00051 Only apply to latest time step. 00052 00053 @param -case <dir>\n 00054 Case directory. 00055 00056 @param -parallel \n 00057 Run in parallel. 00058 00059 @param -help \n 00060 Display help message. 00061 00062 @param -doc \n 00063 Display Doxygen API documentation page for this application. 00064 00065 @param -srcDoc \n 00066 Display Doxygen source documentation page for this application. 00067 00068 \*---------------------------------------------------------------------------*/ 00069 00070 #include <postCalc/calc.H> 00071 #include <finiteVolume/fvc.H> 00072 00073 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00074 00075 void Foam::calc(const argList& args, const Time& runTime, const fvMesh& mesh) 00076 { 00077 IOobject Uheader 00078 ( 00079 "U", 00080 runTime.timeName(), 00081 mesh, 00082 IOobject::MUST_READ 00083 ); 00084 00085 if (Uheader.headerOk()) 00086 { 00087 Info<< " Reading U" << endl; 00088 volVectorField U(Uheader, mesh); 00089 00090 const volTensorField gradU(fvc::grad(U)); 00091 00092 volTensorField SSplusWW = 00093 (symm(gradU) & symm(gradU)) + (skew(gradU) & skew(gradU)); 00094 00095 volScalarField Lambda2 00096 ( 00097 IOobject 00098 ( 00099 "Lambda2", 00100 runTime.timeName(), 00101 mesh, 00102 IOobject::NO_READ, 00103 IOobject::NO_WRITE 00104 ), 00105 -eigenValues(SSplusWW)().component(vector::Y) 00106 ); 00107 00108 Info << " Writing -Lambda2" << endl; 00109 Lambda2.write(); 00110 } 00111 else 00112 { 00113 Info<< " No U" << endl; 00114 } 00115 00116 Info<< "\nEnd\n" << endl; 00117 } 00118 00119 00120 // ************************ vim: set sw=4 sts=4 et: ************************ //