isop.h

#pragma once
 
#include "source.h"
namespace ISOP
{
 namespace detail
 {
  template<int N> double shape_fun(double xi, double eta);
  
  template<> inline double shape_fun<0>(double xi, double eta) { return -0.25*(1-xi)*(1-eta)*(1+xi+eta); }
  template<> inline double shape_fun<1>(double xi, double eta) { return 0.5*(1-xi)*(1+xi)*(1-eta); }
  template<> inline double shape_fun<2>(double xi, double eta) { return -0.25*(1+xi)*(1-eta)*(1-xi+eta); }
  template<> inline double shape_fun<3>(double xi, double eta) { return 0.5*(1+xi)*(1+eta)*(1-eta); }
  template<> inline double shape_fun<4>(double xi, double eta) { return -0.25*(1+xi)*(1+eta)*(1-xi-eta); }
  template<> inline double shape_fun<5>(double xi, double eta) { return 0.5*(1-xi)*(1+xi)*(1+eta); }
  template<> inline double shape_fun<6>(double xi, double eta) { return -0.25*(1-xi)*(1+eta)*(1+xi-eta); }
  template<> inline double shape_fun<7>(double xi, double eta) { return 0.5*(1-xi)*(1+eta)*(1-eta); }
  
  template<int N = 8>
  inline double isop(const double nodes[], double xi, double eta)
  {
   return isop<N-1>(nodes, xi, eta) + nodes[N-1]*shape_fun<N-1>(xi, eta);
  }
  
  template<>
  inline double isop<0>(const double[], double, double)
  {
   return 0;
  }
 }
 
 inline double isop(
        const double nodes[]  
        , double xi           
        , double eta          
                       )
 {
  return detail::isop(nodes, xi, eta);
 }
 
 template<typename Function>
 struct isop_map
 {
  isop_map(Function f, const double* x, const double* y) : f(f), x(x), y(y) {}
  
  double operator()(double xi, double eta)
  {
   return f(isop(x, xi, eta), isop(y, xi, eta));
  }
  
  Function f;
  const double* x;
  const double* y;
 };
 
 template<typename Function>
 isop_map<Function> make_isop_map(Function f, const double* x, const double* y)
 {
  return isop_map<Function>(f, x, y);
 }
 
 double isop_render(
  Source& source       
  ,const double nodx[] 
  ,const double nody[] 
  ,double a_xi         
  ,double b_xi         
  ,double a_eta        
  ,double b_eta        
 );
}