function [err, unew, u_ex] = paraprog(N, M, tau)
    %
    %  Schema d'Euler progressif pour un probleme parabolique unidimensionnel
    %
    %  parametres
    %
    %  N        : nombre de points interieurs dans l'intervalle [0,1]
    %  h        : pas d'espace
    %  M        : nombre de pas de temps
    %  tau      : pas de temps
    %  t        : temps courant
    %  uold     : N-vecteur, uold(i) est une approximation de u(x_i,t_n)
    %  unew     : N-vecteur, unew(i) est une approximation de u(x_i,t_n+1)
    %
    h = 1./(N + 1);
    t = 0.;
    
    alpha = 1.

    %
    % condition initiale (fct w definie ci-dessous)
    %
    for i = 1 : N
      uold(i) = w(i * h);
    end

    %
    % schema d'Euler progressif 
    %
    for n = 1 : M
      unew(1) = ??? ;
      for i = 2 : N - 1
        unew(i) = ??? ;
      end
      unew(N) = ??? ;  
      t = t + tau;

      %
      % reactualiser la solution et imprimer la norme euclidienne de u
      %
      norm2 = 0.;
      for i = 1 : N
        uold(i) = unew(i);
        norm2 = norm2 + uold(i) * uold(i);
      end
      if (mod(n, 50) == 0) 
        fprintf(' pas de temps %d temps %f norm2 %e \n', n, t, norm2);
      end
    end

    %
    % imprimer l'erreur maximum au temps final
    %
    err = 0;
    for i = 1 : N
      x(i) = i * h;
      u_ex(i) = uex(i * h, t);
      erri = abs(unew(i) - uex(i * h, t));
      if (erri > err) 
        err = erri;
      end
    end
    fprintf(' erreur maximum au temps final  %e \n', err);
    plot(x, unew, '*-b', x, u_ex, '-r');
    legend('solution approchee', 'solution exacte');
end

%
% definition de la condition initiale w(x)
%
function init = w(x)
    init = sin(2 * pi * x) ;
end

%
% definition de f(x)
%
function f = f(x, t)
    f = x^3 * sin(2 * pi * x) * exp(-4 * pi^2 * t) ;
%    f= 0;
end

%
% solution u du probleme de la chaleur (utilise pour le calcul de l'erreur)
%
function u_exact = uex(x, t)
    u_exact = sin(2 * pi * x) * exp(-4 * pi^2 * t) ;
end


function [err, u, u_ex] = pararetro(N, M, tau)
    %
    %  Schema d'Euler retrograde pour un probleme parabolique unidimensionnel.
    %  A chaque pas de temps, il s'agit de resoudre le systeme lineaire :
    %            A u^n+1 = u^n + tau f^{n+1}
    %
    %  parametres
    %
    %  N        : nombre d'inconnues du systeme lineaire 
    %  h        : pas d'espace
    %  M        : nombre de pas de temps
    %  tau      : pas de temps
    %  t        : temps courant
    %  u        : N-vecteur, a chaque pas de temps, 
    %             u est second membre du systeme lineaire,
    %             puis solution du systeme lineaire
    %  a        : N-vecteur, diagonale de la matrice A,
    %             puis diagonale de L tq  A=LL^T
    %  c        : (N-1)-vecteur, sur-diagonale et sous-diagonale
    %             de la matrice A, puis sous-diagonale de L tq  A=LL^T
    %
    h = 1/(N + 1);
    t = 0;
    alpha = 1.;

    %
    % condition initiale
    %
    for i = 1 : N
      u(i) = w(i * h);
    end

    %
    % remplissage de la matrice  A
    %
    for i = 1 : N
      a(i) = ???;
    end
    for i = 1 : N-1
      c(i) = ???;
    end

    %
    % decomposition LL^T A = LL^T
    %
    a(1) = sqrt(a(1));
    for i = 1 : N - 1
      c(i) = c(i) / ???; 
      a(i + 1) = sqrt(a(i + 1) - ???);
    end

    %
    % schema d'Euler retrograde: A u^n+1 = u^n + tau f^{n+1}
    %
    for n = 1 : M
      t = t + tau;

      %
      % second membre du systeme lineaire 
      %
      for i = 1 : N
        u(i) = u(i) + ???;
      end
      
      %
      % resolution du systeme lineaire Ly = u^n+ tau f^{n+1}
      %
      u(1) = u(1) / a(1);
      for i = 1 : N - 1
          u(i + 1) = (u(i + 1) - ???) / ???;
      end

    %
    % resolution du systeme lineaire L^T u^n+1 = y
    %
      u(N) = u(N) / a(N);
      for i = N - 1 : -1 : 1
        u(i) = (u(i) - ???) / ???;
      end
    end

    err = 0;
    for i = 1 : N
      x(i) = i * h;
      u_ex(i) = uex(i * h, t);
      erri = abs(u(i) - uex(i * h, t));
      if (erri > err) 
        err = erri;
      end
    end
    fprintf(' erreur maximum au temps final  %e \n', err)
end
