Go to the documentation of this file.
40 optMeshMovementVolumetricBSplinesExternalMotionSolver,
46 optMeshMovementVolumetricBSplinesExternalMotionSolver,
60 const label nCPs(volBSplinesBase_.getTotalControlPointsNumber());
61 dx_.primitiveFieldRef() = vector::zero;
62 cpMovement_ = vector::zero;
64 for (label iCP = 0; iCP < nCPs; iCP++)
72 volBSplinesBase_.boundControlPointMovement(cpMovement_);
79 const label nb = boxes[iNURB].getControlPoints().size();
80 for (label cpI = 0; cpI < nb; ++cpI)
82 const label globalCP = passedCPs + cpI;
85 const label patchI = patchIDs_[pI];
88 boxes[iNURB].patchDxDb(patchI, cpI)
89 & cpMovement_[globalCP]
91 dx_.boundaryField()[patchI].addToInternalField
93 dx_.primitiveFieldRef(),
107 Foam::optMeshMovementVolumetricBSplinesExternalMotionSolver::
108 optMeshMovementVolumetricBSplinesExternalMotionSolver
125 mesh.time().timeName(),
134 cpMovement_(volBSplinesBase_.getTotalControlPointsNumber(),
Zero)
170 scalar maxDisplacement =
gMax(
mag(dx_.primitiveField()));
173 Info<<
"maxAllowedDisplacement/maxDisplacement \t"
174 << getMaxAllowedDisplacement() <<
"/" << maxDisplacement <<
endl;
175 scalar eta = getMaxAllowedDisplacement() / maxDisplacement;
176 Info<<
"Setting eta value to " << eta <<
endl;
186 return volBSplinesBase_.getActiveDesignVariables();
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
Defines the attributes of an object for which implicit objectRegistry management is supported,...
autoPtr< displacementMethod > displMethodPtr_
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
Class constructing a number of volumetric B-Splines boxes, read from dynamicMeshDict....
static constexpr const zero Zero
Global zero (0)
Ostream & endl(Ostream &os)
Add newline and flush stream.
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
vectorField cpMovement_
Movement of control points.
void moveControlPoints(const vectorField &controlPointsMovement)
Move control points. No effect on mesh.
#define forAll(list, i)
Loop across all elements in list.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
virtual scalar computeEta(const scalarField &correction)
Compute eta value based on max displacement.
virtual labelList getActiveDesignVariables() const
Return active design variables.
messageStream Info
Information stream (uses stdout - output is on the master only)
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
volBSplinesBase & volBSplinesBase_
Reference to underlaying volumetric B-Splines morpher.
Mesh data needed to do the Finite Volume discretisation.
void moveMesh()
Calculates surface mesh movement.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
pointVectorField dx_
Boundary movement due to change of NURBS control points.
void computeBoundaryMovement(const scalarField &correction)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
scalarField correction_
Correction of design variables.
Abstract base class for translating an update of the design variables into mesh movement.
defineTypeNameAndDebug(combustionModel, 0)
Type gMax(const FieldField< Field, Type > &f)
void writeControlPoints() const
Write control points to constant and optimisation folders.