Select Git revision
test_qr.c 22.84 KiB
/*******************************************************************************
* This file is part of QuickSched.
* Coypright (c) 2013 Pedro Gonnet (pedro.gonnet@durham.ac.uk)
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*
* *****************************************************************************/
/* Config parameters. */
#include "../config.h"
/* Standard includes. */
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <unistd.h>
#include <math.h>
#include <omp.h>
#include <pthread.h>
#include <cblas.h>
/* Local includes. */
#include "quicksched.h"
/**
* Takes a column major matrix, NOT tile major. size is length of a side of the
* matrix. Only works for square matrices.
* This function is simply for validation and is implemented naively as we know
* of no implementation to retrieve Q from the tiled QR.
*/
double* computeQ(double* HR, int size, int tilesize, double* tau, int tauNum) {
double* Q = malloc(sizeof(double) * size * size);
double* Qtemp = malloc(sizeof(double) * size * size);
double* w = malloc(sizeof(double) * size);
double* ww = malloc(sizeof(double) * size * size);
double* temp = malloc(sizeof(double) * size * size);
int i, k, l, j, n;
bzero(Q, sizeof(double) * size * size);
bzero(Qtemp, sizeof(double) * size * size);
bzero(ww, sizeof(double) * size * size);
for (i = 0; i < size; i++) {
Q[i * size + i] = 1.0;
}
int numcoltile = size / tilesize;
int numrowtile = size / tilesize;
for (k = 0; k < numrowtile; k++) {
for (l = 0; l < tilesize; l++) {
bzero(Qtemp, sizeof(double) * size * size);
for (i = 0; i < size; i++) {
Qtemp[i * size + i] = 1.0;
}
for (i = k; i < numcoltile; i++) {
bzero(w, sizeof(double) * size);
for (j = 0; j < tilesize; j++) {
w[i * tilesize + j] =
HR[(k * tilesize + l) * size + i * tilesize + j]; }
w[k * tilesize + l] = 1.0;
if (k * tilesize + l > i * tilesize) {
for (j = 0; j < k * tilesize + l; j++) w[j] = 0.0;
}
/* Compute (I - tau*w*w')' */
for (j = 0; j < size; j++) {
for (n = 0; n < size; n++) {
if (j != n)
ww[n * size + j] =
-tau[(k * tilesize + l) * tauNum + i] * w[j] * w[n];
else
ww[n * size + j] =
1.0 - tau[(k * tilesize + l) * tauNum + i] * w[j] * w[n];
}
}
/* Qtemp = Qtemp * (I-tau*w*w')' */
cblas_dgemm(CblasColMajor, CblasNoTrans, CblasNoTrans, size, size, size,
1.0, Qtemp, size, ww, size, 0.0, temp, size);
double* b = Qtemp;
Qtemp = temp;
temp = b;
}
cblas_dgemm(CblasColMajor, CblasNoTrans, CblasNoTrans, size, size, size,
1.0, Q, size, Qtemp, size, 0.0, temp, size);
double* b = Q;
Q = temp;
temp = b;
}
}
free(Qtemp);
free(w);
free(ww);
free(temp);
return Q;
}
double* getR(double* HR, int size) {
double* R = malloc(sizeof(double) * size * size);
int i, j;
bzero(R, sizeof(double) * size * size);
for (i = 0; i < size; i++) {
for (j = 0; j <= i; j++) {
R[i * size + j] = HR[i * size + j];
}
}
return R;
}
void printMatrix(double* Matrix, int m, int n, int tilesize) {
int i, j;
for (i = 0; i < m * tilesize; i++) {
for (j = 0; j < n * tilesize; j++) {
printf(" %.3f ", Matrix[j * m * tilesize + i]);
}
printf("\n");
}
}
double* columnToTile(double* columnMatrix, int size, int m, int n,
int tilesize) {
double* TileMatrix;
TileMatrix = malloc(sizeof(double) * size);
if (TileMatrix == NULL) error("failed to allocate TileMatrix");
int i, j, k, l;
for (i = 0; i < n; i++) {
for (j = 0; j < m; j++) {
double* tileStart =
&columnMatrix[i * m * tilesize * tilesize + j * tilesize];
double* tilePos =
&TileMatrix[i * m * tilesize * tilesize + j * tilesize * tilesize];
for (k = 0; k < tilesize; k++) {
tileStart = &columnMatrix[i * m * tilesize * tilesize +
k * m * tilesize + j * tilesize];
for (l = 0; l < tilesize; l++) {
tilePos[k * tilesize + l] = tileStart[l];
}
}
}
}
return TileMatrix;
}
double* tileToColumn(double* tileMatrix, int size, int m, int n, int tilesize) {
double* ColumnMatrix;
ColumnMatrix = (double*)malloc(sizeof(double) * size);
if (ColumnMatrix == NULL) error("failed to allocate ColumnMatrix");
int i, j, k, l;
for (i = 0; i < n; i++) {
for (j = 0; j < m; j++) {
/* Tile on ith column is at i*m*32*32.*/
/* Tile on jth is at j*32*32 */
double* tile =
&tileMatrix[i * m * tilesize * tilesize + j * tilesize * tilesize];
/* Column starts at same position as tile. */
/* Row j*32.*/
double* tilePos =
&ColumnMatrix[i * m * tilesize * tilesize + j * tilesize];
for (k = 0; k < tilesize; k++) {
for (l = 0; l < tilesize; l++) {
tilePos[l] = tile[l];
}
/* Next 32 elements are the position of the tile in the next column.*/
tile = &tile[tilesize];
/* Move to the j*32th position in the next column. */
tilePos = &tilePos[tilesize * m];
}
}
}
return ColumnMatrix;
}
/* Routines for the tiled QR decomposition.*/
/**
*
* @brief Computes the QR decomposition of a tile.
*
* @param cornerTile A pointer to the tile for which the decomposition is
* computed.
* @param tilesize The number of elements on a row/column of the tile.
* @param tauMatrix A pointer to the tau Matrix.
* @param k The value of k for the QR decomposition
* @param tauNum The number of tau values stored for each row of the matrix
* (this is equal to the number of tiles on each column).
*
*
*/
void DGEQRF(double* restrict cornerTile, int tileSize,
double* restrict tauMatrix, int k, int tauNum) {
int i, j, n;
double norm = 0.0, sign, u1, tau, z;
double w[tileSize];
/*Find the householder vector for each row. */
for (i = 0; i < tileSize; i++) {
norm = 0.0;
/*Fill w with the vector.*/
for (j = i; j < tileSize; j++) {
/* ith row is i*tileSize, only want elements on diagonal or below.*/
w[j] = cornerTile[i * tileSize + j];
/*Find the norm as well */
norm = norm + w[j] * w[j];
}
if (w[i] >= 0.0)
sign = -1;
else
sign = 1;
norm = sqrt(norm);
u1 = w[i] - sign * norm;
if (u1 != 0.0) {
for (j = i + 1; j < tileSize; j++) w[j] = w[j] / u1;
} else {
for (j = i + 1; j < tileSize; j++) w[j] = 0.0;
}
if (norm != 0.0)
tau = -sign * u1 / norm;
else
tau = 0.0;
/*Store the below diagonal vector */
for (j = i + 1; j < tileSize; j++) {
cornerTile[i * tileSize + j] = w[j];
}
cornerTile[i * tileSize + i] = sign * norm;
w[i] = 1.0;
/* Apply the householder transformation to the rest of the tile, for
* everything to the right of the diagonal. */
for (j = i + 1; j < tileSize; j++) {
/*Compute w'*A_j*/
z = cornerTile[j * tileSize + i];
for (n = i + 1; n < tileSize; n++) {
z = z + cornerTile[j * tileSize + n] * w[n];
}
/* Tile(m,n) = Tile(m,n) - tau*w[n]* w'A_j */
for (n = i; n < tileSize; n++) {
cornerTile[j * tileSize + n] =
cornerTile[j * tileSize + n] - tau * w[n] * z;
}
}
/* Store tau. We're on k*tileSize+ith row. kth column.*/
tauMatrix[(k * tileSize + i) * tauNum + k] = tau;
}
}
/**
*
* @brief Applies the householder factorisation of the corner to the row tile.
*
* @param cornerTile A pointer to the tile for which the householder is stored.
* @param rowTiles The tile to which the householder is applied.
* @param tilesize The number of elements on a row/column of the tile.
* @param tauMatrix A pointer to the tau Matrix.
* @param jj The value of j for the QR decomposition.
* @param kk The value of k for the QR decomposition.
* @param tauNum The number of tau values stored for each row of the matrix
* (this is equal to the number of tiles on each column).
*
* */
void DLARFT(double* restrict cornerTile, double* restrict rowTile, int tileSize,
int jj, int kk, double* restrict tauMatrix, int tauNum) {
int i, j, n;
double z = 0.0;
double w[tileSize];
/* For each row in the corner Tile*/
for (i = 0; i < tileSize; i++) {
/*Get w for row i */
for (j = i; j < tileSize; j++) {
w[j] = cornerTile[i * tileSize + j];
}
w[i] = 1.0;
/* Apply to the row Tile */
for (j = 0; j < tileSize; j++) {
z = 0.0;
/* Below Diagonal!*/
/*Compute w'*A_j*/
for (n = i; n < tileSize; n++) {
z = z + w[n] * rowTile[j * tileSize + n];
}
for (n = i; n < tileSize; n++) {
rowTile[j * tileSize + n] =
rowTile[j * tileSize + n] -
tauMatrix[(kk * tileSize + i) * tauNum + kk] * w[n] * z;
}
}
}
}
/**
*
* @brief Applies the householder factorisation of the corner to the row tile.
*
* @param cornerTile The corner tile for this value of k.
* @param columnTile The tile for which householders are computed.
* @param tilesize The number of elements on a row/column of the tile.
* @param tauMatrix A pointer to the tau Matrix.
* @param ii The value of i for the QR decomposition.
* @param kk The value of k for the QR decomposition.
* @param tauNum The number of tau values stored for each row of the matrix
* (this is equal to the number of tiles on each column).
*
*
*/
void DTSQRF(double* restrict cornerTile, double* restrict columnTile,
int tilesize, int ii, int kk, double* restrict tauMatrix,
int tauNum) {
int i, j, n;
double norm = 0.0, sign, u1, tau, z;
double w[2*tilesize];
/* For each column compute the householder vector. */
for (i = 0; i < tilesize; i++) {
norm = 0.0;
w[i] = cornerTile[i * tilesize + i];
norm = norm + w[i] * w[i];
for (j = i + 1; j < tilesize; j++) {
w[j] = 0.0;
}
for (j = 0; j < tilesize; j++) {
w[tilesize + j] = columnTile[i * tilesize + j];
norm = norm + w[tilesize + j] * w[tilesize + j];
}
norm = sqrt(norm);
if (w[i] >= 0.0)
sign = -1; else
sign = 1;
u1 = w[i] - sign * norm;
if (u1 != 0.0) {
for (j = i + 1; j < 2 * tilesize; j++) {
w[j] = w[j] / u1;
}
} else {
for (j = i + 1; j < 2 * tilesize; j++) w[j] = 0.0;
}
if (norm != 0)
tau = -sign * u1 / norm;
else
tau = 0.0;
/* Apply to each row to the right.*/
for (j = i; j < tilesize; j++) {
/* Find w'*A_j, w is 0s except for first value with upper tile.*/
z = 1.0 * cornerTile[j * tilesize + i];
for (n = 0; n < tilesize; n++) {
z = z + w[tilesize + n] * columnTile[j * tilesize + n];
}
/* Apply to upper tile.*/
cornerTile[j * tilesize + i] =
cornerTile[j * tilesize + i] - tau * 1.0 * z;
for (n = i + 1; n < tilesize; n++) {
cornerTile[j * tilesize + n] =
cornerTile[j * tilesize + n] - tau * w[n] * z;
}
/* Apply to lower tile.*/
for (n = 0; n < tilesize; n++) {
columnTile[j * tilesize + n] =
columnTile[j * tilesize + n] - tau * w[tilesize + n] * z;
}
}
/* Store w*/
for (j = 0; j < tilesize; j++) {
columnTile[i * tilesize + j] = w[tilesize + j];
}
tauMatrix[(kk * tilesize + i) * tauNum + ii] = tau;
}
}
/**
*
* @brief Applies the householder factorisation of the corner to the row tile.
*
* @param cornerTile A pointer to the tile to have the householder applied.
* @param columnTile The tile in which the householders are stored.
* @param rowTile The upper tile to have the householders applied.
* @param tilesize The number of elements on a row/column of the tile.
* @param tauMatrix A pointer to the tau Matrix.
* @param ii The value of i for the QR decomposition.
* @param jj The value of j for the QR decomposition.
* @param kk The value of k for the QR decomposition.
* @param tauNum The number of tau values stored for each row of the matrix
* (this is equal to the number of tiles on each column).
*
*
*/
void DSSRFT(double* restrict cornerTile, double* restrict columnTile,
double* restrict rowTile, int tilesize, int ii, int jj, int kk,
double* restrict tauMatrix, int tauNum) {
int i, j, n;
double z;
double w[2*tilesize];
for (i = 0; i < tilesize; i++) { for (j = 0; j < i; j++) w[j] = 0.0;
w[i] = 1.0;
for (j = i + 1; j < tilesize; j++) w[j] = 0.0;
for (j = 0; j < tilesize; j++)
w[j + tilesize] = columnTile[i * tilesize + j];
/* Apply householder vector (w) to the tiles.*/
for (j = 0; j < tilesize; j++) {
z = 0.0;
/* Compute w' * A_j */
for (n = 0; n < tilesize; n++) {
z += w[n] * rowTile[j * tilesize + n];
z += w[n + tilesize] * cornerTile[j * tilesize + n];
}
for (n = 0; n < tilesize; n++) {
rowTile[j * tilesize + n] =
rowTile[j * tilesize + n] -
tauMatrix[(kk * tilesize + i) * tauNum + ii] * w[n] * z;
cornerTile[j * tilesize + n] =
cornerTile[j * tilesize + n] -
tauMatrix[(kk * tilesize + i) * tauNum + ii] * w[tilesize + n] * z;
}
}
}
}
/**
* @brief Computed a tiled QR factorization using QuickSched.
*
* @param m Number of tile rows.
* @param n Number of tile columns.
* @param nr_threads Number of threads to use.
*/
void test_qr(int m, int n, int K, int nr_threads, int runs, double* matrix) {
int k, j, i;
double* A, *A_orig, *tau;
struct qsched s;
qsched_task_t* tid, tid_new;
qsched_res_t* rid;
int data[3];
ticks tic, toc_run, tot_setup, tot_run = 0;
enum task_types {
task_DGEQRF,
task_DLARFT,
task_DTSQRF,
task_DSSRFT
};
/* Runner function to pass to the scheduler. */
void runner(int type, void * data) {
/* Decode the task data. */
int* idata = (int*)data;
int i = idata[0], j = idata[1], k = idata[2];
/* Decode and execute the task. */
switch (type) {
case task_DGEQRF:
DGEQRF(&A[(k * m + k) * K * K], K, tau, k, m);
break;
case task_DLARFT:
DLARFT(&A[(k * m + k) * K * K], &A[(j * m + k) * K * K], K, j, k, tau,
m);
break;
case task_DTSQRF:
DTSQRF(&A[(k * m + k) * K * K], &A[(k * m + i) * K * K], K, i, k, tau,
m); break;
case task_DSSRFT:
DSSRFT(&A[(j * m + i) * K * K], &A[(k * m + i) * K * K],
&A[(j * m + k) * K * K], K, i, j, k, tau, m);
break;
default:
error("Unknown task type.");
}
}
/* Allocate and fill the original matrix. */
if ((A = (double*)malloc(sizeof(double)* m* n* K* K)) == NULL ||
(tau = (double*)malloc(sizeof(double)* m* n* K)) == NULL ||
(A_orig = (double*)malloc(sizeof(double) * m * n * K * K)) == NULL)
error("Failed to allocate matrices.");
for (k = 0; k < m * n * K * K; k++) {
if (matrix == NULL)
A_orig[k] = 2 * ((double)rand()) / RAND_MAX - 1.0;
else
A_orig[k] = matrix[k];
}
memcpy(A, A_orig, sizeof(double) * m * n * K * K);
bzero(tau, sizeof(double) * m * n * K);
/* Dump A_orig. */
/* message( "A_orig = [" );
for ( k = 0 ; k < m*K ; k++ ) {
for ( j = 0 ; j < n*K ; j++ )
printf( "%.3f " , A_orig[ j*m*K + k ] );
printf( "\n" );
}
printf( "];\n" ); */
/* Initialize the scheduler. */
qsched_init(&s, nr_threads, qsched_flag_none);
/* Allocate and init the task ID and resource ID matrix. */
tic = getticks();
if ((tid = (qsched_task_t*)malloc(sizeof(qsched_task_t)* m* n)) == NULL ||
(rid = (qsched_res_t*)malloc(sizeof(qsched_res_t) * m * n)) == NULL)
error("Failed to allocate tid/rid matrix.");
for (k = 0; k < m * n; k++) {
tid[k] = qsched_task_none;
rid[k] = qsched_addres(&s, qsched_owner_none, qsched_res_none);
}
/* Build the tasks. */
for (k = 0; k < m && k < n; k++) {
/* Add kth corner task. */
data[0] = k;
data[1] = k;
data[2] = k;
tid_new = qsched_addtask(&s, task_DGEQRF, task_flag_none, data,
sizeof(int) * 3, 2);
qsched_addlock(&s, tid_new, rid[k * m + k]);
if (tid[k * m + k] != -1) qsched_addunlock(&s, tid[k * m + k], tid_new);
tid[k * m + k] = tid_new;
/* Add column tasks on kth row. */
for (j = k + 1; j < n; j++) {
data[0] = k;
data[1] = j;
data[2] = k;
tid_new = qsched_addtask(&s, task_DLARFT, task_flag_none, data,
sizeof(int) * 3, 3);
qsched_addlock(&s, tid_new, rid[j * m + k]);
qsched_adduse(&s, tid_new, rid[k * m + k]);
qsched_addunlock(&s, tid[k * m + k], tid_new);
if (tid[j * m + k] != -1) qsched_addunlock(&s, tid[j * m + k], tid_new); tid[j * m + k] = tid_new;
}
/* For each following row... */
for (i = k + 1; i < m; i++) {
/* Add the row taks for the kth column. */
data[0] = i;
data[1] = k;
data[2] = k;
tid_new = qsched_addtask(&s, task_DTSQRF, task_flag_none, data,
sizeof(int) * 3, 3);
qsched_addlock(&s, tid_new, rid[k * m + i]);
qsched_adduse(&s, tid_new, rid[k * m + k]);
qsched_addunlock(&s, tid[k * m + (i - 1)], tid_new);
if (tid[k * m + i] != -1) qsched_addunlock(&s, tid[k * m + i], tid_new);
tid[k * m + i] = tid_new;
/* Add the inner tasks. */
for (j = k + 1; j < n; j++) {
data[0] = i;
data[1] = j;
data[2] = k;
tid_new = qsched_addtask(&s, task_DSSRFT, task_flag_none, data,
sizeof(int) * 3, 5);
qsched_addlock(&s, tid_new, rid[j * m + i]);
qsched_adduse(&s, tid_new, rid[k * m + i]);
qsched_adduse(&s, tid_new, rid[j * m + k]);
// qsched_addunlock(&s, tid[k * m + i], tid_new);
qsched_addunlock(&s, tid[j * m + i - 1], tid_new);
if (tid[j * m + i] != -1) qsched_addunlock(&s, tid[j * m + i], tid_new);
tid[j * m + i] = tid_new;
}
}
} /* build the tasks. */
tot_setup = getticks() - tic;
/* Dump the number of tasks. */
message("total nr of tasks: %i.", s.count);
message("total nr of deps: %i.", s.count_deps);
message("total nr of res: %i.", s.count_res);
message("total nr of locks: %i.", s.count_locks);
message("total nr of uses: %i.", s.count_uses);
/* Loop over the number of runs. */
for (k = 0; k < runs; k++) {
/* Execute the the tasks. */
tic = getticks();
qsched_run(&s, nr_threads, runner);
toc_run = getticks();
message("%ith run took %lli ticks...", k, toc_run - tic);
tot_run += toc_run - tic;
}
/* Dump A. */
/* message( "A = [" );
for ( k = 0 ; k < m*K ; k++ ) {
for ( j = 0 ; j < n*K ; j++ )
printf( "%.3f " , A[ j*m*K + k ] );
printf( "\n" );
}
printf( "];\n" ); */
/* Dump tau. */
/* message( "tau = [" );
for ( k = 0 ; k < m*K ; k++ ) {
for ( j = 0 ; j < n ; j++ ) printf( "%.3f " , tau[ j*m*K + k ] );
printf( "\n" );
}
printf( "];\n" ); */
/* Dump the tasks. */
/* for ( k = 0 ; k < s.count ; k++ ) {
int *d = (int *)&s.data[ s.tasks[k].data ];
printf( " %i %i %i %i %lli %lli\n" , s.tasks[k].type , s.tasks[k].qid ,
d[0] , d[1] , s.tasks[k].tic , s.tasks[k].toc );
} */
/* Dump the costs. */
message("costs: setup=%lli ticks, run=%lli ticks.", tot_setup,
tot_run / runs);
/* Dump the timers. */
for (k = 0; k < qsched_timer_count; k++)
message("timer %s is %lli ticks.", qsched_timer_names[k],
s.timers[k] / runs);
if (matrix != NULL) {
for (k = 0; k < m * n * K * K; k++) matrix[k] = A[k];
}
/* Test if the decomposition was correct.*/
/*double *tempMatrix = tileToColumn(A, m*n*K*K, m, n, K);
double *Q = computeQ(tempMatrix, m*K, K, tau, m);
double *R = getR(tempMatrix, m*K);
cblas_dgemm(CblasColMajor, CblasNoTrans, CblasNoTrans, m*K, m*K, m*K, 1.0, Q,
m*K, R, m*K, 0.0, tempMatrix, m*K);
free(Q);
Q = tileToColumn(A_orig, m*n*K*K, m, n, K);
for(i = 0; i < m * n * K * K; i++)
{
if(Q[i] != 0 && (Q[i] / tempMatrix[i] > 1.005 || Q[i] / tempMatrix[i] <
0.995))
printf("Not correct at value %i %.3f %.3e %.3e\n", i, A[i], Q[i],
tempMatrix[i]);
}
free(tempMatrix);
free(Q);
free(R);*/
/* Clean up. */
free(A);
free(A_orig);
free(tau);
free(tid);
free(rid);
qsched_free(&s);
}
/* Generates a random matrix. */
double* generateColumnMatrix(int size) {
double* matrix = malloc(sizeof(double) * size * size);
if (matrix == NULL) error("Failed to allocate matrix");
unsigned long int m_z = 35532;
int i;
for (i = 0; i < size * size; i++) {
m_z = (1664525 * m_z + 1013904223) % 4294967296;
matrix[i] = m_z % 100;
if (matrix[i] < 0) matrix[i] += 100;
}
return matrix;
}
/**
* @brief Main function. */
int main(int argc, char* argv[]) {
int c, nr_threads;
int M = 4, N = 4, runs = 1, K = 32;
/* Get the number of threads. */
#pragma omp parallel shared(nr_threads)
{
if (omp_get_thread_num() == 0) nr_threads = omp_get_num_threads();
}
/* Parse the options */
while ((c = getopt(argc, argv, "m:n:k:r:t:")) != -1) switch (c) {
case 'm':
if (sscanf(optarg, "%d", &M) != 1) error("Error parsing dimension M.");
break;
case 'n':
if (sscanf(optarg, "%d", &N) != 1) error("Error parsing dimension M.");
break;
case 'k':
if (sscanf(optarg, "%d", &K) != 1) error("Error parsing tile size.");
break;
case 'r':
if (sscanf(optarg, "%d", &runs) != 1)
error("Error parsing number of runs.");
break;
case 't':
if (sscanf(optarg, "%d", &nr_threads) != 1)
error("Error parsing number of threads.");
omp_set_num_threads(nr_threads);
break;
case '?':
fprintf(stderr, "Usage: %s [-t nr_threads] [-m M] [-n N] [-k K]\n",
argv[0]);
fprintf(stderr, "Computes the tiled QR decomposition of an MxN tiled\n"
"matrix using nr_threads threads.\n");
exit(EXIT_FAILURE);
}
/* Dump arguments. */
message("Computing the tiled QR decomposition of a %ix%i matrix using %i "
"threads (%i runs).",
32 * M, 32 * N, nr_threads, runs);
test_qr(M, N, K, nr_threads, runs, NULL);
}