Go to the documentation of this file.
54 auto tresGrad = tmp<Field<GradType>>
::New(patch_.size(),
Zero);
55 auto& resGrad = tresGrad.ref();
57 const labelList& faceCells = patch_.faceCells();
63 const GeometricField<Type2, fvPatchField, volMesh>&
field =
76 tmp<surfaceInterpolationScheme<Type2>> tinterpScheme
85 GeometricField<Type2, fvsPatchField, surfaceMesh> surfField
92 tmp<vectorField> tnf = patch_.nf();
100 const label cI = faceCells[fI];
101 const cell& cellI =
cells[cI];
102 for (
const label faceI : cellI)
107 const label own = owner[faceI];
123 const label boundaryFaceI = faceI - patchForFlux.start();
127 *surfField.boundaryField()[
patchID][boundaryFaceI];
131 resGrad[fI] /= V[cI];
136 const fvPatchField<Type2>& bField =
field.boundaryField()[patch_.index()];
137 resGrad = nf*bField.snGrad() + (resGrad - nf*(nf & resGrad));
146 if (!addATCUaGradUTerm_)
148 addATCUaGradUTerm_.reset(
new bool(isA<ATCUaGradU>(getATC())));
150 return addATCUaGradUTerm_();
187 const word& solverName
191 managerName_(
"objectiveManager" + solverName),
192 adjointSolverName_(solverName),
193 simulationType_(
"incompressible"),
194 boundaryContrPtr_(
nullptr),
195 addATCUaGradUTerm_(
nullptr)
198 setBoundaryContributionPtr();
214 return adjointSolverName_;
221 return simulationType_;
235 boundaryContrPtr_.
reset
249 <<
"No objectiveManager " << managerName_ <<
" available." <<
nl
250 <<
"Setting boundaryAdjointContributionPtr to nullptr. " <<
nl
251 <<
"OK for decomposePar."
261 return boundaryContrPtr_();
269 patch_.boundaryMesh().mesh().template
270 lookupObject<ATCModel>(
"ATCModel" + adjointSolverName_);
List< label > labelList
A List of labels.
boundaryAdjointContribution & getBoundaryAdjContribution()
Get boundaryContribution.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static tmp< surfaceInterpolationScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return new tmp interpolation scheme.
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
adjointBoundaryCondition(const fvPatch &p, const DimensionedField< Type, volMesh > &iF, const word &solverName)
Construct from field and base name.
A class for handling words, derived from Foam::string.
const fvPatch & patch_
Reference to patch.
const word & simulationType() const
Return the simulationType.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
typeOfRank< typename pTraits< arg1 >::cmptType, direction(pTraits< arg1 >::rank)+direction(pTraits< arg2 >::rank) >::type type
const cellList & cells() const
Base class for solution control classes.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Base class for selecting the adjoint transpose convection model. Inherits from regIOobject to add loo...
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
const ATCModel & getATC() const
ATC type might be useful for a number of BCs. Return here.
#define forAll(list, i)
Loop across all elements in list.
Field< vector > vectorField
Specialisation of Field<T> for vector.
void setBoundaryContributionPtr()
Set the ptr to the correct boundaryAdjointContribution.
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
const polyMesh & mesh() const noexcept
Return the mesh reference.
Generic templated field type.
ITstream & interpolationScheme(const word &name) const
Get interpolation scheme for given name, or default.
word adjointSolverName_
adjointSolver name corresponding to field
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
List< cell > cellList
A List of cells.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
const Type & lookupObject(const word &name, const bool recursive=false) const
word simulationType_
simulationType corresponding to field.
bool addATCUaGradUTerm()
Whether to add the extra term from the UaGradU formulation.
tmp< Field< typename Foam::outerProduct< Foam::vector, Type2 >::type > > computePatchGrad(word name)
Get gradient of field on a specific boundary.
word managerName_
objectiveManager name corresponding to field
Mesh data needed to do the Finite Volume discretisation.
const labelUList & owner() const
Internal face owner. Note bypassing virtual mechanism so.
GeometricField< vector, fvPatchField, volMesh > volVectorField
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Abstract base class for computing contributions of the objective functions to the adjoint boundary co...
virtual tmp< Field< typename Foam::outerProduct< Foam::vector, Type >::type > > dxdbMult() const
Return contribution to sensitivity derivatives.
const word & objectiveManagerName() const
Return objectiveManager name.
void reset(const label nPoints, const label nInternalFaces, const label nFaces, const label nCells)
Reset this primitiveMesh given the primitive array sizes.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
const word & adjointSolverName() const
Return adjointSolverName.
UList< label > labelUList
A UList of labels.
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
#define WarningInFunction
Report a warning using Foam::Warning.
Tensor< scalar > tensor
Tensor of scalars, i.e. Tensor<scalar>.
static autoPtr< boundaryAdjointContribution > New(const word &managerName, const word &adjointSolverName, const word &simulationType, const fvPatch &patch)
Return a reference to the selected turbulence model.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
autoPtr< bool > addATCUaGradUTerm_
Whether to add the extra term from the UaGradU formulation.
const surfaceVectorField & Sf() const
Return cell face area vectors.