Design notes: erosion processes

All (almost all?) modeled erosion processes take as input topographic elevation and compute as output erosion rates integrated over a finite time step. A C++ idiomatic API might be here to use the same function and overload it for each process. The parameters of erosion laws would be defined in their own structures.

Examples for bedrock channel (stream power) erosion and hillslope (diffusion) erosion processes:

struct stream_power_params
{
    double k_coef = 1e-7;
    double m_exp = 0.5;
    double n_exp = 1.;
};

template<class ER, class EL, class FG>
void compute_erosion(ER& erosion,
                     const EL& elevation,
                     const stream_power_params& params,
                     const FG& fgraph);

struct diffusion_params
{
    double k_coef = 1e-3;
};

template<class ER, class EL, class G>
void compute_erosion(ER& erosion,
                     const EL& elevation,
                     const diffusion_params& params,
                     const G& grid);

Some implementations (like implicit-ADI scheme for diffusion) only support a given grid type. We could use SFINAE here:

template<class ER, class EL, class G,
         std::enable_if_t<std::is_same<G, raster_grid_xt<G::xt_selector>>::value, int> = 0>
void compute_erosion(ER& erosion,
                     const EL& elevation,
                     const diffusion_params& params,
                     const G& grid);