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

MRFZone.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 "MRFZone.H"
00027 #include <finiteVolume/fvMesh.H>
00028 #include <finiteVolume/volFields.H>
00029 #include <finiteVolume/surfaceFields.H>
00030 #include <finiteVolume/fvMatrices.H>
00031 #include <OpenFOAM/syncTools.H>
00032 #include <meshTools/faceSet.H>
00033 #include <OpenFOAM/geometricOneField.H>
00034 
00035 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
00036 
00037 defineTypeNameAndDebug(Foam::MRFZone, 0);
00038 
00039 
00040 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
00041 
00042 void Foam::MRFZone::setMRFFaces()
00043 {
00044     const polyBoundaryMesh& patches = mesh_.boundaryMesh();
00045 
00046     // Type per face:
00047     //  0:not in zone
00048     //  1:moving with frame
00049     //  2:other
00050     labelList faceType(mesh_.nFaces(), 0);
00051 
00052     // Determine faces in cell zone
00053     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00054     // (without constructing cells)
00055 
00056     const labelList& own = mesh_.faceOwner();
00057     const labelList& nei = mesh_.faceNeighbour();
00058 
00059     // Cells in zone
00060     boolList zoneCell(mesh_.nCells(), false);
00061 
00062     if (cellZoneID_ != -1)
00063     {
00064         const labelList& cellLabels = mesh_.cellZones()[cellZoneID_];
00065         forAll(cellLabels, i)
00066         {
00067             zoneCell[cellLabels[i]] = true;
00068         }
00069     }
00070 
00071 
00072     label nZoneFaces = 0;
00073 
00074     for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
00075     {
00076         if (zoneCell[own[faceI]] || zoneCell[nei[faceI]])
00077         {
00078             faceType[faceI] = 1;
00079             nZoneFaces++;
00080         }
00081     }
00082 
00083 
00084     labelHashSet excludedPatches(excludedPatchLabels_);
00085 
00086     forAll(patches, patchI)
00087     {
00088         const polyPatch& pp = patches[patchI];
00089 
00090         if (pp.coupled() || excludedPatches.found(patchI))
00091         {
00092             forAll(pp, i)
00093             {
00094                 label faceI = pp.start()+i;
00095 
00096                 if (zoneCell[own[faceI]])
00097                 {
00098                     faceType[faceI] = 2;
00099                     nZoneFaces++;
00100                 }
00101             }
00102         }
00103         else if (!isA<emptyPolyPatch>(pp))
00104         {
00105             forAll(pp, i)
00106             {
00107                 label faceI = pp.start()+i;
00108 
00109                 if (zoneCell[own[faceI]])
00110                 {
00111                     faceType[faceI] = 1;
00112                     nZoneFaces++;
00113                 }
00114             }
00115         }
00116     }
00117 
00118     // Now we have for faceType:
00119     //  0   : face not in cellZone
00120     //  1   : internal face or normal patch face
00121     //  2   : coupled patch face or excluded patch face
00122 
00123     // Sort into lists per patch.
00124 
00125     internalFaces_.setSize(mesh_.nFaces());
00126     label nInternal = 0;
00127 
00128     for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
00129     {
00130         if (faceType[faceI] == 1)
00131         {
00132             internalFaces_[nInternal++] = faceI;
00133         }
00134     }
00135     internalFaces_.setSize(nInternal);
00136 
00137     labelList nIncludedFaces(patches.size(), 0);
00138     labelList nExcludedFaces(patches.size(), 0);
00139 
00140     forAll(patches, patchi)
00141     {
00142         const polyPatch& pp = patches[patchi];
00143 
00144         forAll(pp, patchFacei)
00145         {
00146             label faceI = pp.start()+patchFacei;
00147 
00148             if (faceType[faceI] == 1)
00149             {
00150                 nIncludedFaces[patchi]++;
00151             }
00152             else if (faceType[faceI] == 2)
00153             {
00154                 nExcludedFaces[patchi]++;
00155             }
00156         }
00157     }
00158 
00159     includedFaces_.setSize(patches.size());
00160     excludedFaces_.setSize(patches.size());
00161     forAll(nIncludedFaces, patchi)
00162     {
00163         includedFaces_[patchi].setSize(nIncludedFaces[patchi]);
00164         excludedFaces_[patchi].setSize(nExcludedFaces[patchi]);
00165     }
00166     nIncludedFaces = 0;
00167     nExcludedFaces = 0;
00168 
00169     forAll(patches, patchi)
00170     {
00171         const polyPatch& pp = patches[patchi];
00172 
00173         forAll(pp, patchFacei)
00174         {
00175             label faceI = pp.start()+patchFacei;
00176 
00177             if (faceType[faceI] == 1)
00178             {
00179                 includedFaces_[patchi][nIncludedFaces[patchi]++] = patchFacei;
00180             }
00181             else if (faceType[faceI] == 2)
00182             {
00183                 excludedFaces_[patchi][nExcludedFaces[patchi]++] = patchFacei;
00184             }
00185         }
00186     }
00187 
00188 
00189     if (debug)
00190     {
00191         faceSet internalFaces(mesh_, "internalFaces", internalFaces_);
00192         Pout<< "Writing " << internalFaces.size()
00193             << " internal faces in MRF zone to faceSet "
00194             << internalFaces.name() << endl;
00195         internalFaces.write();
00196 
00197         faceSet MRFFaces(mesh_, "includedFaces", 100);
00198         forAll(includedFaces_, patchi)
00199         {
00200             forAll(includedFaces_[patchi], i)
00201             {
00202                 label patchFacei = includedFaces_[patchi][i];
00203                 MRFFaces.insert(patches[patchi].start()+patchFacei);
00204             }
00205         }
00206         Pout<< "Writing " << MRFFaces.size()
00207             << " patch faces in MRF zone to faceSet "
00208             << MRFFaces.name() << endl;
00209         MRFFaces.write();
00210 
00211         faceSet excludedFaces(mesh_, "excludedFaces", 100);
00212         forAll(excludedFaces_, patchi)
00213         {
00214             forAll(excludedFaces_[patchi], i)
00215             {
00216                 label patchFacei = excludedFaces_[patchi][i];
00217                 excludedFaces.insert(patches[patchi].start()+patchFacei);
00218             }
00219         }
00220         Pout<< "Writing " << excludedFaces.size()
00221             << " faces in MRF zone with special handling to faceSet "
00222             << excludedFaces.name() << endl;
00223         excludedFaces.write();
00224     }
00225 }
00226 
00227 
00228 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00229 
00230 Foam::MRFZone::MRFZone(const fvMesh& mesh, Istream& is)
00231 :
00232     mesh_(mesh),
00233     name_(is),
00234     dict_(is),
00235     cellZoneID_(mesh_.cellZones().findZoneID(name_)),
00236     excludedPatchNames_
00237     (
00238         dict_.lookupOrDefault("nonRotatingPatches", wordList(0))
00239     ),
00240     origin_(dict_.lookup("origin")),
00241     axis_(dict_.lookup("axis")),
00242     omega_(dict_.lookup("omega")),
00243     Omega_("Omega", omega_*axis_)
00244 {
00245     if (dict_.found("patches"))
00246     {
00247         WarningIn("MRFZone(const fvMesh&, Istream&)")
00248             << "Ignoring entry 'patches'\n"
00249             << "    By default all patches within the rotating region rotate.\n"
00250             << "    Optionally supply excluded patches using 'nonRotatingPatches'."
00251             << endl;
00252     }
00253 
00254     const polyBoundaryMesh& patches = mesh_.boundaryMesh();
00255 
00256     axis_ = axis_/mag(axis_);
00257     Omega_ = omega_*axis_;
00258 
00259     excludedPatchLabels_.setSize(excludedPatchNames_.size());
00260 
00261     forAll(excludedPatchNames_, i)
00262     {
00263         excludedPatchLabels_[i] = patches.findPatchID(excludedPatchNames_[i]);
00264 
00265         if (excludedPatchLabels_[i] == -1)
00266         {
00267             FatalErrorIn
00268             (
00269                 "Foam::MRFZone::MRFZone(const fvMesh&, Istream&)"
00270             )   << "cannot find MRF patch " << excludedPatchNames_[i]
00271                 << exit(FatalError);
00272         }
00273     }
00274 
00275 
00276     bool cellZoneFound = (cellZoneID_ != -1);
00277     reduce(cellZoneFound, orOp<bool>());
00278 
00279     if (!cellZoneFound)
00280     {
00281         FatalErrorIn
00282         (
00283             "Foam::MRFZone::MRFZone(const fvMesh&, Istream&)"
00284         )   << "cannot find MRF cellZone " << name_
00285             << exit(FatalError);
00286     }
00287 
00288     setMRFFaces();
00289 }
00290 
00291 
00292 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00293 
00294 void Foam::MRFZone::addCoriolis(fvVectorMatrix& UEqn) const
00295 {
00296     if (cellZoneID_ == -1)
00297     {
00298         return;
00299     }
00300 
00301     const labelList& cells = mesh_.cellZones()[cellZoneID_];
00302     const scalarField& V = mesh_.V();
00303     vectorField& Usource = UEqn.source();
00304     const vectorField& U = UEqn.psi();
00305     const vector& Omega = Omega_.value();
00306 
00307     forAll(cells, i)
00308     {
00309         label celli = cells[i];
00310         Usource[celli] -= V[celli]*(Omega ^ U[celli]);
00311     }
00312 }
00313 
00314 
00315 void Foam::MRFZone::addCoriolis
00316 (
00317     const volScalarField& rho,
00318     fvVectorMatrix& UEqn
00319 ) const
00320 {
00321     if (cellZoneID_ == -1)
00322     {
00323         return;
00324     }
00325 
00326     const labelList& cells = mesh_.cellZones()[cellZoneID_];
00327     const scalarField& V = mesh_.V();
00328     vectorField& Usource = UEqn.source();
00329     const vectorField& U = UEqn.psi();
00330     const vector& Omega = Omega_.value();
00331 
00332     forAll(cells, i)
00333     {
00334         label celli = cells[i];
00335         Usource[celli] -= V[celli]*rho[celli]*(Omega ^ U[celli]);
00336     }
00337 }
00338 
00339 
00340 void Foam::MRFZone::relativeVelocity(volVectorField& U) const
00341 {
00342     const volVectorField& C = mesh_.C();
00343 
00344     const vector& origin = origin_.value();
00345     const vector& Omega = Omega_.value();
00346 
00347     const labelList& cells = mesh_.cellZones()[cellZoneID_];
00348 
00349     forAll(cells, i)
00350     {
00351         label celli = cells[i];
00352         U[celli] -= (Omega ^ (C[celli] - origin));
00353     }
00354 
00355     // Included patches
00356     forAll(includedFaces_, patchi)
00357     {
00358         forAll(includedFaces_[patchi], i)
00359         {
00360             label patchFacei = includedFaces_[patchi][i];
00361             U.boundaryField()[patchi][patchFacei] = vector::zero;
00362         }
00363     }
00364 
00365     // Excluded patches
00366     forAll(excludedFaces_, patchi)
00367     {
00368         forAll(excludedFaces_[patchi], i)
00369         {
00370             label patchFacei = excludedFaces_[patchi][i];
00371             U.boundaryField()[patchi][patchFacei] -=
00372                 (Omega ^ (C.boundaryField()[patchi][patchFacei] - origin));
00373         }
00374     }
00375 }
00376 
00377 
00378 void Foam::MRFZone::absoluteVelocity(volVectorField& U) const
00379 {
00380     const volVectorField& C = mesh_.C();
00381 
00382     const vector& origin = origin_.value();
00383     const vector& Omega = Omega_.value();
00384 
00385     const labelList& cells = mesh_.cellZones()[cellZoneID_];
00386 
00387     forAll(cells, i)
00388     {
00389         label celli = cells[i];
00390         U[celli] += (Omega ^ (C[celli] - origin));
00391     }
00392 
00393     // Included patches
00394     forAll(includedFaces_, patchi)
00395     {
00396         forAll(includedFaces_[patchi], i)
00397         {
00398             label patchFacei = includedFaces_[patchi][i];
00399             U.boundaryField()[patchi][patchFacei] =
00400                 (Omega ^ (C.boundaryField()[patchi][patchFacei] - origin));
00401         }
00402     }
00403 
00404     // Excluded patches
00405     forAll(excludedFaces_, patchi)
00406     {
00407         forAll(excludedFaces_[patchi], i)
00408         {
00409             label patchFacei = excludedFaces_[patchi][i];
00410             U.boundaryField()[patchi][patchFacei] +=
00411                 (Omega ^ (C.boundaryField()[patchi][patchFacei] - origin));
00412         }
00413     }
00414 }
00415 
00416 
00417 void Foam::MRFZone::relativeFlux(surfaceScalarField& phi) const
00418 {
00419     relativeRhoFlux(geometricOneField(), phi);
00420 }
00421 
00422 
00423 void Foam::MRFZone::relativeFlux
00424 (
00425     const surfaceScalarField& rho,
00426     surfaceScalarField& phi
00427 ) const
00428 {
00429     relativeRhoFlux(rho, phi);
00430 }
00431 
00432 
00433 void Foam::MRFZone::absoluteFlux(surfaceScalarField& phi) const
00434 {
00435     absoluteRhoFlux(geometricOneField(), phi);
00436 }
00437 
00438 
00439 void Foam::MRFZone::absoluteFlux
00440 (
00441     const surfaceScalarField& rho,
00442     surfaceScalarField& phi
00443 ) const
00444 {
00445     absoluteRhoFlux(rho, phi);
00446 }
00447 
00448 
00449 void Foam::MRFZone::correctBoundaryVelocity(volVectorField& U) const
00450 {
00451     const vector& origin = origin_.value();
00452     const vector& Omega = Omega_.value();
00453 
00454     // Included patches
00455     forAll(includedFaces_, patchi)
00456     {
00457         const vectorField& patchC = mesh_.Cf().boundaryField()[patchi];
00458 
00459         vectorField pfld(U.boundaryField()[patchi]);
00460 
00461         forAll(includedFaces_[patchi], i)
00462         {
00463             label patchFacei = includedFaces_[patchi][i];
00464 
00465             pfld[patchFacei] = (Omega ^ (patchC[patchFacei] - origin));
00466         }
00467 
00468         U.boundaryField()[patchi] == pfld;
00469     }
00470 }
00471 
00472 
00473 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines