Discussion :

[allolbouchra]

Bonjour, Je demande vos aide s'il vous plait à propos d'un problème dans lequel j'ai eu des problèmes dans le code mais je ne peux pas fixer l'erreur
voilà le code Matlab en pièce jointe et Merci
le code a coriigé

Code :

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
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
%%%%%%%%%%%
%   EF1D  %
%%%%%%%%%%%
% programme Elements finis P1 1D
% pour l'équation -\epsi^2 u"+\p(x) u=\lambda u avec des conditions de Dirichlet aux bords
%
clear all; % effacer toutes les variables
clc; % effacer l'ecran matlab
close all; % fermer toutes les figures
format long; % double precision
% la valeur de p(x)
p=1;
K=8;
%la valeur de epsilon
epsi=0.8;
% on recupere le numero et l'abscisse de chaque noeud (sommet)
coor=load('coor.dat');
% on recupere le numero de l'element et les numeros 
% des sommets de l'élément
elem=load('elem.dat');
[ndlt nul]=size(coor);% ndlt : Nombre de degres de liberte total
dim=nul-1;
[nelem nul]=size(elem);% nelem Nombre d'eléments
ndll=nul-2; % Nombre de degrés de liberté local (par élément)
%%%%% Construction des matrices élémentaires
% Ki matrice de regidite locale
% Ki=(phi1',phi2')^t . (phi1',phi2')
% avec phi1(x)=1-x et phi2(x)=x pour x\in [0,1]
Ki=zeros(ndll,ndll);
Ki(1,1)=1;
Ki(1,2)=-1;
Ki(2,1)=Ki(1,2);% par symétrie
Ki(2,2)=1;
% Mi matrice de masse locale
% Mi=(phi1,phi2)^t . (phi1,phi2)
Mi=zeros(ndll,ndll);
Mi(1,1)=1/3;
Mi(1,2)=1/6;
Mi(2,1)=Mi(1,2);% par symétrie
Mi(2,2)=1/3;
% Matrice globale
K=sparse(ndlt,ndlt);% Matrice de regidite globale ; sparse pour matrice creuse
M=sparse(ndlt,ndlt); % Matrice de masse globale ; sparse pour matrice creuse
X=zeros(ndlt,1); %Le vecteur des abscisses des sommets
X(1)=coor(elem(1,2),2);
hM=-1;%initialisation (car on doit determiner h=max(h_i))
%%%%%%%%%%%%%%%
% ASSEMBLAGE
%%%%%%%%%%%%%%%
for i=1:nelem
    %On parcourt les sommets de l'élément
    I=elem(i,2);
    J=elem(i,3);
    % Les abscisses des sommets
    Xs1=coor(I,2);
    Xs2=coor(J,2);
    % Le pas d'espace
    h=Xs2-Xs1;
    hM=max(h,hM);
    X(i+1)=X(i)+h;
    %Assemblage :
    % Matrice de regidite
    % A FAIRE
    K(I,I)=K(I,I)+(1/h)*Ki(1,1);
    K(I,J)=K(I,J)+(1/h)*Ki(1,2);
    K(J,I)=K(I,J); % par symétrie
    K(J,J)=K(J,J)+(1/h)*Ki(2,2);
    % Matrice de masse
    % A FAIRE
    %*p(ref)
    M(I,I)=M(I,I)+h*Mi(1,1);
    M(I,J)=M(I,J)+h*Mi(1,2);
    M(J,I)=M(I,J); % par symétrie
    M(J,J)=M(J,J)+h*Mi(2,2);
end % Fin de l'assemblage
A=epsi^2.*K+p.*M;
% CONDITION AUX LIMITES DE DIRICHLET
% A FAIRE
% RESOUDRE LE SYSTEME LINEAIRE
A(1,:)=0;
A(:,1)=0;
A(1,1)=1;
A(ndlt,:)=0;
A(:,ndlt)=0;
A(ndlt,ndlt)=1;
% Second Membre
R=chol(M);
T=A*R;
[U lambda]=eigPower(T,X,ndlt);
clear Ki Mi % effacer de la mémoire les matrice elementaires
 
figure(1)% Figure numéro 1
% Solution exacte
lambdap=p+epsi^2*K^2*pi^2;
Up=sin(K*pi*X);
plot(X,U,X,Up)
% A FAIRE
% TRACER SUR LA MEME FIGURE LA SOLUTION EXACTE ET LA SOLUTION APPROCHEE
 
legend('solution numérique','solution exacte')
%figure(2) % Figure numéro 2
% A FAIRE
% TRACER LA VALEUR ABSOLUE DE L'ERREUR ENTRE LA SOLUTION 
% EXACTE ET LA SOLUTION APPROCHEE
%E=U-Up;
%plot(X,E)
 
%legend('erreur')
% CALCULER LA NORME INFINIE DE L4ERREUR

[le fab]

bonjour

c'est quoi l'erreur ?

Fabien