BSpline.h

Go to the documentation of this file.

#ifndef BSPLINE_H
#define BSPLINE_H

#include <vector>

namespace sierra {
namespace nalu {

// Forward declarations
class H5IO;

//====================================================================
//====================================================================

class BSpline{
 public:

  BSpline( const int order,
       const int dimension,
       const bool enableValueClipping=true );

  virtual ~BSpline();

  virtual BSpline * clone() const = 0;

  int get_order() const{return order_;}

  int get_dimension() const{return dim_;}

  virtual double value( const double* ) const = 0;

  double value( std::vector<double> & x ) const{ return value( &x[0] ); }

  virtual void read_hdf5( H5IO & io ) = 0;

  virtual void write_hdf5( H5IO & io ) const = 0;

 protected:

  int order_;
  const int dim_;
  const bool enableValueClipping_;

 private:

  BSpline( const BSpline & );             // no copying
  BSpline & operator=( const BSpline & ); // no assignment

};

//====================================================================
//====================================================================

class BSpline1D : public BSpline{
 public:

  BSpline1D( const int order,
         const std::vector<double> & indepVars,
         const std::vector<double> & depVars,
         const bool allowClipping = true );

  BSpline1D( const bool allowClipping = true );

  BSpline1D( const BSpline1D& );

  BSpline * clone() const{ return new BSpline1D(*this); }

  ~BSpline1D();

  double value( const double* indepVar ) const;
  inline double value( const double & x ) const{ return value(&x); }

  inline const std::vector<double> & get_control_pts() const{ return controlPts_; }
  inline       std::vector<double> & get_control_pts()      { return controlPts_; }
  inline const std::vector<double> & get_knot_vector() const{ return knots_; };

  inline int get_npts() const{ return npts_; }

  inline double get_maxval() const{ return maxIndepVarVal_; }
  inline double get_minval() const{ return minIndepVarVal_; }

  inline const std::vector<double> & get_basisfun() const{ return basisFun_; }

  void sort_inputs( const std::vector<double> & indepVars,
                    const std::vector<double> & depVars,
                    std::vector<double> & sortedIndepVars,
                    std::vector<double> & sortedDepVars ) const;

  void dump();

  void write_hdf5( H5IO & io ) const;
  void  read_hdf5( H5IO & io );

  inline bool operator == (const BSpline1D& a ) const{
    return ( a.npts_ == npts_ &&
         a.maxIndepVarVal_ == maxIndepVarVal_ &&
         a.minIndepVarVal_ == minIndepVarVal_ &&
         a.knots_ == knots_ &&
         a.controlPts_ == controlPts_ );
  }

  inline bool operator != ( const BSpline1D& a ) const{ return !( *this==a ); }

 private:

  void compute_control_pts( const std::vector<double> & indepVars,
                const std::vector<double> & depVars );

  int npts_;
  double maxIndepVarVal_, minIndepVarVal_;
  std::vector<double> knots_, controlPts_;
  mutable std::vector<double> basisFun_;

  BSpline1D& operator=(const BSpline1D&); // no assignment
};

//====================================================================
//====================================================================

class BSpline2D : public BSpline{
 public:

  BSpline2D( const int order,
         const std::vector<double> & indepVars1,
         const std::vector<double> & indepVars2,
         const std::vector<double> & depVars,
         const bool allowClipping = true );

  BSpline2D( const bool allowClipping = true );

  BSpline2D( const BSpline2D& );

  BSpline * clone() const{ return new BSpline2D(*this); }

  ~BSpline2D();

  double value( const double* indepVar ) const;

  void write_hdf5( H5IO & io ) const;
  void  read_hdf5( H5IO & io );

  inline bool operator == (const BSpline2D& a) const{
    bool isEqual = true;
    std::vector<const BSpline1D*>::const_iterator isp  =   dim2Splines_.begin();
    std::vector<const BSpline1D*>::const_iterator ispa = a.dim2Splines_.begin();
    for( ; isp!=dim2Splines_.end(); isp++, ispa++ )
      isEqual = ( (**isp) == (**ispa) ) ? isEqual : false;
    isEqual = (*a.sp1_ == *sp1_) ? isEqual : false;
    return isEqual;
  }

  inline bool operator != ( const BSpline2D& a ) const{ return !( *this==a ); }

 private:

  std::vector<const BSpline1D*> dim2Splines_;
  BSpline1D * sp1_;

  void compute_control_pts( const std::vector<double> & indepVars1,
                const std::vector<double> & indepVars2,
                const std::vector<double> & depVars );

  BSpline2D& operator=(const BSpline2D&); // no assignment
};

//====================================================================
//====================================================================

class BSpline3D : public BSpline{

 public:

  BSpline3D( const int order,
         const std::vector<double> & x1,
         const std::vector<double> & x2,
         const std::vector<double> & x3,
         const std::vector<double> & phi,
         const bool allowClipping = true );

  BSpline3D( const bool allowClipping = true );

  BSpline3D( const BSpline3D& );

  BSpline * clone() const{ return new BSpline3D(*this); }

  ~BSpline3D();

  double value( const double* ) const;

  void write_hdf5( H5IO & io ) const;
  void  read_hdf5( H5IO & io );

  inline bool operator == (const BSpline3D& a) const{
    bool isEqual = true;
    std::vector<const BSpline2D*>::const_iterator isp  = sp2d_.begin();
    std::vector<const BSpline2D*>::const_iterator ispa = a.sp2d_.begin();
    for( ; isp!=sp2d_.end(); isp++, ispa++ )
      isEqual = ( **isp == **ispa ) ? isEqual : false;
    isEqual = (*a.sp1_ == *sp1_) ? isEqual : false;
    return isEqual;
  }

  inline bool operator != ( const BSpline3D& a ) const{ return !( *this==a ); }

 private:

  const int n1_, n2_, n3_;  // number of control points for each dimension
  std::vector<const BSpline2D*> sp2d_;
  BSpline1D * sp1_;

  void compute_control_pts( const std::vector<double> & x1,
                const std::vector<double> & x2,
                const std::vector<double> & x3,
                const std::vector<double> & phi );

  BSpline3D& operator=(const BSpline3D&); // no assignment
};

//====================================================================
//====================================================================

class BSpline4D : public BSpline{

 public:

  BSpline4D( const int order,
         const std::vector<double> & x1,
         const std::vector<double> & x2,
         const std::vector<double> & x3,
         const std::vector<double> & x4,
         const std::vector<double> & phi,
         const bool allowClipping = true );

  BSpline4D( const bool allowClipping = true );

  BSpline4D( const BSpline4D& );

  BSpline * clone() const{ return new BSpline4D(*this); }

  ~BSpline4D();

  double value( const double* x ) const;

  void write_hdf5( H5IO & io ) const;
  void  read_hdf5( H5IO & io );

  inline bool operator == (const BSpline4D& a) const{
    bool isEqual = true;
    std::vector<const BSpline3D*>::const_iterator isp  = sp3d_.begin();
    std::vector<const BSpline3D*>::const_iterator ispa = a.sp3d_.begin();
    for( ; isp!=sp3d_.end(); isp++, ispa++ )
      isEqual = ( **isp == **ispa ) ? isEqual : false;
    isEqual = (*a.sp1_ == *sp1_) ? isEqual : false;
    return isEqual;
  }

  inline bool operator != ( const BSpline4D& a ) const{ return !( *this==a ); }

 private:

  const int n1_, n2_, n3_, n4_;
  std::vector<const BSpline3D*> sp3d_;
  BSpline1D * sp1_;

  void compute_control_pts( const std::vector<double> & x1,
                const std::vector<double> & x2,
                const std::vector<double> & x3,
                const std::vector<double> & x4,
                const std::vector<double> & phi );

  BSpline4D& operator=(const BSpline4D&); // no assignment
};

//====================================================================
//====================================================================

class BSpline5D : public BSpline{

 public:

  BSpline5D( const int order,
         const std::vector<double> & x1,
         const std::vector<double> & x2,
         const std::vector<double> & x3,
         const std::vector<double> & x4,
         const std::vector<double> & x5,
         const std::vector<double> & phi,
         const bool allowClipping = true );

  BSpline5D( const bool allowClipping = true );

  BSpline5D( const BSpline5D& );

  BSpline * clone() const{ return new BSpline5D(*this); }

  ~BSpline5D();

  double value( const double* x ) const;

  void write_hdf5( H5IO & io ) const;
  void  read_hdf5( H5IO & io );

  inline bool operator == (const BSpline5D& a) const{
    bool isEqual = true;
    std::vector<const BSpline4D*>::const_iterator isp  = sp4d_.begin();
    std::vector<const BSpline4D*>::const_iterator ispa = a.sp4d_.begin();
    for( ; isp!=sp4d_.end(); isp++, ispa++ )
      isEqual = ( **isp == **ispa ) ? isEqual : false;
    isEqual = (*a.sp1_ == *sp1_) ? isEqual : false;
    return isEqual;
  }

  inline bool operator != ( const BSpline5D& a ) const{ return !( *this==a ); }

 private:

  const int n1_, n2_, n3_, n4_, n5_;
  std::vector<const BSpline4D*> sp4d_;
  BSpline1D * sp1_;

  void compute_control_pts( const std::vector<double> & x1,
                const std::vector<double> & x2,
                const std::vector<double> & x3,
                const std::vector<double> & x4,
                const std::vector<double> & x5,
                const std::vector<double> & phi );

  BSpline5D& operator=(const BSpline5D&); // no assignment
};

//====================================================================
//====================================================================

} // end nalu namespace
} // end sierra namespace

#endif