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

MixedDiffuseSpecular.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) 2009-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 \*---------------------------------------------------------------------------*/
00025 
00026 #include "MixedDiffuseSpecular.H"
00027 
00028 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00029 
00030 template <class CloudType>
00031 Foam::MixedDiffuseSpecular<CloudType>::MixedDiffuseSpecular
00032 (
00033     const dictionary& dict,
00034     CloudType& cloud
00035 )
00036 :
00037     WallInteractionModel<CloudType>(dict, cloud, typeName),
00038     diffuseFraction_(readScalar(this->coeffDict().lookup("diffuseFraction")))
00039 {}
00040 
00041 
00042 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
00043 
00044 template <class CloudType>
00045 Foam::MixedDiffuseSpecular<CloudType>::~MixedDiffuseSpecular()
00046 {}
00047 
00048 
00049 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00050 
00051 template <class CloudType>
00052 void Foam::MixedDiffuseSpecular<CloudType>::correct
00053 (
00054     const wallPolyPatch& wpp,
00055     const label faceId,
00056     vector& U,
00057     scalar& Ei,
00058     label typeId
00059 )
00060 {
00061     label wppIndex = wpp.index();
00062 
00063     label wppLocalFace = wpp.whichFace(faceId);
00064 
00065     vector nw = wpp.faceAreas()[wppLocalFace];
00066 
00067     // Normal unit vector
00068     nw /= mag(nw);
00069 
00070     // Normal velocity magnitude
00071     scalar U_dot_nw = U & nw;
00072 
00073     CloudType& cloud(this->owner());
00074 
00075     Random& rndGen(cloud.rndGen());
00076 
00077     if (diffuseFraction_ > rndGen.scalar01())
00078     {
00079         // Diffuse reflection
00080 
00081         // Wall tangential velocity (flow direction)
00082         vector Ut = U - U_dot_nw*nw;
00083 
00084         while (mag(Ut) < SMALL)
00085         {
00086             // If the incident velocity is parallel to the face normal, no
00087             // tangential direction can be chosen.  Add a perturbation to the
00088             // incoming velocity and recalculate.
00089 
00090             U = vector
00091             (
00092                 U.x()*(0.8 + 0.2*rndGen.scalar01()),
00093                 U.y()*(0.8 + 0.2*rndGen.scalar01()),
00094                 U.z()*(0.8 + 0.2*rndGen.scalar01())
00095             );
00096 
00097             U_dot_nw = U & nw;
00098 
00099             Ut = U - U_dot_nw*nw;
00100         }
00101 
00102         // Wall tangential unit vector
00103         vector tw1 = Ut/mag(Ut);
00104 
00105         // Other tangential unit vector
00106         vector tw2 = nw^tw1;
00107 
00108         scalar T = cloud.boundaryT().boundaryField()[wppIndex][wppLocalFace];
00109 
00110         scalar mass = cloud.constProps(typeId).mass();
00111 
00112         scalar iDof = cloud.constProps(typeId).internalDegreesOfFreedom();
00113 
00114         U =
00115             sqrt(CloudType::kb*T/mass)
00116            *(
00117                 rndGen.GaussNormal()*tw1
00118               + rndGen.GaussNormal()*tw2
00119               - sqrt(-2.0*log(max(1 - rndGen.scalar01(), VSMALL)))*nw
00120             );
00121 
00122         U += cloud.boundaryU().boundaryField()[wppIndex][wppLocalFace];
00123 
00124         Ei = cloud.equipartitionInternalEnergy(T, iDof);
00125     }
00126     else
00127     {
00128         // Specular reflection
00129 
00130         if (U_dot_nw > 0.0)
00131         {
00132             U -= 2.0*U_dot_nw*nw;
00133         }
00134     }
00135 
00136 }
00137 
00138 
00139 // ************************************************************************* //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines