Solver of the adjoint to the eikonal PDE. More...
Public Member Functions | |
TypeName ("adjointEikonalSolver") | |
Runtime type information. More... | |
adjointEikonalSolver (const fvMesh &mesh, const dictionary &dict, adjointSolver &adjointSolver, 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... | |
boundaryVectorField & | distanceSensitivities () |
Return the sensitivity term depending on da. More... | |
tmp< volTensorField > | getFISensitivityTerm () const |
Return the volume-based sensitivity term depending on da. More... | |
tmp< scalarField > | topologySensitivities (const word &designVarsName) const |
Return sensitivity contribution to topology optimisation. More... | |
const volScalarField & | da () |
Return the adjoint distance field. More... | |
const volScalarField & | d () |
Return the distance field. More... | |
tmp< volVectorField > | gradEikonal () |
Return the gradient of the eikonal equation. More... | |
Protected Member Functions | |
wordList | patchTypes () const |
Return the boundary condition types for da. More... | |
tmp< surfaceScalarField > | computeYPhi () |
Compute convecting velocity. More... | |
void | read () |
Read options each time a new solution is found. More... | |
Protected Attributes | |
const fvMesh & | mesh_ |
dictionary | dict_ |
adjointSolver & | adjointSolver_ |
const labelHashSet & | sensitivityPatchIDs_ |
label | nEikonalIters_ |
scalar | tolerance_ |
scalar | epsilon_ |
labelHashSet | wallPatchIDs_ |
volScalarField | da_ |
volScalarField | source_ |
autoPtr< boundaryVectorField > | distanceSensPtr_ |
Wall face sens w.r.t. (x,y.z) More... | |
Solver of the adjoint to the eikonal PDE.
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; . . } }
Definition at line 143 of file adjointEikonalSolver.H.
adjointEikonalSolver | ( | const fvMesh & | mesh, |
const dictionary & | dict, | ||
adjointSolver & | adjointSolver, | ||
const labelHashSet & | sensitivityPatchIDs | ||
) |
Construct from components.
Definition at line 117 of file adjointEikonalSolver.C.
References Foam::read().
|
virtualdefault |
|
protected |
Return the boundary condition types for da.
Definition at line 47 of file adjointEikonalSolver.C.
References fvMesh::boundary(), adjointEikonalSolver::mesh_, UPtrList< T >::size(), adjointEikonalSolver::wallPatchIDs_, and fvPatchFieldBase::zeroGradientType().
|
protected |
Compute convecting velocity.
Definition at line 76 of file adjointEikonalSolver.C.
References adjointEikonalSolver::adjointSolver_, fvMesh::boundary(), GeometricField< Type, PatchField, GeoMesh >::boundaryFieldRef(), adjointEikonalSolver::d(), Foam::dimless, Foam::fvc::grad(), Foam::fvc::interpolate(), adjointEikonalSolver::mesh_, tmp< T >::New(), IOobjectOption::NO_READ, IOobjectOption::NO_REGISTER, IOobjectOption::NO_WRITE, patches, fvMesh::Sf(), td(), fvMesh::time(), Time::timeName(), adjointEikonalSolver::wallPatchIDs_, adjointSolver::yWall(), and Foam::Zero.
Referenced by adjointEikonalSolver::solve().
|
protected |
Read options each time a new solution is found.
Definition at line 64 of file adjointEikonalSolver.C.
References adjointEikonalSolver::dict_, Foam::constant::electromagnetic::e, adjointEikonalSolver::epsilon_, dictionary::getOrDefault(), adjointEikonalSolver::mesh_, adjointEikonalSolver::nEikonalIters_, schemesLookup::schemesDict(), dictionary::subDict(), and adjointEikonalSolver::tolerance_.
Referenced by adjointEikonalSolver::readDict(), and adjointEikonalSolver::solve().
TypeName | ( | "adjointEikonalSolver" | ) |
Runtime type information.
|
virtual |
Read dict if changed.
Definition at line 174 of file adjointEikonalSolver.C.
References dict, adjointEikonalSolver::dict_, and adjointEikonalSolver::read().
void accumulateIntegrand | ( | const scalar | dt | ) |
Accumulate source term.
Definition at line 183 of file adjointEikonalSolver.C.
References adjointSolver::adjointEikonalSource(), adjointEikonalSolver::adjointSolver_, and adjointEikonalSolver::source_.
void solve | ( | ) |
Calculate the adjoint distance field.
Definition at line 190 of file adjointEikonalSolver.C.
References adjointEikonalSolver::adjointSolver_, adjointEikonalSolver::computeYPhi(), adjointEikonalSolver::d(), adjointEikonalSolver::da_, Foam::ensightOutput::debug, Foam::dimLength, Foam::dimTime, Foam::fvm::div(), Foam::endl(), adjointEikonalSolver::epsilon_, fvOptions, Foam::gMax(), Foam::Info, Foam::fvc::laplacian(), Foam::fvm::laplacian(), Foam::mag(), adjointEikonalSolver::mesh_, adjointEikonalSolver::nEikonalIters_, options::New(), IOobjectOption::NO_READ, IOobjectOption::NO_REGISTER, IOobjectOption::NO_WRITE, Time::printExecutionTime(), adjointEikonalSolver::read(), fvMatrix< Type >::relax(), fvMatrix< Type >::solve(), adjointEikonalSolver::source_, Foam::fvm::SuSp(), td(), fvMesh::time(), Time::timeName(), adjointEikonalSolver::tolerance_, regIOobject::write(), and adjointSolver::yWall().
void reset | ( | ) |
Reset source term.
Definition at line 260 of file adjointEikonalSolver.C.
References DimensionedField< Type, GeoMesh >::dimensions(), adjointEikonalSolver::distanceSensPtr_, adjointEikonalSolver::source_, and Foam::Zero.
boundaryVectorField & distanceSensitivities | ( | ) |
Return the sensitivity term depending on da.
Definition at line 267 of file adjointEikonalSolver.C.
References adjointEikonalSolver::adjointSolver_, fvMesh::boundary(), GeometricField< Type, PatchField, GeoMesh >::boundaryField(), adjointEikonalSolver::d(), adjointEikonalSolver::da_, adjointEikonalSolver::distanceSensPtr_, Foam::endl(), Foam::Info, adjointEikonalSolver::mesh_, adjointEikonalSolver::sensitivityPatchIDs_, td(), and adjointSolver::yWall().
tmp< volTensorField > getFISensitivityTerm | ( | ) | const |
Return the volume-based sensitivity term depending on da.
Definition at line 290 of file adjointEikonalSolver.C.
References adjointEikonalSolver::adjointSolver_, IOobjectOption::AUTO_WRITE, GeometricField< Type, PatchField, GeoMesh >::correctBoundaryConditions(), adjointEikonalSolver::d(), adjointEikonalSolver::da_, DimensionedField< Type, GeoMesh >::dimensions(), Foam::dimLength, Foam::fvc::div(), Foam::endl(), adjointEikonalSolver::epsilon_, Foam::fvc::grad(), Foam::Info, Foam::fvc::interpolate(), adjointEikonalSolver::mesh_, IOobjectOption::NO_READ, GeometricField< Type, PatchField, GeoMesh >::ref(), fvMesh::Sf(), td(), fvMesh::time(), Time::timeName(), adjointEikonalSolver::wallPatchIDs_, adjointSolver::yWall(), Foam::Zero, and fvPatchFieldBase::zeroGradientType().
Referenced by shapeDesignVariables::assembleSensitivities().
tmp< scalarField > topologySensitivities | ( | const word & | designVarsName | ) | const |
Return sensitivity contribution to topology optimisation.
Definition at line 351 of file adjointEikonalSolver.C.
References fvOptions, options::New(), sensitivityTopO::postProcessSens(), td(), and Foam::Zero.
const volScalarField & da | ( | ) |
Return the adjoint distance field.
Definition at line 371 of file adjointEikonalSolver.C.
References adjointEikonalSolver::da_.
const volScalarField& d | ( | ) |
Return the distance field.
Referenced by adjointEikonalSolver::computeYPhi(), adjointEikonalSolver::distanceSensitivities(), adjointEikonalSolver::getFISensitivityTerm(), adjointEikonalSolver::gradEikonal(), and adjointEikonalSolver::solve().
tmp< volVectorField > gradEikonal | ( | ) |
Return the gradient of the eikonal equation.
Definition at line 377 of file adjointEikonalSolver.C.
References adjointEikonalSolver::adjointSolver_, adjointEikonalSolver::d(), Foam::fvc::grad(), tmp< T >::New(), td(), and adjointSolver::yWall().
|
protected |
|
protected |
Definition at line 166 of file adjointEikonalSolver.H.
Referenced by adjointEikonalSolver::read(), and adjointEikonalSolver::readDict().
|
protected |
Definition at line 168 of file adjointEikonalSolver.H.
Referenced by adjointEikonalSolver::accumulateIntegrand(), adjointEikonalSolver::computeYPhi(), adjointEikonalSolver::distanceSensitivities(), adjointEikonalSolver::getFISensitivityTerm(), adjointEikonalSolver::gradEikonal(), and adjointEikonalSolver::solve().
|
protected |
Definition at line 170 of file adjointEikonalSolver.H.
Referenced by adjointEikonalSolver::distanceSensitivities().
|
protected |
Definition at line 172 of file adjointEikonalSolver.H.
Referenced by adjointEikonalSolver::read(), and adjointEikonalSolver::solve().
|
protected |
Definition at line 174 of file adjointEikonalSolver.H.
Referenced by adjointEikonalSolver::read(), and adjointEikonalSolver::solve().
|
protected |
Definition at line 176 of file adjointEikonalSolver.H.
Referenced by adjointEikonalSolver::getFISensitivityTerm(), adjointEikonalSolver::read(), and adjointEikonalSolver::solve().
|
protected |
Definition at line 178 of file adjointEikonalSolver.H.
Referenced by adjointEikonalSolver::computeYPhi(), adjointEikonalSolver::getFISensitivityTerm(), and adjointEikonalSolver::patchTypes().
|
protected |
Definition at line 180 of file adjointEikonalSolver.H.
Referenced by adjointEikonalSolver::da(), adjointEikonalSolver::distanceSensitivities(), adjointEikonalSolver::getFISensitivityTerm(), and adjointEikonalSolver::solve().
|
protected |
Definition at line 182 of file adjointEikonalSolver.H.
Referenced by adjointEikonalSolver::accumulateIntegrand(), adjointEikonalSolver::reset(), and adjointEikonalSolver::solve().
|
protected |
Wall face sens w.r.t. (x,y.z)
Definition at line 187 of file adjointEikonalSolver.H.
Referenced by adjointEikonalSolver::distanceSensitivities(), and adjointEikonalSolver::reset().