Calcul matriciel dans un algo

**thtghgh** · 20/12/2009, 16h07

Bonjour,

Dans un algorithme, je dois calculer cela

for i=1..m
solve H(1:i,1:i)y=s(1:i)
x=x+V(:,1:i)y
stop

Alors je vous donne mon code complet (il y a d'autres calculs avant) car je ne comprends pas d'ou vient mon erreur. Pouvez vous m'aider?
Merci

Annexe - Code Matlab :

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
%%TP Analyse numérique 2
%
% Algorithme du GMRES (avec m redémarrages) préconditionné par la gauche
%
% Entrée : matrice A, vecteur b, vecteur initial x0,
% tolérance
%
% Sortie : vecteur solution x, nombre d'itérations
 
function [solution,it] = gmres_m_pre(A,b,x0,tol,m)
 
n=size(A,2); % Taille de la matrice
 
 
sn=zeros(n,1);
cs=zeros(n,1);
x=x0;
r=b-A*x; % Vecteur 
r0e=r; % Vecteur
p=r; %Vecteur
 
j=1;
it=0;
 
while ((norm(b-A*x)/norm(b)>tol) & (j<=n))
 
    it=it+1;
    M=ssor(A);
    z0=zeros(n,1);
    tol=0.00001;
    rs=inv(M)*r;
 
    v=rs/norm(rs);
    e=zeros(n,1);
    e(1)=1;
    e1=e;
    s=norm(r)*e1;
 
    for i=1:m
 
    z0=zeros(n,1);
    tol=0.00001;
    w=inv(M)*A*v
 
    for k=1:i
        h(k,i)=dot(w,v);
        w=w-h(k,i)*v;
    end
    h(i+1,i)=norm(w);
    v=w/h(i+1,i);
 
    for k=1:i-1
        delta=cs*h(k,i)+sn*h(k+1,i);
        h(k+1,i)=-sn*h(k,i)+cs*h(k+1,i);
        h(k,i)=delta;
    end
 
    kro=sqrt(h(i,i)^2+h(i+1,i)^2);
    cs=h(i,i)/kro;
    sn=h(i+1,i)/kro;
    h(i,i)=cs*h(i,i)+sn*h(i+1,i);
    h(i+1,i)=0;
    temp=s;
    s=-sn*temp;
    temp=cs*temp;
    error=abs(s)/norm(b);
 
 
 
    if (error<=tol)
 
        y= inv(h(1:i,1:i))*s(1:i); % On résoud H(1:i,1:i)y=s(1:i)
 
       x=x+v(:,1:i)*y;
       stop
    end
    end
    if (error<=tol)
        stop
    end
 
    y=inv(h(1:m,1:m))*s(1:m);
    x=x+v(:,1:m)*y;
 
    r=inv(M)*(b-A*x);
    s(i+1)=norm(r);
    error=s(i+1)/norm(b);
end
end

Si vous avez besoin de la fonction ssor.m je vous donnerai le code, ou de l'algo papier je pourrais le mettre à disposition!
J'ai travaillé cet algo mais je ne vois pas d'ou vient mon erreur en fait

Merci

**duf42** · 20/12/2009, 18h01

Bonjour,

Pourrais-tu nous préciser ton erreur?

Merci,
Duf

**thtghgh** · 20/12/2009, 18h34

Excuse moi!
Alors lorsque je teste ma fonction, ca maffiche cette erreur

Code :

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

1
2
3
4
??? Index exceeds matrix dimensions.
 
Error in ==> gmres_m_pre at 85
    x=x+v(:,1:m)*y;

gmres_m_pre étant la fonction que j'ai mise plus haut.

L'algo test est celui ci :

Code :

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

1
2
3
4
5
6
7
8
9
10
11
12
13
 
A=une matrice quelconqie
n=size(A,1)
b=ones(n,1);
x0=zeros(n,1);
tol=0.0000001;
 
[sol3,it3]=gmres_m_pre(A,b,x0,tol,2)
 
th=A\b;
 
err3=norm(sol3-th)
ite3(i,1)=it3