40namespace incompressible
105 sensitivityPatchIDs_,
121 distanceSensPtr.
reset(createZeroBoundaryPtr<vector>(
mesh_));
124 for (
const label patchI : sensitivityPatchIDs_)
126 distanceSensPtr()[patchI] = sens[patchI];
134 meshMovementSensPtr.
reset(createZeroBoundaryPtr<vector>(
mesh_));
137 for (
const label patchI : sensitivityPatchIDs_)
139 meshMovementSensPtr()[patchI] = sens[patchI];
145 for (
const label patchI : sensitivityPatchIDs_)
178 for (
const label patchI : sensitivityPatchIDs_)
184 vectorField& pointPatchSens = wallPointSensVecPtr_()[patchI];
194 const labelList& meshPoints = patch.patch().meshPoints();
201 const labelListList& patchPointFaces = patch.patch().pointFaces();
204 const label patchStartIndex = patch.start();
213 const labelList& pointFaces = patchPointFaces[ppI];
216 label localFaceIndex = pointFaces[pfI];
217 label globalFaceIndex = patchStartIndex + localFaceIndex;
218 const face& faceI = faces[globalFaceIndex];
225 if (faceI[facePointI] == meshPoints[ppI])
235 const tensor& deltaCf = deltaNormals[0];
236 pointPatchSens[ppI] += facePatchSens[localFaceIndex] & deltaCf;
244 const tensor& deltaSf = deltaNormals[1];
245 pointPatchSens[ppI] +=
246 dSfdbMultPatch[localFaceIndex] & deltaSf;
249 const tensor& deltaNf = deltaNormals[2];
250 pointPatchSens[ppI] +=
251 dnfdbMultPatch[localFaceIndex] & deltaNf;
265 for (
const label patchI : sensitivityPatchIDs_)
272 const labelList& meshPoints = patch.patch().meshPoints();
276 const labelListList& patchPointFaces = patch.patch().pointFaces();
281 const labelList& pointFaces = patchPointFaces[ppI];
284 const label localFaceIndex = pointFaces[pfI];
287 pointNormals[meshPoints[ppI]] += nf[localFaceIndex];
288 pointMagSf[meshPoints[ppI]] += magSf[localFaceIndex];
342 includeSurfaceArea_(false),
343 includePressureTerm_(false),
344 includeGradStressTerm_(false),
345 includeTransposeStresses_(false),
346 useSnGradInTranposeStresses_(false),
347 includeDivTerm_(false),
348 includeDistance_(false),
349 includeMeshMovement_(false),
350 includeObjective_(false),
351 eikonalSolver_(nullptr),
352 meshMovementSolver_(nullptr),
361 wallPointSensVecPtr_.reset(createZeroBoundaryPointFieldPtr<vector>(
mesh_));
362 wallPointSensNormalPtr_.reset
364 createZeroBoundaryPointFieldPtr<scalar>(
mesh_)
366 wallPointSensNormalVecPtr_.reset
368 createZeroBoundaryPointFieldPtr<vector>(
mesh_)
372 label nTotalPoints(0);
373 for (
const label patchI : sensitivityPatchIDs_)
388 if (sensitivity::readDict(
dict))
431 adjointTurbulence->wallShapeSensitivities();
437 <<
" Calculating adjoint sensitivity. " <<
endl;
454 for (
const label patchI : sensitivityPatchIDs_)
480 if (isA<wallFvPatch>(patch))
483 gradUbf[patchI] = tnf*
U.boundaryField()[patchI].snGrad();
492 for (label idir = 0; idir < pTraits<vector>::nComponents; ++idir)
496 for (
const label patchI : sensitivityPatchIDs_)
516 tgradUa.
cref().boundaryField();
517 for (
const label patchI : sensitivityPatchIDs_)
526 : (gradUabf[patchI] & nf)
530 *(gradUaNf &
U.boundaryField()[patchI].snGrad())*tnf;
536 for (
const label patchI : sensitivityPatchIDs_)
554 &
U.boundaryField()[patchI].snGrad()
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);
609 + adjointTMsensitivities[patchI]
641 for (
const label patchI : sensitivityPatchIDs_)
646 const label globaPointI = meshPoints[ppI];
647 pointSensGlobal[globaPointI] +=
648 wallPointSensVecPtr_()[patchI][ppI];
662 for (
const label patchI : sensitivityPatchIDs_)
666 wallPointSensVecPtr_()[patchI].map(pointSensGlobal, meshPoints);
671 for (
const label patchI : sensitivityPatchIDs_)
677 const labelList& meshPoints = patch.meshPoints();
682 vectorField patchPointNormals(pointNormals, meshPoints);
683 patchPointNormals /=
mag(patchPointNormals) + VSMALL;
686 wallPointSensVecPtr_()[patchI] /=
689 wallPointSensNormalPtr_()[patchI] =
690 wallPointSensVecPtr_()[patchI]
692 wallPointSensNormalVecPtr_()[patchI] =
693 wallPointSensNormalPtr_()[patchI]
709 3*wallPointSensNormalVecPtr_()[patchI].size()
712 forAll(wallPointSensNormalVecPtr_()[patchI], ptI)
714 patchScalarSens[3*ptI] =
715 wallPointSensNormalVecPtr_()[patchI][ptI].x();
716 patchScalarSens[3*ptI + 1] =
717 wallPointSensNormalVecPtr_()[patchI][ptI].y();
718 patchScalarSens[3*ptI + 2] =
719 wallPointSensNormalVecPtr_()[patchI][ptI].z();
723 forAll(procPatchSens, procI)
725 const scalarField& procSens = procPatchSens[procI];
730 nPassedDVs += procSens.
size();
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
Generic GeometricBoundaryField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
void setSize(const label n)
Alias for resize()
label nProcs() const noexcept
Number of ranks associated with PstreamBuffers.
static void allGatherList(const List< commsStruct > &comms, List< T > &values, const int tag, const label comm)
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
void size(const label n)
Older name for setAddressableSize.
Base class for adjoint solvers.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
T & ref()
Return reference to the managed object without nullptr checking.
void reset(autoPtr< T > &&other) noexcept
Delete managed object and set to new given pointer.
Differentiation of the mesh data structure.
vectorField makeFaceCentresAndAreas_d(const pointField &p, const pointField &p_d)
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
A face is a list of labels corresponding to mesh vertices.
pointField points(const UList< point > &pts) const
Return the points corresponding to this face.
virtual bool write()
Write the output fields.
Mesh data needed to do the Finite Volume discretisation.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
const volScalarField & pa() const
Return const reference to pressure.
const volVectorField & Ua() const
Return const reference to velocity.
Base class for incompressibleAdjoint solvers.
virtual void additionalSensitivityMapTerms(boundaryVectorField &sensitivityMap, const labelHashSet &patchIDs, const scalar dt)
const autoPtr< incompressibleAdjoint::adjointRASModel > & adjointTurbulence() const
Return const reference to the adjointRASModel.
const volVectorField & U() const
Return const reference to velocity.
const volScalarField & p() const
Return const reference to pressure.
const autoPtr< incompressible::RASModelVariables > & RASModelVariables() const
Return const reference to the turbulence model variables.
Solver of the adjoint to the eikonal PDE.
Solver of the adjoint to the Laplace grid displacement equation.
Abstract base class for adjoint-based sensitivities in incompressible flows.
virtual void clearSensitivities()
Zero sensitivity fields and their constituents.
incompressibleAdjointSolver & adjointSolver_
objectiveManager & objectiveManager_
incompressibleAdjointVars & adjointVars_
const incompressibleVars & primalVars_
Calculation of adjoint based sensitivities at wall points.
bool includeTransposeStresses_
Include the transpose part of the adjoint stresses.
virtual void clearSensitivities()
Zero sensitivity fields and their constituents.
bool includeDistance_
Include distance variation in sens computation.
bool useSnGradInTranposeStresses_
Use snGrad in the transpose part of the adjoint stresses.
autoPtr< adjointEikonalSolver > eikonalSolver_
virtual void assembleSensitivities()
Assemble sensitivities.
void constructGlobalPointNormalsAndAreas(vectorField &pointNormals, scalarField &pointMagSf)
Construct globally correct point normals and point areas.
void setSuffixName()
Set suffix name for sensitivity fields.
autoPtr< boundaryVectorField > wallFaceSens_
The face-based part of the sensitivities.
bool includePressureTerm_
Include the adjoint pressure term in sens computation.
bool includeGradStressTerm_
Include the term containing the grad of the stress at the boundary.
void finaliseFaceMultiplier()
bool includeSurfaceArea_
Include surface area in sens computation.
void read()
Read controls and update solver pointers if necessary.
bool includeMeshMovement_
Include mesh movement variation in sens computation.
bool includeDivTerm_
Include the term from the deviatoric part of the stresses.
virtual void accumulateIntegrand(const scalar dt)
Accumulate sensitivity integrands.
autoPtr< boundaryVectorField > dnfdbMult_
autoPtr< adjointMeshMovementSolver > meshMovementSolver_
void finalisePointSensitivities()
autoPtr< boundaryVectorField > dSfdbMult_
Multipliers of d(Sf)/db and d(nf)/db.
bool includeObjective_
Include terms directly emerging from the objective function.
PtrList< objective > & getObjectiveFunctions()
Return reference to objective functions.
virtual const faceList & faces() const
Return raw faces.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const pointField & points() const
Return raw points.
A patch is a list of labels that address the faces in the global face list.
label nPoints() const noexcept
Number of mesh points.
int myProcNo() const noexcept
Return processor number.
const dictionary & dict() const
Return the construction dictionary.
void clearSensitivities()
Zero sensitivity fields and their constituents.
void setSuffix(const word &suffix)
Set suffix.
void write()
Write sensitivity fields.
A class for managing temporary objects.
void clear() const noexcept
const word & solverName() const
Return solver name.
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
#define DebugInfo
Report an information message using Foam::Info.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
const dimensionSet dimVelocity
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)
void unzipRow(const FieldField< Field, SymmTensor< Cmpt > > &input, const direction idx, FieldField< Field, Vector< Cmpt > > &result)
Extract a symmTensor field field row (x,y,z) == (0,1,2)
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
autoPtr< typename GeometricField< Type, fvPatchField, volMesh >::Boundary > createZeroBoundaryPtr(const fvMesh &mesh, bool printAllocation=false)
#define forAll(list, i)
Loop across all elements in list.
A non-counting (dummy) refCount.