volBSplinesBase.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-2020 PCOpt/NTUA
9 Copyright (C) 2013-2020 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 "volBSplinesBase.H"
31
32// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33
34namespace Foam
35{
36
37// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38
40
41// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42
44(
45 const fvMesh& mesh
46)
47:
49 volume_(0),
50 activeDesignVariables_(0)
51{
52 const dictionary NURBSdict
53 (
55 (
57 (
58 "dynamicMeshDict",
59 mesh.time().constant(),
60 mesh,
63 false
64 )
65 ).subDict("volumetricBSplinesMotionSolverCoeffs")
66 );
67 // Read box names and allocate size
68 wordList controlBoxes(NURBSdict.toc());
69 volume_.setSize(controlBoxes.size());
70
71 // Populate NURBS volumes
72 label iBox(0);
73 for (const word& boxName : controlBoxes)
74 {
75 if (NURBSdict.isDict(boxName))
76 {
77 volume_.set
78 (
79 iBox,
81 (
82 NURBSdict.subDict(boxName),
83 mesh,
84 true
85 )
86 );
87 volume_[iBox].writeParamCoordinates();
88 iBox++;
89 }
90 }
91 volume_.setSize(iBox);
92
93 // Determine active design variables
95 label iActive(0);
96 const labelList startCpID(getStartCpID());
97 forAll(volume_, boxI)
98 {
99 const label start(3*startCpID[boxI]);
100 const boolList& isActiveVar = volume_[boxI].getActiveDesignVariables();
101 forAll(isActiveVar, varI)
102 {
103 if (isActiveVar[varI])
104 {
105 activeDesignVariables_[iActive++] = start + varI;
106 }
107 }
108 }
110}
111
112
113// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
114
116{
117 return volume_;
118}
119
120
122{
123 return volume_;
124}
125
126
127const NURBS3DVolume& volBSplinesBase::box(const label boxI) const
128{
129 return volume_[boxI];
130}
131
132
134{
135 return volume_[boxI];
136}
137
138
139const vectorField& volBSplinesBase::getControlPoints(const label& iNURB) const
140{
141 return volume_[iNURB].getControlPoints();
142}
143
144
146{
147 vectorField totalCPs(0);
148 forAll(volume_, iNURB)
149 {
150 totalCPs.append(volume_[iNURB].getControlPoints());
151 }
152
153 return totalCPs;
154}
155
156
158{
159 label nCPs(0);
160 forAll(volume_, iNURB)
161 {
162 nCPs += volume_[iNURB].getControlPoints().size();
163 }
164
165 return nCPs;
166}
167
168
170{
171 return volume_.size();
172}
173
174
176{
177 // Allocate an extra entry to track in which box a CP might be
178 labelList startID(getNumberOfBoxes() + 1);
179 startID[0] = 0;
180 forAll(volume_, iNURB)
181 {
182 startID[iNURB+1] =
183 startID[iNURB] + volume_[iNURB].getControlPoints().size();
184 }
185
186 return startID;
187}
188
189
190label volBSplinesBase::findBoxID(const label cpI) const
191{
192 const labelList startCPID(getStartCpID());
193 for (label iBox = 0; iBox < startCPID.size() - 1 ; ++iBox)
194 {
195 if (cpI >= startCPID[iBox] || cpI < startCPID[iBox + 1])
196 {
197 return iBox;
198 }
199 }
200
202 << "Invalid control point ID " << cpI << endl
203 << exit(FatalError);
204 return -1;
205}
206
207
209{
211}
212
213
215(
216 const vectorField& controlPointsMovement,
217 const labelList& patchesToBeMoved
218)
219{
220 scalar maxDisplacement(0);
221 label pastControlPoints(0);
222 forAll(volume_, iNURB)
223 {
224 const label nb(volume_[iNURB].getControlPoints().size());
225 vectorField localControlPointsMovement(nb, Zero);
226
227 // Set localControlPointsMovement
228 forAll(localControlPointsMovement, iCPM)
229 {
230 localControlPointsMovement[iCPM] =
231 controlPointsMovement[pastControlPoints + iCPM];
232 }
233
234 maxDisplacement = max
235 (
236 maxDisplacement,
237 volume_[iNURB].computeMaxBoundaryDisplacement
238 (
239 localControlPointsMovement,
240 patchesToBeMoved
241 )
242 );
243
244 pastControlPoints += nb;
245 }
246
247 return maxDisplacement;
248}
249
250
252(
253 vectorField& controlPointsMovement
254) const
255{
256 label pastControlPoints(0);
257 forAll(volume_, iNURB)
258 {
259 const label nb(volume_[iNURB].getControlPoints().size());
260 vectorField localControlPointsMovement(nb, Zero);
261
262 // Set localControlPointsMovement
263 forAll(localControlPointsMovement, iCPM)
264 {
265 localControlPointsMovement[iCPM] =
266 controlPointsMovement[pastControlPoints + iCPM];
267 }
268
269 volume_[iNURB].boundControlPointMovement(localControlPointsMovement);
270
271 // Transfer bounding back to controlPointMovement
272 forAll(localControlPointsMovement, iCPM)
273 {
274 controlPointsMovement[pastControlPoints + iCPM] =
275 localControlPointsMovement[iCPM];
276 }
277
278 pastControlPoints += nb;
279 }
280}
281
282
284(
285 const vectorField& controlPointsMovement
286)
287{
288 label pastControlPoints(0);
289 forAll(volume_, iNURB)
290 {
291 const label nb(volume_[iNURB].getControlPoints().size());
292 vectorField localControlPointsMovement(nb, Zero);
293
294 // Set localControlPointsMovement
295 forAll(localControlPointsMovement, iCPM)
296 {
297 localControlPointsMovement[iCPM] =
298 controlPointsMovement[pastControlPoints + iCPM];
299 }
300
301 const vectorField newCps
302 (
303 volume_[iNURB].getControlPoints()
304 + localControlPointsMovement
305 );
306
307 volume_[iNURB].setControlPoints(newCps);
308
309 pastControlPoints += nb;
310 }
311}
312
313
315{
316 for (const NURBS3DVolume& box : volume_)
317 {
318 box.writeCps("cpsBsplines" + mesh_.time().timeName());
319 }
320}
321
322
324{
325 // Does nothing
326 return true;
327}
328
329
331{
332 // Does nothing
333}
334
335
336// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
337
338} // End namespace Foam
339
340// ************************************************************************* //
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
@ MUST_READ_IF_MODIFIED
Definition: IOobject.H:180
void setSize(const label n)
Alias for resize()
Definition: List.H:218
void append(const T &val)
Append an element at the end of the list.
Definition: ListI.H:175
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
Definition: MeshObject.H:91
NURBS3DVolume morpher. Includes support functions for gradient computations Base class providing supp...
Definition: NURBS3DVolume.H:76
const vectorField & getControlPoints() const
Get control points.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
const word & constant() const
Return constant name.
Definition: TimePathsI.H:96
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
Definition: dictionary.C:460
bool isDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Check if entry is found and is a sub-dictionary.
Definition: dictionaryI.H:147
wordList toc() const
Return the table of contents.
Definition: dictionary.C:602
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:162
void updateMesh()
Update for new mesh topology.
Creates the parallel distribution map by describing the source and target objects using box shapes.
Definition: box.H:62
Class constructing a number of volumetric B-Splines boxes, read from dynamicMeshDict....
virtual bool movePoints()
Dummy function required by MeshObject.
void moveControlPoints(const vectorField &controlPointsMovement)
Move control points. No effect on mesh.
const vectorField & getControlPoints(const label &iNURB) const
Get reference to control points.
const PtrList< NURBS3DVolume > & boxes() const
Get const reference to the vol. B-splines boxes.
vectorField getAllControlPoints() const
Get control points from all boxes.
const labelList & getActiveDesignVariables() const
Get active design variables.
void writeControlPoints() const
Write control points to constant and optimisation folders.
scalar computeMaxBoundaryDisplacement(const vectorField &controlPointsMovement, const labelList &patchesToBeMoved)
labelList getStartCpID() const
Get start CP ID for each box.
label findBoxID(const label cpI) const
Find box of certain control point.
NURBS3DVolume & boxRef(const label boxI)
Get non-const reference to a specific box.
PtrList< NURBS3DVolume > & boxesRef()
Get non-const reference to the vol. B-splines boxes.
labelList activeDesignVariables_
Active design variables numbering for all boxes.
label getTotalControlPointsNumber() const
Get cumulative number of control points from all boxes.
void boundControlPointMovement(vectorField &controlPointsMovement) const
Bound control points movement.
PtrList< NURBS3DVolume > volume_
List with volumetric B-splines boxes.
label getNumberOfBoxes() const
Get number of boxes.
A class for handling words, derived from Foam::string.
Definition: word.H:68
#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.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333