|
@@ -68,6 +68,281 @@ static struct starpu_codelet cl22 =
|
|
|
.color = 0x00ff00,
|
|
|
};
|
|
|
|
|
|
+static void run_cholesky(starpu_data_handle_t **data_handles, int rank, int nodes)
|
|
|
+{
|
|
|
+ unsigned k, m, n;
|
|
|
+ unsigned unbound_prio = STARPU_MAX_PRIO == INT_MAX && STARPU_MIN_PRIO == INT_MIN;
|
|
|
+
|
|
|
+ for (k = 0; k < nblocks; k++)
|
|
|
+ {
|
|
|
+ starpu_iteration_push(k);
|
|
|
+
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl11,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k) : STARPU_MAX_PRIO,
|
|
|
+ STARPU_RW, data_handles[k][k],
|
|
|
+ 0);
|
|
|
+
|
|
|
+ for (m = k+1; m<nblocks; m++)
|
|
|
+ {
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl21,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k - m) : (m == k+1)?STARPU_MAX_PRIO:STARPU_DEFAULT_PRIO,
|
|
|
+ STARPU_R, data_handles[k][k],
|
|
|
+ STARPU_RW, data_handles[m][k],
|
|
|
+ 0);
|
|
|
+
|
|
|
+ starpu_mpi_cache_flush(MPI_COMM_WORLD, data_handles[k][k]);
|
|
|
+ if (my_distrib(k, k, nodes) == rank)
|
|
|
+ starpu_data_wont_use(data_handles[k][k]);
|
|
|
+
|
|
|
+ for (n = k+1; n<nblocks; n++)
|
|
|
+ {
|
|
|
+ if (n <= m)
|
|
|
+ {
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl22,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k - m - n) : ((n == k+1) && (m == k+1))?STARPU_MAX_PRIO:STARPU_DEFAULT_PRIO,
|
|
|
+ STARPU_R, data_handles[n][k],
|
|
|
+ STARPU_R, data_handles[m][k],
|
|
|
+ STARPU_RW | STARPU_COMMUTE, data_handles[m][n],
|
|
|
+ 0);
|
|
|
+ }
|
|
|
+ }
|
|
|
+
|
|
|
+ starpu_mpi_cache_flush(MPI_COMM_WORLD, data_handles[m][k]);
|
|
|
+ if (my_distrib(m, k, nodes) == rank)
|
|
|
+ starpu_data_wont_use(data_handles[m][k]);
|
|
|
+ }
|
|
|
+ starpu_iteration_pop();
|
|
|
+ }
|
|
|
+}
|
|
|
+
|
|
|
+/* TODO: generate from compiler polyhedral analysis of classical algorithm */
|
|
|
+static void run_cholesky_column(starpu_data_handle_t **data_handles, int rank, int nodes)
|
|
|
+{
|
|
|
+ unsigned k, m, n;
|
|
|
+ unsigned unbound_prio = STARPU_MAX_PRIO == INT_MAX && STARPU_MIN_PRIO == INT_MIN;
|
|
|
+
|
|
|
+ /* Column */
|
|
|
+ for (n = 0; n<nblocks; n++)
|
|
|
+ {
|
|
|
+ starpu_iteration_push(n);
|
|
|
+
|
|
|
+ /* Row */
|
|
|
+ for (m = n; m<nblocks; m++)
|
|
|
+ {
|
|
|
+ for (k = 0; k < n; k++)
|
|
|
+ {
|
|
|
+ /* Accumulate updates from TRSMs */
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl22,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k - m - n) : ((n == k+1) && (m == k+1))?STARPU_MAX_PRIO:STARPU_DEFAULT_PRIO,
|
|
|
+ STARPU_R, data_handles[n][k],
|
|
|
+ STARPU_R, data_handles[m][k],
|
|
|
+ STARPU_RW | STARPU_COMMUTE, data_handles[m][n],
|
|
|
+ 0);
|
|
|
+ }
|
|
|
+ k = n;
|
|
|
+ if (m > n)
|
|
|
+ {
|
|
|
+ /* non-diagonal block, solve */
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl21,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k - m) : (m == k+1)?STARPU_MAX_PRIO:STARPU_DEFAULT_PRIO,
|
|
|
+ STARPU_R, data_handles[k][k],
|
|
|
+ STARPU_RW, data_handles[m][k],
|
|
|
+ 0);
|
|
|
+ }
|
|
|
+ else
|
|
|
+ {
|
|
|
+ /* diagonal block, factorize */
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl11,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k) : STARPU_MAX_PRIO,
|
|
|
+ STARPU_RW, data_handles[k][k],
|
|
|
+ 0);
|
|
|
+ }
|
|
|
+ }
|
|
|
+
|
|
|
+ starpu_iteration_pop();
|
|
|
+ }
|
|
|
+
|
|
|
+ /* Submit flushes, StarPU will fit them according to the progress */
|
|
|
+ starpu_mpi_cache_flush_all_data(MPI_COMM_WORLD);
|
|
|
+ for (m = 0; m < nblocks; m++)
|
|
|
+ for (n = 0; n < nblocks ; n++)
|
|
|
+ starpu_data_wont_use(data_handles[m][n]);
|
|
|
+}
|
|
|
+
|
|
|
+/* TODO: generate from compiler polyhedral analysis of classical algorithm */
|
|
|
+static void run_cholesky_antidiagonal(starpu_data_handle_t **data_handles, int rank, int nodes)
|
|
|
+{
|
|
|
+ unsigned a, c;
|
|
|
+ unsigned k, m, n;
|
|
|
+ unsigned unbound_prio = STARPU_MAX_PRIO == INT_MAX && STARPU_MIN_PRIO == INT_MIN;
|
|
|
+
|
|
|
+ /* double-antidiagonal number:
|
|
|
+ * - a=0 contains (0,0) plus (1,0)
|
|
|
+ * - a=1 contains (2,0), (1,1) plus (3,0), (2, 1)
|
|
|
+ * - etc.
|
|
|
+ */
|
|
|
+ for (a = 0; a < nblocks; a++)
|
|
|
+ {
|
|
|
+ starpu_iteration_push(a);
|
|
|
+
|
|
|
+ unsigned nfirst;
|
|
|
+ if (2*a < nblocks)
|
|
|
+ nfirst = 0;
|
|
|
+ else
|
|
|
+ nfirst = 2*a - (nblocks-1);
|
|
|
+
|
|
|
+ /* column within first antidiagonal for a */
|
|
|
+ for (n = nfirst; n <= a; n++)
|
|
|
+ {
|
|
|
+ /* row */
|
|
|
+ m = 2*a-n;
|
|
|
+
|
|
|
+ /* Accumulate updates from TRSMs */
|
|
|
+ for (k = 0; k < n; k++)
|
|
|
+ {
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl22,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k - m - n) : ((n == k+1) && (m == k+1))?STARPU_MAX_PRIO:STARPU_DEFAULT_PRIO,
|
|
|
+ STARPU_R, data_handles[n][k],
|
|
|
+ STARPU_R, data_handles[m][k],
|
|
|
+ STARPU_RW | STARPU_COMMUTE, data_handles[m][n],
|
|
|
+ 0);
|
|
|
+ }
|
|
|
+
|
|
|
+ /* k = n */
|
|
|
+ if (n < a)
|
|
|
+ {
|
|
|
+ /* non-diagonal block, solve */
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl21,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k - m) : (m == k+1)?STARPU_MAX_PRIO:STARPU_DEFAULT_PRIO,
|
|
|
+ STARPU_R, data_handles[k][k],
|
|
|
+ STARPU_RW, data_handles[m][k],
|
|
|
+ 0);
|
|
|
+ }
|
|
|
+ else
|
|
|
+ {
|
|
|
+ /* diagonal block, factorize */
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl11,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k) : STARPU_MAX_PRIO,
|
|
|
+ STARPU_RW, data_handles[k][k],
|
|
|
+ 0);
|
|
|
+ }
|
|
|
+ }
|
|
|
+
|
|
|
+ /* column within second antidiagonal for a */
|
|
|
+ for (n = nfirst; n <= a; n++)
|
|
|
+ {
|
|
|
+ /* row */
|
|
|
+ m = 2*a-n + 1;
|
|
|
+
|
|
|
+ if (m >= nblocks)
|
|
|
+ /* Skip first item when even number of tiles */
|
|
|
+ continue;
|
|
|
+
|
|
|
+ /* Accumulate updates from TRSMs */
|
|
|
+ for (k = 0; k < n; k++)
|
|
|
+ {
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl22,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k - m - n) : ((n == k+1) && (m == k+1))?STARPU_MAX_PRIO:STARPU_DEFAULT_PRIO,
|
|
|
+ STARPU_R, data_handles[n][k],
|
|
|
+ STARPU_R, data_handles[m][k],
|
|
|
+ STARPU_RW | STARPU_COMMUTE, data_handles[m][n],
|
|
|
+ 0);
|
|
|
+ }
|
|
|
+ /* non-diagonal block, solve */
|
|
|
+ k = n;
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl21,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k - m) : (m == k+1)?STARPU_MAX_PRIO:STARPU_DEFAULT_PRIO,
|
|
|
+ STARPU_R, data_handles[k][k],
|
|
|
+ STARPU_RW, data_handles[m][k],
|
|
|
+ 0);
|
|
|
+ }
|
|
|
+
|
|
|
+ starpu_iteration_pop();
|
|
|
+ }
|
|
|
+
|
|
|
+ /* Submit flushes, StarPU will fit them according to the progress */
|
|
|
+ starpu_mpi_cache_flush_all_data(MPI_COMM_WORLD);
|
|
|
+ for (m = 0; m < nblocks; m++)
|
|
|
+ for (n = 0; n < nblocks ; n++)
|
|
|
+ starpu_data_wont_use(data_handles[m][n]);
|
|
|
+}
|
|
|
+
|
|
|
+/* TODO: generate from compiler polyhedral analysis of classical algorithm */
|
|
|
+static void run_cholesky_prio(starpu_data_handle_t **data_handles, int rank, int nodes)
|
|
|
+{
|
|
|
+ unsigned a;
|
|
|
+ int k, m, n;
|
|
|
+ unsigned unbound_prio = STARPU_MAX_PRIO == INT_MAX && STARPU_MIN_PRIO == INT_MIN;
|
|
|
+
|
|
|
+ /*
|
|
|
+ * This is basically similar to above, except that we shift k according to the priorities set in the algorithm, so that prio ~ 2*a or 2*a+1
|
|
|
+ * double-antidiagonal number:
|
|
|
+ * - a=0 contains (0,0) plus (1,0)
|
|
|
+ * - a=1 contains (2,0), (1,1) plus (3,0), (2, 1)
|
|
|
+ * - etc.
|
|
|
+ */
|
|
|
+ for (a = 0; a < 4*nblocks; a++)
|
|
|
+ {
|
|
|
+ starpu_iteration_push(a);
|
|
|
+
|
|
|
+ for (k = 0; k < nblocks; k++)
|
|
|
+ {
|
|
|
+ n = k;
|
|
|
+ /* Should be m = a-k-n; for potrf and trsm to respect
|
|
|
+ priorities, but needs to be this for dependencies */
|
|
|
+ m = a-2*k-n;
|
|
|
+
|
|
|
+ if (m < 0 || m >= nblocks)
|
|
|
+ continue;
|
|
|
+
|
|
|
+ if (m == n)
|
|
|
+ {
|
|
|
+ /* diagonal block, factorize */
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl11,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k) : STARPU_MAX_PRIO,
|
|
|
+ STARPU_RW, data_handles[k][k],
|
|
|
+ 0);
|
|
|
+ }
|
|
|
+ else
|
|
|
+ {
|
|
|
+ /* non-diagonal block, solve */
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl21,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k - m) : (m == k+1)?STARPU_MAX_PRIO:STARPU_DEFAULT_PRIO,
|
|
|
+ STARPU_R, data_handles[k][k],
|
|
|
+ STARPU_RW, data_handles[m][k],
|
|
|
+ 0);
|
|
|
+ }
|
|
|
+
|
|
|
+ /* column within antidiagonal for a */
|
|
|
+ for (n = k + 1; n < nblocks; n++)
|
|
|
+ {
|
|
|
+ /* row */
|
|
|
+ m = a-2*k-n;
|
|
|
+
|
|
|
+ if (m >= n && m < nblocks)
|
|
|
+ {
|
|
|
+ /* Update */
|
|
|
+ starpu_mpi_task_insert(MPI_COMM_WORLD, &cl22,
|
|
|
+ STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k - m - n) : ((n == k+1) && (m == k+1))?STARPU_MAX_PRIO:STARPU_DEFAULT_PRIO,
|
|
|
+ STARPU_R, data_handles[n][k],
|
|
|
+ STARPU_R, data_handles[m][k],
|
|
|
+ STARPU_RW | STARPU_COMMUTE, data_handles[m][n],
|
|
|
+ 0);
|
|
|
+ }
|
|
|
+ }
|
|
|
+
|
|
|
+ }
|
|
|
+
|
|
|
+ starpu_iteration_pop();
|
|
|
+ }
|
|
|
+
|
|
|
+ /* Submit flushes, StarPU will fit them according to the progress */
|
|
|
+ starpu_mpi_cache_flush_all_data(MPI_COMM_WORLD);
|
|
|
+ for (m = 0; m < nblocks; m++)
|
|
|
+ for (n = 0; n < nblocks ; n++)
|
|
|
+ starpu_data_wont_use(data_handles[m][n]);
|
|
|
+}
|
|
|
+
|
|
|
/*
|
|
|
* code to bootstrap the factorization
|
|
|
* and construct the DAG
|
|
@@ -79,8 +354,6 @@ void dw_cholesky(float ***matA, unsigned ld, int rank, int nodes, double *timing
|
|
|
starpu_data_handle_t **data_handles;
|
|
|
unsigned k, m, n;
|
|
|
|
|
|
- unsigned unbound_prio = STARPU_MAX_PRIO == INT_MAX && STARPU_MIN_PRIO == INT_MIN;
|
|
|
-
|
|
|
/* create all the DAG nodes */
|
|
|
|
|
|
data_handles = malloc(nblocks*sizeof(starpu_data_handle_t *));
|
|
@@ -91,7 +364,7 @@ void dw_cholesky(float ***matA, unsigned ld, int rank, int nodes, double *timing
|
|
|
for(n = 0; n < nblocks ; n++)
|
|
|
{
|
|
|
int mpi_rank = my_distrib(m, n, nodes);
|
|
|
- if (mpi_rank == rank)
|
|
|
+ if (mpi_rank == rank || (check && rank == 0))
|
|
|
{
|
|
|
//fprintf(stderr, "[%d] Owning data[%d][%d]\n", rank, n, m);
|
|
|
starpu_matrix_data_register(&data_handles[m][n], STARPU_MAIN_RAM, (uintptr_t)matA[m][n],
|
|
@@ -119,50 +392,17 @@ void dw_cholesky(float ***matA, unsigned ld, int rank, int nodes, double *timing
|
|
|
starpu_mpi_barrier(MPI_COMM_WORLD);
|
|
|
start = starpu_timing_now();
|
|
|
|
|
|
- for (k = 0; k < nblocks; k++)
|
|
|
+ switch (submission)
|
|
|
{
|
|
|
- starpu_iteration_push(k);
|
|
|
-
|
|
|
- starpu_mpi_task_insert(MPI_COMM_WORLD, &cl11,
|
|
|
- STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k) : STARPU_MAX_PRIO,
|
|
|
- STARPU_RW, data_handles[k][k],
|
|
|
- 0);
|
|
|
-
|
|
|
- for (m = k+1; m<nblocks; m++)
|
|
|
- {
|
|
|
- starpu_mpi_task_insert(MPI_COMM_WORLD, &cl21,
|
|
|
- STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k - m) : (m == k+1)?STARPU_MAX_PRIO:STARPU_DEFAULT_PRIO,
|
|
|
- STARPU_R, data_handles[k][k],
|
|
|
- STARPU_RW, data_handles[m][k],
|
|
|
- 0);
|
|
|
-
|
|
|
- starpu_mpi_cache_flush(MPI_COMM_WORLD, data_handles[k][k]);
|
|
|
- if (my_distrib(k, k, nodes) == rank)
|
|
|
- starpu_data_wont_use(data_handles[k][k]);
|
|
|
-
|
|
|
- for (n = k+1; n<nblocks; n++)
|
|
|
- {
|
|
|
- if (n <= m)
|
|
|
- {
|
|
|
- starpu_mpi_task_insert(MPI_COMM_WORLD, &cl22,
|
|
|
- STARPU_PRIORITY, noprio ? STARPU_DEFAULT_PRIO : unbound_prio ? (int)(2*nblocks - 2*k - m - n) : ((n == k+1) && (m == k+1))?STARPU_MAX_PRIO:STARPU_DEFAULT_PRIO,
|
|
|
- STARPU_R, data_handles[n][k],
|
|
|
- STARPU_R, data_handles[m][k],
|
|
|
- STARPU_RW | STARPU_COMMUTE, data_handles[m][n],
|
|
|
- 0);
|
|
|
- }
|
|
|
- }
|
|
|
-
|
|
|
- starpu_mpi_cache_flush(MPI_COMM_WORLD, data_handles[m][k]);
|
|
|
- if (my_distrib(m, k, nodes) == rank)
|
|
|
- starpu_data_wont_use(data_handles[m][k]);
|
|
|
- }
|
|
|
- starpu_iteration_pop();
|
|
|
+ case TRIANGLES: run_cholesky(data_handles, rank, nodes); break;
|
|
|
+ case COLUMNS: run_cholesky_column(data_handles, rank, nodes); break;
|
|
|
+ case ANTIDIAGONALS: run_cholesky_antidiagonal(data_handles, rank, nodes); break;
|
|
|
+ case PRIOS: run_cholesky_prio(data_handles, rank, nodes); break;
|
|
|
+ default: STARPU_ABORT();
|
|
|
}
|
|
|
|
|
|
starpu_mpi_wait_for_all(MPI_COMM_WORLD);
|
|
|
starpu_mpi_barrier(MPI_COMM_WORLD);
|
|
|
-
|
|
|
end = starpu_timing_now();
|
|
|
|
|
|
for (m = 0; m < nblocks; m++)
|
|
@@ -170,7 +410,7 @@ void dw_cholesky(float ***matA, unsigned ld, int rank, int nodes, double *timing
|
|
|
for(n = 0; n < nblocks ; n++)
|
|
|
{
|
|
|
/* Get back data on node 0 for the check */
|
|
|
- if (check)
|
|
|
+ if (check && data_handles[m][n])
|
|
|
starpu_mpi_get_data_on_node(MPI_COMM_WORLD, data_handles[m][n], 0);
|
|
|
|
|
|
if (data_handles[m][n])
|
|
@@ -248,24 +488,20 @@ void dw_cholesky_check_computation(float ***matA, int rank, int nodes, int *corr
|
|
|
{
|
|
|
for (m = 0; m < nblocks; m++)
|
|
|
{
|
|
|
- int mpi_rank = my_distrib(m, n, nodes);
|
|
|
- if (mpi_rank == rank)
|
|
|
+ for (nn = BLOCKSIZE*n ; nn < BLOCKSIZE*(n+1); nn++)
|
|
|
{
|
|
|
- for (nn = (size/nblocks)*n ; nn < (size/nblocks)*n+(size/nblocks); nn++)
|
|
|
+ for (mm = BLOCKSIZE*m ; mm < BLOCKSIZE*(m+1); mm++)
|
|
|
{
|
|
|
- for (mm = (size/nblocks)*m ; mm < (size/nblocks)*m+(size/nblocks); mm++)
|
|
|
+ if (nn <= mm)
|
|
|
{
|
|
|
- if (nn <= mm)
|
|
|
+ float orig = (1.0f/(1.0f+nn+mm)) + ((nn == mm)?1.0f*size:0.0f);
|
|
|
+ float err = fabsf(test_mat[mm +nn*size] - orig) / orig;
|
|
|
+ if (err > epsilon)
|
|
|
{
|
|
|
- float orig = (1.0f/(1.0f+nn+mm)) + ((nn == mm)?1.0f*size:0.0f);
|
|
|
- float err = fabsf(test_mat[mm +nn*size] - orig) / orig;
|
|
|
- if (err > epsilon)
|
|
|
- {
|
|
|
- FPRINTF(stderr, "[%d] Error[%u, %u] --> %2.20f != %2.20f (err %2.20f)\n", rank, nn, mm, test_mat[mm +nn*size], orig, err);
|
|
|
- *correctness = 0;
|
|
|
- *flops = 0;
|
|
|
- break;
|
|
|
- }
|
|
|
+ FPRINTF(stderr, "[%d] Error[%u, %u] --> %2.20f != %2.20f (err %2.20f)\n", rank, nn, mm, test_mat[mm +nn*size], orig, err);
|
|
|
+ *correctness = 0;
|
|
|
+ *flops = 0;
|
|
|
+ break;
|
|
|
}
|
|
|
}
|
|
|
}
|