Discussion :

[seanp223]

Bonjour,

J'essaie avec mes modestes connaissances en programmation LUA d'écrire un script pour la résolution d'un problème particulier par la méthode des différences finies (Finite Differences Method - FDM).

L'application de cette méthode passe par la discrétisation d'un espace, puis la génération d'un système de n équations et n inconnues, donc solvable sous certaines conditions.

J'ai donc le système Ax = b (avec A une matrice nxn, x un vecteur colonne de taille n, et b un vecteur de taille n).

Il existe certes plusieurs outils pour résoudre un tel système, mais je souhaiterais implémenter de moi-même un algorithme de résolution. Après quelques recherches, je suis tombé sur la Méthode de la surrelaxation, SOR.
J'ai produit le script LUA suivant, qui a l'air de marcher, mais...

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
 
--Fonction évaluant la norme d'un vecteur
function Norm(m)
	local n = 0
	for i,v in ipairs(m) do
	n = n + v^2
	end
	return math.sqrt(n)
end
 
--Système à résoudre
local A = {
 {4,2,-4,0,2,4,0,0,0},
 {2,2,-1,-2,1,3,2,0,-6},
 {-4,-1,14,1,-8,-3,5,6,20},
 {0,-2,1,6,-1,-4,-3,3,23},
 {2,1,-8,-1,22,4,-10,-3,9},
 {4,3,-3,-4,4,11,1,-4,-22},
 {0,2,5,-3,-10,1,14,2,-15},
 {0,0,6,3,-3,-4,2,19,45}
}
 
local p,q,ACCURRACY=0,0,1E-4
local w=1.86 --paramètre de relaxation
local n = table.getn(A) --taille du système
local MAX = 10000 --nombre maximal d'itérations
local t = MAX
local x = {}  
	for k = 1,n do x[k]=0 end  --vecteur solution.
 
--Processus itératif
repeat
t=t-1
	q = Norm(x)
	for i=1,n do
	local sigma = 0
		for j=1,i-1 do
		sigma = sigma+A[i][j]*x[j]
		end
		for j=i+1,n do
		sigma = sigma+A[i][j]*x[j]
		end
	x[i] = ((1-w)*x[i])+((w/A[i][i])*(A[i][n+1]-sigma))
	end
	p = Norm(x)
until ((math.abs(p-q)<ACCURRACY) or (t==0))
 
--Affichage du vecteur X obtenu en sortie
for k,v in ipairs(x) do
print(k,v)
end

Ci-joint au message, j'ai composé une image qui synthétise tout ce que je pose comme soucis par la suite...Cela permettra sans nul doute de voir plus clair dans mes élucubrations...

Alors...Voila..Avec le système défini dans le code, le résultat est vite obtenu...Avec d'autres systèmes, ça a l'air de toujours marcher...Sauf que pour un même système d'équations, je remarque que l'algorithme ne marche plus dès que je permute deux lignes...Pourtant , c'est bien le même système...!

Et aussi, chose étrange...Je remarque que lorsque les solutions du système sont des valeurs entières, l'algorithme n'a pas l'air de marcher.Pourtant, SOR est définie comme une méthode de résolution d'un système, sans restriction..Alors ça devrait fonctionner, quelque soit le cas...

Je suis donc persuadé que le problème vient de moi (mauvais code, ou mauvaise interprétation), et je m'en remets aux experts pour me guider vers une solution...

Merci!

[FR119492]

Salut!

l'algorithme ne marche plus dès que je permute deux lignes

L'algorithme ne marche que si la matrice est symétrique à diagonale dominante, ce qui n'est évidemment plus le cas si tu permutes deux lignes.
Jean-Marc Blanc

[seanp223]

Ah...

J'ai vu par ici.
Effectivement, lorsque je passe par la méthode FDM, il est assez probable que le système produit soit ordonné de manière à respecter le critère "diagonale dominante".

Mais que dire de l'autre aspect du problème ..? Je remarque que lorsque le vecteur solution n'a que des coefficients entiers, généralement l'aglgorithme n'aboutit pas...

[Aleph69]

Bonjour,

dans le cas 1, la matrice n'est ni symétrique ni à diagonale dominante.
J'ai cru de prime abord qu'elle n'était pas carrée mais tu sembles avoir placé le vecteur second membre en dernière colonne.

EDIT : au passage, les différences finies n'assurent pas nécessairement la dominance diagonale; cela dépend avant tout des équations qu'on discrétise.

[seanp223]

Je m'en suis rendu compte peu de temps après avoir posté...Merci de la précision..

Ceci dit, existe-t-il une autre méthode efficace, que l'on puisse implémenter pour la résolution des systèmes d'équations linéaires ? J'ai déjà pensé à la Méthode de Cramer...Elle me semble assez simple, il suffit d'implémenter un algorithme de calcul de déterminant de matrice nxn...mais la plupart des articles qui l'évoquent la présente comme une méthode peu efficace...Qu'en est-il exactement ? Trop de ressources CPU mobilisées ?

[Aleph69]

Le calcul du déterminant d'une matrice se fait traditionnellement à l'aide d'une factorisation LU (simplifiée). Pour résoudre un système linéaire, mieux vaut donc directement calculer cette factorisation et résoudre les deux systèmes triangulaires associés.

On peut aussi résoudre les systèmes linéaires avec des méthodes itératives appelées "méthodes de Krylov". Les choix ne manquent pas : gradient conjugué, minres, bicg, bicgstab, qmr, tfqmr, gmres, ...

Les algorithmes itératifs pré-cités sont tous décrits dans le livre de référence "Iterative methods for sparse linear systems", dont une version gratuite peut-être téléchargée à partir de ce lien :
http://www-users.cs.umn.edu/~saad/books.html

[FR119492]

Salut!

J'ai déjà pensé à la Méthode de Cramer

C'est certainement le plus mauvais choix que tu puisses faire: pour la plupart des méthodes directes (LU, Cholesky, etc.), le temps de calcul est plus ou moins proportionnel au cube de la taille de la matrice, alors que pour la méthode de Cramer, c'est sa factorielle.

D'autre part, as-tu lu le tutoriel que j'ai écrit à ce sujet?
Jean-Marc Blanc

[seanp223]

@FR119492 : Un lien peut -être vers ce tutoriel ?
@Aleph69 : Merci bien pour le livre!

[FR119492]

Salut!

Un lien peut -être vers ce tutoriel ?

Je constate maintenant que, pour une obscure raison, il a, une fois de plus, disparu de la liste des cours d'algorithmique. Malheureusement, je ne peux rien faire d'autre qu'espérer qu'un responsable du site lira ce message et corrigera cette lacune.
Jean-Marc Blanc

[seanp223]

Entendu donc.
Merci, je considère que mon problème est résolu à présent.