From ad76df43c352745d070ffdc499dc1fa38673b3c2 Mon Sep 17 00:00:00 2001
From: Vaggelis Papoutsis <vaggelisp@gmail.com>
Date: Wed, 12 Feb 2020 19:51:23 +0200
Subject: [PATCH] ENH: refactoring and cleaning of optimisationType

Moved part common to all derived classes (e.g. update) to the base
class to avoid code duplication. Practically, only the protected
updateDesignVariables has to be overwritten in each derived class now.
steadyOptimisation was also affected in a minor way.
---
 .../optimisationManager/singleRun/singleRun.H |   7 +-
 .../steadyOptimisation/steadyOptimisation.C   |  35 ++--
 .../steadyOptimisation/steadyOptimisation.H   |   8 +-
 .../optimisationTypeIncompressible.C          |  85 +++++++--
 .../optimisationTypeIncompressible.H          |  23 ++-
 .../shapeOptimisationIncompressible.C         | 175 ++++--------------
 .../shapeOptimisationIncompressible.H         |  26 ++-
 7 files changed, 155 insertions(+), 204 deletions(-)

diff --git a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/singleRun/singleRun.H b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/singleRun/singleRun.H
index 19c3ff91234..9f5f6708ad5 100644
--- a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/singleRun/singleRun.H
+++ b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/singleRun/singleRun.H
@@ -5,8 +5,8 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2007-2019 PCOpt/NTUA
-    Copyright (C) 2013-2019 FOSS GP
+    Copyright (C) 2007-2020 PCOpt/NTUA
+    Copyright (C) 2013-2020 FOSS GP
     Copyright (C) 2019 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
@@ -107,8 +107,7 @@ public:
         virtual bool update();
 
         //- Update design variables.
-        //  Might employ a line search to find a correction satisfying the step
-        //  convergence criteria
+        //  Does nothing
         virtual void updateDesignVariables();
 };
 
diff --git a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/steadyOptimisation/steadyOptimisation.C b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/steadyOptimisation/steadyOptimisation.C
index 7d7a8dba615..470756d7f15 100644
--- a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/steadyOptimisation/steadyOptimisation.C
+++ b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/steadyOptimisation/steadyOptimisation.C
@@ -5,8 +5,8 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2007-2019 PCOpt/NTUA
-    Copyright (C) 2013-2019 FOSS GP
+    Copyright (C) 2007-2020 PCOpt/NTUA
+    Copyright (C) 2013-2020 FOSS GP
     Copyright (C) 2019 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
@@ -68,8 +68,13 @@ void Foam::steadyOptimisation::updateOptTypeSource()
 
 // * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
 
-void Foam::steadyOptimisation::lineSearchUpdate(scalarField& direction)
+void Foam::steadyOptimisation::lineSearchUpdate()
 {
+    // Compute direction of update
+    tmp<scalarField> tdirection = optType_->computeDirection();
+    scalarField& direction = tdirection.ref();
+
+    // Grab reference to line search
     autoPtr<lineSearch>& lineSrch = optType_->getLineSearch();
 
     // Store starting point
@@ -121,7 +126,7 @@ void Foam::steadyOptimisation::lineSearchUpdate(scalarField& direction)
         else
         {
             // If maximum number of iteration has been reached, continue
-            if (iter == lineSrch->maxIters()-1)
+            if (iter == lineSrch->maxIters() - 1)
             {
                 Info<< "Line search reached max. number of iterations.\n"
                     << "Proceeding to the next optimisation cycle" << endl;
@@ -141,15 +146,10 @@ void Foam::steadyOptimisation::lineSearchUpdate(scalarField& direction)
 }
 
 
-void Foam::steadyOptimisation::fixedStepUpdate(scalarField& direction)
+void Foam::steadyOptimisation::fixedStepUpdate()
 {
-    // Update based on fixed step
-    optType_->update(direction);
-
-    // If direction has been scaled (say by setting the initial eta), the
-    // old correction has to be updated
-    optType_->updateOldCorrection(direction);
-    optType_->write();
+    // Update design variables
+    optType_->update();
 
     // Solve primal equations
     solvePrimalEquations();
@@ -223,20 +223,15 @@ bool Foam::steadyOptimisation::update()
 
 void Foam::steadyOptimisation::updateDesignVariables()
 {
-    // Compute direction of update
-    tmp<scalarField> tdirection = optType_->computeDirection();
-    scalarField& direction = tdirection.ref();
-    autoPtr<lineSearch>& lineSrch = optType_->getLineSearch();
-
     // Update design variables using either a line-search scheme or
     // a fixed-step update
-    if (lineSrch.valid())
+    if (optType_->getLineSearch().valid())
     {
-        lineSearchUpdate(direction);
+        lineSearchUpdate();
     }
     else
     {
-        fixedStepUpdate(direction);
+        fixedStepUpdate();
     }
 
     // Reset adjoint sensitivities in all adjoint solver managers
diff --git a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/steadyOptimisation/steadyOptimisation.H b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/steadyOptimisation/steadyOptimisation.H
index 549b0036878..0fad239c9c9 100644
--- a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/steadyOptimisation/steadyOptimisation.H
+++ b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationManager/steadyOptimisation/steadyOptimisation.H
@@ -5,8 +5,8 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2007-2019 PCOpt/NTUA
-    Copyright (C) 2013-2019 FOSS GP
+    Copyright (C) 2007-2020 PCOpt/NTUA
+    Copyright (C) 2013-2020 FOSS GP
     Copyright (C) 2019 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
@@ -75,10 +75,10 @@ protected:
     // Protected Member Functions
 
         //- Update design variables using a line-search
-        void lineSearchUpdate(scalarField& direction);
+        void lineSearchUpdate();
 
         //- Update design variables using a fixed step
-        void fixedStepUpdate(scalarField& direction);
+        void fixedStepUpdate();
 
 
 public:
diff --git a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/optimisationType/optimisationTypeIncompressible.C b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/optimisationType/optimisationTypeIncompressible.C
index 0cfd8af5de4..06b2de5394b 100644
--- a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/optimisationType/optimisationTypeIncompressible.C
+++ b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/optimisationType/optimisationTypeIncompressible.C
@@ -5,8 +5,8 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2007-2019 PCOpt/NTUA
-    Copyright (C) 2013-2019 FOSS GP
+    Copyright (C) 2007-2020 PCOpt/NTUA
+    Copyright (C) 2013-2020 FOSS GP
     Copyright (C) 2019 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
@@ -157,14 +157,79 @@ autoPtr<optimisationType> optimisationType::New
 
 // * * * * * * * * * * * * * * *  Member Functions   * * * * * * * * * * * * //
 
+void optimisationType::update()
+{
+    // Compute update of the design variables
+    tmp<scalarField> tcorrection(computeDirection());
+    scalarField& correction = tcorrection.ref();
+
+    // Update design variables given the correction
+    update(correction);
+
+    // If direction has been scaled (say by setting the initial eta), the
+    // old correction has to be updated
+    updateOldCorrection(correction);
+    write();
+}
+
+
+void optimisationType::update(scalarField& direction)
+{
+    // Compute eta if needed
+    computeEta(direction);
+
+    // Multiply with line search step, if necessary
+    scalarField correction(direction);
+    if (lineSearch_.valid())
+    {
+        correction *= lineSearch_->step();
+    }
+
+    // Update the design variables
+    updateDesignVariables(correction);
+}
+
+
 tmp<scalarField> optimisationType::computeDirection()
 {
-    // Sum contributions
+    // Sum contributions for sensitivities and objective/constraint values
     scalarField objectiveSens;
     PtrList<scalarField> constraintSens;
     scalar objectiveValue(Zero);
     scalarField constraintValues;
 
+    updateGradientsAndValues
+    (
+        objectiveSens,
+        constraintSens,
+        objectiveValue,
+        constraintValues
+    );
+
+    // Based on the sensitivities, return design variables correction
+    updateMethod_->setObjectiveDeriv(objectiveSens);
+    updateMethod_->setConstraintDeriv(constraintSens);
+    updateMethod_->setObjectiveValue(objectiveValue);
+    updateMethod_->setConstraintValues(constraintValues);
+    tmp<scalarField> tcorrection
+    (
+        new scalarField(objectiveSens.size(), Zero)
+    );
+    scalarField& correction = tcorrection.ref();
+    correction = updateMethod_->returnCorrection();
+
+    return tcorrection;
+}
+
+
+void optimisationType::updateGradientsAndValues
+(
+    scalarField& objectiveSens,
+    PtrList<scalarField>& constraintSens,
+    scalar& objectiveValue,
+    scalarField& constraintValues
+)
+{
     for (adjointSolverManager& adjSolvManager : adjointSolvManagers_)
     {
         const scalar opWeight = adjSolvManager.operatingPointWeight();
@@ -212,20 +277,6 @@ tmp<scalarField> optimisationType::computeDirection()
         }
         constraintValues += opWeight*cValues();
     }
-
-    // Based on the sensitivities, return design variables correction
-    updateMethod_->setObjectiveDeriv(objectiveSens);
-    updateMethod_->setConstraintDeriv(constraintSens);
-    updateMethod_->setObjectiveValue(objectiveValue);
-    updateMethod_->setConstraintValues(constraintValues);
-    tmp<scalarField> tcorrection
-    (
-        new scalarField(objectiveSens.size(), Zero)
-    );
-    scalarField& correction = tcorrection.ref();
-    correction = updateMethod_->returnCorrection();
-
-    return tcorrection;
 }
 
 
diff --git a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/optimisationType/optimisationTypeIncompressible.H b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/optimisationType/optimisationTypeIncompressible.H
index 8dd392c0e7e..1c1a75619fb 100644
--- a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/optimisationType/optimisationTypeIncompressible.H
+++ b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/optimisationType/optimisationTypeIncompressible.H
@@ -5,8 +5,8 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2007-2019 PCOpt/NTUA
-    Copyright (C) 2013-2019 FOSS GP
+    Copyright (C) 2007-2020 PCOpt/NTUA
+    Copyright (C) 2013-2020 FOSS GP
     Copyright (C) 2019 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
@@ -69,7 +69,11 @@ protected:
         autoPtr<volScalarField> sourcePtr_;
         autoPtr<lineSearch> lineSearch_;
 
-        virtual void computeEta(scalarField& correction)=0;
+        //- Update the design variables given their correction
+        virtual void updateDesignVariables(scalarField& correction) = 0;
+
+        //- Compute eta if not set in the first step
+        virtual void computeEta(scalarField& correction) = 0;
 
 
 private:
@@ -132,10 +136,10 @@ public:
         virtual ~optimisationType() = default;
 
         //- Update design variables
-        virtual void update() = 0;
+        virtual void update();
 
         //- Update design variables based on a given correction
-        virtual void update(scalarField& correction) = 0;
+        virtual void update(scalarField& correction);
 
         //- Store design variables, as the starting point for line search
         virtual void storeDesignVariables() = 0;
@@ -146,6 +150,15 @@ public:
         //- Compute update direction
         virtual tmp<scalarField> computeDirection();
 
+        //- Compute cumulative objective and constraint gradients
+        virtual void updateGradientsAndValues
+        (
+            scalarField& objectiveSens,
+            PtrList<scalarField>& constraintSens,
+            scalar& objectiveValue,
+            scalarField& constraintValues
+        );
+
         //- Compute the merit function of the optimisation problem.
         //  Could be different than the objective function in case of
         //  constraint optimisation
diff --git a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/shapeOptimisation/shapeOptimisationIncompressible.C b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/shapeOptimisation/shapeOptimisationIncompressible.C
index d183c40af4b..6af773aca9c 100644
--- a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/shapeOptimisation/shapeOptimisationIncompressible.C
+++ b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/shapeOptimisation/shapeOptimisationIncompressible.C
@@ -5,8 +5,8 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2007-2019 PCOpt/NTUA
-    Copyright (C) 2013-2019 FOSS GP
+    Copyright (C) 2007-2020 PCOpt/NTUA
+    Copyright (C) 2013-2020 FOSS GP
     Copyright (C) 2019-2020 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
@@ -49,6 +49,39 @@ addToRunTimeSelectionTable
 
 // * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
 
+void shapeOptimisation::updateDesignVariables(scalarField& correction)
+{
+    // Communicate the movement to optMeshMovement
+    optMeshMovement_->setCorrection(correction);
+
+    if (updateGeometry_)
+    {
+        // Update the mesh
+        optMeshMovement_->moveMesh();
+
+        if (writeEachMesh_)
+        {
+            Info<< "  Writing new mesh points " << endl;
+            pointIOField points
+            (
+                IOobject
+                (
+                   "points",
+                    mesh_.pointsInstance(),
+                    mesh_.meshSubDir,
+                    mesh_,
+                    IOobject::NO_READ,
+                    IOobject::NO_WRITE,
+                    false
+                ),
+                mesh_.points()
+            );
+            points.write();
+        }
+    }
+}
+
+
 void shapeOptimisation::computeEta
 (
     scalarField& correction
@@ -146,144 +179,6 @@ shapeOptimisation::shapeOptimisation
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-void shapeOptimisation::update()
-{
-    Info<< nl << "Moving Mesh..." << nl << endl;
-
-    // Sum contributions
-    scalarField objectiveSens(0);
-    PtrList<scalarField> constraintSens(0);
-    scalar objectiveValue(Zero);
-    scalarField constraintValues(0);
-    forAll(adjointSolvManagers_, amI)
-    {
-        adjointSolverManager& adjSolvManager(adjointSolvManagers_[amI]);
-        const scalar opWeight(adjSolvManager.operatingPointWeight());
-
-        // Allocate objective sens size if necessary
-        tmp<scalarField> tadjointSolverManagerSens =
-            adjSolvManager.aggregateSensitivities();
-
-        if (objectiveSens.empty())
-        {
-            objectiveSens.setSize(tadjointSolverManagerSens().size(), Zero);
-        }
-        objectiveSens += opWeight*tadjointSolverManagerSens();
-        objectiveValue += opWeight*adjSolvManager.objectiveValue();
-
-        // Allocate constraint sens size if necessary
-        PtrList<scalarField> adjointSolverManagerConstSens
-            = adjSolvManager.constraintSensitivities();
-        tmp<scalarField> cValues = adjSolvManager.constraintValues();
-        if (constraintSens.empty())
-        {
-            constraintSens.setSize(adjointSolverManagerConstSens.size());
-            forAll(constraintSens, cI)
-            {
-                constraintSens.set
-                (
-                    cI,
-                    new scalarField
-                    (
-                        adjointSolverManagerConstSens[cI].size(),
-                        Zero
-                    )
-                );
-                constraintValues.setSize(cValues().size());
-                constraintValues = Zero;
-            }
-        }
-
-        forAll(constraintSens, cI)
-        {
-            constraintSens[cI] += opWeight*adjointSolverManagerConstSens[cI];
-        }
-        constraintValues += opWeight*cValues();
-    }
-
-    // Based on the sensitivities, return design variables correction
-    updateMethod_->setObjectiveDeriv(objectiveSens);
-    updateMethod_->setConstraintDeriv(constraintSens);
-    updateMethod_->setObjectiveValue(objectiveValue);
-    updateMethod_->setConstraintValues(constraintValues);
-    scalarField& correction = updateMethod_->returnCorrection();
-
-    // Computed  eta if needed
-    computeEta(correction);
-    updateMethod_->writeCorrection();
-
-    // Communicate the movement to optMeshMovement
-    optMeshMovement_->setCorrection(correction);
-    if (updateGeometry_)
-    {
-        optMeshMovement_->moveMesh();
-
-        if (writeEachMesh_)
-        {
-            Info<< "  Writing new mesh points " << endl;
-            pointIOField points
-            (
-                IOobject
-                (
-                    "points",
-                     mesh_.pointsInstance(),
-                     mesh_.meshSubDir,
-                     mesh_,
-                     IOobject::NO_READ,
-                     IOobject::NO_WRITE,
-                     false
-                ),
-                mesh_.points()
-            );
-            points.write();
-        }
-    }
-}
-
-
-void shapeOptimisation::update(scalarField& direction)
-{
-    // Computed eta if needed
-    computeEta(direction);
-
-    // Multiply with line search step, if necessary
-    scalarField correction = direction;
-    if (lineSearch_.valid())
-    {
-        correction *= lineSearch_->step();
-    }
-
-    // Communicate the movement to optMeshMovement
-    optMeshMovement_->setCorrection(correction);
-
-    if (updateGeometry_)
-    {
-        // Update the mesh
-        optMeshMovement_->moveMesh();
-
-        if (writeEachMesh_)
-        {
-            Info<< "  Writing new mesh points " << endl;
-            pointIOField points
-            (
-                IOobject
-                (
-                   "points",
-                    mesh_.pointsInstance(),
-                    mesh_.meshSubDir,
-                    mesh_,
-                    IOobject::NO_READ,
-                    IOobject::NO_WRITE,
-                    false
-                ),
-                mesh_.points()
-            );
-            points.write();
-        }
-    }
-}
-
-
 void shapeOptimisation::storeDesignVariables()
 {
     optMeshMovement_->storeDesignVariables();
diff --git a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/shapeOptimisation/shapeOptimisationIncompressible.H b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/shapeOptimisation/shapeOptimisationIncompressible.H
index 5152941fb14..421feaba671 100644
--- a/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/shapeOptimisation/shapeOptimisationIncompressible.H
+++ b/src/optimisation/adjointOptimisation/adjoint/optimisation/optimisationType/incompressible/shapeOptimisation/shapeOptimisationIncompressible.H
@@ -5,8 +5,8 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2007-2019 PCOpt/NTUA
-    Copyright (C) 2013-2019 FOSS GP
+    Copyright (C) 2007-2020 PCOpt/NTUA
+    Copyright (C) 2013-2020 FOSS GP
     Copyright (C) 2019 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
@@ -76,6 +76,10 @@ protected:
 
     // Protected Member Functions
 
+        //- Update the design variables given their correction
+        virtual void updateDesignVariables(scalarField& correction);
+
+        //- Compute eta if not set in the first step
         virtual void computeEta(scalarField& correction);
 
 
@@ -113,20 +117,14 @@ public:
 
     // Member Functions
 
-       //- Master function. Calls all the others
-       void update();
-
-       //- Update design variables based on a given correction
-       virtual void update(scalarField& correction);
-
-       //- Store design variables, as the starting point for line search
-       virtual void storeDesignVariables();
+        //- Store design variables, as the starting point for line search
+        virtual void storeDesignVariables();
 
-       //- Store design variables, as the starting point for line search
-       virtual void resetDesignVariables();
+        //- Store design variables, as the starting point for line search
+        virtual void resetDesignVariables();
 
-       //- Write useful quantities to files
-       virtual void write();
+        //- Write useful quantities to files
+        virtual void write();
 };
 
 
-- 
GitLab