Go to the documentation of this file.
48 objectiveIncompressible,
60 const word& adjointSolverName,
61 const word& primalSolverName
67 mesh_.boundaryMesh().patchSet
72 momentDirection_(
dict.get<
vector>(
"direction")),
73 rotationCentre_(
dict.get<
vector>(
"rotationCenter")),
74 Aref_(
dict.get<scalar>(
"Aref")),
75 lRef_(
dict.get<scalar>(
"lRef")),
76 rhoInf_(
dict.get<scalar>(
"rhoInf")),
77 UInf_(
dict.get<scalar>(
"UInf")),
78 invDenom_(2./(rhoInf_*UInf_*UInf_*Aref_*lRef_)),
81 Foam::createZeroFieldPtr<vector>
88 Foam::createZeroFieldPtr<vector>
95 Foam::createZeroFieldPtr<vector>
100 devReff_(vars_.turbulence()->devReff()())
103 if (momentPatches_.empty())
106 <<
"No valid patch name on which to minimize " <<
type() <<
endl
111 Info<<
"Minimizing " <<
type() <<
" in patches:" <<
endl;
112 for (
const label patchI : momentPatches_)
114 Info<<
"\t " << mesh_.boundary()[patchI].name() <<
endl;
119 bdJdpPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
120 bdSdbMultPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
121 bdxdbMultPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
122 bdxdbDirectMultPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
138 for (
const label patchI : momentPatches_)
143 pressureMoment +=
gSum
145 rhoInf_*(dx ^ Sf)*
p.boundaryField()[patchI]
149 viscousMoment +=
gSum
155 cumulativeMoment = pressureMoment + viscousMoment;
157 scalar moment = cumulativeMoment & momentDirection_;
158 scalar Cm = moment*invDenom_;
160 "Moment|Coeff " << moment <<
"|" << Cm <<
endl;
175 devReff_ = turbVars->devReff(lamTransp,
U)();
182 for (
const label patchI : momentPatches_)
186 bdJdpPtr_()[patchI] = (momentDirection_ ^ dx)*invDenom_*rhoInf_;
195 for (
const label patchI : momentPatches_)
199 bdSdbMultPtr_()[patchI] =
204 (momentDirection_^dx) &
210 + rhoInf_ * (momentDirection_^dx) *
p.boundaryField()[patchI]
233 if (isA<wallFvPatch>(
patch))
237 gradU.boundaryFieldRef()[patchI] =
238 nf *
U.boundaryField()[patchI].snGrad();
244 stressXPtr_().replace(0, stress.
component(0));
245 stressXPtr_().replace(1, stress.
component(1));
246 stressXPtr_().replace(2, stress.
component(2));
248 stressYPtr_().replace(0, stress.
component(3));
249 stressYPtr_().replace(1, stress.
component(4));
250 stressYPtr_().replace(2, stress.
component(5));
252 stressZPtr_().replace(0, stress.
component(6));
253 stressZPtr_().replace(1, stress.
component(7));
254 stressZPtr_().replace(2, stress.
component(8));
262 for (
const label patchI : momentPatches_)
269 bdxdbMultPtr_()[patchI] =
289 for (
const label patchI : momentPatches_)
299 ((
p.boundaryField()[patchI]*nf)
303 bdxdbDirectMultPtr_()[patchI] =
304 (force^momentDirection_)*invDenom_*rhoInf_;
int debug
Static debugging option.
tmp< GeometricField< cmptType, PatchField, GeoMesh > > component(const direction) const
Return a component of the field.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
A class for handling words, derived from Foam::string.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
const incompressibleVars & vars_
addToRunTimeSelectionTable(objectiveIncompressible, objectiveForce, dictionary)
Ostream & endl(Ostream &os)
Add newline and flush stream.
void update_dxdbMultiplier()
Update delta(x)/delta b multiplier.
Type gSum(const FieldField< Field, Type > &f)
A simple single-phase transport model based on viscosityModel.
const volScalarField & p() const
Return const reference to pressure.
#define forAll(list, i)
Loop across all elements in list.
const autoPtr< incompressible::RASModelVariables > & RASModelVariables() const
Return const reference to the turbulence model variables.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
defineTypeNameAndDebug(objectiveForce, 0)
scalar J_
Objective function value and weight.
messageStream Info
Information stream (stdout output on master, null elsewhere)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
const singlePhaseTransportModel & laminarTransport() const
Return const reference to transport model.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const volVectorField & U() const
Return const reference to velocity.
Abstract base class for objective functions in incompressible flows.
void update_dSdbMultiplier()
Update delta(n dS)/delta b multiplier.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
objectiveMoment(const fvMesh &mesh, const dictionary &dict, const word &adjointSolverName, const word &primalSolverName)
from components
Macros for easy insertion into run-time selection tables.
Mesh data needed to do the Finite Volume discretisation.
const autoPtr< incompressible::turbulenceModel > & turbulence() const
Return const reference to the turbulence model.
void update_meanValues()
Update mean drag and lift values.
errorManipArg< error, int > exit(error &err, const int errNo=1)
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
void update_dxdbDirectMultiplier()
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
#define DebugInfo
Report an information message using Foam::Info.
virtual tmp< volScalarField > nu() const
Return the laminar viscosity.
const std::string patch
OpenFOAM patch number as a std::string.
scalar J()
Return the objective function value.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
void update_boundarydJdp()
Update values to be added to the adjoint wall velocity.
A List of wordRe with additional matching capabilities.
autoPtr< boundaryVectorField > bdJdpPtr_
Adjoint (intlet,wall) velocity.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const volScalarField & pInst() const
Return const reference to pressure.