Go to the documentation of this file.
58 const word& adjointSolverName,
59 const word& primalSolverName
65 mesh_.boundaryMesh().patchSet
70 forceDirection_(
dict.get<
vector>(
"direction")),
71 Aref_(
dict.get<scalar>(
"Aref")),
72 rhoInf_(
dict.get<scalar>(
"rhoInf")),
73 UInf_(
dict.get<scalar>(
"UInf")),
76 Foam::createZeroFieldPtr<vector>
83 Foam::createZeroFieldPtr<vector>
90 Foam::createZeroFieldPtr<vector>
97 if (forcePatches_.empty())
100 <<
"No valid patch name on which to minimize " <<
type() <<
endl
105 Info<<
"Minimizing " <<
type() <<
" in patches:" <<
endl;
106 for (
const label patchI : forcePatches_)
108 Info<<
"\t " << mesh_.boundary()[patchI].name() <<
endl;
113 bdJdpPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
114 bdSdbMultPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
115 bdxdbMultPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
116 bdJdStressPtr_.reset(createZeroBoundaryPtr<tensor>(mesh_));
137 pressureForce +=
gSum
144 devReff.boundaryField()[patchI]
149 cumulativeForce = pressureForce + viscousForce;
158 <<
"Force|Coeff " << force <<
"|" << Cforce <<
endl;
189 bdSdbMultPtr_()[patchI] =
215 volTensorField::Boundary& gradUbf = gradU.boundaryFieldRef();
222 if (isA<wallFvPatch>(
patch))
225 gradUbf[patchI] = nf*
U.boundaryField()[patchI].snGrad();
255 bdxdbMultPtr_()[patchI] =
int debug
Static debugging option.
const Cmpt & x() const
Access to the vector x component.
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.
Info<< "Reading field U\n"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);volScalarField rho(IOobject("rho", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), thermo.rho());volVectorField rhoU(IOobject("rhoU", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *U);volScalarField rhoE(IOobject("rhoE", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *(e+0.5 *magSqr(U)));surfaceScalarField pos(IOobject("pos", runTime.timeName(), mesh), mesh, dimensionedScalar("pos", dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar("neg", dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\n"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
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)
void update_dxdbMultiplier()
Update delta(x)/delta b multiplier.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Type gSum(const FieldField< Field, Type > &f)
A simple single-phase transport model based on viscosityModel.
autoPtr< volVectorField > stressZPtr_
const volScalarField & p() const
Return const reference to pressure.
void update_dSdbMultiplier()
Update delta(n dS)/delta b multiplier.
const Cmpt & z() const
Access to the vector z component.
#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 (uses stdout - output is on the master only)
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.
autoPtr< volVectorField > stressXPtr_
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.
labelHashSet forcePatches_
Mesh data needed to do the Finite Volume discretisation.
const autoPtr< incompressible::turbulenceModel > & turbulence() const
Return const reference to the turbulence model.
errorManipArg< error, int > exit(error &err, const int errNo=1)
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
const vector & forceDirection() const
Return force direction.
const Cmpt & y() const
Access to the vector y component.
objectiveForce(const fvMesh &mesh, const dictionary &dict, const word &adjointSolverName, const word &primalSolverName)
Construct from components.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
#define DebugInfo
Report an information message using Foam::Info.
void update_boundarydJdp()
Update values to be added to the adjoint wall velocity.
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.
virtual scalar denom() const
Return denominator, without density.
A List of wordRe with additional matching capabilities.
void update_dJdStressMultiplier()
Update dJ/dStress multiplier.
autoPtr< volVectorField > stressYPtr_
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.
const surfaceVectorField & Sf() const
Return cell face area vectors.