Discussion :

[saad.pizzicato]

Bonsoir les amis,

Dans le cadre de ma thèses je suis entrain de réaliser un programme sur MATLAB: ayant a transformer une base de fonctions W0(i,x) i=1 à 5 en une base Orthonormée de manière à avoir:
Intégrale de 0 à 1 de W(i,x)*W(j,x).dx = Symbole de Chronecker i,j

J'ai programmé ce qui suit:
1/ pour la définition de mes fonctions:

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
 
%Fonctions de bases w0 et ses derivee w1 et w2%
 
 
    for i=1:n
    Alpha1(i)= -cos(nu(i))+ cosh(nu(i));
 
    Alpha2(i)= sinh(nu(i))-sin(nu(i));
 
    for kk=1:length(x)
 
    w0(i,kk)= (cosh(nu(i)*x(kk))-cos(nu(i)*x(kk)))/Alpha1(i)-(sinh(nu(i)*x(kk))-sin(nu(i)*x(kk)))/Alpha2(i);
    w1(i,kk)= nu(i)*((sinh(nu(i)*x(kk))+sin(nu(i)*x(kk)))/Alpha1(i)-(cosh(nu(i)*x(kk))-cos(nu(i)*x(kk)))/Alpha2(i));
    w2(i,kk)= nu(i)*nu(i)*((cosh(nu(i)*x(kk))+cos(nu(i)*x(kk)))/Alpha1(i)-(sinh(nu(i)*x(kk))+sin(nu(i)*x(kk)))/Alpha2(i));
 
    end
    end

Pour l'algorithme d'Orthonormalisation de GRAM SHMIDT:

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
 
% Normalisation des fonctions de base w0  par l'algorithme de Gramm Shmidt%
 
 
    w0carre=w0.^2;
 
 
%Norme w0%
for i=1:n;
w0int(i)=trapz(0:dx:1,w0carre(i,:)); 
Normw0(i)=sqrt(w0int(i));
end
 
%Orthonormalisation%
 
for kk=1:length(x);
    w0n(1,kk)=w0(1,kk)/Normw0(1);    %w0n(1,:) normalisé%
end
 
for kk=1:length(x);
    Provis21n(kk)=w0n(1,kk)*w0(2,kk);
end
 
ALPHA21n=trapz(0:dx:1,Provis21n);
 
for kk=1:length(x);
    w02nn(2,kk)=w0(2,kk)-(ALPHA21n*w0n(1,kk));
    w0n(2,kk)= w02nn(2,kk)/Normw0(2);  %w0n(2,:) normalisé%
end
 
for kk=1:length(x);
    Provis31n(kk)=w0n(1,kk)*w0(3,kk);
    Provis32n(kk)=w0n(2,kk)*w0(3,kk);
end
 
ALPHA31n=trapz(0:dx:1,Provis31n);
ALPHA32n=trapz(0:dx:1,Provis32n);
 
for kk=1:length(x);
    w03nn(3,kk)=w0(3,kk)-(ALPHA31n*w0n(1,kk))-(ALPHA32n*w0n(2,kk));
    w0n(3,kk)= w03nn(3,kk)/Normw0(3);  %w0n(3,:) normalisé%
end
 
for kk=1:length(x);
    Provis41n(kk)=w0n(1,kk)*w0(4,kk);
    Provis42n(kk)=w0n(2,kk)*w0(4,kk);
    Provis43n(kk)=w0n(3,kk)*w0(4,kk);
end
 
ALPHA41n=trapz(0:dx:1,Provis41n);
ALPHA42n=trapz(0:dx:1,Provis42n);
ALPHA43n=trapz(0:dx:1,Provis43n);
 
for kk=1:length(x);
    w04nn(4,kk)=w0(4,kk)-(ALPHA41n*w0n(1,kk))-(ALPHA42n*w0n(2,kk))-(ALPHA43n*w0n(3,kk));
    w0n(4,kk)= w04nn(4,kk)/Normw0(4);  %w0n(4,:) normalisé%
end
 
for kk=1:length(x);
    Provis51n(kk)=w0n(1,kk)*w0(5,kk);
    Provis52n(kk)=w0n(2,kk)*w0(5,kk);
    Provis53n(kk)=w0n(3,kk)*w0(5,kk);
    Provis54n(kk)=w0n(4,kk)*w0(5,kk);
end
 
ALPHA51n=trapz(0:dx:1,Provis51n);
ALPHA52n=trapz(0:dx:1,Provis52n);
ALPHA53n=trapz(0:dx:1,Provis53n);
ALPHA54n=trapz(0:dx:1,Provis54n);
 
for kk=1:length(x);
    w05nn(5,kk)=w0(5,kk)-(ALPHA51n*w0n(1,kk))-(ALPHA52n*w0n(2,kk))-(ALPHA53n*w0n(3,kk))-(ALPHA54n*w0n(4,kk));
    w0n(5,kk)= w04nn(4,kk)/Normw0(5);  %w0n(5,:) normalisé%
end

Mon souci c'est qu'après exécution de l'algorithme pour mes fonction j'obtiens pas des 0 pour Intégrale de 0 à 1 W(i,x).W(j,x).dx lorsque i≠j.

J'obtiens en effet des valeurs comme:

ans =

1.3310e-018


ans =

8.9951e-018


ans =

-1.3925e-010

Je veux savoir si l'erreur provient de ma programmation de l'algorithme? ou que c'est une erreur issue de l'intégration (Par trapèze avec un pas de 1/1000)? ou d'autre chose???


Merci d'avance les amis....

[FLB]

Salut,
quand tu fais des calculs flottants sur Matlab, tout ce qui est de l'ordre de

Code :

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

eps

vaut 0, voir cet article de la faq. Ceci explique au moins les 2 premières valeurs, et j'imagine que la 3ème est aussi le résultat d'un cumul d'approximation numérique.