//#####################################################################
// Copyright 2005-2007, Geoffrey Irving, Andrew Selle, Tamar Shinar.
// This file is part of PhysBAM whose distribution is governed by the license contained in the accompanying file PHYSBAM_COPYRIGHT.txt.
//#####################################################################
// Class BOUNDARY_CONDITION_DOUBLE_FINE
//#####################################################################
#ifndef __BOUNDARY_CONDITION_DOUBLE_FINE__
#define __BOUNDARY_CONDITION_DOUBLE_FINE__

#include <Core/Arrays/ARRAYS_ND.h>
#include <Core/Data_Structures/HASHTABLE.h>
#include <Core/Grids/GRID.h>
#include <Core/Meshes/MESHES_FORWARD.h>
#include <functional>
namespace PhysBAM{

template<class TV> class IMPLICIT_OBJECT;

template<class TV>
class BOUNDARY_CONDITION_DOUBLE_FINE
{
    typedef VECTOR<int,TV::m> IV;typedef typename TV::SCALAR T;
    typedef MESH_OBJECT<TV,TV::m-1> T_SURFACE;
    typedef MESH_OBJECT<TV,TV::m> T_OBJECT;

public:
    enum bc_enum {bc_slip=-3,bc_noslip=-2,bc_free=-1}; // Do not change order
    enum replace_enum {average,replace,skip};
    GRID<TV> mac_grid; // original grid
    GRID<TV> grid; // double fine, node grid.
    int ghost;

    // Indexed by double-fine nodes.
    // Entries between mac-grid faces store boundary condition values.
    // The rest store bits indicating whether flood fill can pass across the entry in a given direction.
    ARRAY<char,IV> bc_type,bc_type_current;

    HASHTABLE<IV,pair<T,int> > bc_p;

    // true=normal direction
    HASHTABLE<FACE_INDEX<TV::m>,tuple<T,int,bool> > bc_u;

    // Sufficient data to solve BC data extending out ghost cells.
    BOUNDARY_CONDITION_DOUBLE_FINE(const GRID<TV>& mac_grid,int ghost);
    BOUNDARY_CONDITION_DOUBLE_FINE(const BOUNDARY_CONDITION_DOUBLE_FINE&) = delete;
    void operator=(const BOUNDARY_CONDITION_DOUBLE_FINE&) = delete;
    ~BOUNDARY_CONDITION_DOUBLE_FINE()=default;

    void Reset();

    struct CB_DATA
    {
        FACE_INDEX<TV::m> face;
        IV cell;
        TV X;
        int e;
    };

    static IV Face_To_Fine_Node(const FACE_INDEX<TV::m>& face)
    {
        return (2*face.index+1).Add_Axis(face.axis,-1);
    }

    static IV Cell_To_Fine_Node(const IV& cell)
    {
        return 2*cell+1;
    }

    static IV Node_To_Fine_Node(const IV& node)
    {
        return 2*node;
    }

    static FACE_INDEX<TV::m> Fine_Node_To_Face(const IV& node, int axis)
    {
        return FACE_INDEX<TV::m>(axis,(node.Add_Axis(axis,1)-1)/2);
    }

    static IV Fine_Node_To_Cell(const IV& node)
    {
        return (node-1)/2;
    }

    static IV Fine_Node_To_Node(const IV& node)
    {
        return node/2;
    }

    // f is only called if boundary condition is used
    void Set(const IMPLICIT_OBJECT<TV>* io,char type,std::function<T(const CB_DATA& data)> f=0,
        bool thin=false,bool invert=false,T contour=0,replace_enum rep=average);
    void Set(const T_SURFACE& surface,char type,std::function<T(const CB_DATA& data)> f=0,
        bool thin=false,replace_enum rep=average);
    void Set(T_OBJECT& object,char type,std::function<T(const CB_DATA& data)> f=0,
        bool thin=false,replace_enum rep=average);
    void Set_Domain_Walls(int side_mask,char type,std::function<T(const CB_DATA& data)> f=0,
        replace_enum rep=average);
    void Set_Region(IV mac_cell_seed,char type);

    void Get_Pressure_Boundary_Conditions(ARRAY<bool,IV>& psi_D,
        ARRAY<bool,FACE_INDEX<TV::m> >& psi_N,ARRAY<T,IV>& p,
        ARRAY<T,FACE_INDEX<TV::m> >& u,int ghost_in) const;
    void Get_Viscosity_Boundary_Conditions(ARRAY<bool,IV>& psi_D,
        ARRAY<bool,FACE_INDEX<TV::m> >& psi_N,ARRAY<T,IV>& u,int axis,int ghost_in) const;
    void Commit_BC(const HASHTABLE<IV,pair<T,int> >& bc_p_cur,
        const HASHTABLE<FACE_INDEX<TV::m>,tuple<T,int,bool> >& bc_u_cur,
        replace_enum rep);
    T* Add_Local_BC(HASHTABLE<IV,pair<T,int> >& bc_p_cur,
        HASHTABLE<FACE_INDEX<TV::m>,tuple<T,int,bool> >& bc_u_cur,
        const CB_DATA& data,char type,bool normal_face,replace_enum rep) const;
//#####################################################################
};
}
#endif
