Bezier.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) 2007-2019, 2022 PCOpt/NTUA
9 Copyright (C) 2013-2019, 2022 FOSS GP
10 Copyright (C) 2019 OpenCFD Ltd.
11-------------------------------------------------------------------------------
12License
13 This file is part of OpenFOAM.
14
15 OpenFOAM is free software: you can redistribute it and/or modify it
16 under the terms of the GNU General Public License as published by
17 the Free Software Foundation, either version 3 of the License, or
18 (at your option) any later version.
19
20 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
21 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
22 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
23 for more details.
24
25 You should have received a copy of the GNU General Public License
26 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
27
28\*---------------------------------------------------------------------------*/
29
30#include "volFields.H"
31#include "surfaceFields.H"
32#include "pointFields.H"
33#include "Bezier.H"
34#include "Pstream.H"
35
36// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
37
38namespace Foam
39{
40
41// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42
44
45// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
46
48:
49 mesh_(mesh),
50 dict_(dict),
51 nBezier_(dict_.subDict("Bezier").get<label>("nBezier")),
52 dxidXj_(nBezier_),
53 confineXmovement_
54 (
55 dict_.subDict("Bezier").get<boolList>("confineXmovement")
56 ),
57 confineYmovement_
58 (
59 dict_.subDict("Bezier").get<boolList>("confineYmovement")
60 ),
61 confineZmovement_
62 (
63 dict_.subDict("Bezier").get<boolList>("confineZmovement")
64 ),
65 confineMovement_(3, boolList(nBezier_, false)),
66 activeDesignVariables_(3*nBezier_)
67{
68 if
69 (
73 )
74 {
76 << "confineMovement lists sizes "
77 << confineXmovement_.size() << " "
78 << confineYmovement_.size() << " "
79 << confineZmovement_.size() << " "
80 << "are incompatible with nBezier " << nBezier_
81 << endl << endl
82 << exit(FatalError);
83 }
87
88 // Determine active design variables
89 label iActive(0);
90 for (label iDir = 0; iDir < 3; ++iDir)
91 {
92 for (label iCP = 0; iCP < nBezier_; ++iCP)
93 {
94 if (confineMovement_[iDir][iCP])
95 {
96 activeDesignVariables_[iActive++] = iDir*nBezier_ + iCP;
97 }
98 }
99 }
101
102 // Read bezier parameterization
103 for (label iCP = 0; iCP < nBezier_; ++iCP)
104 {
105 dxidXj_.set
106 (
107 iCP,
109 (
111 (
112 "dxidXj_"+name(iCP),
113 mesh_.time().timeName(),
114 mesh_,
117 ),
119 )
120 );
121 }
122}
123
124
125// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
126
127label Bezier::nBezier() const
128{
129 return nBezier_;
130}
131
132
134{
135 return dxidXj_;
136}
137
138
140{
141 return confineXmovement_;
142}
143
144
146{
147 return confineYmovement_;
148}
149
150
152{
153 return confineZmovement_;
154}
155
156
158{
159 return confineMovement_;
160}
161
162
164(
165 const label patchI,
166 const label cpI,
167 bool returnDimensionedNormalSens
168) const
169{
170 const fvPatch& patch = mesh_.boundary()[patchI];
171 const polyPatch& ppatch = patch.patch();
172
173 // Return field
174 auto tdndbSens = tmp<tensorField>::New(patch.size(), Zero);
175 auto& dndbSens = tdndbSens.ref();
176
177 // Auxiliary quantities
179 const label patchStart = ppatch.start();
180 const tensorField& dxdbInt = dxidXj_[cpI].primitiveField();
181
182 // Loop over patch faces
183 forAll(patch, fI)
184 {
185 const face& fGlobal = mesh_.faces()[fI + patchStart];
186 const pointField facePoints = fGlobal.points(mesh_.points());
187
188 // Loop over face points
189 tensorField facePointDerivs(facePoints.size(), Zero);
190 forAll(fGlobal, pI)
191 {
192 facePointDerivs[pI] = dxdbInt[fGlobal[pI]];
193 }
194
195 // Determine whether to return variance of dimensioned or unit normal
196 if (returnDimensionedNormalSens)
197 {
198 dndbSens[fI] =
200 (
201 facePoints,
202 facePointDerivs
203 )[1];
204 }
205 else
206 {
207 dndbSens[fI] =
209 (
210 facePoints,
211 facePointDerivs
212 )[2];
213 }
214 }
215 return tdndbSens;
216}
217
218
220(
221 const label patchI,
222 const label cpI,
223 const label idir,
224 bool returnDimensionedNormalSens
225) const
226{
227 const fvPatch& patch = mesh_.boundary()[patchI];
228 const polyPatch& ppatch = patch.patch();
229
230 // Return field
231 auto tdndbSens = tmp<vectorField>::New(patch.size(), Zero);
232 auto& dndbSens = tdndbSens.ref();
233
234 // Auxiliary quantities
236 const label patchStart = ppatch.start();
237 const tensorField& dxdbInt = dxidXj_[cpI].primitiveField();
238 vectorField dxdbDir(dxdbInt.size(), Zero);
239 unzipRow(dxdbInt, vector::components(idir), dxdbDir);
240
241 // Loop over patch faces
242 forAll(patch, fI)
243 {
244 const face& fGlobal = mesh_.faces()[fI + patchStart];
245 const pointField facePoints = fGlobal.points(mesh_.points());
246
247 // Loop over face points
248 vectorField facePointDerivs(facePoints.size(), Zero);
249 forAll(fGlobal, pI)
250 {
251 facePointDerivs[pI] = dxdbDir[fGlobal[pI]];
252 }
253
254 // Determine whether to return variance of dimensioned or unit normal
255 if (returnDimensionedNormalSens)
256 {
257 dndbSens[fI] =
259 (
260 facePoints, facePointDerivs
261 )[1];
262 }
263 else
264 {
265 dndbSens[fI] =
267 (
268 facePoints,
269 facePointDerivs
270 )[2];
271 }
272 }
273
274 return tdndbSens;
275}
276
277
279(
280 const label patchI,
281 const label cpI,
282 bool useChainRule
283) const
284{
285 const polyPatch& patch = mesh_.boundary()[patchI].patch();
286
287 // Return field
288 auto tdxdbFace = tmp<tensorField>::New(patch.size(), Zero);
289 auto& dxdbFace = tdxdbFace.ref();
290 if (useChainRule)
291 {
292 // Auxiliary quantities
294 const label patchStart = patch.start();
295 const tensorField& dxdbInt = dxidXj_[cpI].primitiveField();
296
297 // Loop over patch faces
298 forAll(patch, fI)
299 {
300 const face& fGlobal = mesh_.faces()[fI + patchStart];
301 const pointField facePoints = fGlobal.points(mesh_.points());
302
303 // Loop over face points
304 tensorField facePointDerivs(facePoints.size(), Zero);
305 forAll(fGlobal, pI)
306 {
307 facePointDerivs[pI] = dxdbInt[fGlobal[pI]];
308 }
309 dxdbFace[fI] =
311 (
312 facePoints,
313 facePointDerivs
314 )[0];
315 }
316 }
317 // Simple average of dxdb in points. Older and less accurate
318 else
319 {
322 (
323 dxidXj_[cpI].boundaryField()[patchI].patchInternalField()()
324 )();
325 }
326 return tdxdbFace;
327}
328
329
331(
332 const label patchI,
333 const label cpI,
334 const label idir,
335 bool useChainRule
336) const
337{
338 const polyPatch& patch = mesh_.boundary()[patchI].patch();
339
340 // Return field
341 auto tdxdbFace = tmp<vectorField>::New(patch.size(), Zero);
342 auto& dxdbFace = tdxdbFace.ref();
343
344 if (useChainRule)
345 {
346 // Auxiliary quantities
348 const label patchStart = patch.start();
349 const tensorField& dxdbInt = dxidXj_[cpI].primitiveField();
350 vectorField dxdbDir(dxdbInt.size(), Zero);
351 dxdbDir.replace(0, dxdbInt.component(3*idir));
352 dxdbDir.replace(1, dxdbInt.component(3*idir + 1));
353 dxdbDir.replace(2, dxdbInt.component(3*idir + 2));
354
355 // Loop over patch faces
356 forAll(patch, fI)
357 {
358 const face& fGlobal = mesh_.faces()[fI + patchStart];
359 const pointField facePoints = fGlobal.points(mesh_.points());
360
361 // Loop over face points
362 vectorField facePointDerivs(facePoints.size(), Zero);
363 forAll(fGlobal, pI)
364 {
365 facePointDerivs[pI] = dxdbDir[fGlobal[pI]];
366 }
367 dxdbFace[fI] =
369 (
370 facePoints,
371 facePointDerivs
372 )[0];
373 }
374 }
375 // Simple average of dxdb in points. Older and less accurate
376 else
377 {
379 vectorField dxdb(patch.nPoints(), Zero);
380 dxdb.replace
381 (
382 idir,
383 dxidXj_[cpI].boundaryField()[patchI].patchInternalField()().
384 component(0)
385 );
386 dxdbFace = patchInter.pointToFaceInterpolate(dxdb)();
387 }
388 return tdxdbFace;
389}
390
391
393(
394 const label globalFaceI,
395 const label cpI
396) const
397{
398 const face& faceI(mesh_.faces()[globalFaceI]);
399 tensorField fPoints_d(faceI.size(), Zero);
400 forAll(faceI, fpI)
401 {
402 fPoints_d[fpI] = dxidXj_[cpI].primitiveField()[faceI[fpI]];
403 }
404 return fPoints_d;
405}
406
407
409(
410 const label globalFaceI,
411 const label cpI,
412 const label idir
413) const
414{
415 const face& faceI(mesh_.faces()[globalFaceI]);
416 vectorField fPoints_d(faceI.size(), Zero);
417 forAll(faceI, fpI)
418 {
419 const tensor& dxdbTensor = dxidXj_[cpI].primitiveField()[faceI[fpI]];
420 fPoints_d[fpI].x() = dxdbTensor.component(3*idir);
421 fPoints_d[fpI].y() = dxdbTensor.component(3*idir + 1);
422 fPoints_d[fpI].z() = dxdbTensor.component(3*idir + 2);
423 }
424 return fPoints_d;
425}
426
427
429{
431}
432
433
434// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
435
436} // End namespace Foam
437
438// ************************************************************************* //
Calculation of adjoint based sensitivities for Bezier control points.
Definition: Bezier.H:60
boolList confineYmovement_
Definition: Bezier.H:83
PtrList< pointTensorField > dxidXj_
Definition: Bezier.H:80
const fvMesh & mesh_
Definition: Bezier.H:76
label nBezier() const
Number of Bezier control points.
Definition: Bezier.C:127
const boolList & confineZmovement() const
Confine z movement.
Definition: Bezier.C:151
label nBezier_
Definition: Bezier.H:79
boolListList confineMovement_
Definition: Bezier.H:85
const boolList & confineYmovement() const
Confine y movement.
Definition: Bezier.C:145
const boolListList & confineMovement() const
Info about confining movement in all directions.
Definition: Bezier.C:157
boolList confineXmovement_
Definition: Bezier.H:82
const boolList & confineXmovement() const
Confine x movement.
Definition: Bezier.C:139
labelList activeDesignVariables_
Definition: Bezier.H:86
tmp< tensorField > dxdbFace(const label patchI, const label cpI, bool useChainRule=true) const
dxdb tensor for a Bezier parameterized patch
Definition: Bezier.C:279
tmp< tensorField > dndbBasedSensitivities(const label patchI, const label cpI, bool returnDimensionedNormalSens=true) const
Definition: Bezier.C:164
boolList confineZmovement_
Definition: Bezier.H:84
tensorField facePoints_d(const label globalFaceI, const label cpI) const
For a given (global) face ID, return the change of the face points.
Definition: Bezier.C:393
PtrList< pointTensorField > & dxidXj()
dx/db tensor for all control points
Definition: Bezier.C:133
const labelList & getActiveDesignVariables() const
Return active design variables.
Definition: Bezier.C:428
void replace(const direction, const UList< cmptType > &)
Replace a component field of the field.
Definition: Field.C:557
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
Definition: Field.C:545
Generic GeometricField class.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
void setSize(const label n)
Alias for resize()
Definition: List.H:218
Interpolation class within a primitive patch. Allows interpolation from points to faces and vice vers...
tmp< Field< Type > > pointToFaceInterpolate(const Field< Type > &pf) const
Interpolate from points to faces.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
const Cmpt & component(const direction) const
Definition: VectorSpaceI.H:87
components
Component labeling enumeration.
Definition: Vector.H:81
Differentiation of the mesh data structure.
Definition: deltaBoundary.H:59
vectorField makeFaceCentresAndAreas_d(const pointField &p, const pointField &p_d)
Definition: deltaBoundary.C:72
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
pointField points(const UList< point > &pts) const
Return the points corresponding to this face.
Definition: faceI.H:87
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:712
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1108
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1083
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:75
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:364
A class for managing temporary objects.
Definition: tmp.H:65
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Namespace for OpenFOAM.
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
error FatalError
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
void unzipRow(const FieldField< Field, SymmTensor< Cmpt > > &input, const direction idx, FieldField< Field, Vector< Cmpt > > &result)
Extract a symmTensor field field row (x,y,z) == (0,1,2)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Foam::surfaceFields.