42 volumetricBSplinesDesignVariables,
59 return volBSplinesBase_.getActiveDesignVariables();
66 labelList activeVarsInCPs = volBSplinesBase_.getActiveDesignVariables();
70 activeDesignVariables_ =
71 constraint_().computeActiveDesignVariables(activeVarsInCPs);
80 scalarField cpsScalar(constraint_().designVariablesToControlPoints(dvs));
90 cpI.x() = cpsScalar[varID++];
91 cpI.y() = cpsScalar[varID++];
92 cpI.z() = cpsScalar[varID++];
94 boxI.setControlPoints(cps);
102 scalarField cpsScalar(3*volBSplinesBase_.getTotalControlPointsNumber());
103 const PtrList<NURBS3DVolume>& boxes = volBSplinesBase_.boxes();
105 for (
const NURBS3DVolume& boxI : boxes)
108 for (
const vector& cpI : cps)
110 cpsScalar[varID++] = cpI.x();
111 cpsScalar[varID++] = cpI.y();
112 cpsScalar[varID++] = cpI.z();
118 (constraint_().controlPointsToDesignVariables(cpsScalar));
128 scalarField cpsScalar(3*volBSplinesBase_.getTotalControlPointsNumber());
133 const label nCPs(boxI.getControlPoints().size());
134 for (label cpI = 0; cpI < nCPs; ++cpI)
136 cpsScalar[varID++] = controlPoints[cpI].x();
137 cpsScalar[varID++] = controlPoints[cpI].y();
138 cpsScalar[varID++] = controlPoints[cpI].z();
144 (constraint_().controlPointsToDesignVariables(cpsScalar));
155 readBounds(lowerBounds_,
"lower", -1);
156 readBounds(upperBounds_,
"upper", 1);
159 constraint_().computeBounds(lowerBounds_, upperBounds_);
166 const word& boundsName,
171 if (dict_.found(boundsName +
"CPBounds"))
175 vector CPBounds(dict_.get<
vector>(boundsName +
"CPBounds"));
180 const label nCPs(boxI.getControlPoints().size());
181 for (label iCP = 0; iCP < nCPs; ++iCP)
183 bounds()[varID++] = CPBounds.x();
184 bounds()[varID++] = CPBounds.y();
185 bounds()[varID++] = CPBounds.z();
198 <<
"Reading " << boundsName <<
"Bounds from dict " <<
endl;
200 (
new scalarField(boundsName +
"Bounds", *
this, getVars().size()));
203 else if (nonOverlappingCPs_)
206 <<
"Setting " << boundsName <<
"Bounds from nonOverlappingCPs" 209 const PtrList<NURBS3DVolume>& boxes = volBSplinesBase_.boxesRef();
211 for (
const NURBS3DVolume& boxI : boxes)
214 const Vector<label> nCPsDir = boxI.nCPsPerDirection();
216 for (label idir = 0; idir < 3; ++idir)
219 (
max(cps.component(idir)) -
min(cps.component(idir)))
220 /scalar(nCPsDir[idir] - 1);
222 const label nCPs(boxI.getControlPoints().size());
223 for (label iCP = 0; iCP < nCPs; ++iCP)
226 bounds()[varID++] =
cp.x() +
sign*0.5*dists.x();
227 bounds()[varID++] =
cp.y() +
sign*0.5*dists.y();
228 bounds()[varID++] =
cp.z() +
sign*0.5*dists.z();
243 const PtrList<NURBS3DVolume>& boxes = volBSplinesBase_.boxesRef();
245 for (
const NURBS3DVolume& boxI : boxes)
249 word(
"optimisation")/word(
"controlPoints")/boxI.name()
250 +
name + mesh_.time().timeName() +
".csv" 253 file<<
"\"Points : 0\", \"Points : 1\", \"Points : 2\"," 254 <<
"\"i\", \"j\", \"k\""<<
endl;
257 const label nCPsU = boxI.basisU().nCPs();
258 const label nCPsV = boxI.basisV().nCPs();
261 const label
k = cpI/label(nCPsU*nCPsV);
262 const label j = (cpI -
k*nCPsU*nCPsV)/nCPsU;
263 const label i = (cpI -
k*nCPsU*nCPsV - j*nCPsU);
265 file<< bounds[3*cpI + passed] <<
", " 266 << bounds[3*cpI + 1 + passed] <<
", " 267 << bounds[3*cpI + 2 + passed] <<
", " 272 passed += 3*cps.size();
283 displacementMethod& dm = displMethodPtr_.ref();
285 if (isA<displacementMethodvolumetricBSplinesMotionSolver>(dm))
288 displMethodPtr_->setControlField(cpMovement);
293 tmp<vectorField> tnewPoints =
294 volBSplinesBase_.computeBoundaryDisplacement
297 parametertisedPatches_.toc()
306 mesh_.time().timeName(),
316 for (
const label pI : parametertisedPatches_)
320 dx.primitiveFieldRef(),
321 vectorField(newPoints, mesh_.boundaryMesh()[pI].meshPoints())
326 displMethodPtr_->setMotionField(dx);
333 Foam::volumetricBSplinesDesignVariables::volumetricBSplinesDesignVariables
336 const dictionary&
dict 344 "volumetricBSplinesDesignVariables",
348 IOobject::READ_IF_PRESENT,
353 volBSplinesBase_(const_cast<volBSplinesBase&>(volBSplinesBase::
New(
mesh))),
354 nonOverlappingCPs_(dict_.getOrDefault<bool>(
"nonOverlappingCPs", false)),
355 updateBounds_(dict_.getOrDefault<bool>(
"updateBounds", true)),
356 constraint_(morphingBoxConstraint::
New(
mesh,
dict, *this))
389 volBSplinesBase_.getTotalControlPointsNumber(),
402 cpMovement[iCP].x() = correctionCPs[3*iCP];
403 cpMovement[iCP].y() = correctionCPs[3*iCP + 1];
404 cpMovement[iCP].z() = correctionCPs[3*iCP + 2];
406 volBSplinesBase_.boundControlPointMovement(cpMovement);
415 tmp<vectorField> tcpMovement = controlPointMovement(
correction);
419 setDisplacement(cpMovement);
433 designVariablesToControlPoints();
440 return constraint_().computeEta(
correction, maxInitChange_());
457 constraint_().postProcessSens
467 constraint_().updateBounds(lowerBounds_, upperBounds_);
476 lowerBounds_().writeEntry(
"lowerBounds",
os);
477 writeBounds(lowerBounds_(),
"lowerBounds");
481 upperBounds_().writeEntry(
"upperBounds",
os);
482 writeBounds(upperBounds_(),
"upperBounds");
484 return constraint_().writeData(
os);
493 const displacementMethod& dm = displMethodPtr_();
494 if (isA<displacementMethodvolumetricBSplinesMotionSolver>(dm))
496 Vector<label> decomposed = volBSplinesBase_.decomposeDV(varID);
497 const label boxI = decomposed.x();
498 const label cpILocal = decomposed.y();
499 const label dir = decomposed.z();
514 Vector<label> decomposed = volBSplinesBase_.decomposeDV(varID);
515 const label boxI = decomposed.x();
516 const label cpILocal = decomposed.y();
517 const label dir = decomposed.z();
520 (volBSplinesBase_.boxRef(boxI).patchDxDbFace(patchI, cpILocal));
531 Vector<label> decomposed = volBSplinesBase_.decomposeDV(varID);
532 const label boxI = decomposed.
x();
533 const label cpILocal = decomposed.
y();
534 const label dir = decomposed.
z();
538 volBSplinesBase_.boxRef(boxI).
539 dndbBasedSensitivities(patchI, cpILocal,
false)
551 Vector<label> decomposed = volBSplinesBase_.decomposeDV(varID);
552 const label boxI = decomposed.
x();
553 const label cpILocal = decomposed.
y();
554 const label dir = decomposed.
z();
558 volBSplinesBase_.boxRef(boxI).dndbBasedSensitivities(patchI, cpILocal)
569 Vector<label> decomposed = volBSplinesBase_.decomposeDV(varID);
570 const label boxI = decomposed.
x();
571 const label cpILocal = decomposed.
y();
572 const label dir = decomposed.
z();
586 mesh_.time().timeName(),
dimensionedScalar sign(const dimensionedScalar &ds)
void readBounds(autoPtr< scalar > lowerBoundPtr=nullptr, autoPtr< scalar > upperBoundPtr=nullptr)
Read bounds for design variables, if present.
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix by subtracting the matrix multiplied by the current fi...
autoPtr< morphingBoxConstraint > constraint_
Constraint imposed on the movement of the design variables.
virtual bool globalSum() const
Whether to use global sum when computing matrix-vector products in update methods.
virtual void evolveNumber()
For design variables with a dynamic character (i.e. changing.
A list of keyword definitions, which are a keyword followed by a number of values (eg...
Abstract base class for adjoint-based sensitivities.
static const pointMesh & New(const polyMesh &mesh, Args &&... args)
Get existing or create MeshObject registered with typeName.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
virtual void resetDesignVariables()
Reset to starting point of line search.
virtual scalar computeEta(scalarField &correction)
Compute eta if not set in the first step.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy the source to the destination (recursively if necessary).
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
volBSplinesBase & volBSplinesBase_
Reference to underlaying volumetric B-Splines morpher.
void designVariablesToControlPoints()
Set control points based on current design variables values.
virtual tmp< scalarField > assembleSensitivities(adjointSensitivity &adjointSens)
Add part of sensitivity derivatives related to geometry variations.
GeometricField< vector, pointPatchField, pointMesh > pointVectorField
void writeEntry(const word &keyword, Ostream &os) const
Write the field as a dictionary entry.
const Cmpt & y() const noexcept
Access to the vector y component.
label k
Boltzmann constant.
Ignore writing from objectRegistry::writeObject()
const dimensionSet dimless
Dimensionless.
const adjointSolver & getAdjointSolver() const
Const access to adjoint solver.
NURBS3DVolume morpher. Includes support functions for gradient computations Base class providing supp...
Macros for easy insertion into run-time selection tables.
GeometricField< tensor, pointPatchField, pointMesh > pointTensorField
virtual tmp< scalarField > assembleSensitivities(adjointSensitivity &adjointSens)
Assemble the sensitivity derivatives, by also applying possible constraints.
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find an entry (const access) with the given keyword.
#define forAll(list, i)
Loop across all elements in list.
void reset(T *p=nullptr) noexcept
Delete managed object and set to new given pointer.
virtual void update(scalarField &correction)
Update design variables based on a given correction.
void writeBounds(const scalarField &bounds, const word &name) const
Write current bounds to file.
friend Ostream & operator(Ostream &, const Field< scalar > &)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
virtual tmp< volVectorField > dCdb(const label varID) const
Get dCdb for given design variable.
A class for handling words, derived from Foam::string.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
void controlPointsToDesignVariables()
Set the design variables based on the current control points.
static tmp< T > New(Args &&... args)
Construct tmp with forwarding arguments.
void readBounds(autoPtr< scalar > lowerBoundPtr=nullptr, autoPtr< scalar > upperBoundPtr=nullptr)
Read bounds for design variables, if present.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define DebugInfo
Report an information message using Foam::Info.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
label size() const noexcept
The number of elements in the container.
const Cmpt & x() const noexcept
Access to the vector x component.
void setActiveDesignVariables()
Set IDs of active design variables.
OBJstream os(runTime.globalPath()/outputName)
defineTypeNameAndDebug(combustionModel, 0)
virtual tmp< vectorField > dxdbFace(const label patchI, const label varID) const
Get dxdb for given design variable and patch.
tmp< vectorField > controlPointMovement(const scalarField &correction)
Transform correction of design variables to control points movement.
void unzipCol(const FieldField< Field, SymmTensor< Cmpt >> &input, const direction idx, FieldField< Field, Vector< Cmpt >> &result)
Extract a symmTensor field field column (x,y,z) == (0,1,2)
virtual tmp< vectorField > dSdb(const label patchI, const label varID) const
Get dSdb for given design variable and patch.
label getTotalControlPointsNumber() const
Get cumulative number of control points from all boxes.
const Cmpt & z() const noexcept
Access to the vector z component.
void setDisplacement(const vectorField &cpMovement)
Set the field driving the displacement method.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers...
const word & solverName() const
Return the solver name.
static bool master(const label communicator=worldComm)
True if process corresponds to the master rank in the communicator.
virtual label sensSize() const
Size of the active control points.
Internal & ref(const bool updateAccessTime=true)
Same as internalFieldRef()
components
Component labeling enumeration.
Field< vector > vectorField
Specialisation of Field<T> for vector.
virtual void resetDesignVariables()
Reset to starting point of line search.
void operator+=(const UList< scalar > &)
tmp< pointTensorField > getDxDb(const label cpI)
Get dxCartesiandb for a certain control point.
virtual tmp< vectorField > dxdbVol(const label varID) const
Get dxdb for all mesh points.
A class for managing temporary objects.
virtual bool writeData(Ostream &os) const
Write fields to support continuation.
virtual const labelList & activeSensitivities() const
Components of the active control points.
Defines the attributes of an object for which implicit objectRegistry management is supported...
const Boundary & boundaryField() const noexcept
Return const-reference to the boundary field.
Do not request registration (bool: false)
virtual tmp< vectorField > dndb(const label patchI, const label varID) const
Get dndb for given design variable and patch.
addToRunTimeSelectionTable(functionObject, pointHistory, dictionary)
static constexpr const zero Zero
Global zero (0)