Merge branch 'feature_turbmod_variants' of https://github.com/su2code/SU2 into feature_turbmod_variants

Author: pcarruscag

Common/include/option_structure.hpp

@@ -2254,6 +2254,8 @@ struct StreamwisePeriodicValues {
   su2double Streamwise_Periodic_MassFlow;           /*!< \brief Value of current massflow [kg/s] which results in a delta p and therefore an artificial body force vector. */
   su2double Streamwise_Periodic_IntegratedHeatFlow; /*!< \brief Value of of the net sum of heatflow [W] into the domain. */
   su2double Streamwise_Periodic_InletTemperature;   /*!< \brief Area avg static Temp [K] at the periodic inlet. Used for adaptive outlet heatsink. */
+  su2double Streamwise_Periodic_BoundaryArea;       /*!< \brief Global Surface area of the streamwise periodic interface. */
+  su2double Streamwise_Periodic_AvgDensity;         /*!< \brief Area avg density on the periodic interface. */
 };
 
 /*!

Common/src/CConfig.cpp

@@ -4803,8 +4803,6 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i
       SU2_MPI::Error("A MARKER_PERIODIC pair has to be set with KIND_STREAMWISE_PERIODIC != NONE.", CURRENT_FUNCTION);
     if (Energy_Equation && Streamwise_Periodic_Temperature && nMarker_Isothermal != 0)
       SU2_MPI::Error("No MARKER_ISOTHERMAL marker allowed with STREAMWISE_PERIODIC_TEMPERATURE= YES, only MARKER_HEATFLUX & MARKER_SYM.", CURRENT_FUNCTION);
-    if (DiscreteAdjoint && Kind_Streamwise_Periodic == ENUM_STREAMWISE_PERIODIC::MASSFLOW)
-      SU2_MPI::Error("Discrete Adjoint currently not validated for prescribed MASSFLOW.", CURRENT_FUNCTION);
     if (Ref_Inc_NonDim != DIMENSIONAL)
       SU2_MPI::Error("Streamwise Periodicity only works with \"INC_NONDIM= DIMENSIONAL\", the nondimensionalization with source terms doesn;t work in general.", CURRENT_FUNCTION);
     if (Axisymmetric)

SU2_CFD/include/numerics/turbulent/turb_sources.hpp

@@ -35,9 +35,18 @@
  */
 struct CSAVariables {
   /*--- List of constants ---*/
-  const su2double cv1_3 = pow(7.1, 3), k2 = pow(0.41, 2), cb1 = 0.1355, cw2 = 0.3, ct3 = 1.2, ct4 = 0.5,
-                  cw3_6 = pow(2, 6), sigma = 2.0 / 3.0, cb2 = 0.622, cb2_sigma = cb2 / sigma,
-                  cw1 = cb1 / k2 + (1 + cb2) / sigma, cr1 = 0.5;
+  const su2double cv1_3 = pow(7.1, 3);
+  const su2double k2 = pow(0.41, 2);
+  const su2double cb1 = 0.1355;
+  const su2double cw2 = 0.3;
+  const su2double ct3 = 1.2;
+  const su2double ct4 = 0.5;
+  const su2double cw3_6 = pow(2, 6);
+  const su2double sigma = 2.0 / 3.0;
+  const su2double cb2 = 0.622;
+  const su2double cb2_sigma = cb2 / sigma;
+  const su2double cw1 = cb1 / k2 + (1 + cb2) / sigma;
+  const su2double cr1 = 0.5;
 
   /*--- List of auxiliary functions ---*/
   su2double ft2, d_ft2, r, d_r, g, d_g, glim, fw, d_fw, Ji, d_Ji, S, Shat, d_Shat, fv1, d_fv1, fv2, d_fv2;
@@ -176,7 +185,6 @@ class CSourceBase_TurbSA : public CNumerics {
         const su2double re_v = density * var.dist_i_2 / laminar_viscosity * var.Omega;
         const su2double re_theta = re_v / 2.193;
         const su2double re_theta_t = 803.73 * pow(tu + 0.6067, -1.027);  // MENTER correlation
-        // re_theta_t = 163.0 + exp(6.91-tu); //ABU-GHANNAM & SHAW correlation
 
         const su2double term1 = sqrt(max(re_theta - re_theta_t, 0.0) / (chi_1 * re_theta_t));
         const su2double term2 = sqrt(max((nu_t * chi_2) / nu, 0.0));

SU2_CFD/include/solvers/CIncEulerSolver.hpp

@@ -39,7 +39,7 @@
 class CIncEulerSolver : public CFVMFlowSolverBase<CIncEulerVariable, ENUM_REGIME::INCOMPRESSIBLE> {
 protected:
   vector<CFluidModel*> FluidModel;   /*!< \brief fluid model used in the solver. */
-  StreamwisePeriodicValues SPvals;
+  StreamwisePeriodicValues SPvals, SPvalsUpdated;
 
   /*!
    * \brief Preprocessing actions common to the Euler and NS solvers.
@@ -395,8 +395,30 @@ class CIncEulerSolver : public CFVMFlowSolverBase<CIncEulerVariable, ENUM_REGIME
   inline bool GetHasHybridParallel() const final { return true; }
 
   /*!
-   * \brief Get values for streamwise periodic flow: delta P, m_dot, inlet T, integrated heat.
-   * \return Struct holding 4 su2doubles.
+   * \brief Get values for streamwise periodic flow: delta P, m_dot, inlet T, integrated heat, etc.
+   * \return Struct holding streamwise periodic values.
    */
   StreamwisePeriodicValues GetStreamwisePeriodicValues() const final { return SPvals; }
+
+  /*!
+   * \brief Register In- or Output.
+   * \param[in] input - Boolean whether In- or Output should be registered.
+   * \param[in] config - The particular config.
+   * \returns The number of extra variables.
+   */
+  unsigned long RegisterSolutionExtra(bool input, const CConfig* config) final;
+
+  /*!
+   * \brief Seed the adjoint of the extra solution at the output.
+   * \param[in] adj_sol - Vector containing the adjoint solution to seed.
+   * \param[in] config - The particular config.
+   */
+  void SetAdjoint_SolutionExtra(const su2activevector& adj_sol, const CConfig* config) final;
+
+  /*!
+   * \brief Extract the adjoint of the extra solution at the input.
+   * \param[out] adj_sol - Vector to store the adjoint into.
+   * \param[in] config - The particular config.
+   */
+  void ExtractAdjoint_SolutionExtra(su2activevector& adj_sol, const CConfig* config) final;
 };

SU2_CFD/include/solvers/CSolver.hpp

@@ -3543,6 +3543,28 @@ class CSolver {
    */
   inline virtual void ExtractAdjoint_Solution(CGeometry *geometry, CConfig *config, bool CrossTerm){}
 
+  /*!
+   * \brief Register In- or Output.
+   * \param[in] input - Boolean whether In- or Output should be registered.
+   * \param[in] config - The particular config.
+   * \returns The number of extra variables.
+   */
+  virtual unsigned long RegisterSolutionExtra(bool input, const CConfig* config) { return 0; }
+
+  /*!
+   * \brief Seed the adjoint of the extra solution at the output.
+   * \param[in] adj_sol - Vector containing the adjoint solution to seed.
+   * \param[in] config - The particular config.
+   */
+  virtual void SetAdjoint_SolutionExtra(const su2activevector& adj_sol, const CConfig* config) {}
+
+  /*!
+   * \brief Extract the adjoint of the extra solution at the input.
+   * \param[out] adj_sol - Vector to store the adjoint into.
+   * \param[in] config - The particular config.
+   */
+  virtual void ExtractAdjoint_SolutionExtra(su2activevector& adj_sol, const CConfig* config) {}
+
   /*!
    * \brief  A virtual member.
    * \param[in] geometry - Geometrical definition of the problem.
@@ -4281,8 +4303,8 @@ class CSolver {
   inline virtual bool GetHasHybridParallel() const { return false; }
 
   /*!
-   * \brief Get values for streamwise periodc flow: delta P, m_dot, inlet T, integrated heat.
-   * \return Struct holding 4 su2doubles.
+   * \brief Get values for streamwise periodic flow: delta P, m_dot, inlet T, integrated heat, etc.
+   * \return Struct holding streamwise periodic values.
    */
   virtual StreamwisePeriodicValues GetStreamwisePeriodicValues() const { return StreamwisePeriodicValues(); }
 

SU2_CFD/include/variables/CDiscAdjVariable.hpp

@@ -54,9 +54,10 @@ class CDiscAdjVariable : public CVariable {
   CDiscAdjVariable(const su2double* sol, unsigned long npoint, unsigned long ndim, unsigned long nvar, CConfig *config);
 
   /*!
-   * \brief Destructor of the class.
+   * \brief Allocate extra adjoint variables.
+   * \param[in] nVarExtra - Number of extra variables.
    */
-  ~CDiscAdjVariable() override = default;
+  void AllocateAdjointSolutionExtra(unsigned long nVarExtra);
 
   /*!
    * \brief Set the sensitivity at the node

SU2_CFD/include/variables/CHeatVariable.hpp

@@ -37,8 +37,6 @@
  */
 class CHeatVariable final : public CVariable {
 protected:
-  MatrixType Solution_Direct;  /*!< \brief Direct solution container for use in the adjoint Heat solver. */
-
   CVectorOfMatrix& Gradient_Reconstruction;  /*!< \brief Reference to the gradient of the primitive variables for MUSCL reconstruction for the convective term */
   CVectorOfMatrix Gradient_Aux;              /*!< \brief Auxiliary structure to store a second gradient for reconstruction, if required. */
 

SU2_CFD/include/variables/CVariable.hpp

@@ -94,6 +94,13 @@ class CVariable {
   su2matrix<int> AD_InputIndex;    /*!< \brief Indices of Solution variables in the adjoint vector. */
   su2matrix<int> AD_OutputIndex;   /*!< \brief Indices of Solution variables in the adjoint vector after having been updated. */
 
+  VectorType SolutionExtra; /*!< \brief Stores adjoint solution for extra solution variables.
+                                        Currently only streamwise periodic pressure-drop for massflow prescribed flows. */
+  VectorType ExternalExtra; /*!< \brief External storage for the adjoint value (i.e. for the OF mainly */
+
+  VectorType SolutionExtra_BGS_k; /*!< \brief Intermediate storage, enables cross term extraction as that is also pushed to Solution. */
+
+ protected:
   unsigned long nPoint = 0;  /*!< \brief Number of points in the domain. */
   unsigned long nDim = 0;      /*!< \brief Number of dimension of the problem. */
   unsigned long nVar = 0;        /*!< \brief Number of variables of the problem. */
@@ -406,6 +413,31 @@ class CVariable {
     for(unsigned long iVar = 0; iVar < nVar; iVar++) External(iPoint,iVar) += val_sol[iVar];
   }
 
+  /*!
+   * \brief Store the adjoint solution of the extra adjoint into the external container.
+   */
+  void Set_ExternalExtra_To_SolutionExtra() {
+    assert(SolutionExtra.size() == ExternalExtra.size());
+    for (auto iEntry = 0ul; iEntry < SolutionExtra.size(); iEntry++)
+      ExternalExtra[iEntry] = SolutionExtra[iEntry];
+  }
+
+  /*!
+   * \brief Add the external contribution to the solution for the extra adjoint solutions.
+   */
+  void Add_ExternalExtra_To_SolutionExtra() {
+    assert(SolutionExtra.size() == ExternalExtra.size());
+    for (auto iEntry = 0ul; iEntry < SolutionExtra.size(); iEntry++)
+      SolutionExtra[iEntry] += ExternalExtra[iEntry];
+  }
+
+  /*!
+   * \brief Return the extra adjoint solution.
+   * \return Reference to extra adjoint solution.
+   */
+  inline VectorType& GetSolutionExtra() { return SolutionExtra; }
+  inline const VectorType& GetSolutionExtra() const { return SolutionExtra; }
+
   /*!
    * \brief Update the variables using a conservative format.
    * \param[in] iPoint - Point index.

SU2_CFD/src/drivers/CDiscAdjMultizoneDriver.cpp

@@ -430,7 +430,7 @@ void CDiscAdjMultizoneDriver::Run() {
       /*--- Compute residual from Solution and Solution_BGS_k and update the latter. ---*/
 
       SetResidual_BGS(iZone);
-
+      Set_BGSSolution_k_To_Solution(iZone);
     }
 
     /*--- Set the multizone output. ---*/

SU2_CFD/src/drivers/CMultizoneDriver.cpp

@@ -403,8 +403,10 @@ bool CMultizoneDriver::OuterConvergence(unsigned long OuterIter) {
     auto solvers = solver_container[iZone][INST_0][MESH_0];
 
     for (unsigned short iSol = 0; iSol < MAX_SOLS; iSol++){
-      if (solvers[iSol] != nullptr)
+      if (solvers[iSol] != nullptr) {
         solvers[iSol]->ComputeResidual_Multizone(geometry_container[iZone][INST_0][MESH_0], config_container[iZone]);
+        solvers[iSol]->GetNodes()->Set_BGSSolution_k();
+      }
     }
 
     /*--- Make sure that everything is loaded into the output container. ---*/