//##################################################################### 
// Copyright 2021, Craig Schroeder.
// This file is part of PhysBAM whose distribution is governed by the license contained in the accompanying file PHYSBAM_COPYRIGHT.txt.
//#####################################################################
#ifndef __DIV_FREE_INTERPOLATION_INTERP__
#define __DIV_FREE_INTERPOLATION_INTERP__

#include <Core/Interpolation/GENERAL_INTERPOLATION_AXIS_FACE.h>

namespace PhysBAM {

template<int d> struct DIV_FREE_INTERPOLATION_INTERP_STENCIL;

// interpolating polynomial
template<> 
struct DIV_FREE_INTERPOLATION_INTERP_STENCIL<2>
{
    constexpr static int stencil_size=2;
    template<int a,int d,class T> static VECTOR<VECTOR<T,2+(a==d)>,3> Value(T t);
    template<int a,int d,class T> static VECTOR<VECTOR<T,2+(a==d)>,3> Diff(T t);
    template<int a,int d,class T> static VECTOR<VECTOR<T,2+(a==d)>,3> Hess(T t);
};

// interpolating polynomial
template<> 
struct DIV_FREE_INTERPOLATION_INTERP_STENCIL<3>
{
    constexpr static int stencil_size=2;
    template<int a,int d,class T> static VECTOR<VECTOR<T,2+(a==d)>,4> Value(T t);
    template<int a,int d,class T> static VECTOR<VECTOR<T,2+(a==d)>,4> Diff(T t);
    template<int a,int d,class T> static VECTOR<VECTOR<T,2+(a==d)>,4> Hess(T t);
};

template<class TV,class T2,class T_FACE_ARRAY=ARRAY<T2,FACE_INDEX<TV::m> > >
using DIV_FREE_INTERPOLATION_INTERP=GENERAL_INTERPOLATION_AXIS_FACE<TV,T2,DIV_FREE_INTERPOLATION_INTERP_STENCIL<TV::m>,false,T_FACE_ARRAY>;
}

#endif
