Go to the documentation of this file.
42 namespace incompressible
52 adjointSensitivity::adjointSensitivity
63 primalVars_(primalVars),
64 adjointVars_(adjointVars),
82 Info<<
"adjointSensitivity type : " << modelType <<
endl;
84 auto* ctorPtr = dictionaryConstructorTable(modelType);
93 *dictionaryConstructorTablePtr_
187 if (isA<wallFvPatch>(
patch))
192 nf*
U.boundaryField()[patchI].snGrad();
201 createZeroFieldPtr<vector>(
mesh_,
"stressX", stress.dimensions())
205 createZeroFieldPtr<vector>(
mesh_,
"stressY", stress.dimensions())
209 createZeroFieldPtr<vector>(
mesh_,
"stressZ", stress.dimensions())
212 stressXPtr().replace(0, stress.
component(0));
213 stressXPtr().replace(1, stress.
component(1));
214 stressXPtr().replace(2, stress.
component(2));
216 stressYPtr().replace(0, stress.
component(3));
217 stressYPtr().replace(1, stress.
component(4));
218 stressYPtr().replace(2, stress.
component(5));
220 stressZPtr().replace(0, stress.
component(6));
221 stressZPtr().replace(1, stress.
component(7));
222 stressZPtr().replace(2, stress.
component(8));
233 "objectiveContributions",
245 objectiveContributions +=
246 functions[funcI].weight()
247 *functions[funcI].gradDxDbMultiplier();
260 - nuEff*(gradUa & (gradU +
T(gradU)))
291 - nuEff*(gradU & (gradUa +
T(gradUa)))
299 + objectiveContributions;
301 flowTerm.correctBoundaryConditions();
318 "adjointMeshMovementSource",
339 return (tadjointMeshMovementSource);
class for managing incompressible objective functions.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
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)
virtual const scalarField & calculateSensitivities()
Calculates and returns sensitivity fields.
Base class for incompressibleAdjoint solvers.
A class for handling words, derived from Foam::string.
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)
static options & New(const fvMesh &mesh)
Construct fvOptions and register to database if not present.
static word timeName(const scalar t, const int precision=precision_)
virtual void assembleSensitivities()=0
Assemble sensitivities.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Ostream & endl(Ostream &os)
Add newline and flush stream.
const volScalarField & p() const
Return const reference to pressure.
defineRunTimeSelectionTable(adjointSensitivity, dictionary)
Base class for selecting the adjoint transpose convection model. Inherits from regIOobject to add loo...
Class including all adjoint fields for incompressible flows.
#define forAll(list, i)
Loop across all elements in list.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
dimensionedScalar pow3(const dimensionedScalar &ds)
messageStream Info
Information stream (stdout output on master, null elsewhere)
static autoPtr< adjointSensitivity > New(const fvMesh &mesh, const dictionary &dict, incompressibleVars &primalVars, incompressibleAdjointVars &adjointVars, objectiveManager &objectiveManager)
Return a reference to the selected turbulence model.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual tmp< volTensorField > getFISensitivityTerm() const =0
Get the FI sensitivity derivatives term coming from the ATC.
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.
const Type & lookupObject(const word &name, const bool recursive=false) const
const volVectorField & Ua() const
Return const reference to velocity.
autoPtr< volScalarField > fieldSensPtr_
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
incompressibleAdjointVars & adjointVars_
Mesh data needed to do the Finite Volume discretisation.
void postProcessSens(Field< Type > &sensField, const word &fieldName=word::null, const word &designVariablesName=word::null)
Post process sensitivity field related to the fvOption.
incompressibleVars & primalVars_
errorManipArg< error, int > exit(error &err, const int errNo=1)
tmp< GeometricField< Type, PatchField, GeoMesh > > T() const
Return transpose (only if it is a tensor field)
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const std::string patch
OpenFOAM patch number as a std::string.
Macros to ease declaration of run-time selection tables.
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
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.
tmp< volVectorField > adjointMeshMovementSource()
Compute source term for adjoint mesh movement equation.
const scalarField & getSensitivities() const
Returns the sensitivity fields.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
const Time & time() const
Return the top-level database.
virtual void write(const word &baseName=word::null)
Write sensitivity fields.
const word & adjointSolverName() const
Return name of the adjointSolver.
Abstract base class for adjoint sensitivities.
dimensioned< tensor > dimensionedTensor
Dimensioned tensor obtained from generic dimensioned type.
tmp< volTensorField > computeGradDxDbMultiplier()
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_