Merge remote-tracking branch 'origin/develop' into feature_expliciteuler_scalarsolver

Author: maxaehle

SU2_CFD/include/output/CFlowCompOutput.hpp

@@ -38,11 +38,9 @@ class CVariable;
  */
 class CFlowCompOutput final: public CFlowOutput {
 private:
-
   TURB_MODEL turb_model; //!< Kind of turbulence model
 
 public:
-
   /*!
    * \brief Constructor of the class
    * \param[in] config - Definition of the particular problem.
@@ -55,18 +53,6 @@ class CFlowCompOutput final: public CFlowOutput {
    */
   void LoadHistoryData(CConfig *config, CGeometry *geometry, CSolver **solver) override;
 
-  /*!
-   * \brief Set the values of the volume output fields for a surface point.
-   * \param[in] config - Definition of the particular problem.
-   * \param[in] geometry - Geometrical definition of the problem.
-   * \param[in] solver - The container holding all solution data.
-   * \param[in] iPoint - Index of the point.
-   * \param[in] iMarker - Index of the surface marker.
-   * \param[in] iVertex - Index of the vertex on the marker.
-   */
-  void LoadSurfaceData(CConfig *config, CGeometry *geometry, CSolver **solver,
-                       unsigned long iPoint, unsigned short iMarker, unsigned long iVertex) override;
-
   /*!
    * \brief Set the available volume output fields
    * \param[in] config - Definition of the particular problem.

SU2_CFD/include/output/CFlowIncOutput.hpp

@@ -38,44 +38,25 @@ class CVariable;
  */
 class CFlowIncOutput final: public CFlowOutput {
 private:
-
   TURB_MODEL turb_model;     /*!< \brief The kind of turbulence model*/
   bool heat;                 /*!< \brief Boolean indicating whether have a heat problem*/
   bool weakly_coupled_heat;  /*!< \brief Boolean indicating whether have a weakly coupled heat equation*/
   unsigned short streamwisePeriodic;   /*!< \brief Boolean indicating whether it is a streamwise periodic simulation. */
   bool streamwisePeriodic_temperature; /*!< \brief Boolean indicating streamwise periodic temperature is used. */
 
 public:
-
   /*!
    * \brief Constructor of the class
    * \param[in] config - Definition of the particular problem.
    */
   CFlowIncOutput(CConfig *config, unsigned short nDim);
 
-  /*!
-   * \brief Destructor of the class.
-   */
-  ~CFlowIncOutput(void) override;
-
   /*!
    * \brief Load the history output field values
    * \param[in] config - Definition of the particular problem.
    */
   void LoadHistoryData(CConfig *config, CGeometry *geometry, CSolver **solver) override;
 
-  /*!
-   * \brief Set the values of the volume output fields for a surface point.
-   * \param[in] config - Definition of the particular problem.
-   * \param[in] geometry - Geometrical definition of the problem.
-   * \param[in] solver - The container holding all solution data.
-   * \param[in] iPoint - Index of the point.
-   * \param[in] iMarker - Index of the surface marker.
-   * \param[in] iVertex - Index of the vertex on the marker.
-   */
-  void LoadSurfaceData(CConfig *config, CGeometry *geometry, CSolver **solver,
-                       unsigned long iPoint, unsigned short iMarker, unsigned long iVertex) override;
-
   /*!
    * \brief Set the available volume output fields
    * \param[in] config - Definition of the particular problem.

SU2_CFD/include/output/CFlowOutput.hpp

@@ -40,6 +40,18 @@ class CFlowOutput : public CFVMOutput{
    */
   CFlowOutput(const CConfig *config, unsigned short nDim, bool femOutput);
 
+  /*!
+   * \brief Set the values of the volume output fields for a surface point.
+   * \param[in] config - Definition of the particular problem.
+   * \param[in] geometry - Geometrical definition of the problem.
+   * \param[in] solver - The container holding all solution data.
+   * \param[in] iPoint - Index of the point.
+   * \param[in] iMarker - Index of the surface marker.
+   * \param[in] iVertex - Index of the vertex on the marker.
+   */
+  void LoadSurfaceData(CConfig *config, CGeometry *geometry, CSolver **solver,
+                       unsigned long iPoint, unsigned short iMarker, unsigned long iVertex) override;
+
   /*!
    * \brief Add flow surface output fields
    * \param[in] config - Definition of the particular problem.
@@ -55,6 +67,61 @@ class CFlowOutput : public CFVMOutput{
    */
   void SetAnalyzeSurface(const CSolver *solver, const CGeometry *geometry, CConfig *config, bool output);
 
+  /*!
+   * \brief Add turbulence history fields for the linear solver (FVMComp, FVMInc, FVMNEMO).
+   */
+  void AddHistoryOutputFields_TurbRMS_RES(const CConfig* config);
+
+  /*!
+   * \brief Add turbulence history fields for the linear solver (FVMComp, FVMInc, FVMNEMO).
+   */
+  void AddHistoryOutputFields_TurbMAX_RES(const CConfig* config);
+
+  /*!
+   * \brief Add turbulence history fields for the linear solver (FVMComp, FVMInc, FVMNEMO).
+   */
+  void AddHistoryOutputFields_TurbBGS_RES(const CConfig* config);
+
+  /*!
+   * \brief Add turbulence history fields for the linear solver (FVMComp, FVMInc, FVMNEMO).
+   */
+  void AddHistoryOutputFields_TurbLinsol(const CConfig* config);
+
+  /*!
+   * \brief Set all turbulence history field values.
+   */
+  void LoadHistoryData_Turb(const CConfig* config, const CSolver* const* solver);
+
+  /*!
+   * \brief Add turbulence volume solution fields for a point (FVMComp, FVMInc, FVMNEMO).
+   * \note The order of fields in restart files is fixed. Therefore the split-up.
+   * \param[in] config - Definition of the particular problem.
+   */
+  void SetVolumeOutputFields_TurbSolution(const CConfig* config);
+
+  /*!
+   * \brief Add turbulence volume solution fields for a point (FVMComp, FVMInc, FVMNEMO).
+   * \note The order of fields in restart files is fixed. Therefore the split-up.
+   * \param[in] config - Definition of the particular problem.
+   */
+  void SetVolumeOutputFields_TurbResidual(const CConfig* config);
+
+  /*!
+   * \brief Add turbulence volume limiter fields (and more) for a point (FVMComp, FVMInc, FVMNEMO).
+   * \param[in] config - Definition of the particular problem.
+   */
+  void SetVolumeOutputFields_TurbLimiter(const CConfig* config);
+
+  /*!
+   * \brief Set all turbulence volume field values for a point.
+   * \param[in] config - Definition of the particular problem.
+   * \param[in] solver - The container holding all solution data.
+   * \param[in] geometry - Geometrical definition of the problem.
+   * \param[in] iPoint - Index of the point.
+   */
+  void LoadVolumeData_Turb(const CConfig* config, const CSolver* const* solver, const CGeometry* geometry,
+                           const unsigned long iPoint);
+
   /*!
    * \brief Add aerodynamic coefficients as output fields
    * \param[in] config - Definition of the particular problem.
@@ -179,7 +246,7 @@ class CFlowOutput : public CFVMOutput{
    * \param iPoint
    * \param node_flow
    */
-  void LoadTimeAveragedData(unsigned long iPoint, CVariable *node_flow);
+  void LoadTimeAveragedData(unsigned long iPoint, const CVariable *node_flow);
 
   /*!
    * \brief Write additional output for fixed CL mode.

SU2_CFD/include/output/CNEMOCompOutput.hpp

@@ -37,12 +37,10 @@ class CVariable;
  */
 class CNEMOCompOutput final: public CFlowOutput {
 private:
-
   TURB_MODEL turb_model;   /*!< \brief Kind of turbulence model */
   unsigned short iSpecies, /*!< \brief Species index */
                  nSpecies; /*!< \brief Number of species */
 public:
-
   /*!
    * \brief Constructor of the class
    * \param[in] config - Definition of the particular problem.
@@ -55,18 +53,6 @@ class CNEMOCompOutput final: public CFlowOutput {
    */
   void LoadHistoryData(CConfig *config, CGeometry *geometry, CSolver **solver) override;
 
-  /*!
-   * \brief Set the values of the volume output fields for a surface point.
-   * \param[in] config - Definition of the particular problem.
-   * \param[in] geometry - Geometrical definition of the problem.
-   * \param[in] solver - The container holding all solution data.
-   * \param[in] iPoint - Index of the point.
-   * \param[in] iMarker - Index of the surface marker.
-   * \param[in] iVertex - Index of the vertex on the marker.
-   */
-  void LoadSurfaceData(CConfig *config, CGeometry *geometry, CSolver **solver,
-                       unsigned long iPoint, unsigned short iMarker, unsigned long iVertex) override;
-
   /*!
    * \brief Set the available volume output fields
    * \param[in] config - Definition of the particular problem.

SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp

@@ -972,95 +972,8 @@ class CFVMFlowSolverBase : public CSolver {
 
   /*!
    * \brief Evaluate the vorticity and strain rate magnitude.
-   * \tparam VelocityOffset - Index in the primitive variables where the velocity starts.
    */
-  void ComputeVorticityAndStrainMag(const CConfig& config, unsigned short iMesh,
-                                    const size_t VelocityOffset = 1) {
-
-    auto& StrainMag = nodes->GetStrainMag();
-
-    ompMasterAssignBarrier(StrainMag_Max,0.0, Omega_Max,0.0);
-
-    su2double strainMax = 0.0, omegaMax = 0.0;
-
-    SU2_OMP_FOR_STAT(omp_chunk_size)
-    for (unsigned long iPoint = 0; iPoint < nPoint; ++iPoint) {
-
-      const auto VelocityGradient = nodes->GetGradient_Primitive(iPoint, VelocityOffset);
-      auto Vorticity = nodes->GetVorticity(iPoint);
-
-      /*--- Vorticity ---*/
-
-      Vorticity[0] = 0.0;
-      Vorticity[1] = 0.0;
-      Vorticity[2] = VelocityGradient(1,0)-VelocityGradient(0,1);
-
-      if (nDim == 3) {
-        Vorticity[0] = VelocityGradient(2,1)-VelocityGradient(1,2);
-        Vorticity[1] = -(VelocityGradient(2,0)-VelocityGradient(0,2));
-      }
-
-      /*--- Strain Magnitude ---*/
-
-      AD::StartPreacc();
-      AD::SetPreaccIn(VelocityGradient, nDim, nDim);
-
-      su2double Div = 0.0;
-      for (unsigned long iDim = 0; iDim < nDim; iDim++)
-        Div += VelocityGradient(iDim, iDim);
-      Div /= 3.0;
-
-      StrainMag(iPoint) = 0.0;
-
-      /*--- Add diagonal part ---*/
-
-      for (unsigned long iDim = 0; iDim < nDim; iDim++) {
-        StrainMag(iPoint) += pow(VelocityGradient(iDim, iDim) - Div, 2);
-      }
-      if (nDim == 2) {
-        StrainMag(iPoint) += pow(Div, 2);
-      }
-
-      /*--- Add off diagonals ---*/
-
-      StrainMag(iPoint) += 2.0*pow(0.5*(VelocityGradient(0,1) + VelocityGradient(1,0)), 2);
-
-      if (nDim == 3) {
-        StrainMag(iPoint) += 2.0*pow(0.5*(VelocityGradient(0,2) + VelocityGradient(2,0)), 2);
-        StrainMag(iPoint) += 2.0*pow(0.5*(VelocityGradient(1,2) + VelocityGradient(2,1)), 2);
-      }
-
-      StrainMag(iPoint) = sqrt(2.0*StrainMag(iPoint));
-      AD::SetPreaccOut(StrainMag(iPoint));
-
-      /*--- Max is not differentiable, so we not register them for preacc. ---*/
-      strainMax = max(strainMax, StrainMag(iPoint));
-      omegaMax = max(omegaMax, GeometryToolbox::Norm(3, Vorticity));
-
-      AD::EndPreacc();
-    }
-    END_SU2_OMP_FOR
-
-    if ((iMesh == MESH_0) && (config.GetComm_Level() == COMM_FULL)) {
-      SU2_OMP_CRITICAL {
-        StrainMag_Max = max(StrainMag_Max, strainMax);
-        Omega_Max = max(Omega_Max, omegaMax);
-      }
-      END_SU2_OMP_CRITICAL
-
-      SU2_OMP_BARRIER
-      SU2_OMP_MASTER {
-        su2double MyOmega_Max = Omega_Max;
-        su2double MyStrainMag_Max = StrainMag_Max;
-
-        SU2_MPI::Allreduce(&MyStrainMag_Max, &StrainMag_Max, 1, MPI_DOUBLE, MPI_MAX, SU2_MPI::GetComm());
-        SU2_MPI::Allreduce(&MyOmega_Max, &Omega_Max, 1, MPI_DOUBLE, MPI_MAX, SU2_MPI::GetComm());
-      }
-      END_SU2_OMP_MASTER
-      SU2_OMP_BARRIER
-    }
-
-  }
+  void ComputeVorticityAndStrainMag(const CConfig& config, unsigned short iMesh);
 
   /*!
    * \brief Destructor.

SU2_CFD/include/solvers/CFVMFlowSolverBase.inl

@@ -592,6 +592,94 @@ void CFVMFlowSolverBase<V, R>::ImplicitEuler_Iteration(CGeometry *geometry, CSol
   CompleteImplicitIteration(geometry, nullptr, config);
 }
 
+template <class V, ENUM_REGIME R>
+void CFVMFlowSolverBase<V, R>::ComputeVorticityAndStrainMag(const CConfig& config, unsigned short iMesh) {
+
+  auto& StrainMag = nodes->GetStrainMag();
+
+  ompMasterAssignBarrier(StrainMag_Max,0.0, Omega_Max,0.0);
+
+  su2double strainMax = 0.0, omegaMax = 0.0;
+
+  SU2_OMP_FOR_STAT(omp_chunk_size)
+  for (unsigned long iPoint = 0; iPoint < nPoint; ++iPoint) {
+
+    const auto VelocityGradient = nodes->GetVelocityGradient(iPoint);
+    auto Vorticity = nodes->GetVorticity(iPoint);
+
+    /*--- Vorticity ---*/
+
+    Vorticity[0] = 0.0;
+    Vorticity[1] = 0.0;
+    Vorticity[2] = VelocityGradient(1,0)-VelocityGradient(0,1);
+
+    if (nDim == 3) {
+      Vorticity[0] = VelocityGradient(2,1)-VelocityGradient(1,2);
+      Vorticity[1] = -(VelocityGradient(2,0)-VelocityGradient(0,2));
+    }
+
+    /*--- Strain Magnitude ---*/
+
+    AD::StartPreacc();
+    AD::SetPreaccIn(VelocityGradient, nDim, nDim);
+
+    su2double Div = 0.0;
+    for (unsigned long iDim = 0; iDim < nDim; iDim++)
+      Div += VelocityGradient(iDim, iDim);
+    Div /= 3.0;
+
+    StrainMag(iPoint) = 0.0;
+
+    /*--- Add diagonal part ---*/
+
+    for (unsigned long iDim = 0; iDim < nDim; iDim++) {
+      StrainMag(iPoint) += pow(VelocityGradient(iDim, iDim) - Div, 2);
+    }
+    if (nDim == 2) {
+      StrainMag(iPoint) += pow(Div, 2);
+    }
+
+    /*--- Add off diagonals ---*/
+
+    StrainMag(iPoint) += 2.0*pow(0.5*(VelocityGradient(0,1) + VelocityGradient(1,0)), 2);
+
+    if (nDim == 3) {
+      StrainMag(iPoint) += 2.0*pow(0.5*(VelocityGradient(0,2) + VelocityGradient(2,0)), 2);
+      StrainMag(iPoint) += 2.0*pow(0.5*(VelocityGradient(1,2) + VelocityGradient(2,1)), 2);
+    }
+
+    StrainMag(iPoint) = sqrt(2.0*StrainMag(iPoint));
+    AD::SetPreaccOut(StrainMag(iPoint));
+
+    /*--- Max is not differentiable, so we not register them for preacc. ---*/
+    strainMax = max(strainMax, StrainMag(iPoint));
+    omegaMax = max(omegaMax, GeometryToolbox::Norm(3, Vorticity));
+
+    AD::EndPreacc();
+  }
+  END_SU2_OMP_FOR
+
+  if ((iMesh == MESH_0) && (config.GetComm_Level() == COMM_FULL)) {
+    SU2_OMP_CRITICAL {
+      StrainMag_Max = max(StrainMag_Max, strainMax);
+      Omega_Max = max(Omega_Max, omegaMax);
+    }
+    END_SU2_OMP_CRITICAL
+
+    SU2_OMP_BARRIER
+    SU2_OMP_MASTER {
+      su2double MyOmega_Max = Omega_Max;
+      su2double MyStrainMag_Max = StrainMag_Max;
+
+      SU2_MPI::Allreduce(&MyStrainMag_Max, &StrainMag_Max, 1, MPI_DOUBLE, MPI_MAX, SU2_MPI::GetComm());
+      SU2_MPI::Allreduce(&MyOmega_Max, &Omega_Max, 1, MPI_DOUBLE, MPI_MAX, SU2_MPI::GetComm());
+    }
+    END_SU2_OMP_MASTER
+    SU2_OMP_BARRIER
+  }
+
+}
+
 template <class V, ENUM_REGIME R>
 void CFVMFlowSolverBase<V, R>::SetInletAtVertex(const su2double* val_inlet, unsigned short iMarker,
                                                 unsigned long iVertex) {

SU2_CFD/include/variables/CFlowVariable.hpp

@@ -236,6 +236,7 @@ class CFlowVariable : public CVariable {
    * \return Vorticity array.
    */
   inline su2double* GetVorticity(unsigned long iPoint) final { return Vorticity[iPoint]; }
+  inline const su2double* GetVorticity(unsigned long iPoint) const final { return Vorticity[iPoint]; }
 
   /*!
    * \brief Get the magnitude of rate of strain.