////////////////////////////////////////////////////////////////
// Copyright G. Allaire, Octobre 2017
//
// Sur divers schemas de differences finies
// pour l'equation de transport:
// mise en evidence d'instabilite, décentrement aval
//
////////////////////////////////////////////////////////////////
//
// equation de transport
// u,t + v u,x = 0
//
////////////////////////////////////////////////////////////////
//
lg = 10. ; // lg = demie longueur du domaine
dx = 0.1 ; // dx = pas d'espace
nx = lg/dx ; // nx = demi nombre de mailles
v = -1. ; // v = vitesse
cfl = 0.9 ; // cfl
dt = dx*cfl ; // dt = pas de temps
tfinal = 1. ; // tfinal = temps de simulation
nt = int(tfinal/dt) ; // nt = nombre de pas de temps effectues
//
// initialisation
//
x=zeros(1,2*nx+1) ;
u0=zeros(1,2*nx+1) ;
for i=1:2*nx+1
  x(i) = (i-nx-1)*dx  ;
  u0(i) = max(0.,1.-x(i)**2) ;
end
u=u0 ;
up=u0 ;
um=u0 ;
uexacte=u0 ;

clf()

a=get("current_axes"); //get the handle of the newly created axes
a.axes_visible="on"; // makes the axes visible
a.font_size=6; //set the tics label font size
a.x_location="bottom"; //set the x axis position
a.data_bounds=[-lg,-0.1;lg,1.3]; //set the boundary values for the x, y and z coordinates.
a.sub_tics=[5,0];
a.thickness=2;
//a.labels_font_color=5;

plot(x,u,'b-d')
legends(["solution exacte";"schéma explicite décentré aval"],[[id;1],[-5;1]],opt="ur",font_size=6) ;
title ('ADVECTION, SCHEMA EXPLICITE DECENTRE AVAL' ,'fontsize',5) ;

halt() ;
////////////////////////////////////////////////////////////////
// schema explicite decentre aval
////////////////////////////////////////////////////////////////
tfinal = 1.2 ; // tfinal = temps de simulation
nt = int(tfinal/dt) ; // nt = nombre de pas de temps effectues
//
u=u0 ;
//
for n=1:nt
//
up = shift('+1',u) ;
um = shift('-1',u) ;

u=u - v*dt/(dx)*(u-um) ;

if modulo(n,1) == 0
        uexacte = zeros(1,2*nx+1) ;
        for i=1:2*nx+1
                 j = i - n*dt*v/dx ;
                 j = min(j,2*nx+1) ;
                 uexacte(i) = u0(j) ;
        end
        clf()
        a=get("current_axes");
        a.data_bounds=[-lg,-0.1;lg,1.3];
        a.font_size=6;
        a.thickness=2;
        plot(x,u,'b-d',x,uexacte,'r-')
        legends(["solution exacte";"schéma explicite décentré aval"],[[id;1],[-5;1]],opt="ur",font_size=6) ;
        title ('ADVECTION, SCHEMA EXPLICITE DECENTRE AVAL' ,'fontsize',5) ;
        sleep(100) ;
end
//
end


// comparaison au temps final
        uexacte = zeros(1,2*nx+1) ;
        for i=1:2*nx+1
                 j = i - nt*dt*v/dx ;
                 j = min(j,2*nx+1) ;
                 uexacte(i) = u0(j) ;
        end
        clf()
        a=get("current_axes");
        a.data_bounds=[-lg,-0.1;lg,1.3];
        a.font_size=6;
        a.thickness=2;
        plot(x,u,'b-d',x,uexacte,'r-')
        legends(["solution exacte";"schéma explicite décentré aval"],[[id;1],[-5;1]],opt="ur",font_size=6) ;
        title ('ADVECTION, SCHEMA EXPLICITE DECENTRE AVAL' ,'fontsize',5) ;