- avoids redundant dictionary searching
STYLE: remove dictionary lookupOrDefaultCompat wrapper
- deprecated and replaced by getOrDefaultCompat (2019-05).
The function is usually specific to internal keyword upgrading
(version compatibility) and unlikely to exist in any user code.
STYLE: qualify format/version/compression with IOstreamOption not IOstream
STYLE: reduce number of lookups when scanning {fa,fv}Solution
STYLE: call IOobject::writeEndDivider as static
in RASModelVariables were doing this by checking whether the
corresponding pointer was allocated. In some cases, however, even if the
field does not exist, the pointer is not null, leading to the wrong
output. Made the correspding functions virtual and overwritten their
return values in the derived classes. Kept the initial implementation in
base to facilitate the clone function.
in cases with more than one primal or adjoint solvers
TUT: removed all occurances of useSolverNameForFields
from the optimisation tutorials since it is now set
automatically.
in the sensitivity patches, symmetry::evaluate() needs access to the
internalField which does exist, leading to wrong memory access.
Fixed by specifying a calculated type fvPatchField for all patches when
creating a boundaryField<Type>
Using a symmetry(Plane) as a sensitivity patch is quite rare and
borderline wrong, but this provides a fix nonetheless.
The multiplier of grad(dxdb) is a volTensorField which, by itself, is
memory consuming. The function computing it though was sloppy in terms
of memory management, constituting the peak memory consumption during an
adjoint optimisation. Initial changes to remedy the problem include the
deallocation of some of the volTensorFields included in the computation
of grad(dxdb) once unneeded, the utilisation of volSymmTensorFields
instead of volTensorFields where possible and avoiding allocating some
unnecessary intermediate fields.
Actions to further reduce memory consumption:
- For historical reasons, the code computes/stores the transpose of
grad(dxdb), which is then transposed when used in the computation of
the FI or the ESI sensitivity derivatives. This redundant
transposition can be avoid, saving the allocation of an additional
volTensorField, but the changes need to permeate a number of places in
the code that contribute to grad(dxdb) (e.g. ATC, adjoint turbulence
models, adjoint MRF, etc).
- Allocation of unnecessary pointers in the objective class should be
avoided.
- ATCstandard, ATCUaGradU:
the ATC is now added as a dimensioned field and not as an fvMatrix
to UaEqn. This get rid of many unnecessary allocations.
- ATCstandard:
gradU is cached within the class to avoid its re-computation in
every adjoint iteration of the steady state solver.
- Inlined a number of functions within the primal and adjoint solvers.
This probably has a negligible effect since they likely were inlined
by the compiler either way.
- The momentum diffusivity at the boundary, used by the adjoint boundary
conditions, was computed for the entire field and, then, only the
boundary field of each adjoint boundary condition was used. If many
outlet boundaries exist, the entire nuEff field would be computed as
many times as the number of boundaries, leading to an unnecessary
computational overhead.
- Outlet boundary conditions (both pressure and velocity) use the local
patch gradient to compute their fluxes. This patch gradient requires
the computation of the adjacent cell gradient, which is done on the
fly, on a per patch basis. To compute this patch adjacent gradient
however, the field under the grad sign is interpolated on the entire
mesh. If many outlets exist, this leads to a huge computational
overhead. Solved by caching the interpolated field to the database and
re-using it, in a way similar to the caching of gradient fields (see
fvc::grad).
WIP: functions returning references to primal and adjoint boundary
fields within boundaryAdjointContributions seem to have a non-negligible
overhead for cases with many patches. No easy work-around here since
these are virtual and cannot be inlined.
WIP: introduced the code structure for caching the contributions to
the adjoint boundary conditions that depend only on the primal fields
and reusing. The process needs to be completed and evaluated, to make
sure that the extra code complexity is justified by gains in
performance.
is now appended by the name of the adjoint solver, if more than one
exist. This was necessary for an accurate continuation since, before
these changes, only the ma field of the last solver was written. As a
result, when restarting the first adjoint solver was reading the ma
field of the last one. No changes are needed in fvSolution and fvSchemes
w.r.t. the previous code version.
as a step towards machine-accuracy continuation of the optimisation
loop.
Additionally, control points are now written under the time/uniform
folder, to be in-line with rest of the code structure for continuation.
As a side-effect, the controlPointsDefinition in
constant/dynamicMeshDict does not need to be changed to 'fromFile'
anymore in order to perform the continuation. The 'fromFile' option is
still valid if the user wants to supply the control points manually but,
as with all other controlPointsDefinitions, it will be disregarded if the
proper file exists under the time/uniform/volumetricBSplines folder.
Before the commit, the sensitivity classes were receiving references of
the (incompressible) primal and adjoint variables. However, if
additional physics was added (energy equation, multiphase, etc), the
infrastructure wasn't convenient for accommodating (new terms in the FI
and E-SI formulations, new terms in the sensitivity map, etc).
Now, the sensitivity classes receive a reference to an
incompressibleAdjointSolver and receive the terms for the FI and
sensitivity maps through there. The latter is still WIP.
Modified adjointSimple to incorporate these changes as well.
Each solver now writes its sensitivity derivatives to its dictionary,
enabling also a binary format. If present, the sensitivities are then
re-read from the dictionary, avoiding thus possible loss of information
due to re-computation.
As a side-effect, sensitivities are computed after the completion of
each adjoint solver, instead of being computed after all adjoint solvers
have been completed.
for incompressible flows. The typical convention of appending the primal
field name with 'a' to form the adjoint field is followed for the
adjoint turbulent kinetic energy (i.e. 'ka') but since this would produce
an ugly variable name for the adjoint to omega (i.e. omegaa), the latter
is abbreviated to 'wa'.
The work is based on
\verbatim
Kavvadias, I., Papoutsis-Kiachagias, E.,
Dimitrakopoulos, G., & Giannakoglou, K. (2014).
The continuous adjoint approach to the k–$omega$ SST turbulence model with
applications in shape optimization
Engineering Optimization, 47(11), 1523-1542.
https://doi.org/10.1080/0305215X.2014.979816
\endverbatim
with changes in the discretisation of
a number of differential operators and the formulation of the adjoint to
the wall functions employed by the primal model.
Regarding the latter, the code assumes (and differentiates) the default
behaviour of nutkWallFunction (i.e. nutWallFunction::blendingType::STEPWISE)
and omegaWallFunction (i.e. omegaWallFunction::blendingType::BINOMIAL2).
Due to the availability of a number of terms required for the
formulation of the wall function for ka, the latter is implemented
within adjointkOmegaSST itself, with contributions from objective functions
implemented within kaqRWallFunction. Wall functions for wa are
implemented within waWallFunction.
The initial implementation of the above-mentioned reference was
performed by Dr. Ioannis Kavvadias
the Jacobian of an objective function, defined at the boundary, wrt nut
and gradU. Also modified the current objectives that include such
contributions
- bundles frequently used 'gather/scatter' patterns more consistently.
- combineAllGather -> combineGather + broadcast
- listCombineAllGather -> listCombineGather + broadcast
- mapCombineAllGather -> mapCombineGather + broadcast
- allGatherList -> gatherList + scatterList
- reduce -> gather + broadcast (ie, allreduce)
- The allGatherList currently wraps gatherList/scatterList, but may be
replaced with a different algorithm in the future.
STYLE: PstreamCombineReduceOps.H is mostly unneeded now
A Helmholtz-like filter is applied to the original field of sensitivity
derivatives. The corresponding PDE is solved on the sensitivity patches,
using the finite area infrastructure. A smoothing radius is needed,
which is computed based on the average 'length' of the boundary faces,
if not provided by the user explicitly.
If an faMesh is provided, it will be used; otherwise it will be created
on the fly based on either an faMeshDefinition dictionary in system or
one constructed internally based on the sensitivity patches.
- simplifies local toggling.
- centralize fileModification static variables into IOobject.
They were previously scattered between IOobject and regIOobject
Affected only the first optimisation cycle, if line search was enabled
If eta was not set explicitly, it was computed after evaluating the
directional derivative of the merit function, which was computed
wrongly, leading to an erroneous value of the extrapolated merit
function value.
fvOptionsAdjoint was needlessly duplicating a lot of the functionality
of fvOptions in order to add an interface for computing sensitivity
contributions emerging from fvOptions. To reduce this code duplication:
- fvOptionsAdjoint was removed
- the corresponding sensitivity contributions have moved to fvOptions through
virtual functions (returning a zero contribution in the base so
backwards compatibility is retained)
- all sensitivity classes that were using fvOptionsAdjoint have been
modified appropriately
- all adjoint solvers are now grabbing a reference to an fvOptionList
from the database instead of constructing an fvOptionsAdjointList
Hence, all fvOptions contributions to the adjoint equations
or the sensitivity derivatives can be given through system/fvOptions,
removing the need for separate sub-dictionaries within optimisationDict.
- Expanded the write function in the base class so that it can manage
input coming from the derived ones. This reduces a lot of code
duplication in the latter but keeps the functionality.
- Added a default width for all entries in the objective files.
- If a normalisation factor or a target is set, they are written on the
header of the objective file.
- Cosmetic/code consistency changes in various files.
- Added preLoop, loop and postLoop functions
- Added preIter, mainIter and postIter functions for each SIMPLE
iteration
- Added addMomentumSource and addPressureSource virtual functions, to
allow for additions by derived classes
fvOptions are no longer a member of incompressiblePrimalSolver but are
looked up from the registry in each iteration of each primal solver.
This means that the main system/fvOptions dictionary is read by ALL
instances of the primal solvers and the latter no longer have their
own fvOptions dict in optimisationDict. This is safe since each fvOption
is applied to a specific field and in case of many primal solvers, the
primal fields are named differently for each of them.
In addition, simple is now split in preLoop, loop and postLoop phase.
Furthermore, each SIMPLE iteration is broken down to
a preIter, mainIter and postIter phase, to allow for different behaviour
by derived classes.
Part of the (E)SI shape sensitivities depends of grad(Ua) & nf computed
on the boundary. Up to now, the code was only computing the normal part
of grad(Ua), to avoid the potentially spurious tangential component
which is computed on the cell center and extrapolated to the boundary
faces. However, for some objectives that are strongly related to the
stresses (e.g. moment, stresses), including also the tangential part of
grad(Ua) is necessary for E-SI to replicate the outcome of FI.
Extensive testing on a number of objectives/cases showed
- No regression when including the tangential part
- Improved behaviour in some rare cases (moment, stresses)
Hence, the tangential part is now included by default. The previous code
behaviour can be replicated by setting the useSnGradInTranposeStresses
flag to true.
- controlPointsDefinition is now controled by a class with
runTimeSelection.
- Added a new controlPointsDefinition option that translates, rotates
and scales a given box. The required entries have the same meaning as
in the Paraview 'Transform' filter, facilitating the transition between the
visual placement of control boxes (e.g. in Paraview) and their setup
in the code.
- Improved performance during the parameterization, sensitivity
computation and grid displacement phases by re-using already computed
basis functions.
