primitiveMeshCellCentresAndVols.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 Description
27  Efficient cell-centre calculation using face-addressing, face-centres and
28  face-areas.
29 
30 \*---------------------------------------------------------------------------*/
31 
32 #include "primitiveMesh.H"
33 #include "primitiveMeshTools.H"
34 
35 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
36 
38 {
39  if (debug)
40  {
41  Pout<< "primitiveMesh::calcCellCentresAndVols() : "
42  << "Calculating cell centres and cell volumes"
43  << endl;
44  }
45 
46  // It is an error to attempt to recalculate cellCentres
47  // if the pointer is already set
48  if (cellCentresPtr_ || cellVolumesPtr_)
49  {
51  << "Cell centres or cell volumes already calculated"
52  << abort(FatalError);
53  }
54 
55  // set the accumulated cell centre to zero vector
56  cellCentresPtr_ = new vectorField(nCells());
57  vectorField& cellCtrs = *cellCentresPtr_;
58 
59  // Initialise cell volumes to 0
60  cellVolumesPtr_ = new scalarField(nCells());
61  scalarField& cellVols = *cellVolumesPtr_;
62 
63  // Make centres and volumes
65  (
66  *this,
67  faceCentres(),
68  faceAreas(),
69  cellCtrs,
70  cellVols
71  );
72 
73  if (debug)
74  {
75  Pout<< "primitiveMesh::calcCellCentresAndVols() : "
76  << "Finished calculating cell centres and cell volumes"
77  << endl;
78  }
79 }
80 
81 
82 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
83 
85 {
86  if (!cellCentresPtr_)
87  {
88  //calcCellCentresAndVols();
89  const_cast<primitiveMesh&>(*this).updateGeom();
90  }
91 
92  return *cellCentresPtr_;
93 }
94 
95 
97 {
98  if (!cellVolumesPtr_)
99  {
100  //calcCellCentresAndVols();
101  const_cast<primitiveMesh&>(*this).updateGeom();
102  }
103 
104  return *cellVolumesPtr_;
105 }
106 
107 
108 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
primitiveMeshTools.H
primitiveMesh.H
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::Pout
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
Foam::vectorField
Field< vector > vectorField
Specialisation of Field<T> for vector.
Definition: primitiveFieldsFwd.H:54
Foam::primitiveMesh::nCells
label nCells() const noexcept
Number of mesh cells.
Definition: primitiveMeshI.H:96
Foam::primitiveMeshTools::makeCellCentresAndVols
static void makeCellCentresAndVols(const primitiveMesh &mesh, const vectorField &fCtrs, const vectorField &fAreas, vectorField &cellCtrs, scalarField &cellVols)
Calculate cell centres and volumes from face properties.
Definition: primitiveMeshTools.C:107
Foam::Field< vector >
Foam::primitiveMesh::calcCellCentresAndVols
void calcCellCentresAndVols() const
Calculate cell centres and volumes.
Definition: primitiveMeshCellCentresAndVols.C:37
cellVols
const scalarField & cellVols
Definition: temperatureAndPressureVariables.H:51
Foam::FatalError
error FatalError
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::primitiveMesh::cellVolumes
const scalarField & cellVolumes() const
Definition: primitiveMeshCellCentresAndVols.C:96
Foam::primitiveMesh::cellCentres
const vectorField & cellCentres() const
Definition: primitiveMeshCellCentresAndVols.C:84
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::primitiveMesh::faceCentres
const vectorField & faceCentres() const
Definition: primitiveMeshFaceCentresAndAreas.C:77
Foam::primitiveMesh::updateGeom
virtual void updateGeom()
Update all geometric data.
Definition: primitiveMesh.C:365
Foam::primitiveMesh::faceAreas
const vectorField & faceAreas() const
Definition: primitiveMeshFaceCentresAndAreas.C:89
Foam::primitiveMesh
Cell-face mesh analysis engine.
Definition: primitiveMesh.H:78