fvcSmooth.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2022 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 "fvcSmooth.H"
27 #include "volFields.H"
28 #include "FvFaceCellWave.H"
29 #include "smoothData.H"
30 #include "sweepData.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
35 (
36  volScalarField& field,
37  const scalar coeff
38 )
39 {
40  const fvMesh& mesh = field.mesh();
41  scalar maxRatio = 1 + coeff;
42 
43  DynamicList<labelPair> changedPatchAndFaces(mesh.nFaces()/100 + 100);
44  DynamicList<smoothData> changedFacesInfo(changedPatchAndFaces.size());
45 
46  const labelUList& owner = mesh.owner();
47  const labelUList& neighbour = mesh.neighbour();
48 
49  forAll(owner, facei)
50  {
51  const label own = owner[facei];
52  const label nbr = neighbour[facei];
53 
54  // Check if owner value much larger than neighbour value or vice versa
55  if (field[own] > maxRatio*field[nbr])
56  {
57  changedPatchAndFaces.append(labelPair(-1, facei));
58  changedFacesInfo.append(smoothData(field[own]));
59  }
60  else if (field[nbr] > maxRatio*field[own])
61  {
62  changedPatchAndFaces.append(labelPair(-1, facei));
63  changedFacesInfo.append(smoothData(field[nbr]));
64  }
65  }
66 
67  // Insert all faces of coupled patches - FvFaceCellWave will correct them
68  forAll(mesh.boundary(), patchi)
69  {
70  const fvPatch& patch = mesh.boundary()[patchi];
71 
72  if (patch.coupled())
73  {
74  forAll(patch, patchFacei)
75  {
76  const label own = patch.faceCells()[patchFacei];
77 
78  changedPatchAndFaces.append(labelPair(patchi, patchFacei));
79  changedFacesInfo.append(smoothData(field[own]));
80  }
81  }
82  }
83 
84  changedPatchAndFaces.shrink();
85  changedFacesInfo.shrink();
86 
87  // Set initial field on cells
88  List<smoothData> cellData(mesh.nCells());
89  forAll(field, celli)
90  {
91  cellData[celli] = field[celli];
92  }
93 
94  // Set initial field on faces
95  List<smoothData> internalFaceData(mesh.nInternalFaces());
96  List<List<smoothData>> patchFaceData
97  (
99  sizesListList<List<List<smoothData>>>
100  (
102  listListSizes(mesh.boundary()),
103  smoothData()
104  )
105  );
106 
107  // Create track data
109  td.maxRatio = maxRatio;
110 
111  // Propagate information over whole domain
113  (
114  mesh,
115  changedPatchAndFaces,
116  changedFacesInfo,
117  internalFaceData,
118  patchFaceData,
119  cellData,
120  mesh.globalData().nTotalCells() + 1, // max iterations
121  td
122  );
123 
124  forAll(field, celli)
125  {
126  field[celli] = cellData[celli].value();
127  }
128 
130 }
131 
132 
134 (
135  volScalarField& field,
136  const volScalarField& alpha,
137  const label nLayers,
138  const scalar alphaDiff,
139  const scalar alphaMax,
140  const scalar alphaMin
141 )
142 {
143  const fvMesh& mesh = field.mesh();
144 
145  DynamicList<labelPair> changedPatchAndFaces(mesh.nFaces()/100 + 100);
146  DynamicList<smoothData> changedFacesInfo(changedPatchAndFaces.size());
147 
148  const labelUList& owner = mesh.owner();
149  const labelUList& neighbour = mesh.neighbour();
150 
151  forAll(owner, facei)
152  {
153  const label own = owner[facei];
154  const label nbr = neighbour[facei];
155 
156  // Check if owner value much larger than neighbour value or vice versa
157  if (mag(alpha[own] - alpha[nbr]) > alphaDiff)
158  {
159  changedPatchAndFaces.append(labelPair(-1, facei));
160  changedFacesInfo.append(smoothData(max(field[own], field[nbr])));
161  }
162  }
163 
164  // Insert all faces of coupled patches - FvFaceCellWave will correct them
165  forAll(mesh.boundary(), patchi)
166  {
167  const fvPatch& patch = mesh.boundary()[patchi];
168 
169  if (patch.coupled())
170  {
171  forAll(patch, patchFacei)
172  {
173  const label own = patch.faceCells()[patchFacei];
174 
175  const scalarField alphapn
176  (
177  alpha.boundaryField()[patchi].patchNeighbourField()
178  );
179 
180  if (mag(alpha[own] - alphapn[patchFacei]) > alphaDiff)
181  {
182  changedPatchAndFaces.append(labelPair(patchi, patchFacei));
183  changedFacesInfo.append(smoothData(field[own]));
184  }
185  }
186  }
187  }
188 
189  changedPatchAndFaces.shrink();
190  changedFacesInfo.shrink();
191 
192  // Set initial field on cells
193  List<smoothData> cellData(mesh.nCells());
194  forAll(field, celli)
195  {
196  cellData[celli] = field[celli];
197  }
198 
199  // Set initial field on faces
200  List<smoothData> internalFaceData(mesh.nInternalFaces());
201  List<List<smoothData>> patchFaceData
202  (
204  sizesListList<List<List<smoothData>>>
205  (
207  listListSizes(mesh.boundary()),
208  smoothData()
209  )
210  );
211 
212  // Create track data
214  td.maxRatio = 1.0;
215 
216  // Propagate information over whole domain
218  (
219  mesh,
220  internalFaceData,
221  patchFaceData,
222  cellData,
223  td
224  );
225 
226  smoothData.setFaceInfo(changedPatchAndFaces, changedFacesInfo);
227 
228  smoothData.iterate(nLayers);
229 
230  forAll(field, celli)
231  {
232  field[celli] = cellData[celli].value();
233  }
234 
236 }
237 
238 
239 void Foam::fvc::sweep
240 (
241  volScalarField& field,
242  const volScalarField& alpha,
243  const label nLayers,
244  const scalar alphaDiff
245 )
246 {
247  const fvMesh& mesh = field.mesh();
248 
249  DynamicList<labelPair> changedPatchAndFaces(mesh.nFaces()/100 + 100);
250  DynamicList<sweepData> changedFacesInfo(changedPatchAndFaces.size());
251 
252  const labelUList& owner = mesh.owner();
253  const labelUList& neighbour = mesh.neighbour();
254  const surfaceVectorField& Cf = mesh.Cf();
255  const surfaceVectorField::Boundary& CfBf = mesh.Cf().boundaryField();
256 
257  forAll(owner, facei)
258  {
259  const label own = owner[facei];
260  const label nbr = neighbour[facei];
261 
262  // Check if owner value much larger than neighbour value or vice versa
263  if (mag(alpha[own] - alpha[nbr]) > alphaDiff)
264  {
265  changedPatchAndFaces.append(labelPair(-1, facei));
266  changedFacesInfo.append
267  (
268  sweepData(max(field[own], field[nbr]), Cf[facei])
269  );
270  }
271  }
272 
273  // Insert all faces of coupled patches - FvFaceCellWave will correct them
274  forAll(mesh.boundary(), patchi)
275  {
276  const fvPatch& patch = mesh.boundary()[patchi];
277 
278  if (patch.coupled())
279  {
280  forAll(patch, patchFacei)
281  {
282  const label own = patch.faceCells()[patchFacei];
283 
284  const scalarField alphapn
285  (
286  alpha.boundaryField()[patchi].patchNeighbourField()
287  );
288 
289  if (mag(alpha[own] - alphapn[patchFacei]) > alphaDiff)
290  {
291  changedPatchAndFaces.append(labelPair(patchi, patchFacei));
292  changedFacesInfo.append
293  (
294  sweepData(field[own], CfBf[patchi][patchFacei])
295  );
296  }
297  }
298  }
299  }
300 
301  changedPatchAndFaces.shrink();
302  changedFacesInfo.shrink();
303 
304  // Set initial field on cells
305  List<sweepData> cellData(mesh.nCells());
306 
307  // Set initial field on faces
308  List<sweepData> internalFaceData(mesh.nInternalFaces());
309  List<List<sweepData>> patchFaceData
310  (
312  sizesListList<List<List<sweepData>>>
313  (
315  listListSizes(mesh.boundary()),
316  sweepData()
317  )
318  );
319 
320  // Propagate information over whole domain
322  (
323  mesh,
324  internalFaceData,
325  patchFaceData,
326  cellData
327  );
328 
329  sweepData.setFaceInfo(changedPatchAndFaces, changedFacesInfo);
330 
331  sweepData.iterate(nLayers);
332 
333  forAll(field, celli)
334  {
335  if (cellData[celli].valid(sweepData.data()))
336  {
337  field[celli] = max(field[celli], cellData[celli].value());
338  }
339  }
340 
342 }
343 
344 
345 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Wave propagation of information through grid. Every iteration information goes through one layer of c...
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
void spread(volScalarField &field, const volScalarField &alpha, const label nLayers, const scalar alphaDiff=0.2, const scalar alphaMax=0.99, const scalar alphaMin=0.01)
Definition: fvcSmooth.C:134
Class used to pass additional data in.
Definition: smoothData.H:61
const surfaceVectorField & Cf() const
Return face centres.
label nInternalFaces() const
const Boundary & boundaryField() const
Return const-reference to the boundary field.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
label nFaces() const
const TrackingData & data() const
Additional data to be passed into container.
label nCells() const
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:63
label nTotalCells() const
Return total number of cells in decomposed mesh.
fvMesh & mesh
void smooth(volScalarField &field, const scalar coeff)
Definition: fvcSmooth.C:35
Helper class used by fvc::sweep function.
Definition: sweepData.H:56
const labelUList & neighbour() const
Internal face neighbour.
Definition: fvMesh.H:423
Helper class used by the fvc::smooth and fvc::spread functions.
Definition: smoothData.H:55
void sweep(volScalarField &field, const volScalarField &alpha, const label nLayers, const scalar alphaDiff=0.2)
Definition: fvcSmooth.C:240
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
virtual const labelUList & faceCells() const
Return faceCells.
Definition: fvPatch.C:99
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:178
Pair< label > labelPair
Label pair.
Definition: labelPair.H:48
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1435
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
void setFaceInfo(const List< labelPair > &changedPatchAndFaces, const List< Type > &changedFacesInfo)
Set initial changed faces.
const Mesh & mesh() const
Return mesh.
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:417
label patchi
scalar maxRatio
Cut off distance.
Definition: smoothData.H:66
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
void correctBoundaryConditions()
Correct boundary field.
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual label iterate(const label maxIter)
Iterate until no changes or maxIter reached. Returns actual.
virtual bool coupled() const
Return true if this patch is coupled.
Definition: fvPatch.H:163
Provides functions smooth spread and sweep which use the FvFaceCellWave algorithm to smooth and redis...
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:800