//#####################################################################
// Copyright 2006, Andrew Selle.
// This file is part of PhysBAM whose distribution is governed by the license contained in the accompanying file PHYSBAM_COPYRIGHT.txt.
//#####################################################################
// Class ADVECTION_MACCORMACK_UNIFORM
//#####################################################################
#ifndef __ADVECTION_MACCORMACK_UNIFORM__
#define __ADVECTION_MACCORMACK_UNIFORM__

#include <Core/Advection/ADVECTION.h>
namespace PhysBAM{

template<class TV,class T2> class BOUNDARY;

template<class TV,class T2,class T_NESTED_ADVECTION>
class ADVECTION_MACCORMACK_UNIFORM:public ADVECTION<TV,T2>
{
    typedef typename TV::SCALAR T;
    typedef VECTOR<int,TV::m> IV;
private:
    // false to use 1st order advection
    const ARRAY<bool,IV>* node_mask;
    const ARRAY<bool,IV>* cell_mask;
    const ARRAY<bool,FACE_INDEX<TV::m> >* face_mask;
    T_NESTED_ADVECTION& nested_advection;
public:
    bool clamp_extrema; // otherwise fall back to 1st order
    bool ensure_second_order; // sacrifice stability for accuracy.
    const BOUNDARY<TV,T>* boundary=0;
    
    ADVECTION_MACCORMACK_UNIFORM(T_NESTED_ADVECTION& nested_advection_input,const BOUNDARY<TV,T>& boundary,const ARRAY<bool,IV>* node_mask_input,const ARRAY<bool,IV>* cell_mask_input,const ARRAY<bool,FACE_INDEX<TV::m> >* face_mask_input);

    void Update_Advection_Equation_Node(const GRID<TV>& grid,ARRAY<T2,IV>& Z,const ARRAY<T2,IV>& Z_ghost,const ARRAY<TV,IV>& V,const T dt,const T time,
        const ARRAY<T2,IV>* Z_min_ghost_input=0,const ARRAY<T2,IV>* Z_max_ghost_input=0,ARRAY<T2,IV>* Z_min_input=0,ARRAY<T2,IV>* Z_max_input=0,
        ARRAY<T2,FACE_INDEX<TV::m> >* flux=0);
    void Update_Advection_Equation_Cell_Lookup(const GRID<TV>& grid,ARRAY<T2,IV>& Z,const ARRAY<T2,IV>& Z_ghost,const ARRAY<typename TV::SCALAR,FACE_INDEX<TV::m> >& face_velocities,
        const T dt,const T time,const ARRAY<T2,IV>* Z_min_ghost_input,const ARRAY<T2,IV>* Z_max_ghost_input,ARRAY<T2,IV>* Z_min_input,
        ARRAY<T2,IV>* Z_max_input,
        ARRAY<T2,FACE_INDEX<TV::m> >* flux);
    void Update_Advection_Equation_Face_Lookup(const GRID<TV>& grid,ARRAY<T,FACE_INDEX<TV::m> >& Z,const ARRAY<typename TV::SCALAR,FACE_INDEX<TV::m> >& Z_ghost,
        const ARRAY<typename TV::SCALAR,FACE_INDEX<TV::m> >& face_velocities,const T dt,const T time,const ARRAY<typename TV::SCALAR,FACE_INDEX<TV::m> >* Z_min_ghost_input,
        const ARRAY<typename TV::SCALAR,FACE_INDEX<TV::m> >* Z_max_ghost_input,ARRAY<T,FACE_INDEX<TV::m> >* Z_min_input,ARRAY<T,FACE_INDEX<TV::m> >* Z_max_input,
        VECTOR<ARRAY<T2,FACE_INDEX<TV::m> >,TV::m>* flux);
//#####################################################################
};
}
#endif
