adjointEikonalSolver Class Reference

Solver of the adjoint to the eikonal PDE. More...

Collaboration diagram for adjointEikonalSolver:
[legend]

Public Member Functions

 TypeName ("adjointEikonalSolver")
 Runtime type information. More...
 
 adjointEikonalSolver (const fvMesh &mesh, const dictionary &dict, const autoPtr< incompressible::RASModelVariables > &RASModelVars, incompressibleAdjointVars &adjointVars, const labelHashSet &sensitivityPatchIDs)
 Construct from components. More...
 
virtual ~adjointEikonalSolver ()=default
 
virtual bool readDict (const dictionary &dict)
 Read dict if changed. More...
 
void accumulateIntegrand (const scalar dt)
 Accumulate source term. More...
 
void solve ()
 Calculate the adjoint distance field. More...
 
void reset ()
 Reset source term. More...
 
boundaryVectorFielddistanceSensitivities ()
 Return the sensitivity term depending on da. More...
 
tmp< volTensorFieldgetFISensitivityTerm () const
 Return the volume-based sensitivity term depending on da. More...
 
const volScalarFieldda ()
 Return the adjoint distance field. More...
 
const volScalarFieldd ()
 Return the distance field. More...
 
tmp< volVectorFieldgradEikonal ()
 Return the gradient of the eikonal equation. More...
 

Protected Member Functions

wordList patchTypes () const
 Return the boundary condition types for da. More...
 
tmp< surfaceScalarFieldcomputeYPhi ()
 Compute convecting velocity. More...
 
void read ()
 Read options each time a new solution is found. More...
 

Protected Attributes

const fvMeshmesh_
 
dictionary dict_
 
const autoPtr< incompressible::RASModelVariables > & RASModelVars_
 
autoPtr< Foam::incompressibleAdjoint::adjointRASModel > & adjointTurbulence_
 
const labelHashSetsensitivityPatchIDs_
 
label nEikonalIters_
 
scalar tolerance_
 
scalar epsilon_
 
labelHashSet wallPatchIDs_
 
volScalarField da_
 
volScalarField source_
 
autoPtr< boundaryVectorFielddistanceSensPtr_
 Wall face sens w.r.t. (x,y.z) More...
 

Detailed Description

Solver of the adjoint to the eikonal PDE.

Reference:

    For the development of the adjoint eikonal PDE and its boundary
    conditions

        Papoutsis-Kiachagias, E. M., & Giannakoglou, K. C. (2014).
        Continuous Adjoint Methods for Turbulent Flows, Applied to Shape
        and Topology Optimization: Industrial Applications.
        Archives of Computational Methods in Engineering, 23(2), 255-299.
        http://doi.org/10.1007/s11831-014-9141-9

To be as consistent as possible, it is recommended to use the advectionDiffusion wallDist method in fvSchemes, instead of the more widely used meshWave

Example of the adjointEikonalSolver specification in optimisationDict:

    optimisation
    {
        sensitivities
        {
            includeDistance true;
            adjointEikonalSolver
            {
                // epsilon should be the same as the one used
                // in fvSchemes/wallDist/advectionDiffusionCoeffs
                epsilon   0.1;
                iters     1000;
                tolerance 1e-6;
            }
        }
    }

Example of the entries in fvSchemes:

    divSchemes
    {
        .
        .
        // avoid bounded schemes since yPhi is not conservative
        div(-yPhi,da) Gauss linearUpwind grad(da);
        .
        .
    }
    laplacianSchemes
    {
        .
        .
        laplacian(yWall,da) Gauss linear corrected;
        .
        .
    }

Also, the solver specification and a relaxation factor for da are required in fvSolution

    da
    {
        solver          PBiCGStab;
        preconditioner  DILU;
        tolerance       1e-9;
        relTol          0.1;
    }

    relaxationFactors
    {
        equations
        {
            .
            .
            da           0.5;
            .
            .
        }
    }
See also
Foam::patchDistMethod::advectionDiffusion Foam::wallDist
Source files

Definition at line 146 of file adjointEikonalSolverIncompressible.H.

Constructor & Destructor Documentation

◆ adjointEikonalSolver()

adjointEikonalSolver ( const fvMesh mesh,
const dictionary dict,
const autoPtr< incompressible::RASModelVariables > &  RASModelVars,
incompressibleAdjointVars adjointVars,
const labelHashSet sensitivityPatchIDs 
)

Construct from components.

Definition at line 117 of file adjointEikonalSolverIncompressible.C.

References adjointEikonalSolver::read().

Here is the call graph for this function:

◆ ~adjointEikonalSolver()

virtual ~adjointEikonalSolver ( )
virtualdefault

Member Function Documentation

◆ patchTypes()

wordList patchTypes ( ) const
protected

Return the boundary condition types for da.

Definition at line 48 of file adjointEikonalSolverIncompressible.C.

References fvMesh::boundary(), adjointEikonalSolver::mesh_, UPtrList< T >::size(), and adjointEikonalSolver::wallPatchIDs_.

Here is the call graph for this function:

◆ computeYPhi()

tmp< surfaceScalarField > computeYPhi ( )
protected

Compute convecting velocity.

Definition at line 78 of file adjointEikonalSolverIncompressible.C.

References fvMesh::boundary(), GeometricField< Type, PatchField, GeoMesh >::boundaryFieldRef(), adjointEikonalSolver::d(), Foam::dimless, Foam::fvc::grad(), Foam::fvc::interpolate(), adjointEikonalSolver::mesh_, Time::New(), IOobject::NO_READ, IOobject::NO_WRITE, patches, adjointEikonalSolver::RASModelVars_, fvMesh::Sf(), fvMesh::time(), Time::timeName(), adjointEikonalSolver::wallPatchIDs_, and Foam::Zero.

Referenced by adjointEikonalSolver::solve().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ read()

void read ( )
protected

Read options each time a new solution is found.

Definition at line 65 of file adjointEikonalSolverIncompressible.C.

References adjointEikonalSolver::dict_, e, adjointEikonalSolver::epsilon_, dictionary::getOrDefault(), adjointEikonalSolver::mesh_, adjointEikonalSolver::nEikonalIters_, schemesLookup::schemesDict(), dictionary::subDict(), and adjointEikonalSolver::tolerance_.

Referenced by adjointEikonalSolver::adjointEikonalSolver(), and adjointEikonalSolver::solve().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ TypeName()

TypeName ( "adjointEikonalSolver"  )

Runtime type information.

◆ readDict()

bool readDict ( const dictionary dict)
virtual

Read dict if changed.

Definition at line 175 of file adjointEikonalSolverIncompressible.C.

References dict, adjointEikonalSolver::dict_, and dictionary::subOrEmptyDict().

Here is the call graph for this function:

◆ accumulateIntegrand()

void accumulateIntegrand ( const scalar  dt)

Accumulate source term.

Definition at line 183 of file adjointEikonalSolverIncompressible.C.

References adjointEikonalSolver::adjointTurbulence_, adjointRASModel::distanceSensitivities(), and adjointEikonalSolver::source_.

Here is the call graph for this function:

◆ solve()

◆ reset()

void reset ( )

Reset source term.

Definition at line 237 of file adjointEikonalSolverIncompressible.C.

References DimensionedField< Type, GeoMesh >::dimensions(), adjointEikonalSolver::distanceSensPtr_, adjointEikonalSolver::source_, and Foam::Zero.

Here is the call graph for this function:

◆ distanceSensitivities()

boundaryVectorField & distanceSensitivities ( )

Return the sensitivity term depending on da.

Definition at line 244 of file adjointEikonalSolverIncompressible.C.

References fvMesh::boundary(), GeometricField< Type, PatchField, GeoMesh >::boundaryField(), adjointEikonalSolver::d(), adjointEikonalSolver::da_, adjointEikonalSolver::distanceSensPtr_, Foam::endl(), Foam::Info, adjointEikonalSolver::mesh_, adjointEikonalSolver::RASModelVars_, and adjointEikonalSolver::sensitivityPatchIDs_.

Here is the call graph for this function:

◆ getFISensitivityTerm()

tmp< volTensorField > getFISensitivityTerm ( ) const

Return the volume-based sensitivity term depending on da.

Definition at line 265 of file adjointEikonalSolverIncompressible.C.

References IOobject::AUTO_WRITE, adjointEikonalSolver::d(), adjointEikonalSolver::da_, DimensionedField< Type, GeoMesh >::dimensions(), Foam::dimLength, Foam::endl(), adjointEikonalSolver::epsilon_, Foam::fvc::grad(), Foam::Info, adjointEikonalSolver::mesh_, IOobject::NO_READ, adjointEikonalSolver::RASModelVars_, tmp< T >::ref(), fvMesh::time(), Time::timeName(), adjointEikonalSolver::wallPatchIDs_, and Foam::Zero.

Here is the call graph for this function:

◆ da()

const volScalarField & da ( )

Return the adjoint distance field.

Definition at line 315 of file adjointEikonalSolverIncompressible.C.

References adjointEikonalSolver::da_.

◆ d()

const volScalarField & d ( )

Return the distance field.

Referenced by adjointEikonalSolver::computeYPhi(), adjointEikonalSolver::distanceSensitivities(), adjointEikonalSolver::getFISensitivityTerm(), adjointEikonalSolver::gradEikonal(), and adjointEikonalSolver::solve().

Here is the caller graph for this function:

◆ gradEikonal()

tmp< volVectorField > gradEikonal ( )

Return the gradient of the eikonal equation.

Definition at line 321 of file adjointEikonalSolverIncompressible.C.

References adjointEikonalSolver::d(), Foam::fvc::grad(), Time::New(), and adjointEikonalSolver::RASModelVars_.

Here is the call graph for this function:

Member Data Documentation

◆ mesh_

◆ dict_

◆ RASModelVars_

◆ adjointTurbulence_

◆ sensitivityPatchIDs_

const labelHashSet& sensitivityPatchIDs_
protected

◆ nEikonalIters_

label nEikonalIters_
protected

◆ tolerance_

scalar tolerance_
protected

◆ epsilon_

◆ wallPatchIDs_

◆ da_

◆ source_

◆ distanceSensPtr_

autoPtr<boundaryVectorField> distanceSensPtr_
protected

Wall face sens w.r.t. (x,y.z)

Definition at line 187 of file adjointEikonalSolverIncompressible.H.

Referenced by adjointEikonalSolver::distanceSensitivities(), and adjointEikonalSolver::reset().


The documentation for this class was generated from the following files: