![]() |
|
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 tetPolyMeshFaceDecomp 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 calcPolyMeshAddressing.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 tetPolyMeshFaceDecomp.C 00088 calcTetPolyMeshFaceDecompGeometry.C 00089 calcTetPolyMeshFaceDecompAddressing.C 00090 calcTetPolyMeshFaceDecompParPointData.C 00091 addParallelPointPatchFaceDecomp.C 00092 00093 \*---------------------------------------------------------------------------*/ 00094 00095 #ifndef tetPolyMeshFaceDecomp_H 00096 #define tetPolyMeshFaceDecomp_H 00097 00098 #include "tetFemSolution.H" 00099 #include "polyMesh.H" 00100 #include "primitiveMesh.H" 00101 #include "tetPolyBoundaryMeshFaceDecomp.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 tetPolyMeshMapperFaceDecomp; 00116 00117 /*---------------------------------------------------------------------------*\ 00118 Class tetPolyMeshFaceDecomp Declaration 00119 \*---------------------------------------------------------------------------*/ 00120 00121 class tetPolyMeshFaceDecomp 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 tetPolyBoundaryMeshFaceDecomp boundary_; 00134 00135 //- Offset in numbering to first face centre 00136 label faceOffset_; 00137 00138 //- Offset in numbering to first cell centre 00139 label cellOffset_; 00140 00141 // Demand-driven data 00142 00143 //- Number of points 00144 mutable label nPoints_; 00145 00146 //- Number of edges 00147 mutable label nEdges_; 00148 00149 //- Number of tetrahedra 00150 mutable label nTets_; 00151 00152 //- Ldu addressing 00153 mutable lduAddressingStore* lduPtr_; 00154 00155 //- Max number of points per cell 00156 mutable label maxNPointsForCell_; 00157 00158 00159 // Parallel mesh analysis tools 00160 00161 //- List of points shared between processor patches 00162 mutable labelList* parPointsPtr_; 00163 00164 //- List of edges shared between processor patches 00165 mutable edgeList* parEdgesPtr_; 00166 00167 00168 // Private Member Functions 00169 00170 //- Disallow default bitwise copy construct 00171 tetPolyMeshFaceDecomp(const tetPolyMeshFaceDecomp&); 00172 00173 //- Disallow default bitwise assignment 00174 void operator=(const tetPolyMeshFaceDecomp&); 00175 00176 00177 // Private member functions to calculate demand driven data 00178 00179 //- Calculate addressing 00180 void calcAddressing() const; 00181 00182 //- Check for and create a parallel point patch 00183 void addParallelPointPatch(); 00184 00185 //- Calculate points shared between parallel patches 00186 void calcParPointData() const; 00187 void clearOutParPointData() const; 00188 00189 //- Clear out all data 00190 void clearOut() const; 00191 00192 public: 00193 00194 // Declare name of the class and its debug switch 00195 ClassName("tetPolyMesh"); 00196 00197 typedef tetPolyMeshFaceDecomp Mesh; 00198 typedef tetPolyBoundaryMeshFaceDecomp BoundaryMesh; 00199 00200 00201 // Constructors 00202 00203 //- Construct from components 00204 explicit tetPolyMeshFaceDecomp(const polyMesh& pMesh); 00205 00206 // Destructor 00207 00208 ~tetPolyMeshFaceDecomp(); 00209 00210 00211 // Member Functions 00212 00213 // Access 00214 00215 //- Return reference to polyMesh 00216 const polyMesh& operator()() const 00217 { 00218 return mesh_; 00219 } 00220 00221 //- Return number of points in decomposition 00222 label nPoints() const; 00223 00224 //- Return number of edges in decomposition 00225 label nEdges() const; 00226 00227 //- Return number of tetrahedra in decomposition 00228 label nTets() const; 00229 00230 //- Return number of cells in polyhedral mesh 00231 label nCells() const 00232 { 00233 return mesh_.nCells(); 00234 } 00235 00236 //- Return number of tetrahedra in decomposition for cell 00237 label nTetsForCell(const label cellID) const; 00238 00239 //- Return number of edges in decomposition for a face 00240 label nEdgesForFace(const label faceID) const; 00241 00242 //- Return number of edges in decomposition connected to a 00243 // given point 00244 label nEdgesForPoint(const label pointID) const; 00245 00246 //- Return ldu addressing 00247 const lduAddressing& ldu() const; 00248 00249 //- Return face offset 00250 label faceOffset() const 00251 { 00252 return faceOffset_; 00253 } 00254 00255 //- Return cell offset 00256 label cellOffset() const 00257 { 00258 return cellOffset_; 00259 } 00260 00261 //- Return list of edge labels coming out of a point 00262 labelList edgesForPoint(const label pointID) const; 00263 00264 //- Return points 00265 tmp<pointField> points() const; 00266 00267 //- Return complete list of cell shapes. All are tetrahedra 00268 cellShapeList tetCells() const; 00269 00270 //- Return reference to boundary mesh 00271 const tetPolyBoundaryMeshFaceDecomp& boundary() const 00272 { 00273 return boundary_; 00274 } 00275 00276 //- Return tetrahedral decomposition for cell 00277 tetCellList tets(const label cellID) const; 00278 00279 //- Return max number of tets in a cell 00280 label maxNPointsForCell() const; 00281 00282 //- Fill buffer with tet decomposition addressing 00283 // Used for FEM matrix assembly. 00284 // localToGlobalBuffer - sized to max number of vertices per cell 00285 // in the mesh 00286 // globalToLocalBuffer - sized to total number of points in 00287 // the mesh and initialised to -1 00288 label addressing 00289 ( 00290 const label cellID, 00291 labelList& localToGlobalBuffer, 00292 labelList& globalToLocalBuffer 00293 ) const; 00294 00295 //- Clear global to local addressing 00296 void clearAddressing 00297 ( 00298 const label cellID, 00299 const label nCellPoints, 00300 labelList& localToGlobalBuffer, 00301 labelList& globalToLocalBuffer 00302 ) const; 00303 00304 //- Fill buffer with dot-products of shape functions 00305 // Used for FEM matrix assembly 00306 void gradNiDotGradNj 00307 ( 00308 const label cellID, 00309 Matrix<scalar>& denseMatrix, 00310 const labelList& globalToLocalBuffer 00311 ) const; 00312 00313 //- Fill buffer with tensor products of shape functions 00314 // Used for FEM matrix assembly 00315 void gradNiGradNj 00316 ( 00317 const label cellID, 00318 Matrix<tensor>& denseMatrix, 00319 const labelList& globalToLocalBuffer 00320 ) const; 00321 00322 //- Fill buffer with the volume integral distributed into vertices 00323 void volIntegral 00324 ( 00325 const label cellID, 00326 scalarField& buffer, 00327 const labelList& globalToLocalBuffer 00328 ) const; 00329 00330 00331 // Parallel mesh analysis data 00332 00333 //- Return parallel info 00334 const parallelInfo& parallelData() const 00335 { 00336 return mesh_.parallelData(); 00337 } 00338 00339 //- Shared parallel points 00340 const labelList& parallelPoints() const; 00341 00342 //- Shared parallel edges 00343 const edgeList& parallelEdges() const; 00344 00345 00346 // Edit 00347 00348 //- Update mesh topology using the morph engine 00349 void updateMesh(const tetPolyMeshMapperFaceDecomp& mapper); 00350 00351 00352 // Member Operators 00353 00354 bool operator!=(const tetPolyMeshFaceDecomp&) const; 00355 bool operator==(const tetPolyMeshFaceDecomp&) const; 00356 }; 00357 00358 00359 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00360 00361 } // End namespace Foam 00362 00363 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00364 00365 #endif 00366 00367 // ************************************************************************* //