Méthodes numériques et langage C Plan du cours Algebre linéaire Espaces de dimension finie

14 Pages·2015·1.39 MB·French

Plan du cours
Rappels d’alg`ebre lin´eaire 4
VecteursdeRn 4
Applicationlin´eaire 7
Matricesetpropri´et´es 8
M´ethodes num´eriques et langage C
Repr´esentation en langage C 16
Alg`ebre lin´eaire Blasetallocationm´emoire 16
Vecteurs 20
Matrices 23
R.Flamary Matricescreuses 25
Op´erations vectorielles et matricielles 27
Op´erationssurlesvecteurs 27
Op´erationsmatrice/vecteur 30
17octobre2018
Op´erationsmatricielles 33
R´esolutions d’´equations lin´eaires 35
Syst`emed’´equationsetmatrices 35
Matricesparticuli`eres 39
PivotdeGauss 42
Factorisationdematrice 48
2/53
Algebre lin´eaire Espaces de dimension finie
Espace vectoriel
Principe UnespacevectorielsurK estestunensembleE munidedeuxlois:
I Outilsmath´ematiquespermettantdemod´eliserlesespacesvectorielsetdes I L’addition+:E2 E,op´erationcommutativeetassociatived’´el´ementneutre0
→
relationslin´eairesentrevecteurs. appel´evecteurnul.
I Relationslin´eairessimples`ainterpr´eter. I Lamultiplicationparunscalaire :K E E,op´erationdistributiveet
· × →
d’´el´ementneutre1.
I Mod´elisationdesdonn´eesnum´eriquessouslaformedetableau.
Unvecteurdev E estd´efinitparunsymboleengrasdanscecours.
∈
Exemples d’utilisation Famille de vecteurs
I Syst`emesderecommandation(Netflix,Amazon,Criteo). I Soitunefamilledevecteur vi i=1,...,m.
{ }
I Moteursderecherche(Google,Bing,Duckduckgo). I Lafamilleestditelin´eairementind´ependantesilaseulecombinaisonlin´eairedela
forme λ v =0estcelleavecλ =0, i.
I Pr´edictiondevaleursr´eelles(M´et´eo,coursdelabourse). i i i i ∀
I Traitementdusignaletdesimages(syst`emeslinaires,filtrage,reconstruction). I LesouPs-espacevectorielengendr´eparunefamille{vi}i=1,...,m d´enot´e
Vec( v )estl’ensembledescombinaisonslin´eairesdev .
i i=1,...,m i
I Simulationsphysiques(m´ethodesd’´el´ementfinis,´equationsdelachaleur). { }
I Lafamille vi i=1,...,m lin´eairementind´ependanteestunebasedeE si
{ }
Vec( v )=E.
i i=1,...,m
{ }
3/53 4/53
Vecteurs dans Rn Produit scalaire
Notation
LesvecteursdeRn peuventˆetreexprim´essouslaforme: I Leproduitscalaireentredeuxvecteursxety deRn estl’applicationbilin´eaire:
v1 x,y =x>y= xiyi
v h i
2 i
v= ..  ou v>=[v1,v2,...,vn] X
 .  I Lanormeeuclidienneaucarr´ed’unvecteurestsonproduitscalaireaveclui
vn mˆeme:
  v 2= v,v =v>v
Enl’absencedusymbole ,touslesvecteurssontsuppos´esvecteurscolonne. k k h i
>
I Unespacevectorielassoci´e`aunproduitscalaireestappel´eunespace
Base canonique et norme
pr´e-hilbertien.
I LabasecanoniquedeRn estlafamilledevecteurs vi i=1,...,n tellequevi=δi I In´egalit´edeCauchy-Schwartz:
{ }
ou` δ estlevecteurdiracnulsurtoutessescomposante`aparlacomposantei
i
´egale`a1. x,y E x,y x y
∀ ∈ |h i| ≤ k kk k
I Lanormed’unvecteurdeRn estd´efiniepar
I L’angleθ entredeuxvecteurspeutˆetred´etermin´egrˆace`a
v 2= v2
k k i x,y = x y cosθ
i h i k kk k
X
I Leproduitscalairepeutˆetreutilis´epourprojeterunvecteursurunautrevecteur.
I Labasecanoniqueestorthonormale,c’est`adireorthogonaleetchaque´el´ementa
unenormeunit´e.
5/53 6/53
Application lin´eaire Matrice
Notation
D´efinition LesmatricesdeRm×n sontd´efiniespar:
SoientE etF deuxespacesvectoriels.Uneapplicationf deE versF estditelin´eaire
sielleestadditiveetcommute`alamultiplicationparlesscalaires: a1,1 a1,2 a1,n
···
a a a
2,1 2,2 2,n
∀x,y∈E, f(x+y)=f(x)+f(y), A= ... ... ·.·.·. ... , ai,j ∈R,∀i,j
x E, λ K, f(λx)=λf(x).  
∀ ∈ ∀ ∈ am,1 am,2 ··· am,n
 
Autrementdit,f pr´eservelescombinaisonslin´eaires.
Onpourra´egalementutiliserA(:,k)pourd´esignerlaki`emecolonnedeAetA(k,:)
pourd´esignerlaki`emelignedeA.
Noyaux et image
Attention,enConutilisel’indexage`a0,a correspondradonc`al’indexentier(0,0).
Pourtouteapplicationlin´eairef deE dansF, 1,1
Op´erations de base
I LesvecteursxdeE telsquef(x)=0formentunsous-espacevectorieldeE,
appel´elenoyaudef etnot´eKer(f). I Transposition
I lesvecteursf(x)pourxdansE formentunsous-espacevectorieldeF,appel´e C=A> ci,j =aj,i
⇔
l’imagedef etnot´eIm(f).
I Addition
I Lesdimensionsdessous-espacessontli´eesparlet´eor`emedurang: C=A+B ci,j =ai,j+bi,j
⇔
dim(E)=dim Ker(f)+dim Im(f) I Multiplicationparunscalaireα
C=αA c =αa
ij ij
⇔
7/53 8/53
14, 6, 2,
−
1 2 3 1 1 1 2 2 1
−
2 4 6 , 2 2 2 , 2 3 2
    − − 
3 6 9 3 3 3 1 2 2
−
     
1 0 2 2
− 22 − −
2 , 1 , , 3 4
−  −  28 − − 
1 0 (cid:20) (cid:21) 2 2
     
Multiplication matricielle Matrice et application lin´eaire
Multiplication matricielle I Touteapplicationlin´eaireentredeuxespacesmunischacund’unebasefinieest
SoientdeuxmatricesA Rm×p et ncolonnes repr´esentableparunematrice.
B∈Rp×n ∈ nes B I Uneapplicationlin´eairef :Rn→Rm peutdoncˆetremisesouslaforme:
g
li f(x)=Ax
p p
pcolonnes
C=AB c = a b
⇔ i,j i,k k,j
kX=1 I L’imageIm(f)def estlesousespaceVec(A(:,k)k=1,...,n)engendr´epar
lescolonnesdeA.
avecC∈Rm×n es I LenoyauKer(f)def estlesousespaceorthogonalauxlignesA(k,:)deA.
Attentionleproduitmatricieln’estpas gn A C
commutatif(AB6=BA). mli I `aCavsalpeaurrtirc´euellileer.ou`A=v> estunvecteurtranspos´e,f :Rn→Restunefonction
Multiplication matrice vecteur f(x)=v>x=hv,xi= vkxk
SiB=bestunvecteurdetaillepalors Xk
Lavaleurdelafonctionestunesommedescomposantex dexpond´er´eesparv.
k
p
c=Ab c = a b I Casparticulierou` m=netleslignesdeAd´efinissentunebaseorthonormalede
⇔ i kX=1 i,k k Conseil:faireledessin. Rn :l’applicationlin´eaireestunchangementdebase.
9/53 10/53
Matrices particuli`eres Exercice 1 : Multiplication matricielle
Soitlesvecteursetmatricesuivants:
Matrice Diagonale (m=n) Matrice triangulaire
1 1 1 1 0 1 2
−
a1,1 0 ··· 0 a1,1 0 ··· 0 u=2, v=1, A=−1 1 −1, B=3 4
 0... a2...,2 ·.·.·. 0...  a2...,1 a2...,2 ·.·.·. 0...  Calculerlesprodu3itssuivants1:  0 −1 1  5 6
   
 0 0 ··· am,n am,1 am,2 ··· am,n u>u= u>v= u>Av=
I Matricediagonale: I Matricetriangulaireinf´erieure:
i6=j⇒aij =0, ∀i,j i<j⇒aij =0, ∀i,j uu>= uv>= AA=
I Matricescalairesiaii=C, i. I A> esttriangulairesup´erieure.
∀
I Matriceidentit´esiaii=1, i. I Strictementtriangulaireinf´erieuresi
∀ Au= Av= B>u= AB=
i j a =0, i,j
ij
≤ ⇒ ∀
11/53 12/53
Applications de l’alg`ebre lin´eaire Applications de l’alg`ebre lin´eaire
Pr´ediction lin´eaire
I Soitlafonctionlin´eaire
Infographie
f(x)=v>x= vkxk
I L’immensemajorit´edesmod`eles3D Xk
utilis´esdanslesjeuxvid´eosont
I xestunvecteurd´ecrivantunobjetet
repr´esent´espar
v leparam`etredumod`ele.
I VecteursdepositiondansR3.
I Pr´edictiondevaleurcontinue.
I Matriced’adjacence.
I Lesalgorithmesd’illuminationfont
Exemples d’utilisation
souventappelauproduitscalaire
entreunvecteurnormal`alasurface I M´et´eo(temp´erature,humidit´e,
etladirectiond’unesourcelumineuse directionduvent`apartirdemesures).
ponctuelle.
I Moteurderecherche(score`apartir
I Lesrotationsetd´eplacementsdes derequˆetewebetpagerank).
membresdesmodels3Dsont
I D´etectiondeSPAM(score`apartir
mod´elis´espardesproduitsmatriciels.
d’email).
13/53 14/53
Applications de l’alg`ebre lin´eaire Repr´esentation m´emoire en langage C
Factorisation de matrice
I Oncherche`afactoriserunematrice Repr´esentation m´emoire
A∈Rm×n souslaforme: I Ilexistedemultiplesimpl´ementationspossibles(pointeurs,tableaux,...).
A UV> I Danscecoursnousnousconcentreronssurl’utilisationdepointeurs.
≈
avecU Rm×p,V Rn×p et I Pourlesmatricesonveutunstockagem´emoirecontigupourpermettreune
∈ ∈ interfaceavecBLAS(interfacestandardpourlesop´erationsd’alg`ebrelin´eaire).
p m,p n
(cid:28) (cid:28)
I Apeutnepasˆetrecompl`etement
observ´ee. Principe g´en´eral
I D´eclarationdesvecteursetmatrices(pointeurs).
Exemples d’utilisation
I Allocationm´emoiredynamiqueaud´ebutduprogramme(malloc).
I RecommandationsAmazon. I Calculseffectu´essurlesvecteursmatricesallou´es.
I RecommandationsNetflix. I Lib´erationm´emoire`alafinduprogramme(free).
I Clusteringder´eseauxsociaux.
I R´eductiondedimension.
15/53 16/53
BLAS (Basic Linear Algebra Subprograms) BLAS et complexit´e
Lesfonctionsd´efiniesdansBLASsontregroup´eesparniveauxdecomplexit´e
Description
algorithmique.
I BLASestunesp´ecificationd´ecrivantlesdiversesop´erationsdebaseenalg`ebre
Niveau 1 (Complexit´e lin´eaire)
lin´eaire.
I Vecteurtaillen.
I L’originedelasp´ecificationestunebiblioth`equeFortranpropos´eeen1979.
I T(n)=O(n),S(n)=O(n)
I Uneimpl´ementationder´ef´erenceenFortranestdisponiblesur
http://www.netlib.org/blas/. I Exemple:produitscalaire,multiplicationparunscalaire,additiondevecteurs.
I UneinterfaceCappel´eecblasest´egalementdisponible.
Niveau 2 (Complexit´e quadratique)
I Typesdefonctionss´epar´eesparleurcomplexit´e(m´emoire/temporelle).
I Vecteurtaillen,Matricedetaillem n
×
Importance de la compatibilit´e I T(n)=O(n2)ouO(mn),S(n)=O(n2)ouO(mn)
I Denombreusesbiblioth`equesfournissentune I Exemple:produitdyadique,multiplicationmatrice/vecteur.
impl´ementationBLASextrˆemementoptimis´ee.
I Choixdelabiblioth`equeauniveausyst`eme Niveau 3 (Complexit´e cubique)
sousLinux. I Matricedetaillen n
×
I Exemples:IntelMKL,Atlas,OpenBLAS. I T(n)=O(n3),S(n)=O(n2)
I Exemple:produitmatriciel.
Sourcehttp://blog.nguyenvq.com/
17/53 18/53
Convention des fonctions BLAS Vecteurs en C
I Lesfonctionsetproc´eduresBLASsontdesfonctionsFortranenmajuscule.
I Nomtypique:xNOMFONCTIONou` xdonneletypedesvecteurs: x x[0] x[1] x[2] x[3] x[4] x[5] x[6] x[7] x[8] x[9]
S Simplepr´ecision(float)
x x x x x x x x x x
D Doublepr´ecision(double) 1 2 3 4 5 6 7 8 9 10
C Complexesimplepr´ecision
Z Complexedoublepr´ecision Principe
I L’interfacecblasfournitdesfonctionsCappelantdirectementlesfonctions
I Lesvecteurssontrepr´esent´espardespointeursversdesflottants.
BLAS.
I D´eclaration:double *x;
I Nomdesfonctionscblas:cblas xnomfonction
I Initialisationm´emoire:x=malloc(n*sizeof(double));(callocpourinit.`a0)
I Lepassagedeparam`etresefaittoujoursparr´ef´erencepourlesmatriceset
vecteurs(ondonnelespointeurs) I Utilisation:x[i]=10.;ou*(x+i)=10.;
I Lespointeursdoiventadresserdesespacesm´emoireallou´es. I Lib´erationm´emoire:free(x);
I DanslecoursnousdonneronslenomdesfonctionsBLASassoci´eesaux I Attention`alib´ererlam´emoireet`anepasfairedemallocdanslaboucled’une
op´erationsdebaseenalg`ebrelin´eaire. proc´edureit´erative.
19/53 20/53
1 $./ex_vector
2 v1=[ 0.00 0.00 1.00 0.00 0.00
0.00 0.00 0.00 0.00 0.00 ]
3 v2=[ 0.00 1.00 2.00 1.00 4.00
1.00 6.00 1.00 8.00 1.00 ]
(n) 4 Modification de v2:
O 5 v2=[ 0.00 3.00 6.00 3.00 12.00
(n) 3.00 18.00 3.00 24.00 3.00 ]
O
(n)
O
(n)
O
(mn)
O
(mn)
O
(1)
O
(mn)
O
Fonctions utiles pour les vecteurs Exercice 2 : Initialisation et impression
Quelleestlasortieimprim´eedansleterminalpourleprogrammesuivant:
Code Allocation et initialisation Code source Sortie
1 double *vinit(int n,double alpha) 1 int n=10;
2 { I v←α1n 2 double *v1;
3 double * res= (double *)malloc(n* I Allocationetmise`aαdes 3 double *v2;
sizeof(double)); composantesd’unvecteur. 4 double alpha=3.;
4 for (int i=0;i<n;i++) 5
5 res[i]=alpha; I Complexit´eT(): 6 v1=vinit(n,0);
6 return res; 7 v1[2]=1;
I Complexit´eS():
7 } 8
9 v2=vinit(n,1);
Code Impression ´ecran 10 for (int i=0;i<n;i+=2)
1 void vprint(int n, double *v) 11 v2[i]*=i;
2 { I Imprimeunvecteurdansleterminal 12
3 printf("["); (suruneligne). 13 printf("v1=");vprint(n,v1);
4 for (int i=0;i<n;i++) 14 printf("v2=");vprint(n,v2);
5 printf(" %2.2f ",v[i]); I Complexit´eT(): 15
6 printf("]\n"); I Complexit´eS(): 16 printf("Modification de v2:\n");
7 } 17 for (int i=0;i<n;i++)
18 v2[i]*=alpha;
19 printf("v2=");vprint(n,v2);
21/53 22/53
Matrices en C Allocation et lib´eration de m´emoire
Code Allocation et initialisation
a a[0] a[1] a[2]
1 double **minit(int m, int n,double
alpha) I A←α1m1>n
2 { I Allocationetmise`aαdes
3 double ** res= (double **)malloc(m* composantesd’unematrice.
a[0][0]a[0][1]a[0][2]a[0][3]a[1][0]a[1][1]a[1][2]a[1][3]a[2][0]a[2][1]a[2][2]a[2][3] sizeof(double *));
a a a a a a a a a a a a 4 double * data=(double *)malloc(m*n* I Complexit´eT():
1,1 1,2 1,3 1,4 2,1 2,2 2,3 2,4 3,1 3,2 3,3 3,4 sizeof(double));
I Complexit´eS():
5 for (int i=0;i<m*n;i++)
Principe 6 data[i]=alpha;
7 for (int i=0;i<m;i++)
I Lesmatricessontrepr´esent´eesparuntableaudepointeurs(verschaquelignede 8 res[i]=data+i*n;
9 return res;
lamatrice). 10 }
I D´eclaration:double **a;
Code Lib´eration m´emoire
I Initialisationm´emoire:2malloc(unpourletableaua,unpourlesdonn´ees) 1
I Utilisation:a[i][j]=10.; 23 v{oid mfree(double **a) I Ltaibb´eleraautiodnedpeoilnatemu´erms)o.ire(donn´eeset
I Lib´erationm´emoire(2free):free(&a[0][0]);free(a); 4 free(&a[0][0]); I Complexit´eT():
5 free(a);
I Uneespacem´emoirecontiguestn´ecessairepouruneinterfaceavecBLAS. 6 } I Complexit´eS():
23/53 24/53
Matrices creuses Matrice creuse (2)
Exemple sous Matlab
Principe
Op´eration Dense Creux
I n=10000
I Stockagem´emoireclassiqueest dense . v>v 0.0001 0.0005
(cid:28) (cid:29) I nv =100 Av 0.0525 0.0007
I Lesmatricesetvecteurscreux(parcimonieux)ontunpetitnombrenv (cid:28)nde I Tempsdecalcul(ensec). AA 19.2145 0.0397
valeursdiff´erentesde0.
I Ilestplusefficacedestockerenm´emoireunelistedespositionsetlesvaleurs
correspondantes. Exemples d’utilisation
I Vecteurv detaillenavecnv valeursdiff´erentesdez´erosformatcompress´e:
is Nnv etvs Rnv
∈ ∈
I MatriceAdetaillem navecnv valeursdiff´erentesdez´erosformatcompress´e:
×
is Nnv, js Nnv etvs Rnv
∈ ∈ ∈
I Larepr´esentationci-dessuss’appellelarepr´esentationparcoordonn´ee.
I Diff´erentstypesderepr´esentation(Yalesparsematrix,CompressedSparse
I Donn´eesweb(pageweb,profil).
Column)sontplusoumoinsefficacesselonlesop´erations`aeffectuer(op´eration
matricielle,remplissage). I Graphederelations(r´eseauxsociaux).
I Objets3Deninfographie.
25/53 I Simulationsnum´eriques(´el´ements 26/53
finis).
Op´erations de base sur des vecteurs Op´erations sur des vecteurs en C
I Convolutiondiscr`ete.
Multiplication par un scalaire (xSCAL) Code
I Multiplicationparunscalaire(xSCAL) 1 void vscale(int n, double *v, double
I v αv alpha)
←
v←αv I FonctionBLASdeNiveau1: 23 { for (int i=0;i<n;i++)
I Sommeetsoustraction void cblas_dscal(int n,double alpha 4 v[i]=alpha*v[i];
, double *v,int incx) 5 }
v x+y
← Appel:cblas_dscal(n,alpha,v,1);
I Produitscalaire(xDOT)
s x>y Source:http://www.cise.ufl.edu/research/sparse/matrices/
← Produit scalaire (xDOT) Code
I Casg´en´eral(xAXPY) 1 double vdot(int n , double *x,double
y αx+y I s x>y *y)
← ←
I Complexit´edetoutescesop´erations: I FonctionBLASdeNiveau1: 23 { double res=0;
double cblas_ddot(int n,double *x, 4 for (int i=0;i<n;i++)
T(n)=O(n), S(n)=O(n) int incx,double *y,int incy); 5 res+=x[i]*y[i];
Appel: 6 return res;
7 }
double s=cblas_ddot(n,x,1,y,1);
27/53 28/53
1 void vaxpy(int n, double alpha ,
double *x, double *y)
2 {
3 for (int i=0;i<n;i++)
4 y[i]+=alpha*x[i];
5 }
1 void mgemv(int m, int n, double
alpha, double **A, double *x
, double beta, double *y)
2 {
3 for (int i=0;i<m;i++)
4 {
5 y[i]=beta*y[i];
6 for (int j=0;j<n;j++)
7 y[i]+=alpha*A[i][j]*x[j];
8 }
9 }
vaxpy(n,x,1,y);
vaxpy(n,x,-1,y);
Solution
y=vinit(n,0);vaxpy(n,x,alpha,y);
T(n)=O(mn), S(n)=O(mn)
Solution
Exercice 3 : xAXPY Op´erations matrice/vecteur
Code Exercice
Mise `a jour de rang 1 (xGER)
1 void vaxpy(int n ,double alpha,
double *x, double *y) I Compl´eterlafonctionxAXPY`agauche.
2 { I Onsuposend´eclar´eetlesvecteursxety A←αxy>+A
34 initialis´es. I Avecα∈R,y∈Rn,x∈Rm etA∈Rm×n.
5 I Montrerquecettefonctionpermet I Complexit´e:T(n)=O(mn), S(n)=O(mn)
6 d’effectuerdessommes,soustractionset
7 multiplicationparunscalaire: I Casparticulier:Asym´etriqueetx=y (xSYR).
8
9 I Sommey x+y: Produit matrice vecteur g´en´eralis´e (xGEMV)
10 } I Soustractio←ny y x:
← −
y αAx+βy
←
I Mult.scalairey αx:
← I Avecα,β R2,y Rm,x Rn etA Rm×n.
∈ ∈ ∈ ∈
I Complexit´e:
I Casparticuliers:
I Asym´etrique(xSYMV).
I Abande(xSBMV).
I Atriangulaire(xSTRMV).
FonctionBLASdeNiveau1:
void cblas_daxpy(int n,double alpha
, double *x,int incx,double *y,int 29/53 30/53
incy);
Mise `a jour de rang 1 (xGER) Exercice 4 : Produit g´en´eralis´e (xGEMV)
Appel:
Code Exercice
cblas_daxpy(n,alpha,x,1,y,1);
1 void mgemv(int m, int n, double
Impl´ementation Code alpha, double **A, double *x y αAx+βy
, double beta, double * y) ←
A←αxy>+A 1 void,*mgdyeo)ru(bilnet*m*,A,indtounb,ledo*uxb,ledoaulbplhea 23 { I αA,β∈RmR×2,n.y∈Rm,x∈Rn et
2 { 4 ∈
II αAFo∈n∈cRtRi,omny×B∈nL.RAnS,dxe∈NRivmeauet1: 3456 } forf(oiArn[ti(]ii[n=jt0];+ji==<a0ml;;pjih<+an+*;)xj[+i+])*y[j]; 56789 II CFvCoBooiLnmdAcSptc_ilo´ebTntlReaArBsN_SlLaPdAOgfSoSeEmndvce(TtiCrNoBaniLnvAsxeSAGa_,uEORiM1DntEV:Rm`ao,rgdiaenurtc,hne,.
void cblas_dger( CBLAS_ORDER order, 10 double alpha, double *A, int lda,
int m, int n, double alpha, double 11 double *x, int incx, double beta,
*x, int incx, double *y, int incy, 12 double *y, int incy);
double *A, int lda); 13 }
Appel:
Appel:
cblas_dgemv(CblasRowMajor,CblasNoTrans
cblas_dger(CblasRowMajor,m,n,alpha, ,m,n,alpha,&A[0][0],n,x,1,beta,y,1);
x,1,y,1,&A[0][0],n);
31/53 32/53
Op´eration matricielle g´en´eralis´ee (xGEMM) Op´eration matricielle g´en´eralis´ee (xGEMM)
Impl´ementation Code
Op´eration 1 void mgemm(int m, int n, int p,
double alpha, double **A, double
C αAB+βC
C αAB+βC ← **B, double beta, double **C)
← 2 {
I α,β∈R2,C∈Rm×n,A∈Rm×p,B∈Rp×n I αB,β∈∈RpR×2n,C∈Rm×n,A∈Rm×p, 345 ff{oorr ((iinntt ij==00;;ij<<mn;;ij++++))
I Complexit´e:T(n)=O(mnp), S(n)=O(mn) I FonctionBLASdeNiveau1: 6 C[i][j]=beta*C[i][j];
void cblas_dgemm(CBLAS_ORDER Order, 7 for (int k=0;k<p;k++)
8 C[i][j]+=alpha*A[i][k]*B[k][j];
Cas particuliers CBLAS_TRANSPOSE TransA, 9 }
CBLAS_TRANSPOSE TransB, int m, int 10 }
I Asym´etrique(xSYMM) n, int p, double alpha, double *A,
I Ahermitienne(xHEMM) int lda, double *B, int ldb, double
beta, double *C, int ldc);
I Atriangulaire(xTRMM)
Appel:
I B=A> (xSYRK) cblas_dgemm(CblasRowMajor,
I C αAB>+α¯BA>+βC(mise`ajoursym´etriquexSYRK2) CblasNoTrans,CblasNoTrans,m,n,p,
←
alpha,&A[0][0],p,&B[0][0],n,beta,&C
[0][0],n);
33/53 34/53
Syst`eme d’´equations lin´eaires Solutions des syst`emes d’´equations lin´eaires
a x +a x +...+a x =b Existence et unicit´e des solutions Exemple n=3
1,1 1 1,2 2 1,n n 1
a x +a x +...+a x =b D´ependdem,netdurangdeA:
2,1 1 2,2 2 2,n n 2
 .
 .... ⇔ Ax = b II mA<denr,anilge<xisnteilunexeisintefinuint´eedinefisnoitlu´etdioens.olution.
.
avecA Rma×mn,1ext1b+amR,m2x.2+...+am,nxn=bm II mA=denraentgA>dne,rialnng’enxi,stileepxaisstedeunsoelusotilountionunique
∈ ∈ solutiondesmoindrescarr´es.
→
I Oncherchelevecteurxdetaillensatisfaisant
lesm´equationscidessus.
R´esolution du syst`eme
I L’ensembledessolutionsd´ependdespropri´et´es
I Casparticulierou` m=netAderangn
delamatriceA.
I Solutionunique:
I Lar´esolutiondusyst`emepeutˆetre
particuli`erementsimpledanscertainscas x?=A−1b
particuliers(Adiagonale,triangulaire). I E´viterl’inversionmatricielledirected`esquepossible
(matricesparticuli`eres).
35/53 36/53
O(n)
O(n2)ouO(n) 1
1
 
1
 
LAPACK (Linear Algebra PACKage) LAPACK et interface C
Convention de nommage
I NomdesfonctionsLAPACK:XMMAAA
Description
X typededonn´ees(mˆemequepourBLAS).
I LAPACKestunebiblioth`equeFortranqui MM typedematrice(GEg´en´erale,GBbande,SYsym´etrique,TR
fournitdesfonctionsdehautniveauenalg`ebre triangulaire,DIdiagonale...
lin´eaire. AAA op´erationdetaille2-3caract`eres.
I Biblioth`equed’accompagnementdeBLAS I Typed’op´erationsAAAsupport´ees:
utilisantlesroutinesBLASlepluspossible. I R´esolutiondesyst`emeslin´eaires(SV).
I InterfaceCofficielle:LAPACKE. I Moindrescarr´es(LS).
I D´ecompositionenvaleurs/vecteurspropresEIG(TRD).
Fonctions propos´ees (Niveau 3) I D´ecompositionLU(TRF).
I D´ecompositionenvaleursinguli`ereSVD(BRD).
I R´esolutiondesyst`emeslin´eaires(casg´en´eral).
I Estimationdesmoindrescarr´es.
LAPACKE
I Factorisationdematrice(EIG,SVD,LU,QR).
I Biblioth`equeofficielled’interfaceC.
I Nomdesfonctions:LAPACKE xmmaa.
I Lestableauxd´efinisdanslecourssontdetype:LAPACK ROW MAJOR
37/53 38/53
Matrice A diagonale Matrice triangulaire
R´esolution du syst`eme Matrice
R´esolution du syst`eme
l 0 0
I SoitA=Lunematricetriangulaireinf´erieure. 1,1 ···
II nSo´eluqtuioatnison´esviidnedn´etpeesndadxekak,nlkat=xefksoar=dbmekbelkaforme a100...,1 a200...,2 ···.···.···. an00...,n I Lcoemsypsots`eamnetelxipn1a´e=raicroelb1m1,p1peuotsaˆenttree:r´esolu SyLst`e=melln2...,,11 lln2...,,22 ··.··.··. ln0...,n
k,k b l x
I Siai,i=0etbi6=0lesyst`emen’admetpasdesolution. x2= 2−l2,22,1 1 ll1,1xx1+l x ==bb1
2,1 1 2,2 2 2
I Siai,i=0etbi=0lesyst`emeadmetuneinfinit´edesolutions. xk = bk− lki=−11lk,ixi ...
Pk,k l x + +l x =b
k,1 1 ··· k,k k k
Code Discussion
1 void diagsv(int n,double **A,double Exercice 5 : Solution
*b,double *x) I Complexit´eT(n): Trouverlasolutionxpourlesyst`emesuivant
2 { I Complexit´eS(n):
3 for (int i=0;i<n;i++) 2 0 0 2 x=
4 x[i]=b[i]/A[i][i]; L= 4 1 0 b= 3
5 }  6 −1 2  7
− − −
   
39/53 40/53