Go to the documentation of this file.
39 namespace incompressible
47 adjointSensitivity, sensitivityBezierFI, dictionary
79 for (
const label patchI : sensitivityPatchIDs_)
89 for (label iter = 0; iter < meshMovementIters_; iter++)
91 Info<<
"Mesh Movement Propagation(direction, CP), ("
92 << idir <<
", " << iCP <<
"), Iteration : "<< iter <<
endl;
100 scalar residual =
mag(mEqn.
solve().initialResidual());
104 mesh_.time().printExecutionTime(
Info);
107 if (residual < meshMovementResidualLimit_)
109 Info<<
"\n***Reached dxdb convergence limit, iteration " << iter
121 sensitivityBezierFI::sensitivityBezierFI
142 flowSens_(3*Bezier_.nBezier(),
Zero),
143 dSdbSens_(3*Bezier_.nBezier(),
Zero),
144 dndbSens_(3*Bezier_.nBezier(),
Zero),
145 dxdbDirectSens_(3*Bezier_.nBezier(),
Zero),
146 dVdbSens_(3*Bezier_.nBezier(),
Zero),
147 distanceSens_(3*Bezier_.nBezier(),
Zero),
148 optionsSens_(3*Bezier_.nBezier(),
Zero),
149 bcSens_(3*Bezier_.nBezier(),
Zero),
151 derivativesFolder_(
"optimisation"/
type() +
"Derivatives"),
153 meshMovementIters_(-1),
154 meshMovementResidualLimit_(1.
e-7),
169 mkDir(derivativesFolder_);
185 distanceSensPtr.
reset
187 createZeroFieldPtr<tensor>
198 const label nDVs = 3*nBezier;
199 for (label iDV = 0; iDV < nDVs; iDV++)
201 label iCP = iDV%nBezier;
202 label idir = iDV/nBezier;
228 for (
const label patchI : sensitivityPatchIDs_)
277 distanceSensPtr().primitiveField()
320 Info<<
"Writing control point sensitivities to file" <<
endl;
331 <<
setw(widthDV) <<
"#dv" <<
" "
332 <<
setw(width) <<
"total" <<
" "
333 <<
setw(width) <<
"flow" <<
" "
334 <<
setw(width) <<
"dSdb" <<
" "
335 <<
setw(width) <<
"dndb" <<
" "
336 <<
setw(width) <<
"dxdbDirect" <<
" "
337 <<
setw(width) <<
"dVdb" <<
" "
338 <<
setw(width) <<
"distance" <<
" "
339 <<
setw(width) <<
"options" <<
" "
344 label lastActive(-1);
346 for (label iDV = 0; iDV < nDVs; iDV++)
348 const label iCP(iDV%nBezier);
349 const label idir(iDV/nBezier);
350 if (!confineMovement[idir][iCP])
352 if (iDV!=lastActive + 1)
358 <<
setw(widthDV) << iDV <<
" "
vectorField optionsDxDbMult_
dx/db multiplier coming from fvOptions
virtual void clearSensitivities()
Zero sensitivity fields and their constituents.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
void reset(T *p=nullptr) noexcept
Delete managed object and set to new given pointer.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
class for managing incompressible objective functions.
const word & solverName() const
Return solver name.
scalarField optionsSens_
Term depending on fvOptions.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
autoPtr< boundaryVectorField > bcDxDbMult_
const boolList & confineXmovement() const
Confine x movement.
A class for handling words, derived from Foam::string.
virtual void clearSensitivities()
Zero sensitivity fields and their constituents.
defineTypeNameAndDebug(adjointEikonalSolver, 0)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
scalarField dVdbSens_
Term depending on delta(V)/delta b.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
bool read(const char *buf, int32_t &val)
Same as readInt32.
volTensorField gradDxDbMult_
grad(dx/db) multiplier
static word timeName(const scalar t, const int precision=precision_)
scalarField dxdbDirectSens_
void read()
Read options and update solver pointers if necessary.
scalarField distanceSens_
Term depending on distance differentiation.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
scalarField divDxDbMult_
div(dx/db) multiplier
Ostream & endl(Ostream &os)
Add newline and flush stream.
Type gSum(const FieldField< Field, Type > &f)
const boolList & confineZmovement() const
Confine z movement.
Dimension set for the base types.
virtual void assembleSensitivities()
Assemble sensitivities.
Class including all adjoint fields for incompressible flows.
autoPtr< adjointEikonalSolver > eikonalSolver_
Adjoint eikonal equation solver.
scalarField dSdbSens_
Term depending on delta(n dS)/delta b.
label nBezier() const
Number of Bezier control points.
tmp< volVectorField > solveMeshMovementEqn(const label iCP, const label idir)
autoPtr< boundaryVectorField > dSfdbMult_
Fields related to direct sensitivities.
addToRunTimeSelectionTable(adjointSensitivity, sensitivityBezier, dictionary)
scalarField flowSens_
Flow related term.
messageStream Info
Information stream (uses stdout - output is on the master only)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
word name(const complex &c)
Return string representation of complex.
scalarField dndbSens_
Term depending on delta(n)/delta b.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
Istream and Ostream manipulators taking arguments.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
tmp< tensorField > dxdbFace(const label patchI, const label cpI, bool useChainRule=true) const
dxdb tensor for a Bezier parameterized patch
static tmp< volVectorField > autoCreateMeshMovementField(const fvMesh &mesh, const word &name, const dimensionSet &dims)
Auto create variable for mesh movement.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
autoPtr< boundaryVectorField > dxdbDirectMult_
Macros for easy insertion into run-time selection tables.
incompressibleAdjointVars & adjointVars_
Mesh data needed to do the Finite Volume discretisation.
Omanip< int > setw(const int i)
dictionary subOrEmptyDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX, const bool mandatory=false) const
scalar meshMovementResidualLimit_
Output to file stream, using an OSstream.
static bool master(const label communicator=0)
Am I the master process.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Base class for Field Integral-based sensitivity derivatives.
static unsigned int defaultPrecision()
Return the default precision.
bool includeDistance_
Include distance variation in sens computation.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
scalarField bcSens_
Term depending on the differenation of boundary conditions.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
const dimensionedScalar e
Elementary charge.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
SolverPerformance< Type > solve(const dictionary &)
Solve returning the solution statistics.
const Time & time() const
Return the top-level database.
const boolList & confineYmovement() const
Confine y movement.
const boolListList & confineMovement() const
Info about confining movement in all directions.
virtual void write(const word &baseName=word::null)
Write sensitivities to file.
autoPtr< boundaryVectorField > dnfdbMult_
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
fv::IOoptionListAdjoint fvOptionsAdjoint(mesh)
bool mkDir(const fileName &pathName, mode_t mode=0777)
Make a directory and return an error if it could not be created.
Type gMax(const FieldField< Field, Type > &f)
Base class for solution control classes.
fileName derivativesFolder_
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
tmp< tensorField > dndbBasedSensitivities(const label patchI, const label cpI, bool returnDimensionedNormalSens=true) const