FreeFOAM The Cross-Platform CFD Toolkit
fvSurfaceMapper.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) 1991-2010 OpenCFD Ltd.
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 Description
25  FV surface mapper.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "fvSurfaceMapper.H"
30 #include <finiteVolume/fvMesh.H>
31 #include <OpenFOAM/mapPolyMesh.H>
32 #include <OpenFOAM/faceMapper.H>
33 
34 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
35 
36 void Foam::fvSurfaceMapper::calcAddressing() const
37 {
38  if
39  (
40  directAddrPtr_
41  || interpolationAddrPtr_
42  || weightsPtr_
43  || insertedObjectLabelsPtr_
44  )
45  {
46  FatalErrorIn("void fvSurfaceMapper::calcAddressing() const)")
47  << "Addressing already calculated"
48  << abort(FatalError);
49  }
50 
51  // Mapping
52 
53  const label oldNInternal = faceMap_.nOldInternalFaces();
54 
55  // Assemble the maps
56  if (direct())
57  {
58  // Direct mapping - slice to size
59  directAddrPtr_ =
60  new labelList
61  (
63  );
64  labelList& addr = *directAddrPtr_;
65 
66  // Adjust for creation of an internal face from a boundary face
67  forAll (addr, faceI)
68  {
69  if (addr[faceI] > oldNInternal)
70  {
71  addr[faceI] = 0;
72  }
73  }
74  }
75  else
76  {
77  // Interpolative mapping - slice to size
78  interpolationAddrPtr_ =
79  new labelListList
80  (
82  );
83  labelListList& addr = *interpolationAddrPtr_;
84 
85  weightsPtr_ =
86  new scalarListList
87  (
88  scalarListList::subList(faceMap_.weights(), size())
89  );
90  scalarListList& w = *weightsPtr_;
91 
92  // Adjust for creation of an internal face from a boundary face
93  forAll (addr, faceI)
94  {
95  if (max(addr[faceI]) >= oldNInternal)
96  {
97  addr[faceI] = labelList(1, 0);
98  w[faceI] = scalarList(1, 1.0);
99  }
100  }
101  }
102 
103  // Inserted objects
104 
105  // If there are, assemble the labels
106  if (insertedObjects())
107  {
108  const labelList& insFaces = faceMap_.insertedObjectLabels();
109 
110  insertedObjectLabelsPtr_ = new labelList(insFaces.size());
111  labelList& ins = *insertedObjectLabelsPtr_;
112 
113  label nIns = 0;
114 
115  forAll (insFaces, faceI)
116  {
117  // If the face is internal, keep it here
118  if (insFaces[faceI] < size())
119  {
120  ins[nIns] = insFaces[faceI];
121  nIns++;
122  }
123  }
124 
125  ins.setSize(nIns);
126  }
127  else
128  {
129  // No inserted objects
130  insertedObjectLabelsPtr_ = new labelList(0);
131  }
132 }
133 
134 
135 void Foam::fvSurfaceMapper::clearOut()
136 {
137  deleteDemandDrivenData(directAddrPtr_);
138  deleteDemandDrivenData(interpolationAddrPtr_);
139  deleteDemandDrivenData(weightsPtr_);
140 
141  deleteDemandDrivenData(insertedObjectLabelsPtr_);
142 }
143 
144 
145 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
146 
147 // Construct from components
148 Foam::fvSurfaceMapper::fvSurfaceMapper
149 (
150  const fvMesh& mesh,
151  const faceMapper& fMapper
152 )
153 :
154  mesh_(mesh),
155  faceMap_(fMapper),
156  directAddrPtr_(NULL),
157  interpolationAddrPtr_(NULL),
158  weightsPtr_(NULL),
159  insertedObjectLabelsPtr_(NULL)
160 {}
161 
162 
163 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
164 
166 {
167  clearOut();
168 }
169 
170 
171 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
172 
174 {
175  if (!direct())
176  {
178  (
179  "const unallocLabelList& fvSurfaceMapper::"
180  "directAddressing() const"
181  ) << "Requested direct addressing for an interpolative mapper."
182  << abort(FatalError);
183  }
184 
185  if (!directAddrPtr_)
186  {
187  calcAddressing();
188  }
189 
190  return *directAddrPtr_;
191 }
192 
193 
195 {
196  if (direct())
197  {
199  (
200  "const labelListList& fvSurfaceMapper::addressing() const"
201  ) << "Requested interpolative addressing for a direct mapper."
202  << abort(FatalError);
203  }
204 
205  if (!interpolationAddrPtr_)
206  {
207  calcAddressing();
208  }
209 
210  return *interpolationAddrPtr_;
211 }
212 
213 
215 {
216  if (direct())
217  {
219  (
220  "const scalarListList& fvSurfaceMapper::weights() const"
221  ) << "Requested interpolative weights for a direct mapper."
222  << abort(FatalError);
223  }
224 
225  if (!weightsPtr_)
226  {
227  calcAddressing();
228  }
229 
230  return *weightsPtr_;
231 }
232 
233 
235 {
236  if (!insertedObjectLabelsPtr_)
237  {
238  calcAddressing();
239  }
240 
241  return *insertedObjectLabelsPtr_;
242 }
243 
244 
245 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
246 
247 
248 // * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * //
249 
250 
251 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
252 
253 
254 // ************************ vim: set sw=4 sts=4 et: ************************ //