(refactor) Improved output and changed xtr output to xoc

blast/src/blBody.cpp

@@ -33,7 +33,7 @@ double Body::getAvgxtr(size_t iregion) const
     double sum = 0.0;
     for (auto &section : sections)
     {
-        sum += section->regions[iregion]->getxtr();
+        sum += section->regions[iregion]->getxoctr();
     }
     return sum / static_cast<double>(sections.size());
 }

blast/src/blBoundaryLayer.h

@@ -76,6 +76,8 @@ public:
     size_t getnNodes() const { return nodes.size(); }
     size_t getnElms() const { return nodes.size() - 1; }
     double getxtr() const { return xtr; }
+    double getxoctr() const  { return nodes[transitionNode]->xoc; }
+    double getxitr() const { return nodes[transitionNode]->xi; }
     double getxtrF() const { return xtrF; }
     double getMaxMach() const { return *std::max_element(Me.begin(), Me.end()); }
     size_t getnVar() const { return nVar; }

blast/src/blDriver.cpp

@@ -68,13 +68,17 @@ int Driver::run()
     int pointExitCode;      // Output of pseudo time integration (0: converged; -1:
                             // unsuccessful, failed; >0: unsuccessful nb iter).
 
-    int numberFailed = 0;
+    int numberFailedPerBody = 0;
     double maxMach = 0.;
 
     solverExitCode = 0;
 
+    if (verbose >= 1)
+        std::cout << "Solving boundary layer for Re " << std::setprecision(1) << Re/1e6 << "e6, " << bodies.size() << " bodies." << std::endl;
+
     for (auto &body : bodies)
     {
+        numberFailedPerBody = 0;
         for (auto &section : body->sections)
         {
             // Loop over the regions (upper, lower and wake).
@@ -112,12 +116,8 @@ int Driver::run()
                     // Impose wake boundary condition.
                     std::vector<std::vector<double>> teConditions;
                     for (auto &iregion : section->regions)
-                    {
                         if (iregion->getName() == "upper" || iregion->getName() == "lower")
-                        {
                             teConditions.push_back(iregion->getTeConditions());
-                        }
-                    }
                     reg->setWakeBC(Re, teConditions);
                     lockTrans = true;
                 }
@@ -130,20 +130,10 @@ int Driver::run()
                     // Initialize solution at point
                     tSolver->initialCondition(inod, reg);
                     // Solve equations
-                    tms["1-Solver"].start();
                     pointExitCode = tSolver->integration(inod, reg);
-                    if (pointExitCode != 0)
-                    {
-                        // reg->printSolution(inod);
-                        // std::cout << "side" << reg->getName() << " inod " << inod << " exitCode " << pointExitCode << std::endl;
-                        solverExitCode = -1;
-                        numberFailed++;
-                        // std::cout << "solverExitCode " << solverExitCode << std::endl;
-                    }
-                    tms["1-Solver"].stop();
+                    if (pointExitCode != 0){solverExitCode = -1; numberFailedPerBody++;}
 
                     // Transition
-                    tms["2-Transition"].start();
                     if (!lockTrans)
                     {
                         // Forced transition
@@ -151,12 +141,6 @@ int Driver::run()
                         {
                             reg->u[inod * reg->getnVar() + 2] = reg->getnCrit();
                             averageTransition(inod, reg, true);
-                            if (verbose > 1)
-                            {
-                                std::cout << std::fixed;
-                                std::cout << std::setprecision(2);
-                                std::cout << reg->getxtr() << " (" << reg->transitionNode << ") ";
-                            }
                             lockTrans = true;
                         }
                         // Free transtion
@@ -170,33 +154,33 @@ int Driver::run()
                             if (reg->u[inod * reg->getnVar() + 2] >= reg->getnCrit())
                             {
                                 averageTransition(inod, reg, false);
-                                if (verbose > 1)
-                                {
-                                    std::cout << std::fixed;
-                                    std::cout << std::setprecision(2);
-                                    std::cout << reg->getxtr() << " (" << reg->transitionNode << ") ";
-                                }
                                 lockTrans = true;
                             }
                         }
                     }
-                    tms["2-Transition"].stop();
                 }
-                tms["3-Blowing"].start();
                 reg->computeBlowingVelocity();
-                tms["3-Blowing"].stop();
                 iRegion++;
             }
-            if (verbose > 1)
-                std::cout << std::endl;
         }
         body->computeDrag();
         Cdf = body->Cdf;
         Cdt = body->Cdt;
         Cdp = body->Cdp;
+        if (verbose >= 2)
+            std::cout << body->getName() << "; " << numberFailedPerBody << " points failed." << std::endl;
+    }
+    if (verbose >= 1)
+    {
+        std::cout << "Maximum Mach number: " << maxMach << std::endl;
+        if (solverExitCode == 0)
+            std::cout << "\x1b[1;32m" << "Boundary layer solver converged." << "\x1b[0m" << std::endl;
+        else if (solverExitCode == -1)
+            std::cout << "\x1b[1;31m" << "Boundary layer solver diverged!" << "\x1b[0m" << std::endl;
+        else if (solverExitCode != 0)
+            std::cout << "\x1b[1;33m" << "Boundary layer solver not fully converged." << "\x1b[0m" << std::endl;
+        std::cout << std::endl;
     }
-    if (numberFailed > 0)
-        std::cout << numberFailed << " points failed" << std::endl;
     return solverExitCode; // exit code
 }
 

blast/src/blSolver.cpp

@@ -134,11 +134,12 @@ int Solver::integration(size_t inod, std::shared_ptr<BoundaryLayer> &bl)
         innerIters++;
     }
 
-    if (std::isnan(normSysRes) || normSysRes / normSysRes0 > 1e-3)
+    if (std::isnan(normSysRes))
     {
         resetSolution(inod, bl, sln0);
         return -1;
     }
+    resetSolution(inod, bl, sln0);
     return innerIters;
 }
 

blast/tests/dart/t_naca0012_2D.py

@@ -110,14 +110,14 @@ def main():
     coupler, isol, vsol = vutils.initBlast(icfg, vcfg)
     tms['pre'].stop()
 
-    # Save in vtk format
-    vutils.save(vsol.bodies[0].sections)
-
     print(ccolors.ANSI_BLUE + 'PySolving...' + ccolors.ANSI_RESET)
     tms['solver'].start()
     aeroCoeffs = coupler.run()
     tms['solver'].stop()
 
+    # Save in vtk format
+    vutils.save(vsol.bodies[0].sections)
+
     # Display results.
     print(ccolors.ANSI_BLUE + 'PyRes...' + ccolors.ANSI_RESET)
     print('      Re        M    alpha       Cl       Cd      Cdp      Cdf       Cm')