%A1DCEC.m
%Edwin Kite

%A 1-D Catchment Evolution Code

%Purpose: Simulate the evolution of a catchment.
%TO DO: Add flexure.

nx = 1e2; %number of gridpoints
dx = 1e3; %distance between grid points, in meters
eod = nx/2; %end of deposition (by aeolian processes)

rod = 0.0001; %rate of deposition (m/yr); 

nt = 2002; %number of timesteps
dt = 1e2; %years

k_erod = 0.01;

w = 2e3; %Full width of catchment (m);

p = 1/3600; %Precip. rate (m/s)

i = 1; %Precip. intermittency (dimensionless);

m = 1; n = 1; %ms and ns
%initialize:
z(1:nx) = 0;

%main loop:
for t = 1:nt
    z(1:eod) = z(1:eod) + rod*dt;
    
    %S05 = [0 diff(z) 0]./dx; %dimensionless slope, defined at half-points. zero-elevation boundary condition at the catchment terminus
    %S = (S05(1:end-1) + S05(2:end))./2; %slope evaluated at elevation points 
    S = [diff(z) 0]./dx; %dimensionless slope. zero-slope boundary condition at the catchment terminus
    
    Q = [p.*w*dx*[1:eod],ones(1,50).*p.*w*dx*eod]; %discharge (no infiltration, no evaporation)
    
    Supwind ;%Upwind-differencing slope
    
    stc = k_erod.*Q.^m.*abs(Supwind).^n; %sediment transport capacity; move downstream.
    
    z(1:end-1) = z(1:end-1) - stc(1:end-1);
    z(2:end) = z(2:end) + stc(1:end-1);
    
    zs(t,:) = z;
    
    %Call diffusive delivery of sediment from ridgelines/hillslopes into
    %channels
    
    %Call flexure
    
end