Discussion :

[jslpfasc2]

Bonjour à tous, j'ai un projet de mathématiques dans lequel je dois programmer sur Scilab l'algorithme d'Uzawa. J'ai donc écrit mon programme mais il ne fonctionne pas et je ne comprends pas ce que j'ai fait de mal. Je suis débutant en Scilab donc je suis perdu. D'abord je vous donne le contexte de l'exercice :

On est dans l'algorithme d'Uzawa pour contraintes d'inégalités Cx<f
Voici l'algorithme :
-On se donne lambda0 un vecteur de R^q
Itérer pour n>=0
(i) trouver xn t.q. grad(J(xn)) + (C')*lambdan = 0
(ii) lambdan+1 = PF(lambdan + rho*(C*xn - f))

où rho désigne le pas et PF la projection sur l'ensemble des contraintes qui est
K = {x / x>=g}
g est définit dans mon programme, que voici :

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
function [x,i,resi]=U1(A,b,C,f,eta,Imax,rho)
 
x=zeros(b);
r=A*x-b;
nr=norm(r);
i=0;
err=eta;
resi=[];
lambda=zeros(n,1);
while (i<Imax & norm(C*x-f)>eta) 
  xold=x;
  i=i+1;
  x=inv(A)*(b-C'*lambda);
  lambda = lambda+rho*(C*x-f);
  r=A*x-b;
  resi=[resi;err];
  if nr>1e10; disp('Explosion!'); break; end;
end
if (err<eta); disp('Methode terminee normalement'); end
if (i==Imax); disp('Imax atteint! Augmenter le nombre d''iterations'); end;
endfunction
 
  n=29;
  h=1/(n+1);
 n=29;
A=zeros(n,n);
A=zeros(28,28);
    A(n,n)=2;
    A(n,n-1)=-1;
    for i=1:n-1 do
        A(i,i)=2;
        A(i,i+1)=-1;
        A(i+1,i)=-1;
    end
 
  b=-20*ones(n,1);
 
  eta=1e-4;
  Imax=10000;
  S=spec(A);
  rho=1000;
 
t=zeros(n,1);
for i=1:n do
    t(i)=i*h;
end
 
g=zeros(n,1);
for i=1:n do
  g(i)=-1+max(0,0.565-10*(t(i)-0.4)^2);
end
 
  C=-eye(n,n);
  f=-g;
 
      disp('rho :');disp(rho);

Code :

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

1
2
3
4
5
6
     [x,i,resi]=U1(A,b,C,f,eta,Imax,rho);
      clf;
      xgrid(3);
      plot2d((1:i),resi,logflag='nl', style=[5]);
      xtitle('nombre d itérations = '+string(i));
      disp('dernier residu :');disp(resi(i));

Le problème c'est que ça m'affiche des valeurs bizarres, et que je ne trouve pas ce qui doit converger vers quoi...

Si vous pouviez juste jeter un oeil au programme et me dire s'il y a une absurdité ce serait vraiment génial. Merci beaucoup!

[jslpfasc2]

Ah en fait c'est bon j'ai trouvé mon erreur, j'ai oublié de diviser la matrice A par le pas, du coup ça me donnait des valeurs beaucoup trop grandes! -_-
Je classe donc ce sujet en résolu et qu'il tombe dans les tréfonds du forum.
Merci beaucoup!