patchCorrectedInterpolationTemplates.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) 2015 OpenFOAM Foundation
9 Copyright (C) 2016 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
30#include "PointData.H"
31#include "PointEdgeWave.H"
34
35// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36
37template <class Type>
38void Foam::patchCorrectedInterpolation::interpolateType
39(
40 const GeometricField<Type, fvPatchField, volMesh>& cellDisplacement,
41 GeometricField<Type, pointPatchField, pointMesh>& pointDisplacement
42) const
43{
44 // Create an uncorrected field
45 GeometricField<Type, pointPatchField, pointMesh>
46 pointUncorrectedDisplacement
47 (
48 IOobject
49 (
50 "pointUncorrectedDisplacement",
51 mesh().time().timeName(),
52 mesh()
53 ),
54 pointDisplacement.mesh(),
55 pointDisplacement.dimensions(),
57 );
58
59 // Interpolate to the uncorrected field, overwriting the fixed value
60 // boundary conditions
61 pointUncorrectedDisplacement ==
62 volPointInterpolation::New(mesh()).interpolate
63 (
64 cellDisplacement,
66 (
67 pointUncorrectedDisplacement.boundaryField().size(),
69 )
70 );
71
72 // Set the point displacement to the uncorrected result everywhere except
73 // for on the boundary
74 pointDisplacement.primitiveFieldRef() =
75 pointUncorrectedDisplacement.primitiveField();
76 pointDisplacement.correctBoundaryConditions();
77
78 // Set the residual displacement as the difference between the boundary
79 // specification and the uncorrected solution
80 // (this reuses the uncorrected displacement field as the residual)
81 pointUncorrectedDisplacement ==
82 pointDisplacement - pointUncorrectedDisplacement;
83
84 // Interpolate the residual from the boundary into the field
85 interpolateDataFromPatchGroups(pointUncorrectedDisplacement);
86
87 // Add the residual to the point displacement and correct the boundary
88 pointDisplacement += pointUncorrectedDisplacement;
89 pointDisplacement.correctBoundaryConditions();
90}
91
92
93template <class Type>
94void Foam::patchCorrectedInterpolation::interpolateDataFromPatchGroups
95(
96 GeometricField<Type, pointPatchField, pointMesh>& data
97) const
98{
99 // Initialise
100 pointScalarField weight
101 (
102 IOobject
103 (
104 "weight",
105 mesh().time().timeName(),
106 mesh()
107 ),
108 data.mesh(),
111 );
112 data = dimensioned<Type>(data.dimensions(), Zero);
113
114 forAll(patchGroups_, patchGroupI)
115 {
116 // Distance and data for the current group
117 pointScalarField patchDistance
118 (
119 IOobject
120 (
121 "patchDistance",
122 mesh().time().timeName(),
123 mesh()
124 ),
125 data.mesh(),
126 dimensionedScalar(data.dimensions(), Zero),
128 );
129 GeometricField<Type, pointPatchField, pointMesh> patchData(data);
130
131 // Wave the data through the mesh
132 propagateDataFromPatchGroup
133 (
134 patchGroupI,
135 patchDistance,
136 patchData
137 );
138
139 // Calculate the weight and add to weighted sum
140 const scalarField patchWeight
141 (
142 1/max(sqr(patchDistance.primitiveField()), SMALL)
143 );
144 data.primitiveFieldRef() += patchWeight*patchData.primitiveField();
145 weight.primitiveFieldRef() += patchWeight;
146 }
147
148 // Complete the average
149 data /= weight;
150}
151
152
153template <class Type>
154void Foam::patchCorrectedInterpolation::propagateDataFromPatchGroup
155(
156 const label patchGroupi,
158 GeometricField<Type, pointPatchField, pointMesh>& data
159) const
160{
161 const labelList& patchGroup(patchGroups_[patchGroupi]);
162
163 // Get the size of the seed info
164 label nSeedInfo(0);
165 forAll(patchGroup, patchGroupi)
166 {
167 const label patchi(patchGroup[patchGroupi]);
168
169 nSeedInfo += data.boundaryField()[patchi].size();
170 }
171
172 // Generate the seed labels and info
173 labelList seedLabels(nSeedInfo);
174 List<PointData<Type>> seedInfo(nSeedInfo);
175 nSeedInfo = 0;
176 forAll(patchGroup, patchGroupi)
177 {
178 const label patchi(patchGroup[patchGroupi]);
179
180 pointPatchField<Type>& patchDataField(data.boundaryFieldRef()[patchi]);
181
182 patchDataField.updateCoeffs();
183
184 const pointPatch& patch(patchDataField.patch());
185 const Field<Type> patchData(patchDataField.patchInternalField());
186
187 forAll(patch.meshPoints(), patchPointi)
188 {
189 const label pointi(patch.meshPoints()[patchPointi]);
190
191 seedLabels[nSeedInfo] = pointi;
192
193 seedInfo[nSeedInfo] =
194 PointData<Type>
195 (
196 mesh().points()[pointi],
197 0,
198 patchData[patchPointi]
199 );
200
201 nSeedInfo++;
202 }
203 }
204
205 // Wave the data through the mesh
206 List<PointData<Type>> allPointInfo(mesh().nPoints());
207 List<PointData<Type>> allEdgeInfo(mesh().nEdges());
208 PointEdgeWave<PointData<Type>>
209 (
210 mesh(),
211 seedLabels,
212 seedInfo,
213 allPointInfo,
214 allEdgeInfo,
216 );
217
218 // Copy result into the fields
219 forAll(allPointInfo, pointi)
220 {
221 distance[pointi] = sqrt(allPointInfo[pointi].distSqr());
222 data[pointi] = allPointInfo[pointi].data();
223 }
224}
225
226// ************************************************************************* //
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
label nTotalPoints() const noexcept
Return total number of points in decomposed mesh. Not.
const fvMesh & mesh() const
Return const-reference to the mesh.
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1310
dynamicFvMesh & mesh
const labelList nEdges(UPstream::listGatherValues< label >(aMesh.nEdges()))
const pointField & points
label nPoints
word timeName
Definition: getTimeIndex.H:3
const std::string patch
OpenFOAM patch number as a std::string.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
scalar distance(const vector &p1, const vector &p2)
Definition: curveTools.C:12
List< word > wordList
A List of words.
Definition: fileName.H:63
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
const dimensionSet dimless
Dimensionless.
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField
List< label > labelList
A List of labels.
Definition: List.H:66
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
dimensionedScalar sqrt(const dimensionedScalar &ds)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
static const char *const typeName
The type name used in ensight case files.