Esercitazioni di MATLAB: Lezione del 4 Marzo 2008

Un mini riassunto di MATLAB

Esempio di programmazione

Filelinearsistemsolve.m

:open:
 1%
 2% Risolve il sistema lineare A*x = b
 3%
 4function x = linearsistemsolve( A, b )
 5
 6  [nr,nc] = size(A) ;
 7
 8  % controllo dimensioni
 9  if nr ~= nc
10    error( sprintf( 'la matrice A (%d x %d) non e` quadrata ', nr,nc) ) ;
11  end
12
13  neq = length(b) ;
14
15  if nr ~= neq
16    error( sprintf( 'la matrice A (%d x %d) non e` compatibile\ncon la lunghezza del vettore b, length(b) = %d ', nr, nc, neq ) ) ;
17  end
18
19  M = [ A b ] ; % costruisco la matrice aumentata
20  for i=1:neq-1
21    % cerco l'elemento di massimo modulo colonna i-esima
22    [V,k] = max(abs(M(i:end,i))) ;
23
24    % rimetto l'indice k alla giusta posizione
25    k = k + i-1 ;
26
27    % scambio la riga i con la riga k
28    BF     = M(i,:) ;
29    M(i,:) = M(k,:) ;
30    M(k,:) = BF ;
31
32    % controllo singolarità 
33    if  V == 0 
34      error( 'Matrice singolare' ) ;  
35    end
36
37    % costruisco la matrice di Frobenius
38    L             = eye(neq) ;
39    L(i+1:end,i)  = -M(i+1:end,i)./M(i,i) ;
40
41    % azzero gli elementi della colonna i
42    M             = L * M ;
43
44  end
45    
46  % alloco il vettore soluzione
47  x = zeros(neq,1) ;
48
49  % a questo punto la matrice è in forma triangolare
50  % applichiamo le operazioni di ritorno
51  x(neq) = M(neq,neq) / M(neq,neq+1) ;
52  for i=neq-1:-1:1
53    x(i) = ( M(i,neq+1) - M(i,i+1:end-1) * x(i+1:end) ) / M(i,i) ;
54  end