[Défis][2] Pourquoi les matrices et les vecteurs ne s'entendent pas mieux ?

C'est bien dommage :aie:

Prenons le cas d'une matrice A et d'un vecteur B :

Code:

---

1 2 3

A = [1 2 ; 3 4 ; 5 6 ; 7 8]   B = [2 4]

---

Code:

---

1 2 3 4 5 6 7 8 9 10 11

A =   1 2 3 4 5 6 7 8     B =   2 4

---

Il est relativement fréquent que l'on soit amené à effectuer une opération entre la première colonne de A ([1 3 5 7]) et la première valeur de B (2) et la seconde colonne de A ([2 4 6 8]) et la seconde valeur de B (4).

Par exemple, une simple soustraction, A - B qui renverrait :

Code:

---

1 2 3 4 5 6

ans =   -1 -2 1 0 3 2 5 4

---

Malheureusement si l'on fait directement la soustraction entre la matrice et le vecteur, MATLAB ne veut rien savoir :

Code:

---

1 2 3

>> A-B ??? Error using ==> minus Matrix dimensions must agree.

---

Il faut donc utiliser une autre méthode : boucle FOR-END, code vectorisé, fonction spéciale...

------------------------------------------------------

Défi : Trouver la méthode la plus rapide pour effectuer une soustraction entre les colonnes d'une matrice et un vecteur


Contrainte : le code devra être contenu dans un fichier m fonction dont le nom sera votre pseudo (si votre pseudo n'est pas accepté comme nom de fonction, prenez un nom similaire :aie:). Cette fonction prendra en argument d'entrée une matrice A et un vecteur B et retournera en argument de sortie une matrice C résultat.

Par exemple, si je (Dut) propose une solution, la première ligne du fichier sera :

Code:

---

function C = Dut(A,B)

---

Notes :

• les noms de variable A,B,C sont donnés à titre indicatif et peuvent être modifiés
• rien n'empêche de faire appel à d'autres fichiers à l'intérieur de matvec.m

Alors, quelle serait la méthode la plus rapide sous MATLAB selon vous ?

A vous de jouer... donnez vos solutions à la suite de ce message:king:

------------------------------------------------------


Les solutions proposées :

mr_samurai :

Code:

---

1 2 3

function C = mr_samurai(A, B)   C = A - ones(size(A,1),1)*B;

---

vinc-mai :

Code:

---

1 2 3

function C = vinc_mai(A,B)   C = A - repmat(B,size(A,1),1);

---

Caro-line :

Code:

---

1 2 3 4 5 6 7

function C = caroline(A, B)   [m,n] = size(A); C = zeros(size(A)); for i=1:m C(i,:) = A(i,:)-B; end

---

tug83

Code:

---

1 2 3

function C = tug83(A,B)   C = bsxfun(@minus, A, B);

---

Dut (voir le fichier MEX ici) :

Code:

---

1 2 3

function C = Dut(A,B);   C = Dutmex(A,B);

---

------------------------------------------------------


Code utilisé pour la comparaison :

Code:

---

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

function [T,users,N] = microDefi2(num)   if nargin == 0 num = 7; end   % Taille croissante des données N = 10.^(0:num-1);   % Données initiales iA = [1 2 ; 3 4 ; 5 6 ; 7 8]; iB = [2 4];   % Membres participant au défi users = {'mr_samurai' 'vinc-mai' 'caroline' 'tug83'}; nusers = numel(users);   mt = 5; nT = 1;   for n=1:numel(N)   % Augmentation de la taille des données au fur et à mesure A = repmat(iA,1,N(n)); B = repmat(iB,1,N(n));   %% Solution de mr_samurai for k=1:mt tic C = mr_samurai(A,B); t(k) = toc; clear C end T(nT,1) = mean(t);   %% Solution de vinc-mai for k=1:mt tic C = vinc_mai(A,B); t(k) = toc; clear C end T(nT,2) = mean(t);   %% Solution de Caro-line for k=1:mt tic C = caroline(A,B); t(k) = toc; clear C end T(nT,3) = mean(t);   %% Solution de tug83 if exist('bsxfun','builtin')==5 % Cette fonction n'est pas disponible avec toutes les versions   for k=1:mt tic C = tug83(A,B); t(k) = toc; clear C end T(nT,4) = mean(t); else T(nT,4) = nan; end   %% Fin des solutions proposées nT = nT +1;   end

---

------------------------------------------------------


Temps d'exécution :

MATLAB R2008a - Windows XP - 2GHz - 2Go Ram

Code:

---

1 2 3 4 5 6 7 8

N mr_samurai vinc-mai caroline tug83 Dut 1 0.000848 0.000932 0.001002 0.000739 0.000249 10 0.000036 0.000076 0.000020 0.000045 0.000012 100 0.000036 0.000081 0.000018 0.000055 0.000019 1000 0.000107 0.000126 0.000053 0.000164 0.000106 10000 0.001491 0.001802 0.000913 0.001691 0.001289 100000 0.021650 0.033351 0.033862 0.018619 0.017151 1000000 0.242351 0.342116 0.294612 0.189642 0.179889

---

MATLAB R2007a - Windows Vista SP1 - 2GHz - 2Go Ram

Code:

---

1 2 3 4 5 6 7 8

N mr_samurai vinc-mai caroline tug83 1 0.166230 0.000183 0.003227 0.002319 10 0.001539 0.000057 0.000020 0.000032 100 0.000035 0.000065 0.000019 0.000045 1000 0.000185 0.000277 0.000150 0.000202 10000 0.003575 0.003802 0.003620 0.002428 100000 0.027099 0.030994 0.028084 0.020660 1000000 0.264510 0.301678 0.280330 0.207564

---

MATLAB R2008b - Linux 64bits - Processeur 8 cœurs

Code:

---

1 2 3 4 5 6 7 8

N mr_samurai vinc-mai caroline tug83 Dut 1 0.000021 0.000043 0.000014 0.000025 0.000009 10 0.000013 0.000039 0.000009 0.000022 0.000008 100 0.000017 0.000042 0.000011 0.000029 0.000014 1000 0.000080 0.000079 0.000042 0.000098 0.000074 10000 0.001122 0.000519 0.000773 0.001603 0.000893 100000 0.017785 0.023169 0.025011 0.011389 0.012374 1000000 0.185730 0.246647 0.260227 0.115151 0.127354

---

Matlab 2008a - Linux 32 bits - 3.2GHz - 2Go

Code:

---

1 2 3 4 5 6 7 8

N mr_samurai vinc-mai caroline tug83 Dut 1 0.000078 0.000135 0.000062 0.000084 0.000327 10 0.000036 0.000131 0.000026 0.000062 0.000052 100 0.000050 0.000142 0.000039 0.000082 0.000073 1000 0.000222 0.000157 0.000055 0.000193 0.000163 10000 0.002153 0.001704 0.000903 0.001920 0.001592 100000 0.023430 0.027308 0.028055 0.019134 0.018004 1000000 0.247297 0.259067 0.275122 0.186231 0.174071

---

• N est l'ordre de grandeur du nombre de colonne de A (A possède 4*N colonnes)
• les résultats sont donnés en secondes

------------------------------------------------------

Salut,

J'ouvre le bal :D .

Code:

---

1 2 3 4 5 6 7 8 9 10 11

function C = mr_samurai(A, B) % Micro-Défi 2 %   if nargin == 0 A = rand(50000,100); B = rand(1,100); end   C = A - ones(size(A,1),1) * B;

---

Plus rapide que REPMAT sur ma machine 8-) .

++

Salut.
J'allais utiliser repmat.

Citation:

---

Envoyé par mr_samurai

Plus rapide que REPMAT sur ma machine 8-) .

---

Je passe mon tour alors.

Citation:

---

Envoyé par vinc-mai

J'allais utiliser repmat.

---

Pas sûr que REPMAT soit plus lent :aie: .

Citation:

---

Envoyé par vinc-mai

Salut.
J'allais utiliser repmat.


Je passe mon tour alors.

---

Faut pas baisser les bras :calin:
En plus la rapidité d'exécution dépendra aussi de la taille des matrices sur lesquelles on travaille : une méthode pourra être bien pour des matrices de grande taille alors qu'une autre sera plus efficace sur des matrices de petite taille.

L'énoncé ne donnant pas la taille des matrices....:aie:

Donc moi je me lance (en piquant le début du samurai) :

Code:

---

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

function C = caroline(A, B) % Micro-Défi 2 % %doit retourner C tq C(i,j) = A(i,j) - B(j)   %Gestion des arguments if nargin == 0 N=50000; M=100; A = rand(N,M); B = rand(1,M); end   %Vérifier que A et B ont le même nombre de colonnes [m,n] = size(A); [p,q] = size(B); if p~=1 || n~=q error('pas bon'); return end   %Note : sans la pré-allocation %Pour N=500 t=0.1 alors qu'avec t=0.001 %Je n'ai pas eu le courage d'attendre pour N=50000 ! C=zeros(size(A)); for i=1:m C(i,:) = A(i,:)-B; end

---

Et dans ce cas en bouclant 500 fois pour faire une moyenne du temps je trouve :

Citation:

---

mr_samurai : 0.1687
caroline : 0.1418

---

:yaisse2:
Mais je ne suis pas sure que ce soit toujours vrai suivant les entrées qu'on prendra.

Voici mon script de test :

Code:

---

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

clc clear   A = rand(50000,100); B = rand(1,100); nb = 100;   %---- S1 tic for z=1:nb C = A - repmat(B,size(A,1),1); end toc clear C   %---- S2 tic for z=1:nb C = zeros(size(A)); for u=1:size(A,1) C(u,:) = A(u,:) - B; end end toc clear C   %---- S3 tic for z=1:nb C = A - ones(size(A,1),1) * B; end toc clear C

---

Ce qui donne :

Code:

---

1 2 3

Elapsed time is 10.358000 seconds. % REMAPT Elapsed time is 10.407000 seconds. % BOUCLE FOR Elapsed time is 8.429000 seconds. % ONES * B

---

A Dut de trancher :aie: .

++

Ben tu dois avoir une bête de course, avec le même code j'obtiens :

Code:

---

1 2 3

Elapsed time is 15.810041 seconds.%REPMAT Elapsed time is 13.523129 seconds.%BOUCLE FOR Elapsed time is 14.859095 seconds.%ONES

---

Et une 2ème fois :

Code:

---

1 2 3

Elapsed time is 12.789178 seconds. Elapsed time is 13.511811 seconds. Elapsed time is 15.273793 seconds.

---

:roll:

Citation:

---

Envoyé par mr_samurai

A Dut de trancher :aie: .

---

Pfff...Pô juste.

Voici les temps que j'obtiens

Code:

---

1 2 3

Elapsed time is 14.831326 seconds.%Repmat Elapsed time is 29.354350 seconds.%Boucle for Elapsed time is 13.074401 seconds.%ones

---

Code:

---

1 2 3

Elapsed time is 14.202258 seconds. Elapsed time is 28.605043 seconds. Elapsed time is 12.782492 seconds.

---

ones semble être la méthode la plus rapide pour l'instant.
La boucle for est vraiment lente chez moi!

Citation:

---

Envoyé par vinc-mai

ones semble être la méthode la plus rapide pour l'instant.
La boucle for est vraiment lente chez moi!

---

Comme toujours le monde entier est contre moi :cry:

Non par contre il serait intéressant que tu nous donnes :
- ta version de MATLAB (Samurai et moi n'avons déjà pas la même)
- ton OS (là on a le même)

Et puis bon pour le faire proprement il faudrait tout fermer à côté mais personnellement j'ai la flemme donc on va laisser le Responsable juger :aie:

Oups, j'ai oublié de préciser ma version:
R2008.a sous linux (2.6.24).

Allez moi aussi je me lance , je vais utiliser bsxfun (une feature dispo depuis la R2007a):

Code:

---

1 2 3 4 5 6 7 8

A = [1 2 ; 3 4 ; 5 6 ; 7 8]; B = [2 4]; C = bsxfun(@minus, A, B) C = -1 -2 1 0 3 2 5 4

---

dois je mettre un tic / toc? :D

Code:

---

1 2 3 4 5 6 7 8 9 10

%---- S4 C = zeros(size(A)); tic for z=1:nb for u=1:size(A,2) C(:, u) = A(:, u) - B(u); end end toc

---

[Edit = Dut] Cette solution est déjà donnée par mr_samurai :?

solution utilisant les coordinees homogenes... lente, juste postee a titre informatif ;)


Juste un commentaire en ce qui concerne la solution avec repmat, il existe des version mex plus rapides en general que celle de matlab, voir par exemple http://research.microsoft.com/~minka...re/lightspeed/

Code:

---

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16

%---- S5   tic for z=1:nb % A en coordinees homogenes hA = [A ones(size(A,1), 1)]; % B en form de matrice hB = diag(ones(1, numel(B)+1)); hB(end,1:end-1) = -B; % translation hC = hA * hB; % retour aux coordinees usuelles hC = hC(1:end, 1:end-1); end toc

---

Le premier message a été mis à jour avec les temps d'exécution pour des matrices de différentes tailles.

Si certains pouvaient faire tourner la fonction microDefi2.m et me retourner la variable T dans un fichier mat avec la configuration de l'ordinateur (OS, processeur), on pourrait faire un semblant de comparatif par machine ;)

A utiliser comme ceci :

Code:

---

T = microDefi2;

---

Si une erreur mémoire survient, essayez comme ceci :

Code:

---

1 2	num = 6; T = microDefi2(num);

---

En faisant décroitre num (entier) jusqu'à ce que ça marche.

Pensez à nettoyer MATLAB avant d'exécuter le code :

Code:

---

clear all

---

Pour ceux qui utilisent une version récente, pensez à mentionner si vous utilisez le multithreading ;)

Citation:

---

Envoyé par paradize3

il existe des version mex plus rapides en general

---

Voila une idée que personne ne semble avoir creusé... sans parler d'utiliser un MEX de repmat, je rappelle que les structures itératives peuvent aussi être intégrées dans un fichier MEX...

Quelqu'un veut essayer ?

Voilà.

MATLAB R2007a - Windows Vista SP1 - 2GHz - 2Go Ram

[Edit=Dut] :merci:

Citation:

---

Envoyé par paradize3

Code:

---

1 2 3 4 5 6 7 8 9 10

%---- S4 C = zeros(size(A)); tic for z=1:nb for u=1:size(A,2) C(:, u) = A(:, u) - B(u); end end toc

---

[Edit = Dut] Cette solution est déjà donnée par mr_samurai :?

---

pas tout a fait!! compare precisemment, et teste :)
je pense que ma version devrait etre plus rapide (dans le cas ou le nombre de lignes est plus grand que le nombre de colonnes), de plus elle est plus efficace au niveau du cache CPU..

Salutations,

Greg

Citation:

---

Envoyé par paradize3

plus rapide dans le cas ou le nombre de lignes est plus grand que le nombre de colonnes

---

ce qui est l'oppose de ton test ;) mais c'est pas grave lol

Sous Linux 64-bit , 8 cores , sous R2008b
J'ai zippé le MAT-file

C'est vrai que pour les processeurs multicoeurs, il peut être intéressant de faire l'essai avec le multithread aussi... ;)