OpenFOAM logo
Open Source CFD Toolkit

tetPolyMeshCellDecomp.H

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-2005 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 the
00013     Free Software Foundation; either version 2 of the License, or (at your
00014     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, write to the Free Software Foundation,
00023     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
00024 
00025 Class
00026     tetPolyMeshCellDecomp
00027 
00028 Description
00029     The tetrahedral decomposition of the mesh based on the underlying polyMesh.
00030     Uses minimal storage - the tet decomposition is provided on a cell-by-cell
00031     basis.
00032 
00033     This class will be used as the underlying structure for FEM mesh motion
00034     and therefore it needs to be able to create and store tet poing-edge
00035     addressing for the solver.
00036 
00037     Note on addressing:
00038         The problem of addressing for Finite elements is two-fold. The
00039         first issue is addressing in the FEM matrix. The FEM coefficients are
00040         associated with edges in the tet decomposition. In order for the FOAM
00041         solver to work correctly, the coefficients need to be ordered in the
00042         upper triangular order. It is convinient to keep the order of edges
00043         and off-diagonal coefficients the same; therefore, the order of edges
00044         in the tet mesh is fixed by the ordering of points. The construction
00045         of the list from the edge list in polyMesh will be described in
00046         calcTetPolyMeshAddressing.C.
00047 
00048         The second issue is the construction of the FEM matrix.
00049         A fast method by Sanjay Mathur (Fluent Wavre 19/Sep/2003) is
00050         currently used:
00051 
00052         1) take a group of elements.  The discretisation will be done
00053         over a group.  A goup should be such that I can allocate a
00054         dense square matrix for all elements in the group.  Do NOT
00055         initialise arrays to zero.  As only a part of the group will
00056         be addressed, it needs to be zeroed out only for the entries
00057         it will have to use.  As the number of vertices per group will
00058         be variable, the local dense matrix changes for every group.
00059 
00060         2).  Allocate group-to-global and global-to-group arrays.  The
00061         second one is full of -1-s for vertices that are not in the
00062         group.  It will be initialised to -1 and is kept outside of
00063         the group.  This is done by going through all the vertices in
00064         the group and looking up a global-to-local index. If the index
00065         is -1, this global index has not been seen before: allocate
00066         next free local index to the vertex.
00067 
00068         3) For every element in the local-to-global, fidn the matrix
00069         row in ldu addressing.  For each element in a row, using
00070         global-to-local, find the local off-diagonal index and zero it
00071         out.  These are the only ones being used!
00072 
00073         4) Fill in local matrix based on the local vertex indices into
00074         the local dense matrix.
00075 
00076         5) Assemble global matrix by doing the following:
00077             - go through all group-to-global indices
00078             - for each group-to-global grab the row
00079             - for each non-zero element of the row, look up
00080             global-to-group.  This gives the x,y coordinates in the
00081             dense matrix and and insert it.
00082 
00083         6) Clear the global-to-local addressing by looking up the
00084         local-to-global indices and zeroing out only those entries.
00085 
00086 SourceFiles
00087     tetPolyMeshCellDecomp.C
00088     calcTetPolyMeshCellDecompGeometry.C
00089     calcTetPolyMeshCellDecompAddressing.C
00090     calcTetPolyMeshCellDecompParPointData.C
00091     addParallelPointPatchCellDecomp.C
00092 
00093 \*---------------------------------------------------------------------------*/
00094 
00095 #ifndef tetPolyMeshCellDecomp_H
00096 #define tetPolyMeshCellDecomp_H
00097 
00098 #include "tetFemSolution.H"
00099 #include "polyMesh.H"
00100 #include "primitiveMesh.H"
00101 #include "tetPolyBoundaryMeshCellDecomp.H"
00102 #include "tetCellList.H"
00103 #include "lduAddressingStore.H"
00104 #include "cellShapeList.H"
00105 #include "FieldFields.H"
00106 
00107 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00108 
00109 namespace Foam
00110 {
00111 
00112 // Class forward declarations
00113 template<class Type> class Matrix;
00114 
00115 class tetPolyMeshMapperCellDecomp;
00116 
00117 /*---------------------------------------------------------------------------*\
00118                            Class tetPolyMeshCellDecomp Declaration
00119 \*---------------------------------------------------------------------------*/
00120 
00121 class tetPolyMeshCellDecomp
00122 :
00123     public tetFemSolution
00124 {
00125     // Permanent data
00126 
00127         // Minimum mesh data
00128 
00129         //- Reference to underlying polyMesh
00130         const polyMesh& mesh_;
00131 
00132         //- Boundary mesh
00133         tetPolyBoundaryMeshCellDecomp boundary_;
00134 
00135         //- Offset in numbering to first cell centre
00136         label cellOffset_;
00137 
00138     // Demand-driven data
00139 
00140         //- Number of points
00141         mutable label nPoints_;
00142 
00143         //- Number of edges
00144         mutable label nEdges_;
00145 
00146         //- Number of tetrahedra
00147         mutable label nTets_;
00148 
00149         //- Ldu addressing
00150         mutable lduAddressingStore* lduPtr_;
00151 
00152         //- Max number of points per cell
00153         mutable label maxNPointsForCell_;
00154 
00155 
00156         // Parallel mesh analysis tools
00157 
00158             //- List of points shared between processor patches
00159             mutable labelList* parPointsPtr_;
00160 
00161             //- List of edges shared between processor patches
00162             mutable edgeList* parEdgesPtr_;
00163 
00164 
00165     // Private Member Functions
00166 
00167         //- Disallow default bitwise copy construct
00168         tetPolyMeshCellDecomp(const tetPolyMeshCellDecomp&);
00169 
00170         //- Disallow default bitwise assignment
00171         void operator=(const tetPolyMeshCellDecomp&);
00172 
00173 
00174     // Private member functions to calculate demand driven data
00175 
00176         //- Calculate addressing
00177         void calcAddressing() const;
00178 
00179         //- Check for and create a parallel point patch
00180         void addParallelPointPatch();
00181 
00182         //- Calculate points shared between parallel patches
00183         void calcParPointData() const;
00184         void clearOutParPointData() const;
00185 
00186         //- Clear out all data
00187         void clearOut() const;
00188 
00189 public:
00190 
00191     // Declare name of the class and its debug switch
00192     ClassName("tetPolyMesh");
00193 
00194     typedef tetPolyMeshCellDecomp Mesh;
00195     typedef tetPolyBoundaryMeshCellDecomp BoundaryMesh;
00196 
00197 
00198     // Constructors
00199 
00200         //- Construct from components
00201         explicit tetPolyMeshCellDecomp(const polyMesh& pMesh);
00202 
00203     // Destructor
00204 
00205         ~tetPolyMeshCellDecomp();
00206 
00207 
00208     // Member Functions
00209 
00210         //- Return reference to polyMesh
00211         const polyMesh& operator()() const
00212         {
00213             return mesh_;
00214         }
00215 
00216         //- Return number of points in decomposition
00217         label nPoints() const;
00218 
00219         //- Return number of edges in decomposition
00220         label nEdges() const;
00221 
00222         //- Return number of tetrahedra in decomposition
00223         label nTets() const;
00224 
00225         //- Return number of cells in polyhedral mesh
00226         label nCells() const
00227         {
00228             return mesh_.nCells();
00229         }
00230 
00231         //- Return number of tetrahedra in decomposition for a cell
00232         label nTetsForCell(const label cellID) const;
00233 
00234         //- Return number of edges in decomposition for a face
00235         label nEdgesForFace(const label faceID) const;
00236 
00237         //- Return number of edges in decomposition connected to a given point
00238         label nEdgesForPoint(const label pointID) const;
00239 
00240         //- Return ldu addressing
00241         const lduAddressing& ldu() const;
00242 
00243         //- Return cell offset
00244         label cellOffset() const
00245         {
00246             return cellOffset_;
00247         }
00248 
00249         //- Return list of edge labels coming out of a point
00250         labelList edgesForPoint(const label pointID) const;
00251 
00252         //- Return points
00253         tmp<pointField> points() const;
00254 
00255         //- Return complete list of cell shapes. All are tetrahedra
00256         cellShapeList tetCells() const;
00257 
00258         //- Return reference to boundary mesh
00259         const tetPolyBoundaryMeshCellDecomp& boundary() const
00260         {
00261             return boundary_;
00262         }
00263 
00264         //- Return tetrahedral decomposition for cell
00265         tetCellList tets(const label cellID) const;
00266 
00267         //- Return max number of tets in a cell
00268         label maxNPointsForCell() const;
00269 
00270         //- Fill buffer with tet decomposition addressing
00271         //  Used for FEM matrix assembly.
00272         //  localToGlobalBuffer - size to accept max number of vertices per cell
00273         //                        in the mesh
00274         //  globalToLocalBuffer - sized to total number of points in the mesh
00275         //                        and initialised to -1
00276         label addressing
00277         (
00278             const label cellID,
00279             labelList& localToGlobalBuffer,
00280             labelList& globalToLocalBuffer
00281         ) const;
00282 
00283         //- Clear global to local addressing
00284         void clearAddressing
00285         (
00286             const label cellID,
00287             const label nCellPoints,
00288             labelList& localToGlobalBuffer,
00289             labelList& globalToLocalBuffer
00290         ) const;
00291 
00292         //- Fill buffer with dot-products of shape functions
00293         // Used for FEM matrix assembly
00294         void gradNiDotGradNj
00295         (
00296             const label cellID,
00297             Matrix<scalar>& denseMatrix,                                       
00298             const labelList& globalToLocalBuffer
00299         ) const;
00300 
00301         //- Fill buffer with tensor products of shape functions
00302         // Used for FEM matrix assembly
00303         void gradNiGradNj
00304         (
00305             const label cellID,
00306             Matrix<tensor>& denseMatrix,                                       
00307             const labelList& globalToLocalBuffer
00308         ) const;
00309 
00310         //- Fill buffer with the volume integral distributed into vertices
00311         void volIntegral
00312         (
00313             const label cellID,
00314             scalarField& buffer,
00315             const labelList& globalToLocalBuffer
00316         ) const;
00317 
00318 
00319         // Parallel mesh analysis data
00320 
00321             //- Return parallel info
00322             const parallelInfo& parallelData() const
00323             {
00324                 return mesh_.parallelData();
00325             }
00326 
00327             //- Shared parallel points
00328             const labelList& parallelPoints() const;
00329 
00330             //- Shared parallel edges
00331             const edgeList& parallelEdges() const;
00332 
00333 
00334         // Edit
00335 
00336             //- Update mesh topology using the morph engine
00337             void updateMesh(const tetPolyMeshMapperCellDecomp& mapper);
00338 
00339 
00340     // Member Operators
00341 
00342         bool operator!=(const tetPolyMeshCellDecomp&) const;
00343         bool operator==(const tetPolyMeshCellDecomp&) const;
00344 };
00345 
00346 
00347 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00348 
00349 } // End namespace Foam
00350 
00351 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00352 
00353 #endif
00354 
00355 // ************************************************************************* //
For further information go to www.openfoam.org