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

Author: jtlau

Common/include/CConfig.hpp

@@ -1119,8 +1119,8 @@ class CConfig {
 
   unsigned long edgeColorGroupSize; /*!< \brief Size of the edge groups colored for OpenMP parallelization of edge loops. */
 
-  unsigned short Kind_InletInterpolationFunction; /*!brief type of spanwise interpolation function to use for the inlet face. */
-  unsigned short Kind_Inlet_InterpolationType;    /*!brief type of spanwise interpolation data to use for the inlet face. */
+  INLET_SPANWISE_INTERP Kind_InletInterpolationFunction; /*!brief type of spanwise interpolation function to use for the inlet face. */
+  INLET_INTERP_TYPE Kind_Inlet_InterpolationType;    /*!brief type of spanwise interpolation data to use for the inlet face. */
   bool PrintInlet_InterpolatedData;               /*!brief option for printing the interpolated data file. */
 
   /*--- libROM configure options ---*/
@@ -8742,12 +8742,12 @@ class CConfig {
   /*!
    * \brief Get the kind of inlet face interpolation function to use.
    */
-  inline unsigned short GetKindInletInterpolationFunction(void) const { return Kind_InletInterpolationFunction; }
+  inline INLET_SPANWISE_INTERP GetKindInletInterpolationFunction(void) const { return Kind_InletInterpolationFunction; }
 
   /*!
    * \brief Get the kind of inlet face interpolation data type.
    */
-  inline unsigned short GetKindInletInterpolationType (void) const  { return Kind_Inlet_InterpolationType; }
+  inline INLET_INTERP_TYPE GetKindInletInterpolationType (void) const  { return Kind_Inlet_InterpolationType; }
 
   /*!
    * \brief Get whether to print inlet interpolated data or not.

Common/include/geometry/CGeometry.hpp

@@ -1198,26 +1198,6 @@ class CGeometry {
    */
   inline vector<vector<unsigned long> > GetPlanarPoints() const {return Plane_points;}
 
-  /*!
-   * \brief Given arrays x[1..n] and y[1..n] containing a tabulated function, i.e., yi = f(xi), with
-              x1 < x2 < . . . < xN , and given values yp1 and ypn for the first derivative of the interpolating
-              function at points 1 and n, respectively, this routine returns an array y2[1..n] that contains
-              the second derivatives of the interpolating function at the tabulated points xi. If yp1 and/or
-              ypn are equal to 1 × 1030 or larger, the routine is signaled to set the corresponding boundary
-              condition for a natural spline, with zero second derivative on that boundary.
-                        Numerical Recipes: The Art of Scientific Computing, Third Edition in C++.
-   */
-  void SetSpline(vector<su2double> &x, vector<su2double> &y, unsigned long n, su2double yp1, su2double ypn, vector<su2double> &y2);
-
-  /*!
-   * \brief Given the arrays xa[1..n] and ya[1..n], which tabulate a function (with the xai’s in order),
-              and given the array y2a[1..n], which is the output from spline above, and given a value of
-              x, this routine returns a cubic-spline interpolated value y.
-              Numerical Recipes: The Art of Scientific Computing, Third Edition in C++.
-   * \return The interpolated value of for x.
-   */
-  su2double GetSpline(vector<su2double> &xa, vector<su2double> &ya, vector<su2double> &y2a, unsigned long n, su2double x);
-
   /*!
    * \brief Compute the intersection between a segment and a plane.
    * \param[in] Segment_P0 - Definition of the particular problem.

Common/include/linear_algebra/blas_structure.hpp

@@ -482,6 +482,31 @@ class CBlasStructure {
     }
   }
 
+  /*!
+   * \brief Algorithm to solve a linear system with a tridiagonal matrix.
+   * \param[in] lower - lower diagonal
+   * \param[in] main - main diagonal
+   * \param[in,out] upper - upper diagonal (modified on exit)
+   * \param[in,out] rhs - right hand side on entry, solution on exit
+   * \note Same size for all vectors. Use row index for lower and upper vector (e.g. lower[0] does not matter).
+   */
+  template<class Vec, class Scalar = su2double>
+  static void tdma(const Vec& lower, const Vec& main, Vec& upper, Vec& rhs) {
+    const int N = main.size();
+
+    upper[0] /= main[0];
+    rhs[0] /= main[0];
+
+    for (int i=1; i<N; i++) {
+      const Scalar denom = 1.0 / (main[i]-lower[i]*upper[i-1]);
+      upper[i] *= denom;
+      rhs[i] = (rhs[i]-lower[i]*rhs[i-1])*denom;
+    }
+
+    for (int i=N-2; i>=0; i--)
+      rhs[i] -= upper[i]*rhs[i+1];
+  }
+
 private:
 
 #if !(defined(HAVE_LIBXSMM) || defined(HAVE_BLAS) || defined(HAVE_MKL)) || (defined(CODI_REVERSE_TYPE) || defined(CODI_FORWARD_TYPE))