Equation de diffusion convection

Fermé
dianbobo Messages postés 151 Date d'inscription mardi 20 avril 2010 Statut Membre Dernière intervention 15 juin 2014 - 28 nov. 2011 à 02:28
bonjour , tout le monde
je ne sais pas si je poste au bon endroit
voilà je voudrais coder en scilab l'équation de diffusion convection
(c(n+1,i)= a1*c(n,i-1)-a2b*c(n,i)+a3*c(n,i+1)) ,j'ai la forme discrétisé de cette équation
voilà comment j'ai implémenté :
n représente le temps et i l'espace.
//************************************      Initialisation des variables          *****************************************//
    //Longueur du domaine
    L= 1000;
    //Pas d'espace
    dx =10
    //Nombre de subdivisions spatiales
    N = L/dx;
    //Pas de temps
    dt=41;
    //Nombre d'itération en temps (Que je fixe moi meme) 
    M=100;
    
    // Vecteur décrivant la discrétisation spatiale
    x=linspace(0,L,N);
    x=x';//Je transpose x pour avoir la meme taille que c1
    //Les différentes constantes
    D=0.1;
    q=0.3;
    u = 0.05;
    //calcul des variables intervenant dans le schéma
    a1 = (dt*u./(q*dx)+dt*D./(q*dx*dx));
    a2 = (2*dt*D./(q*dx*dx)+dt.*u./(q*dx));
    a3 = (dt.*D./(q*dx*dx));
    //initialisation de la solution initiale en tant que vecteur
    c0=zeros(N,1);
    c0(1)=1;
    //initialisation du vecteur c1	
    c1=zeros(N,1);
//************************************      Fin de l'nitialisation des variables          *****************************************//	
    
    
//****************************************     Calcul de la solution        ***********************************************************//
for k=1:M                                      //La boucle en temps
    for i=2:N-1                              // La boucle en espace
        c1(i)=a1.*c0(i-1)-a2.*c0(i)+a3.*c0(i+1);  //Le schéma
     end                                       //fin de la boucle en espace
     c1(1)=c0(1);	 //condition aux limites (bord 1)
     c1(N)=c0(N);	//condition aux limites (bord2)

     c0=c1  ;                                  // c1 est calculé à  chaque pas de temps puis est echangé avec c0
   end	 //Fin de la boucle en temps
   
//****************************************     Fin du Calcul     ***********************************************************// 
   
//*****************************              Représentation Graphique          ****************************************************// 
   
plot2d(x,c1,style = 5)






le code tourne mais j'ai l'impression d'avoir mal écrit
car je n'ai pas la convergence
merci de votre aide