47 objectiveIncompressible,
48 objectiveFlowRatePartition,
59 const word& adjointSolverName,
60 const word& primalSolverName
78 targetFlowRateFraction_(),
79 currentFlowRateFraction_(outletPatches_.size(),
Zero),
81 flowRateDifference_(outletPatches_.size(),
Zero)
89 "targetFractions", {{
"targetPercentages", 2306}},
90 targetFlowRateFraction_
94 const label nOutPatches = outletPatches_.size();
95 targetFlowRateFraction_.setSize(nOutPatches, 1.0/scalar(nOutPatches));
98 if (targetFlowRateFraction_.size() != outletPatches_.size())
101 <<
"Inconsistent sizes for targetFractions and outletPatches" 107 bdJdvPtr_.reset(createZeroBoundaryPtr<vector>(mesh_));
108 bdJdvnPtr_.reset(createZeroBoundaryPtr<scalar>(mesh_));
111 for (
const label patchI : outletPatches_)
113 const fvPatch&
patch = mesh_.boundary()[patchI];
114 unsigned int wordSize = word(
patch.name() +
"Ratio").size();
115 width_ =
max(width_, wordSize);
122 scalar objectiveFlowRatePartition::J()
129 for (
const label patchI : inletPatches_)
132 inletFlowRate_ +=
gSum(
phi.boundaryField()[patchI]);
136 forAll(outletPatches_, pI)
138 const label patchI = outletPatches_[pI];
139 const scalar outletFlowRate =
gSum(
phi.boundaryField()[patchI]);
140 currentFlowRateFraction_[pI] = -outletFlowRate/inletFlowRate_;
141 flowRateDifference_[pI] =
142 targetFlowRateFraction_[pI] - currentFlowRateFraction_[pI];
143 J_ += 0.5*flowRateDifference_[pI]*flowRateDifference_[pI];
150 void objectiveFlowRatePartition::update_boundarydJdv()
152 forAll(outletPatches_, pI)
154 const label patchI = outletPatches_[pI];
156 bdJdvPtr_()[patchI] = nf*flowRateDifference_[pI]/inletFlowRate_;
161 void objectiveFlowRatePartition::update_boundarydJdvn()
163 forAll(outletPatches_, pI)
165 const label patchI = outletPatches_[pI];
166 bdJdvnPtr_()[patchI] = flowRateDifference_[pI]/inletFlowRate_;
171 void objectiveFlowRatePartition::addHeaderInfo()
const 173 objFunctionFilePtr_()
174 <<
setw(width_) <<
"#inletFlowRate" <<
" " 175 <<
setw(width_) << -inletFlowRate_ <<
endl;
176 forAll(outletPatches_, pI)
178 const label patchI = outletPatches_[pI];
180 objFunctionFilePtr_()
182 <<
setw(width_) << targetFlowRateFraction_[pI] <<
endl;
187 void objectiveFlowRatePartition::addHeaderColumns()
const 189 for (
const label patchI : outletPatches_)
192 objFunctionFilePtr_()
198 void objectiveFlowRatePartition::addColumnValues()
const 200 for (
const scalar flowRate : currentFlowRateFraction_)
202 objFunctionFilePtr_()
203 <<
setw(width_) << flowRate <<
" ";
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
List of values generated by applying the access operation to each list item.
bool readIfPresentCompat(const word &keyword, std::initializer_list< std::pair< const char *, int >> compat, T &val, enum keyType::option matchOpt=keyType::REGEX) const
Find an entry if present, and assign to T val using any compatibility names if needed. FatalIOError if it is found and there are excess tokens.
const dictionary & dict() const
Return objective dictionary.
errorManipArg< error, int > exit(error &err, const int errNo=1)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
objectiveFlowRatePartition(const fvMesh &mesh, const dictionary &dict, const word &adjointSolverName, const word &primalSolverName)
from components
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface...
Ostream & endl(Ostream &os)
Add newline and flush stream.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Macros for easy insertion into run-time selection tables.
#define forAll(list, i)
Loop across all elements in list.
Type gSum(const FieldField< Field, Type > &f)
A class for handling words, derived from Foam::string.
defineTypeNameAndDebug(objectivePartialVolume, 1)
A List of wordRe with additional matching capabilities.
Istream and Ostream manipulators taking arguments.
addToRunTimeSelectionTable(objectiveGeometric, objectivePartialVolume, dictionary)
Mesh data needed to do the Finite Volume discretisation.
const std::string patch
OpenFOAM patch number as a std::string.
List< label > sortedToc(const UList< bool > &bools)
Return the (sorted) values corresponding to 'true' entries.
Omanip< int > setw(const int i)
A class for managing temporary objects.
Abstract base class for objective functions in incompressible flows.
static constexpr const zero Zero
Global zero (0)