Go to the documentation of this file.
39 namespace incompressible
53 fixedValueFvPatchScalarField::typeName
58 daTypes[patchi] = zeroGradientFvPatchScalarField::typeName;
95 volVectorField::Boundary& nybf = ny.boundaryFieldRef();
99 nybf[patchi] == -
patches[patchi].nf();
112 adjointEikonalSolver::adjointEikonalSolver
122 dict_(
dict.subOrEmptyDict(
"adjointEikonalSolver")),
123 RASModelVars_(RASModelVars),
124 adjointTurbulence_(adjointTurbulence),
125 sensitivityPatchIDs_(sensitivityPatchIDs),
129 wallPatchIDs_(mesh_.boundaryMesh().findPatchIDs<
wallPolyPatch>()),
135 mesh_.time().timeName(),
149 mesh_.time().timeName(),
157 distanceSensPtr_(createZeroBoundaryPtr<vector>(mesh_))
167 dict_ =
dict.subOrEmptyDict(
"adjointEikonalSolver");
196 Info<<
"Adjoint Eikonal Iteration : " << iter <<
endl;
207 scalar residual = daEqn.
solve().initialResidual();
217 Info<<
"\n***Reached adjoint eikonal convergence limit, iteration "
218 << iter <<
"***\n\n";
235 Info<<
"Calculating distance sensitivities " <<
endl;
245 distanceSens[patchi] =
256 Info<<
"Calculating distance sensitivities " <<
endl;
300 return tdistanceSens;
Defines the attributes of an object for which implicit objectRegistry management is supported,...
void reset()
Reset source term.
tmp< volTensorField > getFISensitivityTerm() const
Return the volume-based sensitivity term depending on da.
void accumulateIntegrand(const scalar dt)
Accumulate source term.
boundaryVectorField & distanceSensitivities()
Return the sensitivity term depending on da.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
defineTypeNameAndDebug(adjointEikonalSolver, 0)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
A class for managing temporary objects.
static constexpr const zero Zero
Global zero.
labelHashSet wallPatchIDs_
bool read(const char *buf, int32_t &val)
Same as readInt32.
void solve()
Calculate the adjoint distance field.
static word timeName(const scalar t, const int precision=precision_)
Ostream & endl(Ostream &os)
Add newline and flush stream.
T lookupOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
const volScalarField & da()
Return the adjoint distance field.
const labelHashSet & sensitivityPatchIDs_
wordList patchTypes(nPatches)
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Solver of the adjoint to the eikonal PDE.
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
wordList patchTypes() const
Return the boundary condition types for da.
dimensionedScalar pow3(const dimensionedScalar &ds)
messageStream Info
Information stream (uses stdout - output is on the master only)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
const volScalarField & d()
Return the distance field.
void read()
Read options each time a new solution is found.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
tmp< surfaceScalarField > computeYPhi()
Compute convecting velocity.
double elapsedCpuTime() const
Return CPU time (in seconds) from the start.
tmp< volVectorField > gradEikonal()
Return the gradient of the eikonal equation.
virtual tmp< volScalarField > distanceSensitivities()=0
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
autoPtr< Foam::incompressibleAdjoint::adjointRASModel > & adjointTurbulence_
Mesh data needed to do the Finite Volume discretisation.
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
virtual bool readDict(const dictionary &dict)
Read dict if changed.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const dimensionedScalar e
Elementary charge.
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
const autoPtr< incompressible::RASModelVariables > & RASModelVars_
const polyBoundaryMesh & patches
A special matrix type and solver, designed for finite volume solutions of scalar equations....
SolverPerformance< Type > solve(const dictionary &)
Solve returning the solution statistics.
const Time & time() const
Return the top-level database.
static const Vector< scalar > zero
autoPtr< boundaryVectorField > distanceSensPtr_
Wall face sens w.r.t. (x,y.z)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
dimensioned< tensor > dimensionedTensor
Dimensioned tensor obtained from generic dimensioned type.
Type gMax(const FieldField< Field, Type > &f)
double elapsedClockTime() const
Returns wall-clock time from clock instantiation.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const surfaceVectorField & Sf() const
Return cell face area vectors.