Discussion :

[zackaririane]

Code Scilab :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
clear; 
function [y]=uex(x,t)
    y=sin(%pi*x)*cos(%pi*t)
endfunction
// definition parametres calculs
 
//maillage
un0=0;
unpl=0
N=49; 
L=1.;
h=L/(N+1);
 
//temps
T=1; 
M=499
dt=T/(M+1); 
 
// Lambda
Lambda= (dt^2)/(h*h) ;
 
//maillage
x=zeros(N,1) ;
for j=1:N;
    x(j)= j*h;
 
end;
 
// definition matrice
 
A=zeros(N,N) ;
A(1,1) = 2+Lambda ;
A(1,2) =  -Lambda ;
for j=2:N-1;
    A(j,j) =2+Lambda  ;
    A(j,j+1) = -Lambda ;
    A(j,j-1) = -Lambda ;
end;
A(N,N-1) = -Lambda;
A(N,N) =  2+Lambda ;
 
//==================================
//boucle en temps
//==================================
// Solution à l'instant initial
un=  zeros(N,1);
for j=1:N;
        un(j)= sin(%pi*x(j));
end;
 
 
un_1=zeros(N,1)
 
for n=1:M;
    t(n)=n*dt;
    t(n+1)=t(n)+dt
 
     for j=1:N;
        un_1(j)= x(j);
    end;
 
    for j=1:N;
        ue(j)= uex(x(j),t(n));
    end;
    // definition second membre
 
    b=  zeros(N,1) ;  
 
     un(1)=0
     un_1(1)= un(j)+dt*un(1)
    for j=2:N  ;
        b(j)=2*un_1(j)-un(j); 
    end;
    b(1)=sin(%pi*x(j))
    // calcul de la solution à l'instant t(n+1)
    unp1 = A\b ; // resolution du systeme lineaire INV de A
 
    // Graphe solution approchée vs solution exacte toutes les 100 itérations 
 
 
    if(modulo(n,499)==0) then
        plot(x,unp1,'-b',x,ue,'+b');    
        xlabel('x');ylabel('u');    
        legend(   'solution exacte');
    end
    //passage iteration suivante
    un=un_1
    un_1=unp1
end;