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

polyMesh.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 "polyMesh.H"
00027 #include <OpenFOAM/Time.H>
00028 #include <OpenFOAM/cellIOList.H>
00029 #include <OpenFOAM/SubList.H>
00030 #include <OpenFOAM/wedgePolyPatch.H>
00031 #include <OpenFOAM/emptyPolyPatch.H>
00032 #include <OpenFOAM/globalMeshData.H>
00033 #include <OpenFOAM/processorPolyPatch.H>
00034 #include <OpenFOAM/OSspecific.H>
00035 #include <OpenFOAM/demandDrivenData.H>
00036 
00037 #include <OpenFOAM/pointMesh.H>
00038 
00039 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
00040 
00041 defineTypeNameAndDebug(Foam::polyMesh, 0);
00042 
00043 
00044 Foam::word Foam::polyMesh::defaultRegion = "region0";
00045 Foam::word Foam::polyMesh::meshSubDir = "polyMesh";
00046 
00047 
00048 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
00049 
00050 void Foam::polyMesh::calcDirections() const
00051 {
00052     for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
00053     {
00054         solutionD_[cmpt] = 1;
00055     }
00056 
00057     // Knock out empty and wedge directions. Note:they will be present on all
00058     // domains.
00059 
00060     label nEmptyPatches = 0;
00061     label nWedgePatches = 0;
00062 
00063     vector emptyDirVec = vector::zero;
00064     vector wedgeDirVec = vector::zero;
00065 
00066     forAll(boundaryMesh(), patchi)
00067     {
00068         if (boundaryMesh()[patchi].size())
00069         {
00070             if (isA<emptyPolyPatch>(boundaryMesh()[patchi]))
00071             {
00072                 nEmptyPatches++;
00073                 emptyDirVec += sum(cmptMag(boundaryMesh()[patchi].faceAreas()));
00074             }
00075             else if (isA<wedgePolyPatch>(boundaryMesh()[patchi]))
00076             {
00077                 const wedgePolyPatch& wpp = refCast<const wedgePolyPatch>
00078                 (
00079                     boundaryMesh()[patchi]
00080                 );
00081 
00082                 nWedgePatches++;
00083                 wedgeDirVec += cmptMag(wpp.centreNormal());
00084             }
00085         }
00086     }
00087 
00088     reduce(nEmptyPatches, maxOp<label>());
00089     reduce(nWedgePatches, maxOp<label>());
00090 
00091     if (nEmptyPatches)
00092     {
00093         reduce(emptyDirVec, sumOp<vector>());
00094 
00095         emptyDirVec /= mag(emptyDirVec);
00096 
00097         for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
00098         {
00099             if (emptyDirVec[cmpt] > 1e-6)
00100             {
00101                 solutionD_[cmpt] = -1;
00102             }
00103             else
00104             {
00105                 solutionD_[cmpt] = 1;
00106             }
00107         }
00108     }
00109 
00110 
00111     // Knock out wedge directions
00112 
00113     geometricD_ = solutionD_;
00114 
00115     if (nWedgePatches)
00116     {
00117         reduce(wedgeDirVec, sumOp<vector>());
00118 
00119         wedgeDirVec /= mag(wedgeDirVec);
00120 
00121         for (direction cmpt=0; cmpt<vector::nComponents; cmpt++)
00122         {
00123             if (wedgeDirVec[cmpt] > 1e-6)
00124             {
00125                 geometricD_[cmpt] = -1;
00126             }
00127             else
00128             {
00129                 geometricD_[cmpt] = 1;
00130             }
00131         }
00132     }
00133 }
00134 
00135 
00136 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00137 
00138 Foam::polyMesh::polyMesh(const IOobject& io)
00139 :
00140     objectRegistry(io),
00141     primitiveMesh(),
00142     points_
00143     (
00144         IOobject
00145         (
00146             "points",
00147             time().findInstance(meshDir(), "points"),
00148             meshSubDir,
00149             *this,
00150             IOobject::MUST_READ,
00151             IOobject::NO_WRITE
00152         )
00153     ),
00154     faces_
00155     (
00156         IOobject
00157         (
00158             "faces",
00159             time().findInstance(meshDir(), "faces"),
00160             meshSubDir,
00161             *this,
00162             IOobject::MUST_READ,
00163             IOobject::NO_WRITE
00164         )
00165     ),
00166     owner_
00167     (
00168         IOobject
00169         (
00170             "owner",
00171             time().findInstance(meshDir(), "faces"),
00172             meshSubDir,
00173             *this,
00174             IOobject::READ_IF_PRESENT,
00175             IOobject::NO_WRITE
00176         )
00177     ),
00178     neighbour_
00179     (
00180         IOobject
00181         (
00182             "neighbour",
00183             time().findInstance(meshDir(), "faces"),
00184             meshSubDir,
00185             *this,
00186             IOobject::READ_IF_PRESENT,
00187             IOobject::NO_WRITE
00188         )
00189     ),
00190     clearedPrimitives_(false),
00191     boundary_
00192     (
00193         IOobject
00194         (
00195             "boundary",
00196             time().findInstance(meshDir(), "boundary"),
00197             meshSubDir,
00198             *this,
00199             IOobject::MUST_READ,
00200             IOobject::NO_WRITE
00201         ),
00202         *this
00203     ),
00204     bounds_(points_),
00205     geometricD_(Vector<label>::zero),
00206     solutionD_(Vector<label>::zero),
00207     pointZones_
00208     (
00209         IOobject
00210         (
00211             "pointZones",
00212             time().findInstance
00213             (
00214                 meshDir(),
00215                 "pointZones",
00216                 IOobject::READ_IF_PRESENT
00217             ),
00218             meshSubDir,
00219             *this,
00220             IOobject::READ_IF_PRESENT,
00221             IOobject::NO_WRITE
00222         ),
00223         *this
00224     ),
00225     faceZones_
00226     (
00227         IOobject
00228         (
00229             "faceZones",
00230             time().findInstance
00231             (
00232                 meshDir(),
00233                 "faceZones",
00234                 IOobject::READ_IF_PRESENT
00235             ),
00236             meshSubDir,
00237             *this,
00238             IOobject::READ_IF_PRESENT,
00239             IOobject::NO_WRITE
00240         ),
00241         *this
00242     ),
00243     cellZones_
00244     (
00245         IOobject
00246         (
00247             "cellZones",
00248             time().findInstance
00249             (
00250                 meshDir(),
00251                 "cellZones",
00252                 IOobject::READ_IF_PRESENT
00253             ),
00254             meshSubDir,
00255             *this,
00256             IOobject::READ_IF_PRESENT,
00257             IOobject::NO_WRITE
00258         ),
00259         *this
00260     ),
00261     globalMeshDataPtr_(NULL),
00262     moving_(false),
00263     changing_(false),
00264     curMotionTimeIndex_(time().timeIndex()),
00265     oldPointsPtr_(NULL)
00266 {
00267     if (exists(owner_.objectPath()))
00268     {
00269         initMesh();
00270     }
00271     else
00272     {
00273         cellIOList cLst
00274         (
00275             IOobject
00276             (
00277                 "cells",
00278                 time().findInstance(meshDir(), "cells"),
00279                 meshSubDir,
00280                 *this,
00281                 IOobject::MUST_READ,
00282                 IOobject::NO_WRITE
00283             )
00284         );
00285 
00286         // Set the primitive mesh
00287         initMesh(cLst);
00288 
00289         owner_.write();
00290         neighbour_.write();
00291     }
00292 
00293     // Calculate topology for the patches (processor-processor comms etc.)
00294     boundary_.updateMesh();
00295 
00296     // Calculate the geometry for the patches (transformation tensors etc.)
00297     boundary_.calcGeometry();
00298 
00299     // Warn if global empty mesh
00300     if (returnReduce(nPoints(), sumOp<label>()) == 0)
00301     {
00302         WarningIn("polyMesh(const IOobject&)")
00303             << "no points in mesh" << endl;
00304     }
00305     if (returnReduce(nCells(), sumOp<label>()) == 0)
00306     {
00307         WarningIn("polyMesh(const IOobject&)")
00308             << "no cells in mesh" << endl;
00309     }
00310 }
00311 
00312 
00313 Foam::polyMesh::polyMesh
00314 (
00315     const IOobject& io,
00316     const Xfer<pointField>& points,
00317     const Xfer<faceList>& faces,
00318     const Xfer<labelList>& owner,
00319     const Xfer<labelList>& neighbour,
00320     const bool syncPar
00321 )
00322 :
00323     objectRegistry(io),
00324     primitiveMesh(),
00325     points_
00326     (
00327         IOobject
00328         (
00329             "points",
00330             instance(),
00331             meshSubDir,
00332             *this,
00333             IOobject::NO_READ,
00334             IOobject::AUTO_WRITE
00335         ),
00336         points
00337     ),
00338     faces_
00339     (
00340         IOobject
00341         (
00342             "faces",
00343             instance(),
00344             meshSubDir,
00345             *this,
00346             IOobject::NO_READ,
00347             IOobject::AUTO_WRITE
00348         ),
00349         faces
00350     ),
00351     owner_
00352     (
00353         IOobject
00354         (
00355             "owner",
00356             instance(),
00357             meshSubDir,
00358             *this,
00359             IOobject::NO_READ,
00360             IOobject::AUTO_WRITE
00361         ),
00362         owner
00363     ),
00364     neighbour_
00365     (
00366         IOobject
00367         (
00368             "neighbour",
00369             instance(),
00370             meshSubDir,
00371             *this,
00372             IOobject::NO_READ,
00373             IOobject::AUTO_WRITE
00374         ),
00375         neighbour
00376     ),
00377     clearedPrimitives_(false),
00378     boundary_
00379     (
00380         IOobject
00381         (
00382             "boundary",
00383             instance(),
00384             meshSubDir,
00385             *this,
00386             IOobject::NO_READ,
00387             IOobject::AUTO_WRITE
00388         ),
00389         *this,
00390         0
00391     ),
00392     bounds_(points_, syncPar),
00393     geometricD_(Vector<label>::zero),
00394     solutionD_(Vector<label>::zero),
00395     pointZones_
00396     (
00397         IOobject
00398         (
00399             "pointZones",
00400             instance(),
00401             meshSubDir,
00402             *this,
00403             IOobject::NO_READ,
00404             IOobject::NO_WRITE
00405         ),
00406         *this,
00407         0
00408     ),
00409     faceZones_
00410     (
00411         IOobject
00412         (
00413             "faceZones",
00414             instance(),
00415             meshSubDir,
00416             *this,
00417             IOobject::NO_READ,
00418             IOobject::NO_WRITE
00419         ),
00420         *this,
00421         0
00422     ),
00423     cellZones_
00424     (
00425         IOobject
00426         (
00427             "cellZones",
00428             instance(),
00429             meshSubDir,
00430             *this,
00431             IOobject::NO_READ,
00432             IOobject::NO_WRITE
00433         ),
00434         *this,
00435         0
00436     ),
00437     globalMeshDataPtr_(NULL),
00438     moving_(false),
00439     changing_(false),
00440     curMotionTimeIndex_(time().timeIndex()),
00441     oldPointsPtr_(NULL)
00442 {
00443     // Check if the faces and cells are valid
00444     forAll (faces_, faceI)
00445     {
00446         const face& curFace = faces_[faceI];
00447 
00448         if (min(curFace) < 0 || max(curFace) > points_.size())
00449         {
00450             FatalErrorIn
00451             (
00452                 "polyMesh::polyMesh\n"
00453                 "(\n"
00454                 "    const IOobject& io,\n"
00455                 "    const pointField& points,\n"
00456                 "    const faceList& faces,\n"
00457                 "    const cellList& cells\n"
00458                 ")\n"
00459             )   << "Face " << faceI << "contains vertex labels out of range: "
00460                 << curFace << " Max point index = " << points_.size()
00461                 << abort(FatalError);
00462         }
00463     }
00464 
00465     // Set the primitive mesh
00466     initMesh();
00467 }
00468 
00469 
00470 Foam::polyMesh::polyMesh
00471 (
00472     const IOobject& io,
00473     const Xfer<pointField>& points,
00474     const Xfer<faceList>& faces,
00475     const Xfer<cellList>& cells,
00476     const bool syncPar
00477 )
00478 :
00479     objectRegistry(io),
00480     primitiveMesh(),
00481     points_
00482     (
00483         IOobject
00484         (
00485             "points",
00486             instance(),
00487             meshSubDir,
00488             *this,
00489             IOobject::NO_READ,
00490             IOobject::AUTO_WRITE
00491         ),
00492         points
00493     ),
00494     faces_
00495     (
00496         IOobject
00497         (
00498             "faces",
00499             instance(),
00500             meshSubDir,
00501             *this,
00502             IOobject::NO_READ,
00503             IOobject::AUTO_WRITE
00504         ),
00505         faces
00506     ),
00507     owner_
00508     (
00509         IOobject
00510         (
00511             "owner",
00512             instance(),
00513             meshSubDir,
00514             *this,
00515             IOobject::NO_READ,
00516             IOobject::AUTO_WRITE
00517         ),
00518         0
00519     ),
00520     neighbour_
00521     (
00522         IOobject
00523         (
00524             "neighbour",
00525             instance(),
00526             meshSubDir,
00527             *this,
00528             IOobject::NO_READ,
00529             IOobject::AUTO_WRITE
00530         ),
00531         0
00532     ),
00533     clearedPrimitives_(false),
00534     boundary_
00535     (
00536         IOobject
00537         (
00538             "boundary",
00539             instance(),
00540             meshSubDir,
00541             *this,
00542             IOobject::NO_READ,
00543             IOobject::AUTO_WRITE
00544         ),
00545         *this,
00546         0
00547     ),
00548     bounds_(points_, syncPar),
00549     geometricD_(Vector<label>::zero),
00550     solutionD_(Vector<label>::zero),
00551     pointZones_
00552     (
00553         IOobject
00554         (
00555             "pointZones",
00556             instance(),
00557             meshSubDir,
00558             *this,
00559             IOobject::NO_READ,
00560             IOobject::NO_WRITE
00561         ),
00562         *this,
00563         0
00564     ),
00565     faceZones_
00566     (
00567         IOobject
00568         (
00569             "faceZones",
00570             instance(),
00571             meshSubDir,
00572             *this,
00573             IOobject::NO_READ,
00574             IOobject::NO_WRITE
00575         ),
00576         *this,
00577         0
00578     ),
00579     cellZones_
00580     (
00581         IOobject
00582         (
00583             "cellZones",
00584             instance(),
00585             meshSubDir,
00586             *this,
00587             IOobject::NO_READ,
00588             IOobject::NO_WRITE
00589         ),
00590         *this,
00591         0
00592     ),
00593     globalMeshDataPtr_(NULL),
00594     moving_(false),
00595     changing_(false),
00596     curMotionTimeIndex_(time().timeIndex()),
00597     oldPointsPtr_(NULL)
00598 {
00599     // Check if faces are valid
00600     forAll (faces_, faceI)
00601     {
00602         const face& curFace = faces_[faceI];
00603 
00604         if (min(curFace) < 0 || max(curFace) > points_.size())
00605         {
00606             FatalErrorIn
00607             (
00608                 "polyMesh::polyMesh\n"
00609                 "(\n"
00610                 "    const IOobject&,\n"
00611                 "    const Xfer<pointField>&,\n"
00612                 "    const Xfer<faceList>&,\n"
00613                 "    const Xfer<cellList>&\n"
00614                 ")\n"
00615             )   << "Face " << faceI << "contains vertex labels out of range: "
00616                 << curFace << " Max point index = " << points_.size()
00617                 << abort(FatalError);
00618         }
00619     }
00620 
00621     // transfer in cell list
00622     cellList cLst(cells);
00623 
00624     // Check if cells are valid
00625     forAll (cLst, cellI)
00626     {
00627         const cell& curCell = cLst[cellI];
00628 
00629         if (min(curCell) < 0 || max(curCell) > faces_.size())
00630         {
00631             FatalErrorIn
00632             (
00633                 "polyMesh::polyMesh\n"
00634                 "(\n"
00635                 "    const IOobject&,\n"
00636                 "    const Xfer<pointField>&,\n"
00637                 "    const Xfer<faceList>&,\n"
00638                 "    const Xfer<cellList>&\n"
00639                 ")\n"
00640             )   << "Cell " << cellI << "contains face labels out of range: "
00641                 << curCell << " Max face index = " << faces_.size()
00642                 << abort(FatalError);
00643         }
00644     }
00645 
00646     // Set the primitive mesh
00647     initMesh(cLst);
00648 }
00649 
00650 
00651 void Foam::polyMesh::resetPrimitives
00652 (
00653     const Xfer<pointField>& points,
00654     const Xfer<faceList>& faces,
00655     const Xfer<labelList>& owner,
00656     const Xfer<labelList>& neighbour,
00657     const labelList& patchSizes,
00658     const labelList& patchStarts,
00659     const bool validBoundary
00660 )
00661 {
00662     // Clear addressing. Keep geometric props for mapping.
00663     clearAddressing();
00664 
00665     // Take over new primitive data.
00666     // Optimized to avoid overwriting data at all
00667     if (&points)
00668     {
00669         points_.transfer(points());
00670         bounds_ = boundBox(points_, validBoundary);
00671     }
00672 
00673     if (&faces)
00674     {
00675         faces_.transfer(faces());
00676     }
00677 
00678     if (&owner)
00679     {
00680         owner_.transfer(owner());
00681     }
00682 
00683     if (&neighbour)
00684     {
00685         neighbour_.transfer(neighbour());
00686     }
00687 
00688 
00689     // Reset patch sizes and starts
00690     forAll(boundary_, patchI)
00691     {
00692         boundary_[patchI] = polyPatch
00693         (
00694             boundary_[patchI].name(),
00695             patchSizes[patchI],
00696             patchStarts[patchI],
00697             patchI,
00698             boundary_
00699         );
00700     }
00701 
00702 
00703     // Flags the mesh files as being changed
00704     setInstance(time().timeName());
00705 
00706     // Check if the faces and cells are valid
00707     forAll (faces_, faceI)
00708     {
00709         const face& curFace = faces_[faceI];
00710 
00711         if (min(curFace) < 0 || max(curFace) > points_.size())
00712         {
00713             FatalErrorIn
00714             (
00715                 "polyMesh::polyMesh::resetPrimitives\n"
00716                 "(\n"
00717                 "    const Xfer<pointField>&,\n"
00718                 "    const Xfer<faceList>&,\n"
00719                 "    const Xfer<labelList>& owner,\n"
00720                 "    const Xfer<labelList>& neighbour,\n"
00721                 "    const labelList& patchSizes,\n"
00722                 "    const labelList& patchStarts\n"
00723                 ")\n"
00724             )   << "Face " << faceI << " contains vertex labels out of range: "
00725                 << curFace << " Max point index = " << points_.size()
00726                 << abort(FatalError);
00727         }
00728     }
00729 
00730 
00731     // Set the primitive mesh from the owner_, neighbour_.
00732     // Works out from patch end where the active faces stop.
00733     initMesh();
00734 
00735 
00736     if (validBoundary)
00737     {
00738         // Note that we assume that all the patches stay the same and are
00739         // correct etc. so we can already use the patches to do
00740         // processor-processor comms.
00741 
00742         // Calculate topology for the patches (processor-processor comms etc.)
00743         boundary_.updateMesh();
00744 
00745         // Calculate the geometry for the patches (transformation tensors etc.)
00746         boundary_.calcGeometry();
00747 
00748         // Warn if global empty mesh
00749         if
00750         (
00751             (returnReduce(nPoints(), sumOp<label>()) == 0)
00752          || (returnReduce(nCells(), sumOp<label>()) == 0)
00753         )
00754         {
00755             FatalErrorIn
00756             (
00757                 "polyMesh::polyMesh::resetPrimitives\n"
00758                 "(\n"
00759                 "    const Xfer<pointField>&,\n"
00760                 "    const Xfer<faceList>&,\n"
00761                 "    const Xfer<labelList>& owner,\n"
00762                 "    const Xfer<labelList>& neighbour,\n"
00763                 "    const labelList& patchSizes,\n"
00764                 "    const labelList& patchStarts\n"
00765                 "    const bool validBoundary\n"
00766                 ")\n"
00767             )   << "no points or no cells in mesh" << endl;
00768         }
00769     }
00770 }
00771 
00772 
00773 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
00774 
00775 Foam::polyMesh::~polyMesh()
00776 {
00777     clearOut();
00778     resetMotion();
00779 }
00780 
00781 
00782 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00783 
00784 const Foam::fileName& Foam::polyMesh::dbDir() const
00785 {
00786     if (objectRegistry::dbDir() == defaultRegion)
00787     {
00788         return parent().dbDir();
00789     }
00790     else
00791     {
00792         return objectRegistry::dbDir();
00793     }
00794 }
00795 
00796 
00797 Foam::fileName Foam::polyMesh::meshDir() const
00798 {
00799     return dbDir()/meshSubDir;
00800 }
00801 
00802 
00803 const Foam::fileName& Foam::polyMesh::pointsInstance() const
00804 {
00805     return points_.instance();
00806 }
00807 
00808 
00809 const Foam::fileName& Foam::polyMesh::facesInstance() const
00810 {
00811     return faces_.instance();
00812 }
00813 
00814 
00815 const Foam::Vector<Foam::label>& Foam::polyMesh::geometricD() const
00816 {
00817     if (geometricD_.x() == 0)
00818     {
00819         calcDirections();
00820     }
00821 
00822     return geometricD_;
00823 }
00824 
00825 
00826 Foam::label Foam::polyMesh::nGeometricD() const
00827 {
00828     return cmptSum(geometricD() + Vector<label>::one)/2;
00829 }
00830 
00831 
00832 const Foam::Vector<Foam::label>& Foam::polyMesh::solutionD() const
00833 {
00834     if (solutionD_.x() == 0)
00835     {
00836         calcDirections();
00837     }
00838 
00839     return solutionD_;
00840 }
00841 
00842 
00843 Foam::label Foam::polyMesh::nSolutionD() const
00844 {
00845     return cmptSum(solutionD() + Vector<label>::one)/2;
00846 }
00847 
00848 
00849 // Add boundary patches. Constructor helper
00850 void Foam::polyMesh::addPatches
00851 (
00852     const List<polyPatch*>& p,
00853     const bool validBoundary
00854 )
00855 {
00856     if (boundaryMesh().size())
00857     {
00858         FatalErrorIn
00859         (
00860             "void polyMesh::addPatches(const List<polyPatch*>&, const bool)"
00861         )   << "boundary already exists"
00862             << abort(FatalError);
00863     }
00864 
00865     // Reset valid directions
00866     geometricD_ = Vector<label>::zero;
00867     solutionD_ = Vector<label>::zero;
00868 
00869     boundary_.setSize(p.size());
00870 
00871     // Copy the patch pointers
00872     forAll (p, pI)
00873     {
00874         boundary_.set(pI, p[pI]);
00875     }
00876 
00877     // parallelData depends on the processorPatch ordering so force
00878     // recalculation. Problem: should really be done in removeBoundary but
00879     // there is some info in parallelData which might be interesting inbetween
00880     // removeBoundary and addPatches.
00881     deleteDemandDrivenData(globalMeshDataPtr_);
00882 
00883     if (validBoundary)
00884     {
00885         // Calculate topology for the patches (processor-processor comms etc.)
00886         boundary_.updateMesh();
00887 
00888         // Calculate the geometry for the patches (transformation tensors etc.)
00889         boundary_.calcGeometry();
00890 
00891         boundary_.checkDefinition();
00892     }
00893 }
00894 
00895 
00896 // Add mesh zones. Constructor helper
00897 void Foam::polyMesh::addZones
00898 (
00899     const List<pointZone*>& pz,
00900     const List<faceZone*>& fz,
00901     const List<cellZone*>& cz
00902 )
00903 {
00904     if (pointZones().size() || faceZones().size() || cellZones().size())
00905     {
00906         FatalErrorIn
00907         (
00908             "void addZones\n"
00909             "(\n"
00910             "    const List<pointZone*>&,\n"
00911             "    const List<faceZone*>&,\n"
00912             "    const List<cellZone*>&\n"
00913             ")"
00914         )   << "point, face or cell zone already exists"
00915             << abort(FatalError);
00916     }
00917 
00918     // Point zones
00919     if (pz.size())
00920     {
00921         pointZones_.setSize(pz.size());
00922 
00923         // Copy the zone pointers
00924         forAll (pz, pI)
00925         {
00926             pointZones_.set(pI, pz[pI]);
00927         }
00928 
00929         pointZones_.writeOpt() = IOobject::AUTO_WRITE;
00930     }
00931 
00932     // Face zones
00933     if (fz.size())
00934     {
00935         faceZones_.setSize(fz.size());
00936 
00937         // Copy the zone pointers
00938         forAll (fz, fI)
00939         {
00940             faceZones_.set(fI, fz[fI]);
00941         }
00942 
00943         faceZones_.writeOpt() = IOobject::AUTO_WRITE;
00944     }
00945 
00946     // Cell zones
00947     if (cz.size())
00948     {
00949         cellZones_.setSize(cz.size());
00950 
00951         // Copy the zone pointers
00952         forAll (cz, cI)
00953         {
00954             cellZones_.set(cI, cz[cI]);
00955         }
00956 
00957         cellZones_.writeOpt() = IOobject::AUTO_WRITE;
00958     }
00959 }
00960 
00961 
00962 const Foam::pointField& Foam::polyMesh::points() const
00963 {
00964     if (clearedPrimitives_)
00965     {
00966         FatalErrorIn("const pointField& polyMesh::points() const")
00967             << "points deallocated"
00968             << abort(FatalError);
00969     }
00970 
00971     return points_;
00972 }
00973 
00974 
00975 const Foam::faceList& Foam::polyMesh::faces() const
00976 {
00977     if (clearedPrimitives_)
00978     {
00979         FatalErrorIn("const faceList& polyMesh::faces() const")
00980             << "faces deallocated"
00981             << abort(FatalError);
00982     }
00983 
00984     return faces_;
00985 }
00986 
00987 
00988 const Foam::labelList& Foam::polyMesh::faceOwner() const
00989 {
00990     return owner_;
00991 }
00992 
00993 
00994 const Foam::labelList& Foam::polyMesh::faceNeighbour() const
00995 {
00996     return neighbour_;
00997 }
00998 
00999 
01000 // Return old mesh motion points
01001 const Foam::pointField& Foam::polyMesh::oldPoints() const
01002 {
01003     if (!oldPointsPtr_)
01004     {
01005         if (debug)
01006         {
01007             WarningIn("const pointField& polyMesh::oldPoints() const")
01008                 << "Old points not available.  Forcing storage of old points"
01009                 << endl;
01010         }
01011 
01012         oldPointsPtr_ = new pointField(points_);
01013         curMotionTimeIndex_ = time().timeIndex();
01014     }
01015 
01016     return *oldPointsPtr_;
01017 }
01018 
01019 
01020 Foam::tmp<Foam::scalarField> Foam::polyMesh::movePoints
01021 (
01022     const pointField& newPoints
01023 )
01024 {
01025     if (debug)
01026     {
01027         Info<< "tmp<scalarField> polyMesh::movePoints(const pointField&) : "
01028             << " Moving points for time " << time().value()
01029             << " index " << time().timeIndex() << endl;
01030     }
01031 
01032     moving(true);
01033 
01034     // Pick up old points
01035     if (curMotionTimeIndex_ != time().timeIndex())
01036     {
01037         // Mesh motion in the new time step
01038         deleteDemandDrivenData(oldPointsPtr_);
01039         oldPointsPtr_ = new pointField(points_);
01040         curMotionTimeIndex_ = time().timeIndex();
01041     }
01042 
01043     points_ = newPoints;
01044 
01045     if (debug)
01046     {
01047         // Check mesh motion
01048         if (primitiveMesh::checkMeshMotion(points_, true))
01049         {
01050             Info<< "tmp<scalarField> polyMesh::movePoints"
01051                 << "(const pointField&) : "
01052                 << "Moving the mesh with given points will "
01053                 << "invalidate the mesh." << nl
01054                 << "Mesh motion should not be executed." << endl;
01055         }
01056     }
01057 
01058     points_.writeOpt() = IOobject::AUTO_WRITE;
01059     points_.instance() = time().timeName();
01060 
01061 
01062     tmp<scalarField> sweptVols = primitiveMesh::movePoints
01063     (
01064         points_,
01065         oldPoints()
01066     );
01067 
01068     // Adjust parallel shared points
01069     if (globalMeshDataPtr_)
01070     {
01071         globalMeshDataPtr_->movePoints(points_);
01072     }
01073 
01074     // Force recalculation of all geometric data with new points
01075 
01076     bounds_ = boundBox(points_);
01077     boundary_.movePoints(points_);
01078 
01079     pointZones_.movePoints(points_);
01080     faceZones_.movePoints(points_);
01081     cellZones_.movePoints(points_);
01082 
01083     // Reset valid directions (could change with rotation)
01084     geometricD_ = Vector<label>::zero;
01085     solutionD_ = Vector<label>::zero;
01086 
01087 
01088     // Hack until proper callbacks. Below are all the polyMeh MeshObjects with a
01089     // movePoints function.
01090 
01091     // pointMesh
01092     if (thisDb().foundObject<pointMesh>(pointMesh::typeName))
01093     {
01094         const_cast<pointMesh&>
01095         (
01096             thisDb().lookupObject<pointMesh>
01097             (
01098                 pointMesh::typeName
01099             )
01100         ).movePoints(points_);
01101     }
01102 
01103     return sweptVols;
01104 }
01105 
01106 
01107 // Reset motion by deleting old points
01108 void Foam::polyMesh::resetMotion() const
01109 {
01110     curMotionTimeIndex_ = 0;
01111     deleteDemandDrivenData(oldPointsPtr_);
01112 }
01113 
01114 
01115 // Return parallel info
01116 const Foam::globalMeshData& Foam::polyMesh::globalData() const
01117 {
01118     if (!globalMeshDataPtr_)
01119     {
01120         if (debug)
01121         {
01122             Pout<< "polyMesh::globalData() const : "
01123                 << "Constructing parallelData from processor topology" << nl
01124                 << "This needs the patch faces to be correctly matched"
01125                 << endl;
01126         }
01127         // Construct globalMeshData using processorPatch information only.
01128         globalMeshDataPtr_ = new globalMeshData(*this);
01129     }
01130 
01131     return *globalMeshDataPtr_;
01132 }
01133 
01134 
01135 // Remove all files and some subdirs (eg, sets)
01136 void Foam::polyMesh::removeFiles(const fileName& instanceDir) const
01137 {
01138     fileName meshFilesPath = thisDb().path()/instanceDir/meshDir();
01139 
01140     rm(meshFilesPath/"points");
01141     rm(meshFilesPath/"faces");
01142     rm(meshFilesPath/"owner");
01143     rm(meshFilesPath/"neighbour");
01144     rm(meshFilesPath/"cells");
01145     rm(meshFilesPath/"boundary");
01146     rm(meshFilesPath/"pointZones");
01147     rm(meshFilesPath/"faceZones");
01148     rm(meshFilesPath/"cellZones");
01149     rm(meshFilesPath/"meshModifiers");
01150     rm(meshFilesPath/"parallelData");
01151 
01152     // remove subdirectories
01153     if (isDir(meshFilesPath/"sets"))
01154     {
01155         rmDir(meshFilesPath/"sets");
01156     }
01157 }
01158 
01159 void Foam::polyMesh::removeFiles() const
01160 {
01161     removeFiles(instance());
01162 }
01163 
01164 
01165 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines