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

inverseFaceDistanceDiffusivity.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 \*---------------------------------------------------------------------------*/
00025 
00026 #include "inverseFaceDistanceDiffusivity.H"
00027 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00028 #include <OpenFOAM/HashSet.H>
00029 #include <meshTools/wallPoint.H>
00030 #include <OpenFOAM/MeshWave.H>
00031 
00032 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
00033 
00034 namespace Foam
00035 {
00036     defineTypeNameAndDebug(inverseFaceDistanceDiffusivity, 0);
00037 
00038     addToRunTimeSelectionTable
00039     (
00040         motionDiffusivity,
00041         inverseFaceDistanceDiffusivity,
00042         Istream
00043     );
00044 }
00045 
00046 
00047 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00048 
00049 Foam::inverseFaceDistanceDiffusivity::inverseFaceDistanceDiffusivity
00050 (
00051     const fvMotionSolver& mSolver,
00052     Istream& mdData
00053 )
00054 :
00055     uniformDiffusivity(mSolver, mdData),
00056     patchNames_(mdData)
00057 {
00058     correct();
00059 }
00060 
00061 
00062 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
00063 
00064 Foam::inverseFaceDistanceDiffusivity::~inverseFaceDistanceDiffusivity()
00065 {}
00066 
00067 
00068 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00069 
00070 void Foam::inverseFaceDistanceDiffusivity::correct()
00071 {
00072     const polyMesh& mesh = mSolver().mesh();
00073     const polyBoundaryMesh& bdry = mesh.boundaryMesh();
00074 
00075     labelHashSet patchSet(bdry.size());
00076 
00077     label nPatchFaces = 0;
00078 
00079     forAll (patchNames_, i)
00080     {
00081         label pID = bdry.findPatchID(patchNames_[i]);
00082 
00083         if (pID > -1)
00084         {
00085             patchSet.insert(pID);
00086             nPatchFaces += bdry[pID].size();
00087         }
00088     }
00089 
00090     List<wallPoint> faceDist(nPatchFaces);
00091     labelList changedFaces(nPatchFaces);
00092 
00093     nPatchFaces = 0;
00094 
00095     forAllConstIter(labelHashSet, patchSet, iter)
00096     {
00097         const polyPatch& patch = bdry[iter.key()];
00098 
00099         const vectorField::subField fc = patch.faceCentres();
00100 
00101         forAll(fc, patchFaceI)
00102         {
00103             changedFaces[nPatchFaces] = patch.start() + patchFaceI;
00104 
00105             faceDist[nPatchFaces] = wallPoint(fc[patchFaceI], 0);
00106 
00107             nPatchFaces++;
00108         }
00109     }
00110     faceDist.setSize(nPatchFaces);
00111     changedFaces.setSize(nPatchFaces);
00112 
00113     MeshWave<wallPoint> waveInfo
00114     (
00115         mesh,
00116         changedFaces,
00117         faceDist,
00118         mesh.globalData().nTotalCells() // max iterations
00119     );
00120 
00121     const List<wallPoint>& faceInfo = waveInfo.allFaceInfo();
00122     const List<wallPoint>& cellInfo = waveInfo.allCellInfo();
00123 
00124     for (label faceI=0; faceI<mesh.nInternalFaces(); faceI++)
00125     {
00126         scalar dist = faceInfo[faceI].distSqr();
00127 
00128         faceDiffusivity_[faceI] = 1.0/sqrt(dist);
00129     }
00130 
00131     forAll(faceDiffusivity_.boundaryField(), patchI)
00132     {
00133         fvsPatchScalarField& bfld = faceDiffusivity_.boundaryField()[patchI];
00134 
00135         const unallocLabelList& faceCells = bfld.patch().faceCells();
00136 
00137         if (patchSet.found(patchI))
00138         {
00139             forAll(bfld, i)
00140             {
00141                 scalar dist = cellInfo[faceCells[i]].distSqr();
00142                 bfld[i] = 1.0/sqrt(dist);
00143             }
00144         }
00145         else
00146         {
00147             label start = bfld.patch().patch().start();
00148 
00149             forAll(bfld, i)
00150             {
00151                 scalar dist = faceInfo[start+i].distSqr();
00152                 bfld[i] = 1.0/sqrt(dist);
00153             }
00154         }
00155     }
00156 }
00157 
00158 
00159 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines