optMeshMovementVolumetricBSplinesExternalMotionSolver.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
32
33
34// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35
36namespace Foam
37{
39 (
41 0
42 );
44 (
48 );
49}
50
51
52// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
53
56(
57 const scalarField& correction
58)
59{
63
64 for (label iCP = 0; iCP < nCPs; iCP++)
65 {
66 cpMovement_[iCP].x() = correction[3*iCP];
67 cpMovement_[iCP].y() = correction[3*iCP + 1];
68 cpMovement_[iCP].z() = correction[3*iCP + 2];
69 }
70
71 // Bound control point movement for non-active CPs
73
74 // Compute boundary movement
75 label passedCPs(0);
77 forAll(boxes, iNURB)
78 {
79 const label nb = boxes[iNURB].getControlPoints().size();
80 for (label cpI = 0; cpI < nb; ++cpI)
81 {
82 const label globalCP = passedCPs + cpI;
83 forAll(patchIDs_, pI)
84 {
85 const label patchI = patchIDs_[pI];
86 vectorField boundaryMovement
87 (
88 boxes[iNURB].patchDxDb(patchI, cpI)
89 & cpMovement_[globalCP]
90 );
91 dx_.boundaryField()[patchI].addToInternalField
92 (
94 boundaryMovement
95 );
96 }
97 }
98
99 // Increment number of passed sensitivities
100 passedCPs += nb;
101 }
102}
103
104
105// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
106
107Foam::optMeshMovementVolumetricBSplinesExternalMotionSolver::
108optMeshMovementVolumetricBSplinesExternalMotionSolver
109(
110 fvMesh& mesh,
111 const dictionary& dict,
112 const labelList& patchIDs
113)
114:
115 optMeshMovement(mesh, dict, patchIDs),
116 volBSplinesBase_
117 (
119 ),
120 dx_
121 (
123 (
124 "dx",
125 mesh.time().timeName(),
126 mesh,
127 IOobject::NO_READ,
128 IOobject::NO_WRITE,
129 false
130 ),
133 ),
134 cpMovement_(volBSplinesBase_.getTotalControlPointsNumber(), Zero)
135{}
136
137
138// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
139
141{
142 // Compute boundary movement
143 computeBoundaryMovement(correction_);
144
145 // Set boundary movement of motion solver
146 displMethodPtr_->setMotionField(dx_);
147
148 // Positions of control points have not changed since only the boundary dx
149 // has been computed.
150 // Use correction to update them
151 volBSplinesBase_.moveControlPoints(cpMovement_);
152
153 // Write control points to files
154 volBSplinesBase_.writeControlPoints();
155
156 // Move the mesh and check quality
158}
159
160
161Foam::scalar
163(
165)
166{
167 computeBoundaryMovement(correction);
168
169 // Get maximum boundary movement
170 scalar maxDisplacement = gMax(mag(dx_.primitiveField()));
171
172 // Compute eta value
173 Info<< "maxAllowedDisplacement/maxDisplacement \t"
174 << getMaxAllowedDisplacement() << "/" << maxDisplacement << endl;
175 scalar eta = getMaxAllowedDisplacement() / maxDisplacement;
176 Info<< "Setting eta value to " << eta << endl;
177
178 return eta;
179}
180
181
184const
185{
186 return volBSplinesBase_.getActiveDesignVariables();
187}
188
189
190// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
label size() const noexcept
The number of elements in the list.
Definition: UPtrListI.H:106
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Converts NURBS volume control points update to actual mesh movement. Internal points are moved based ...
pointVectorField dx_
Boundary movement due to change of NURBS control points.
virtual scalar computeEta(const scalarField &correction)
Compute eta value based on max displacement.
volBSplinesBase & volBSplinesBase_
Reference to underlaying volumetric B-Splines morpher.
Abstract base class for translating an update of the design variables into mesh movement.
virtual void moveMesh()
labelList patchIDs_
IDs of patches to be moved.
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:55
Class constructing a number of volumetric B-Splines boxes, read from dynamicMeshDict....
PtrList< NURBS3DVolume > & boxesRef()
Get non-const reference to the vol. B-splines boxes.
label getTotalControlPointsNumber() const
Get cumulative number of control points from all boxes.
void boundControlPointMovement(vectorField &controlPointsMovement) const
Bound control points movement.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
dynamicFvMesh & mesh
word timeName
Definition: getTimeIndex.H:3
Namespace for OpenFOAM.
const dimensionSet dimless
Dimensionless.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Type gMax(const FieldField< Field, Type > &f)
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
A non-counting (dummy) refCount.
Definition: refCount.H:59