Go to the documentation of this file.
49 objectiveIncompressible,
61 const word& adjointSolverName,
62 const word& primalSolverName
73 bdJdpPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
74 bdJdvPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
75 bdJdvnPtr_.reset(createZeroBoundaryPtr<scalar>(mesh_));
76 bdJdvtPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
97 <<
"No patches provided to PtLosses. "
98 <<
"Choosing them according to the patch mass flows" <<
nl;
107 const scalar mass =
gSum(phiPatch);
108 if (
mag(mass) > SMALL)
110 objectiveReportPatches.append(patchI);
114 patches_.
transfer(objectiveReportPatches);
116 patchPt_.setSize(patches_.size());
118 if (patches_.empty())
121 <<
"No valid patch name on which to minimize " <<
type() <<
endl
126 Info<<
"Minimizing " <<
type() <<
" in patches:" <<
endl;
146 const label patchI = patches_[oI];
150 (
U.boundaryField()[patchI] & Sf)
152 p.boundaryField()[patchI]
153 + 0.5*
magSqr(
U.boundaryField()[patchI])
156 patchPt_[oI] =
mag(pt);
170 const label patchI = patches_[oI];
175 bdJdpPtr_()[patchI] = -(
U.boundaryField()[patchI] & nf)*nf;
187 const label patchI = patches_[oI];
194 - (
p.boundaryField()[patchI] + 0.5*
magSqr(Ub))*nf
207 const label patchI = patches_[oI];
213 -
p.boundaryField()[patchI]
214 - 0.5*
magSqr(
U.boundaryField()[patchI])
215 -
sqr(
U.boundaryField()[patchI] & nf);
226 const label patchI = patches_[oI];
232 bdJdvtPtr_()[patchI] = -Un*(
U.boundaryField()[patchI] - Un*nf);
245 if (objFunctionFilePtr_.empty())
247 setObjectiveFilePtr();
248 objFunctionFilePtr_() <<
setw(4) <<
"#" <<
" ";
249 objFunctionFilePtr_() <<
setw(width) <<
"ptLosses" <<
" ";
250 objFunctionFilePtr_() <<
setw(width) <<
"JCycle" <<
" ";
253 label patchI = patches_[oI];
254 objFunctionFilePtr_()
257 objFunctionFilePtr_() <<
endl;
261 objFunctionFilePtr_() <<
setw(width) <<
J_ <<
" ";
262 objFunctionFilePtr_() <<
setw(width) << JCycle() <<
" ";
265 objFunctionFilePtr_() <<
setw(width) << patchPt_[pI] <<
" ";
267 objFunctionFilePtr_() <<
endl;
int debug
Static debugging option.
void initialize()
Return the objectiveReportPatches.
const surfaceScalarField & phiInst() const
Return const reference to volume flux.
A class for handling words, derived from Foam::string.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
const incompressibleVars & vars_
addToRunTimeSelectionTable(objectiveIncompressible, objectiveForce, dictionary)
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
autoPtr< boundaryScalarField > bdJdvnPtr_
Adjoint outlet pressure.
void update_boundarydJdv()
Update values to be added to the adjoint outlet velocity.
autoPtr< boundaryVectorField > bdJdvPtr_
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const Type & value() const
Return const reference to value.
Type gSum(const FieldField< Field, Type > &f)
const volVectorField & UInst() const
Return const reference to velocity.
autoPtr< boundaryVectorField > bdJdvtPtr_
Adjoint outlet velocity.
const volScalarField & p() const
Return const reference to pressure.
#define forAll(list, i)
Loop across all elements in list.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
defineTypeNameAndDebug(objectiveForce, 0)
scalar J_
Objective function value and weight.
messageStream Info
Information stream (uses stdout - output is on the master only)
void transfer(List< T > &list)
Istream and Ostream manipulators taking arguments.
const volVectorField & U() const
Return const reference to velocity.
Abstract base class for objective functions in incompressible flows.
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.
Mesh data needed to do the Finite Volume discretisation.
Omanip< int > setw(const int i)
virtual bool write(const bool valid=true) const
Write objective function values and its contrituents.
List< Key > sortedToc() const
The table of contents (the keys) in sorted order.
errorManipArg< error, int > exit(error &err, const int errNo=1)
static bool master(const label communicator=0)
Am I the master process.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
scalar J()
Return the objective function value.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
void update_boundarydJdp()
Update values to be added to the adjoint inlet velocity.
static unsigned int defaultPrecision()
Return the default precision.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
A List of wordRe with additional matching capabilities.
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool useGroups=true) const
Return the set of patch IDs corresponding to the given names.
const Time & time() const
Return the top-level database.
autoPtr< boundaryVectorField > bdJdpPtr_
Adjoint (intlet,wall) velocity.
objectivePtLosses(const fvMesh &mesh, const dictionary &dict, const word &adjointSolverName, const word &primalSolverName)
from components
void update_boundarydJdvn()
Update values to be added to the adjoint outlet pressure.
#define WarningInFunction
Report a warning using Foam::Warning.
void update_boundarydJdvt()
Update values to be added to the adjoint outlet tangential velocity.
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
const volScalarField & pInst() const
Return const reference to pressure.