function [ tv, uv ] = crank_nicolson_modelproblem( lambda, y0, t0, tf, Nh )
% CRANK_NICOLSON_MODELPROBLEM Crank-Nicolson method for the model problem 
% ODE in the form
% y'(t) = lambda y(t),  t \in (t0,tf)
% y(0) = y_0
%
%  [ tv, uv ] = crank_nicolson_modelproblem( lambda, y0, t0, tf, Nh )
%  Inputs: lambda = real parameter (negative)
%          y0     = initial value
%          t0     = initial time
%          tf     = final time
%          Nh     = number of time subintervals
%  Output: tv     = vector of time steps (1 x (Nh+1))
%          uv     = vector of approximate solution at times tv
%

tv = linspace( t0, tf, Nh + 1 );
h = ( tf - t0 ) / Nh;

uv = zeros( 1, Nh + 1 );
uv( 1 ) = y0;

alpha = h * lambda / 2;

for n = 1 : Nh    
    uv( n + 1 ) = uv( n ) * ( 1 + alpha ) / ( 1 - alpha ); 
end

return