Go to the documentation of this file.
39 namespace incompressible
48 sensitivitySurfacePoints,
104 sensitivityPatchIDs_,
120 distanceSensPtr.
reset(createZeroBoundaryPtr<vector>(
mesh_));
123 for (
const label patchI : sensitivityPatchIDs_)
125 distanceSensPtr()[patchI] = sens[patchI];
133 meshMovementSensPtr.
reset(createZeroBoundaryPtr<vector>(
mesh_));
136 for (
const label patchI : sensitivityPatchIDs_)
138 meshMovementSensPtr()[patchI] = sens[patchI];
144 for (
const label patchI : sensitivityPatchIDs_)
177 for (
const label patchI : sensitivityPatchIDs_)
183 vectorField& pointPatchSens = wallPointSensVecPtr_()[patchI];
203 const label patchStartIndex =
patch.start();
212 const labelList& pointFaces = patchPointFaces[ppI];
215 label localFaceIndex = pointFaces[pfI];
216 label globalFaceIndex = patchStartIndex + localFaceIndex;
217 const face& faceI = faces[globalFaceIndex];
224 if (faceI[facePointI] == meshPoints[ppI])
234 const tensor& deltaCf = deltaNormals[0];
235 pointPatchSens[ppI] += facePatchSens[localFaceIndex] & deltaCf;
243 const tensor& deltaSf = deltaNormals[1];
244 pointPatchSens[ppI] +=
245 dSfdbMultPatch[localFaceIndex] & deltaSf;
248 const tensor& deltaNf = deltaNormals[2];
249 pointPatchSens[ppI] +=
250 dnfdbMultPatch[localFaceIndex] & deltaNf;
264 for (
const label patchI : sensitivityPatchIDs_)
280 const labelList& pointFaces = patchPointFaces[ppI];
283 const label localFaceIndex = pointFaces[pfI];
286 pointNormals[meshPoints[ppI]] += nf[localFaceIndex];
287 pointMagSf[meshPoints[ppI]] += magSf[localFaceIndex];
332 sensitivitySurfacePoints::sensitivitySurfacePoints
350 includeSurfaceArea_(
false),
351 includePressureTerm_(
false),
352 includeGradStressTerm_(
false),
353 includeTransposeStresses_(
false),
354 useSnGradInTranposeStresses_(
false),
355 includeDivTerm_(
false),
356 includeDistance_(
false),
357 includeMeshMovement_(
false),
358 includeObjective_(
false),
359 eikonalSolver_(
nullptr),
360 meshMovementSolver_(
nullptr),
361 wallFaceSens_(createZeroBoundaryPtr<vector>(mesh_)),
362 dSfdbMult_(createZeroBoundaryPtr<vector>(mesh_)),
363 dnfdbMult_(createZeroBoundaryPtr<vector>(mesh_))
369 wallPointSensVecPtr_.reset(createZeroBoundaryPointFieldPtr<vector>(mesh_));
370 wallPointSensNormalPtr_.reset
372 createZeroBoundaryPointFieldPtr<scalar>(mesh_)
374 wallPointSensNormalVecPtr_.reset
376 createZeroBoundaryPointFieldPtr<vector>(mesh_)
380 label nTotalPoints(0);
381 for (
const label patchI : sensitivityPatchIDs_)
383 label
nPoints = mesh_.boundaryMesh()[patchI].nPoints();
427 <<
" Calculating auxilary quantities " <<
endl;
439 if (isA<wallFvPatch>(
patch))
444 nf*
U.boundaryField()[patchI].snGrad();
453 createZeroFieldPtr<vector>(
mesh_,
"stressX", stress.dimensions())
457 createZeroFieldPtr<vector>(
mesh_,
"stressY", stress.dimensions())
461 createZeroFieldPtr<vector>(
mesh_,
"stressZ", stress.dimensions())
464 stressXPtr().replace(0, stress.
component(0));
465 stressXPtr().replace(1, stress.
component(1));
466 stressXPtr().replace(2, stress.
component(2));
468 stressYPtr().replace(0, stress.
component(3));
469 stressYPtr().replace(1, stress.
component(4));
470 stressYPtr().replace(2, stress.
component(5));
472 stressZPtr().replace(0, stress.
component(6));
473 stressZPtr().replace(1, stress.
component(7));
474 stressZPtr().replace(2, stress.
component(8));
494 adjointTurbulence->wallShapeSensitivities();
500 <<
" Calculating adjoint sensitivity. " <<
endl;
504 for (
const label patchI : sensitivityPatchIDs_)
522 &
U.boundaryField()[patchI].snGrad()
524 * nuEff.boundaryField()[patchI]
537 gradStressTerm = (-((Uab & nf)*gradp.
boundaryField()[patchI]));
552 (gradUa.boundaryField()[patchI] & nf)
555 nuEff.boundaryField()[patchI]
556 *(gradUaNf &
U.boundaryField()[patchI].snGrad())
563 scalar(1./3.)*nuEff.boundaryField()[patchI]
566 &
U.boundaryField()[patchI].snGrad()
578 &
U.boundaryField()[patchI].snGrad()
589 const scalar wei = functions[funcI].weight();
592 wei*functions[funcI].dxdbDirectMultiplier(patchI);
596 wei*dt*functions[funcI].dSdbMultiplier(patchI);
598 wei*dt*functions[funcI].dndbMultiplier(patchI);
610 + adjointTMsensitivities[patchI]
638 for (
const label patchI : sensitivityPatchIDs_)
643 const label globaPointI = meshPoints[ppI];
644 pointSensGlobal[globaPointI] +=
645 wallPointSensVecPtr_()[patchI][ppI];
659 for (
const label patchI : sensitivityPatchIDs_)
663 wallPointSensVecPtr_()[patchI].map(pointSensGlobal, meshPoints);
668 for (
const label patchI : sensitivityPatchIDs_)
679 vectorField patchPointNormals(pointNormals, meshPoints);
680 patchPointNormals /=
mag(patchPointNormals) + VSMALL;
683 wallPointSensVecPtr_()[patchI] /=
686 wallPointSensNormalPtr_()[patchI] =
687 wallPointSensVecPtr_()[patchI]
689 wallPointSensNormalVecPtr_()[patchI] =
690 wallPointSensNormalPtr_()[patchI]
706 3*wallPointSensNormalVecPtr_()[patchI].size()
709 forAll(wallPointSensNormalVecPtr_()[patchI], ptI)
711 patchScalarSens[3*ptI] =
712 wallPointSensNormalVecPtr_()[patchI][ptI].x();
713 patchScalarSens[3*ptI + 1] =
714 wallPointSensNormalVecPtr_()[patchI][ptI].y();
715 patchScalarSens[3*ptI + 2] =
716 wallPointSensNormalVecPtr_()[patchI][ptI].z();
721 forAll(procPatchSens, procI)
723 const scalarField& procSens = procPatchSens[procI];
728 nPassedDVs += procSens.size();
autoPtr< boundaryVectorField > wallFaceSens_
The face-based part of the sensitivities.
const dictionary & dict() const
Return the construction dictionary.
virtual const pointField & points() const
Return raw points.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
class for managing incompressible objective functions.
bool includeTransposeStresses_
Include the transpose part of the adjoint stresses.
tmp< GeometricField< cmptType, PatchField, GeoMesh > > component(const direction) const
Return a component of the field.
const word & solverName() const
Return solver name.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Differentiation of the mesh data structure.
A class for handling words, derived from Foam::string.
autoPtr< boundaryVectorField > dnfdbMult_
defineTypeNameAndDebug(adjointEikonalSolver, 0)
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
A class for managing temporary objects.
virtual void accumulateIntegrand(const scalar dt)
Accumulate sensitivity integrands.
static constexpr const zero Zero
Global zero (0)
void constructGlobalPointNormalsAndAreas(vectorField &pointNormals, scalarField &pointMagSf)
Construct globally correct point normals and point areas.
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
bool read(const char *buf, int32_t &val)
Same as readInt32.
bool includePressureTerm_
Include the adjoint pressure term in sens computation.
autoPtr< adjointMeshMovementSolver > meshMovementSolver_
bool useSnGradInTranposeStresses_
Use snGrad in the transpose part of the adjoint stresses.
autoPtr< boundaryVectorField > dSfdbMult_
Multipliers of d(Sf)/db and d(nf)/db.
void read()
Read controls and update solver pointers if necessary.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const volScalarField & p() const
Return const reference to pressure.
Class including all adjoint fields for incompressible flows.
label nPoints() const noexcept
Number of mesh points.
#define forAll(list, i)
Loop across all elements in list.
const autoPtr< incompressible::RASModelVariables > & RASModelVariables() const
Return const reference to the turbulence model variables.
bool includeDivTerm_
Include the term from the deviatoric part of the stresses.
Solver of the adjoint to the eikonal PDE.
addToRunTimeSelectionTable(adjointSensitivity, sensitivityBezier, dictionary)
void write()
Write sensitivity fields.
Abstract base class for adjoint-based sensitivities in incompressible flows.
virtual bool readDict(const dictionary &dict)
Read dictionary if changed.
autoPtr< adjointEikonalSolver > eikonalSolver_
A patch is a list of labels that address the faces in the global face list.
void setSize(const label n)
Alias for resize()
virtual bool readDict(const dictionary &dict)
Read dict if changed.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const volScalarField & pa() const
Return const reference to pressure.
const volVectorField & U() const
Return const reference to velocity.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
virtual void clearSensitivities()
Zero sensitivity fields and their constituents.
bool includeMeshMovement_
Include mesh movement variation in sens computation.
const volVectorField & Ua() const
Return const reference to velocity.
bool includeSurfaceArea_
Include surface area in sens computation.
bool includeObjective_
Include terms directly emerging from the objective function.
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.
vectorField makeFaceCentresAndAreas_d(const pointField &p, const pointField &p_d)
incompressibleAdjointVars & adjointVars_
Mesh data needed to do the Finite Volume discretisation.
incompressibleVars & primalVars_
virtual void write(const word &baseName=word::null)
Write sensitivity fields.
void reset(autoPtr< T > &&other) noexcept
Delete managed object and set to new given pointer.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
bool includeDistance_
Include distance variation in sens computation.
virtual const faceList & faces() const
Return raw faces.
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
bool includeGradStressTerm_
Include the term containing the grad of the stress at the boundary.
#define DebugInfo
Report an information message using Foam::Info.
static int myProcNo(const label communicator=worldComm)
Number of this process (starting from masterNo() = 0)
const std::string patch
OpenFOAM patch number as a std::string.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
const autoPtr< incompressibleAdjoint::adjointRASModel > & adjointTurbulence() const
Return const reference to the adjointRASModel.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Solver of the adjoint to the Laplace grid displacement equation.
void clearSensitivities()
Zero sensitivity fields and their constituents.
pointField points(const UList< point > &pts) const
Return the points corresponding to this face.
static const word null
An empty word.
void setSuffixName()
Set suffix name for sensitivity fields.
static const Vector< scalar > zero
A face is a list of labels corresponding to mesh vertices.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
void setSuffix(const word &suffix)
Set suffix.
void finaliseFaceMultiplier()
virtual void clearSensitivities()
Zero sensitivity fields and their constituents.
void finalisePointSensitivities()
static label nProcs(const label communicator=worldComm)
Number of processes in parallel run, and 1 for serial run.
virtual void assembleSensitivities()
Assemble sensitivities.
Base class for solution control classes.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
virtual void write(const word &baseName=word::null)
Write sensitivity fields.
PtrList< objective > & getObjectiveFunctions()
Return reference to objective functions.
objectiveManager & objectiveManager_