singleLayerRegion.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 | Copyright (C) 2011-2016 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "singleLayerRegion.H"
27 #include "fvMesh.H"
28 #include "Time.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace regionModels
36 {
37  defineTypeNameAndDebug(singleLayerRegion, 0);
38 }
39 }
40 
41 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
42 
43 void Foam::regionModels::singleLayerRegion::constructMeshObjects()
44 {
45  // construct patch normal vectors
46  nHatPtr_.reset
47  (
48  new volVectorField
49  (
50  IOobject
51  (
52  "nHat",
53  time_.timeName(),
54  regionMesh(),
56  NO_WRITE
57  ),
58  regionMesh(),
59  dimensionedVector("zero", dimless, Zero),
61  )
62  );
63 
64  // construct patch areas
65  magSfPtr_.reset
66  (
67  new volScalarField
68  (
69  IOobject
70  (
71  "magSf",
72  time_.timeName(),
73  regionMesh(),
75  NO_WRITE
76  ),
77  regionMesh(),
78  dimensionedScalar("zero", dimArea, 0.0),
80  )
81  );
82 }
83 
84 
85 void Foam::regionModels::singleLayerRegion::initialise()
86 {
87  if (debug)
88  {
89  Pout<< "singleLayerRegion::initialise()" << endl;
90  }
91 
92  label nBoundaryFaces = 0;
93  const polyBoundaryMesh& rbm = regionMesh().boundaryMesh();
94  volVectorField& nHat = nHatPtr_();
95  volScalarField& magSf = magSfPtr_();
96  forAll(intCoupledPatchIDs_, i)
97  {
98  const label patchi = intCoupledPatchIDs_[i];
99  const polyPatch& pp = rbm[patchi];
100  const labelList& fCells = pp.faceCells();
101 
102  nBoundaryFaces += fCells.size();
103 
104  UIndirectList<vector>(nHat, fCells) = pp.faceNormals();
105  UIndirectList<scalar>(magSf, fCells) = mag(pp.faceAreas());
106  }
107  nHat.correctBoundaryConditions();
108  magSf.correctBoundaryConditions();
109 
110  if (nBoundaryFaces != regionMesh().nCells())
111  {
113  << "Number of primary region coupled boundary faces not equal to "
114  << "the number of cells in the local region" << nl << nl
115  << "Number of cells = " << regionMesh().nCells() << nl
116  << "Boundary faces = " << nBoundaryFaces << nl
117  << abort(FatalError);
118  }
119 
120  scalarField passiveMagSf(magSf.size(), 0.0);
121  passivePatchIDs_.setSize(intCoupledPatchIDs_.size(), -1);
122  forAll(intCoupledPatchIDs_, i)
123  {
124  const label patchi = intCoupledPatchIDs_[i];
125  const polyPatch& ppIntCoupled = rbm[patchi];
126  if (ppIntCoupled.size() > 0)
127  {
128  label cellId = rbm[patchi].faceCells()[0];
129  const cell& cFaces = regionMesh().cells()[cellId];
130 
131  label facei = ppIntCoupled.start();
132  label faceO = cFaces.opposingFaceLabel(facei, regionMesh().faces());
133 
134  label passivePatchi = rbm.whichPatch(faceO);
135  passivePatchIDs_[i] = passivePatchi;
136  const polyPatch& ppPassive = rbm[passivePatchi];
137  UIndirectList<scalar>(passiveMagSf, ppPassive.faceCells()) =
138  mag(ppPassive.faceAreas());
139  }
140  }
141 
142  Pstream::listCombineGather(passivePatchIDs_, maxEqOp<label>());
143  Pstream::listCombineScatter(passivePatchIDs_);
144 
145  magSf.field() = 0.5*(magSf + passiveMagSf);
146  magSf.correctBoundaryConditions();
147 }
148 
149 
150 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
151 
153 {
154  return regionModel::read();
155 }
156 
157 
158 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
159 
160 Foam::regionModels::singleLayerRegion::singleLayerRegion
161 (
162  const fvMesh& mesh,
163  const word& regionType
164 )
165 :
166  regionModel(mesh, regionType),
167  nHatPtr_(NULL),
168  magSfPtr_(NULL),
169  passivePatchIDs_()
170 {}
171 
172 
173 Foam::regionModels::singleLayerRegion::singleLayerRegion
174 (
175  const fvMesh& mesh,
176  const word& regionType,
177  const word& modelName,
178  bool readFields
179 )
180 :
181  regionModel(mesh, regionType, modelName, false),
182  nHatPtr_(NULL),
183  magSfPtr_(NULL),
184  passivePatchIDs_()
185 {
186  if (active_)
187  {
188  constructMeshObjects();
189  initialise();
190 
191  if (readFields)
192  {
193  read();
194  }
195  }
196 }
197 
198 
199 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
200 
202 {}
203 
204 
205 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
206 
208 {
209  if (!nHatPtr_.valid())
210  {
212  << "Region patch normal vectors not available"
213  << abort(FatalError);
214  }
215 
216  return nHatPtr_();
217 }
218 
219 
221 {
222  if (!magSfPtr_.valid())
223  {
225  << "Region patch areas not available"
226  << abort(FatalError);
227  }
228 
229  return magSfPtr_();
230 }
231 
232 
233 const Foam::labelList&
235 {
236  return passivePatchIDs_;
237 }
238 
239 
240 // ************************************************************************* //
virtual const volScalarField & magSf() const
Return the face area magnitudes / [m2].
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
static const char *const typeName
Definition: Field.H:94
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
defineTypeNameAndDebug(regionModel, 0)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
bool read(const char *, int32_t &)
Definition: int32IO.C:85
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual bool read()
Read control parameters from dictionary.
Definition: regionModel.C:171
List< label > labelList
A List of labels.
Definition: labelList.H:56
static const zero Zero
Definition: zero.H:91
errorManip< error > abort(error &err)
Definition: errorManip.H:131
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
static const char nl
Definition: Ostream.H:262
virtual const volVectorField & nHat() const
Return the patch normal vectors.
label patchi
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
label cellId
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual const labelList & passivePatchIDs() const
Return the list of patch IDs opposite to internally.
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
Namespace for OpenFOAM.
virtual bool read()
Read control parameters from dictionary.