address review comments

Author: pcarruscag

…include/toolboxes/CQuasiNewtonDriver.hpp …oolboxes/CQuasiNewtonInvLeastSquares.hppCommon/include/toolboxes/CQuasiNewtonDriver.hpp renamed to Common/include/toolboxes/CQuasiNewtonInvLeastSquares.hpp Common/include/toolboxes/CQuasiNewtonDriver.hpp renamed to Common/include/toolboxes/CQuasiNewtonInvLeastSquares.hpp

@@ -1,5 +1,5 @@
 /*!
- * \file CQuasiNewtonDriver.hpp
+ * \file CQuasiNewtonInvLeastSquares.hpp
  * \brief Implements a method to accelerate and stabilize the convergence
  * of fixed point iterations, the history of past iterates is used to compute
  * a least squares approximation of the inverse of the Jacobian, which is then
@@ -46,7 +46,7 @@
  * as the new input of the FP, run the FP, etc.
  */
 template<class Scalar_t>
-class CQuasiNewtonDriver {
+class CQuasiNewtonInvLeastSquares {
 public:
   using Scalar = Scalar_t;
   using Index = typename su2matrix<Scalar>::Index;
@@ -57,7 +57,7 @@ class CQuasiNewtonDriver {
 
   enum: size_t {BLOCK_SIZE = 1024};     /*!< \brief Loop tiling parameter. */
   std::vector<su2matrix<Scalar> > X, R; /*!< \brief Input and residual history of the FP. */
-  su2matrix<Scalar> corr;               /*!< \brief Correction to the natural FP result. */
+  su2matrix<Scalar> work;               /*!< \brief Work matrix (FP result, correction, and approx solution). */
   su2vector<Scalar> mat, rhs, sol;      /*!< \brief Matrix, rhs, and solution of the normal equations. */
   Index iSample = 0;                    /*!< \brief Current sample index. */
   Index nPtDomain = 0;                  /*!< \brief Local size of the history, considered in dot products. */
@@ -73,7 +73,7 @@ class CQuasiNewtonDriver {
 
   void computeNormalEquations() {
     /*--- Size for the dot products. ---*/
-    const auto end = std::min<Index>(nPtDomain,corr.rows())*corr.cols();
+    const auto end = std::min<Index>(nPtDomain,work.rows())*work.cols();
 
     mat = Scalar(0); rhs = Scalar(0);
 
@@ -90,7 +90,7 @@ class CQuasiNewtonDriver {
     }
 
     /*--- MPI reduction of the dot products. ---*/
-    if (nPtDomain < corr.rows()) {
+    if (nPtDomain < work.rows()) {
       const auto type = (sizeof(Scalar) < sizeof(double))? MPI_FLOAT : MPI_DOUBLE;
 
       su2vector<Scalar> tmp(mat.size());
@@ -125,6 +125,7 @@ class CQuasiNewtonDriver {
         for (Index k = 0; k < blkSize; ++k) {
           sum += (ri1[k]-ri0[k]) * (rj1[k]-rj0[k]);
         }
+        /*--- 1D index of (i,j) in lower triangular storage. ---*/
         const auto iCoeff = i*(i+1)/2 + j;
         mat(iCoeff) += sum;
       }
@@ -143,10 +144,10 @@ class CQuasiNewtonDriver {
 
 public:
   /*! \brief Default construction without allocation. */
-  CQuasiNewtonDriver() = default;
+  CQuasiNewtonInvLeastSquares() = default;
 
   /*! \brief Construction with allocation, see "resize". */
-  CQuasiNewtonDriver(Index nsample, Index npt, Index nvar, Index nptdomain = 0) {
+  CQuasiNewtonInvLeastSquares(Index nsample, Index npt, Index nvar, Index nptdomain = 0) {
     resize(nsample, npt, nvar, nptdomain);
   }
 
@@ -162,7 +163,7 @@ class CQuasiNewtonDriver {
       SU2_MPI::Error("Invalid quasi-Newton parameters", CURRENT_FUNCTION);
     iSample = 0;
     nPtDomain = nptdomain? nptdomain : npt;
-    corr.resize(npt,nvar);
+    work.resize(npt,nvar);
     X.clear();
     R.clear();
     for (Index i = 0; i < nsample; ++i) {
@@ -186,10 +187,10 @@ class CQuasiNewtonDriver {
    * \brief Access the current fixed-point result.
    * \note Use these to STORE the result of running the FP.
    */
-  su2matrix<Scalar>& FPresult() { return corr; }
-  const su2matrix<Scalar>& FPresult() const { return corr; }
-  Scalar& FPresult(Index iPt, Index iVar) { return corr(iPt,iVar); }
-  const Scalar& FPresult(Index iPt, Index iVar) const { return corr(iPt,iVar); }
+  su2matrix<Scalar>& FPresult() { return work; }
+  const su2matrix<Scalar>& FPresult() const { return work; }
+  Scalar& FPresult(Index iPt, Index iVar) { return work(iPt,iVar); }
+  const Scalar& FPresult(Index iPt, Index iVar) const { return work(iPt,iVar); }
 
   /*!
    * \brief Access the current solution approximation.
@@ -207,9 +208,9 @@ class CQuasiNewtonDriver {
   const su2matrix<Scalar>& compute() {
     /*--- Compute FP residual, clear correction. ---*/
     SU2_OMP_SIMD
-    for (Index i = 0; i < corr.size(); ++i) {
-      R[iSample].data()[i] = corr.data()[i] - X[iSample].data()[i];
-      corr.data()[i] = Scalar(0);
+    for (Index i = 0; i < work.size(); ++i) {
+      R[iSample].data()[i] = work.data()[i] - X[iSample].data()[i];
+      work.data()[i] = Scalar(0);
     }
 
     if (iSample > 0) {
@@ -229,9 +230,9 @@ class CQuasiNewtonDriver {
         const auto x0 = X[k].data();
         const auto r0 = R[k].data();
         SU2_OMP_SIMD
-        for (Index i = 0; i < corr.size(); ++i) {
+        for (Index i = 0; i < work.size(); ++i) {
           Scalar dy = r1[i]-r0[i] + x1[i]-x0[i];
-          corr.data()[i] += sol(k) * dy;
+          work.data()[i] += sol(k) * dy;
         }
       }
     }
@@ -244,9 +245,9 @@ class CQuasiNewtonDriver {
 
     /*--- Set new solution. ---*/
     SU2_OMP_SIMD
-    for (Index i = 0; i < corr.size(); ++i)
-      corr.data()[i] += R[iSample].data()[i] + X[iSample].data()[i];
-    std::swap(X[++iSample], corr);
+    for (Index i = 0; i < work.size(); ++i)
+      work.data()[i] += R[iSample].data()[i] + X[iSample].data()[i];
+    std::swap(X[++iSample], work);
 
     return solution();
   }

SU2_CFD/src/drivers/CDiscAdjMultizoneDriver.cpp

@@ -31,7 +31,7 @@
 #include "../../include/output/COutputLegacy.hpp"
 #include "../../include/output/COutput.hpp"
 #include "../../include/iteration/CIterationFactory.hpp"
-#include "../../../Common/include/toolboxes/CQuasiNewtonDriver.hpp"
+#include "../../../Common/include/toolboxes/CQuasiNewtonInvLeastSquares.hpp"
 
 CDiscAdjMultizoneDriver::CDiscAdjMultizoneDriver(char* confFile,
                                                  unsigned short val_nZone,
@@ -164,7 +164,7 @@ void CDiscAdjMultizoneDriver::Run() {
 
   unsigned long wrt_sol_freq = 9999;
   unsigned long nOuterIter = driver_config->GetnOuter_Iter();
-  vector<CQuasiNewtonDriver<passivedouble> > QNDriver(nZone);
+  vector<CQuasiNewtonInvLeastSquares<passivedouble> > fixPtCorrector(nZone);
 
   for (iZone = 0; iZone < nZone; iZone++) {
 
@@ -182,10 +182,10 @@ void CDiscAdjMultizoneDriver::Run() {
     /*--- Prepare quasi-Newton drivers. ---*/
 
     if (config_container[iZone]->GetnQuasiNewtonSamples() > 1) {
-      QNDriver[iZone].resize(config_container[iZone]->GetnQuasiNewtonSamples(),
-                             geometry_container[iZone][INST_0][MESH_0]->GetnPoint(),
-                             GetTotalNumberOfVariables(iZone, true),
-                             geometry_container[iZone][INST_0][MESH_0]->GetnPointDomain());
+      fixPtCorrector[iZone].resize(config_container[iZone]->GetnQuasiNewtonSamples(),
+                                   geometry_container[iZone][INST_0][MESH_0]->GetnPoint(),
+                                   GetTotalNumberOfVariables(iZone, true),
+                                   geometry_container[iZone][INST_0][MESH_0]->GetnPointDomain());
     }
   }
 
@@ -264,9 +264,9 @@ void CDiscAdjMultizoneDriver::Run() {
 
       /*--- Reset QN driver for new inner iterations. ---*/
 
-      if (QNDriver[iZone].size()) {
-        QNDriver[iZone].reset();
-        if(restart && (iOuterIter==1)) GetAllSolutions(iZone, true, QNDriver[iZone]);
+      if (fixPtCorrector[iZone].size()) {
+        fixPtCorrector[iZone].reset();
+        if(restart && (iOuterIter==1)) GetAllSolutions(iZone, true, fixPtCorrector[iZone]);
       }
 
       for (unsigned long iInnerIter = 0; iInnerIter < nInnerIter[iZone]; iInnerIter++) {
@@ -299,10 +299,10 @@ void CDiscAdjMultizoneDriver::Run() {
 
         /*--- Use QN driver to improve the solution. ---*/
 
-        if (QNDriver[iZone].size()) {
-          GetAllSolutions(iZone, true, QNDriver[iZone].FPresult());
-          QNDriver[iZone].compute();
-          if(iInnerIter) SetAllSolutions(iZone, true, QNDriver[iZone]);
+        if (fixPtCorrector[iZone].size()) {
+          GetAllSolutions(iZone, true, fixPtCorrector[iZone].FPresult());
+          fixPtCorrector[iZone].compute();
+          if(iInnerIter) SetAllSolutions(iZone, true, fixPtCorrector[iZone]);
         }
 
         /*--- This is done explicitly here for multizone cases, only in inner iterations and not when

SU2_CFD/src/drivers/CDiscAdjSinglezoneDriver.cpp

@@ -32,7 +32,7 @@
 #include "../../include/output/COutput.hpp"
 #include "../../include/iteration/CIterationFactory.hpp"
 #include "../../include/iteration/CTurboIteration.hpp"
-#include "../../../Common/include/toolboxes/CQuasiNewtonDriver.hpp"
+#include "../../../Common/include/toolboxes/CQuasiNewtonInvLeastSquares.hpp"
 
 CDiscAdjSinglezoneDriver::CDiscAdjSinglezoneDriver(char* confFile,
                                                    unsigned short val_nZone,
@@ -152,14 +152,14 @@ void CDiscAdjSinglezoneDriver::Run() {
 
   const bool steady = !config->GetTime_Domain();
 
-  CQuasiNewtonDriver<passivedouble> QNDriver;
+  CQuasiNewtonInvLeastSquares<passivedouble> fixPtCorrector;
   if (config->GetnQuasiNewtonSamples() > 1) {
-    QNDriver.resize(config->GetnQuasiNewtonSamples(),
-                    geometry_container[ZONE_0][INST_0][MESH_0]->GetnPoint(),
-                    GetTotalNumberOfVariables(ZONE_0,true),
-                    geometry_container[ZONE_0][INST_0][MESH_0]->GetnPointDomain());
+    fixPtCorrector.resize(config->GetnQuasiNewtonSamples(),
+                          geometry_container[ZONE_0][INST_0][MESH_0]->GetnPoint(),
+                          GetTotalNumberOfVariables(ZONE_0,true),
+                          geometry_container[ZONE_0][INST_0][MESH_0]->GetnPointDomain());
 
-    if (TimeIter != 0) GetAllSolutions(ZONE_0, true, QNDriver);
+    if (TimeIter != 0) GetAllSolutions(ZONE_0, true, fixPtCorrector);
   }
 
   for (auto Adjoint_Iter = 0ul; Adjoint_Iter < nAdjoint_Iter; Adjoint_Iter++) {
@@ -207,9 +207,9 @@ void CDiscAdjSinglezoneDriver::Run() {
 
     /*--- Correct the solution with the quasi-Newton approach. ---*/
 
-    if (QNDriver.size()) {
-      GetAllSolutions(ZONE_0, true, QNDriver.FPresult());
-      SetAllSolutions(ZONE_0, true, QNDriver.compute());
+    if (fixPtCorrector.size()) {
+      GetAllSolutions(ZONE_0, true, fixPtCorrector.FPresult());
+      SetAllSolutions(ZONE_0, true, fixPtCorrector.compute());
     }
 
   }