LASs - Linear Algebra Routines on OmpSs  1.0.0
LASs
kdtrsm.c
Go to the documentation of this file.
1 #include "../include/lass.h"
2 
3 /**
4  *
5  * @file kdtrsm.c
6  *
7  * @brief LASs-DDSs kdtrsm routine.
8  *
9  * LASs-DDSs is a software package provided by:
10  * Barcelona Supercomputing Center - Centro Nacional de Supercomputacion
11  *
12  * @author Pedro Valero-Lara pedro.valero@bsc.es
13  * @author Boro Sofranac boro.sofranac@bsc.es
14  * @date 2017-12-14
15  * @reviewer
16  * @modified
17  *
18  **/
19 
20 /**
21  *
22  * @ingroup DDSS
23  *
24  * Solves one of the matrix equations:
25  *
26  * op( A ) * X = ALPHA * B, or X * op( A ) = ALPHA * B
27  *
28  * where op( A ) is one of:
29  *
30  * op( A ) = A or
31  * op( A ) = A**T
32  *
33  * ALPHA is a scalar, X and B are M by N matrices, A is a unit, or non-unit,
34  * upper or lower triangular matrix. The matrix X is overwritten on B
35  *
36  **/
37 
38 /**
39  *
40  * @param[in]
41  * SIDE enum DDSS_SIDE.
42  * SIDE specifies the position of the triangular A matrix in
43  * the operations:
44  * - Left: op( A ) * X = ALPHA * B.
45  * - Right: X * op( A ) = ALPHA * B.
46  *
47  * @param[in]
48  * UPLO enum DDSS_UPLO.
49  * UPLO specifies the form of A is stored:
50  * - Lower: Lower triangle of A is stored. The upper traingular part is
51  * not referenced.
52  * - Upper: Upper triangle of A is stored. The lower triangular part is
53  * not referenced.
54  *
55  * @param[in]
56  * TRANS_A enum DDSS_TRANS.
57  * TRANS_A specifies the form of op( A ) to be used:
58  * - NoTrans: op( A ) = A.
59  * - Trans: op( A ) = A**T.
60  *
61  * @param[in]
62  * DIAG enum DDSS_DIAG.
63  * DIAG specifies whether or not A is unit triangular as follows:
64  * - Unit: A is assumed to be unit triangular.
65  * - NonUnit: A is not assumed to be unit triangular.
66  *
67  * @param[in]
68  * M int.
69  * M specifies the number of rows of B. M must be at least zero.
70  *
71  * @param[in]
72  * N int.
73  * N specifies the number of columns of B. N must be at least zero.
74  *
75  * @param[in]
76  * ALPHA double.
77  * ALPHA specifies the scalar alpha.
78  *
79  * @param[in]
80  * A double *.
81  * A is a pointer to a matrix of dimension K by K, where K is M when
82  * SIDE = Left and is N otherwise.
83  * When UPLO = Uppper the strictly lower triangular part of A is not
84  * referenced and when UPLO = Lower the strictly upper triangular part
85  * of A is not referenced.
86  * Note that when DIAG = Unit, the diagonal elements of A are not
87  * referenced either, but are assumed to be unity.
88  *
89  * @param[in]
90  * LDA int.
91  * LDA specifies the number of columns of A ( row-major order ).
92  * When SIDE = Left then LDA must be at least max( 1, M ), otherwise
93  * LDA must be at least max( 1, N ).
94  *
95  * @param[in, out]
96  * B double *.
97  * B is a pointer to a matrix of dimension M by N.
98  * On exit the matrix B is overwritten by the solution matrix X.
99  *
100  * @param[in]
101  * LDB int.
102  * LDB specifies the number of columns of B ( row-major order ).
103  * LDB must be at least max( 1, N ).
104  *
105  **/
106 
107 /**
108  *
109  * @retval Success successful exit
110  * @retval NoSuccess unsuccessful exit
111  *
112  **/
113 
114 /**
115  *
116  * @sa ddss_tile
117  * @sa ddss_flat2tiled
118  * @sa ddss_dsymflat2tiled
119  * @sa ddss_tiled2flat
120  *
121  **/
122 
123 enum LASS_RETURN kdtrsm( enum DDSS_SIDE SIDE, enum DDSS_UPLO UPLO,
124  enum DDSS_TRANS TRANS_A, enum DDSS_DIAG DIAG,
125  int M, int N,
126  const double ALPHA, double *A, int LDA,
127  double *B, int LDB )
128 {
129 
130  // Local variables
131  int k;
132  int mt, nt;
133  int ki, ni, mi, nni, mmi;
134  int Am, An;
135  int Bm, Bn;
136  int tile_size_m;
137  int tile_size_mm;
138  int tile_size_n;
139  int tile_size_nn;
140  int tile_size_k;
141  double alpha;
142 
143  // Number of tiles
144  if ( M % TILE_SIZE == 0 )
145  {
146  mt = M / TILE_SIZE;
147  }
148  else
149  {
150  mt = ( M / TILE_SIZE ) + 1;
151  }
152 
153  if ( N % TILE_SIZE == 0 )
154  {
155  nt = N / TILE_SIZE;
156  }
157  else
158  {
159  nt = ( N / TILE_SIZE ) + 1;
160  }
161 
162  /****************************
163  --Tile matrices declaration--
164  ****************************/
165 
166  if ( SIDE == Left )
167  {
168  Am = An = mt;
169  k = M;
170  }
171  else
172  {
173  Am = An = nt;
174  k = N;
175  }
176 
177  Bm = mt;
178  Bn = nt;
179 
180  /***************************
181  --Tile matrices allocation--
182  ***************************/
183 
184  double ( *TILE_A )[An][TILE_SIZE * TILE_SIZE] = malloc(
185  Am * An * TILE_SIZE * TILE_SIZE * sizeof( double ) );
186 
187  if ( TILE_A == NULL )
188  {
189  fprintf( stderr, "Failure in ddss_dtile_alloc for matrix TILE_A\n" );
190  return NoSuccess;
191  }
192 
193  double ( *TILE_B )[Bn][TILE_SIZE * TILE_SIZE] = malloc(
194  Bm * Bn * TILE_SIZE * TILE_SIZE * sizeof( double ) );
195 
196  if ( TILE_B == NULL )
197  {
198  fprintf( stderr, "Failure in ddss_dtile_alloc for matrix TILE_B\n" );
199  return NoSuccess;
200  }
201 
202  /*********************************************
203  --From flat data layout to tiled data layout--
204  *********************************************/
205 
206  // From flat matrix A to tile matrix TILE_A
207  ddss_dsymflat2tiled( k, k, A, LDA, Am, An, TILE_A, UPLO );
208 
209  // From flat matrix B to tile matrix TILE_B
210  ddss_dflat2tiled( M, N, B, LDB, Bm, Bn, TILE_B );
211 
212  /*************
213  --DTRSM tile--
214  *************/
215 
216  if ( SIDE == Left )
217  {
218  if ( UPLO == Upper )
219  {
220  // --SIDE = Left & UPLO = Upper & TRANS_A = NoTrans--
221  if ( TRANS_A == NoTrans )
222  {
223  for ( ki = mt - 1; ki >= 0; ki-- )
224  {
225  tile_size_k = ddss_tile_size( M, ki );
226  if ( ki == mt - 1 )
227  {
228  alpha = ALPHA;
229  }
230  else
231  {
232  alpha = 1.0;
233  }
234 
235  for ( ni = 0; ni < nt; ni++ )
236  {
237  tile_size_n = ddss_tile_size( N, ni );
238 
239  #pragma oss task in( TILE_A[ki][ki] ) \
240  inout( TILE_B[ki][ni] ) \
241  shared( TILE_A, TILE_B ) \
242  firstprivate( ki, ni ) \
243  label( dtrsm )
244  cblas_dtrsm( CblasRowMajor,
245  ( CBLAS_SIDE ) SIDE, ( CBLAS_UPLO ) UPLO,
246  ( CBLAS_TRANSPOSE ) TRANS_A,
247  ( CBLAS_DIAG ) DIAG,
248  tile_size_k,
249  tile_size_n,
250  alpha, TILE_A[ki][ki], tile_size_k,
251  TILE_B[ki][ni], tile_size_n );
252  }
253 
254  for ( mi = ki - 1; mi >= 0; mi-- )
255  {
256  tile_size_m = ddss_tile_size( M, mi );
257  for ( nni = 0; nni < nt; nni++ )
258  {
259  tile_size_nn = ddss_tile_size( N, nni );
260 
261  #pragma oss task in( TILE_A[mi][ki] ) \
262  in( TILE_B[ki][nni] ) \
263  inout( TILE_B[mi][nni] ) \
264  shared( TILE_A, TILE_B ) \
265  firstprivate( ki, mi, nni ) \
266  label( dgemm )
267  cblas_dgemm( CblasRowMajor,
268  CblasNoTrans, CblasNoTrans,
269  tile_size_m,
270  tile_size_nn,
271  tile_size_k,
272  -1.0, TILE_A[mi][ki], tile_size_k,
273  TILE_B[ki][nni], tile_size_nn,
274  alpha, TILE_B[mi][nni], tile_size_nn );
275  }
276  }
277  }
278  }
279 
280  // --SIDE = Left & UPLO = Upper & TRANS_A = Trans--
281  else
282  {
283  for ( ki = 0; ki < mt; ki++ )
284  {
285  tile_size_k = ddss_tile_size( M, ki );
286  if ( ki == 0 )
287  {
288  alpha = ALPHA;
289  }
290  else
291  {
292  alpha = 1.0;
293  }
294 
295  for ( ni = 0; ni < nt; ni++ )
296  {
297  tile_size_n = ddss_tile_size( N, ni );
298 
299  #pragma oss task in( TILE_A[ki][ki] ) \
300  inout( TILE_B[ki][ni] ) \
301  shared( TILE_A, TILE_B ) \
302  firstprivate( ki, ni ) \
303  label( dtrsm )
304  cblas_dtrsm( CblasRowMajor,
305  ( CBLAS_SIDE ) SIDE, ( CBLAS_UPLO ) UPLO,
306  ( CBLAS_TRANSPOSE ) TRANS_A,
307  ( CBLAS_DIAG ) DIAG,
308  tile_size_k,
309  tile_size_n,
310  alpha, TILE_A[ki][ki], tile_size_k,
311  TILE_B[ki][ni], tile_size_n );
312  }
313 
314  for ( mi = ki + 1; mi < mt; mi++ )
315  {
316  tile_size_m = ddss_tile_size( M, mi );
317  for ( nni = 0; nni < nt; nni++ )
318  {
319  tile_size_nn = ddss_tile_size( N, nni );
320 
321  #pragma oss task in( TILE_A[ki][mi]) \
322  in( TILE_B[ki][nni] ) \
323  inout( TILE_B[mi][nni] ) \
324  shared( TILE_A, TILE_B ) \
325  firstprivate( ki, mi, nni ) \
326  label( dgemm )
327  cblas_dgemm( CblasRowMajor,
328  CblasTrans, CblasNoTrans,
329  tile_size_m,
330  tile_size_nn,
331  tile_size_k,
332  -1.0, TILE_A[ki][mi], tile_size_m,
333  TILE_B[ki][nni], tile_size_nn,
334  alpha, TILE_B[mi][nni], tile_size_nn );
335  }
336  }
337  }
338  }
339  }
340 
341  // --SIDE = Left & UPLO = Lower & TRANS_A = NoTrans--
342  else
343  {
344  if ( TRANS_A == NoTrans )
345  {
346  for ( ki = 0; ki < mt; ki++ )
347  {
348  tile_size_k = ddss_tile_size( M, ki );
349  if ( ki == 0 )
350  {
351  alpha = ALPHA;
352  }
353  else
354  {
355  alpha = 1.0;
356  }
357 
358  for ( ni = 0; ni < nt; ni++ )
359  {
360  tile_size_n = ddss_tile_size( N, ni );
361 
362  #pragma oss task in( TILE_A[ki][ki] ) \
363  inout( TILE_B[ki][ni] ) \
364  shared( TILE_A, TILE_B ) \
365  firstprivate( ki, ni ) \
366  label( dtrsm )
367  cblas_dtrsm( CblasRowMajor,
368  ( CBLAS_SIDE ) SIDE, ( CBLAS_UPLO ) UPLO,
369  ( CBLAS_TRANSPOSE ) TRANS_A,
370  ( CBLAS_DIAG ) DIAG,
371  tile_size_k,
372  tile_size_n,
373  alpha, TILE_A[ki][ki], tile_size_k,
374  TILE_B[ki][ni], tile_size_n );
375 
376  }
377  for ( nni = 0; nni < nt; nni++ )
378  {
379  tile_size_nn = ddss_tile_size( N, nni );
380 
381  for ( mi = ki + 1; mi < mt; mi++ )
382  {
383 
384  tile_size_m = ddss_tile_size( M, mi );
385 
386  #pragma oss task in( TILE_A[mi][ki] ) \
387  in( TILE_B[ki][nni] ) \
388  inout( TILE_B[mi][nni] ) \
389  shared( TILE_A, TILE_B ) \
390  firstprivate( ki, mi, nni ) \
391  label( dgemm )
392  cblas_dgemm( CblasRowMajor,
393  CblasNoTrans, CblasNoTrans,
394  tile_size_m,
395  tile_size_nn,
396  tile_size_k,
397  -1.0, TILE_A[mi][ki], tile_size_k,
398  TILE_B[ki][nni], tile_size_nn,
399  alpha, TILE_B[mi][nni], tile_size_nn );
400  }
401  }
402  }
403  }
404 
405  // --SIDE = Left & UPLO = Lower & TRANS_A = Trans--
406  else
407  {
408  for ( ki = mt - 1; ki >= 0; ki-- )
409  {
410  tile_size_k = ddss_tile_size( M, ki );
411  if ( ki == mt - 1 )
412  {
413  alpha = ALPHA;
414  }
415  else
416  {
417  alpha = 1.0;
418  }
419 
420  for ( ni = 0; ni < nt; ni++ )
421  {
422  tile_size_n = ddss_tile_size( N, ni );
423 
424  #pragma oss task in( TILE_A[ki][ki] ) \
425  inout( TILE_B[ki][ni] ) \
426  shared( TILE_A, TILE_B ) \
427  firstprivate( ki, ni ) \
428  label( dtrsm )
429  cblas_dtrsm( CblasRowMajor,
430  ( CBLAS_SIDE ) SIDE, ( CBLAS_UPLO ) UPLO,
431  ( CBLAS_TRANSPOSE ) TRANS_A,
432  ( CBLAS_DIAG ) DIAG,
433  tile_size_k,
434  tile_size_n,
435  alpha, TILE_A[ki][ki], tile_size_k,
436  TILE_B[ki][ni], tile_size_n );
437  }
438 
439  for ( mi = ki - 1; mi >= 0; mi-- )
440  {
441  tile_size_m = ddss_tile_size( M, mi );
442  for ( nni = 0; nni < nt; nni++ )
443  {
444  tile_size_nn = ddss_tile_size( N, nni );
445 
446  #pragma oss task in( TILE_A[ki][mi] ) \
447  in( TILE_B[ki][nni] ) \
448  inout( TILE_B[mi][nni] ) \
449  shared( TILE_A, TILE_B ) \
450  firstprivate( ki, mi, nni ) \
451  label( dgemm )
452  cblas_dgemm( CblasRowMajor,
453  CblasTrans, CblasNoTrans,
454  tile_size_m,
455  tile_size_nn,
456  tile_size_k,
457  -1.0, TILE_A[ki][mi], tile_size_m,
458  TILE_B[ki][nni], tile_size_nn,
459  alpha, TILE_B[mi][nni], tile_size_nn );
460  }
461  }
462  }
463  }
464  }
465  }
466 
467  // --SIDE = Right & UPLO = Upper & TRANS_A = NoTrans--
468  else
469  {
470  if ( UPLO == Upper )
471  {
472  if ( TRANS_A == NoTrans )
473  {
474  for ( ki = 0; ki < nt; ki++ )
475  {
476  tile_size_k = ddss_tile_size( N, ki );
477  if ( ki == 0 )
478  {
479  alpha = ALPHA;
480  }
481  else
482  {
483  alpha = 1.0;
484  }
485 
486  for ( mi = 0; mi < mt; mi++ )
487  {
488  tile_size_m = ddss_tile_size( M, mi );
489 
490  #pragma oss task in( TILE_A[ki][ki] ) \
491  inout( TILE_B[mi][ki] ) \
492  shared( TILE_A, TILE_B ) \
493  firstprivate( ki, mi ) \
494  label( dtrsm )
495  cblas_dtrsm( CblasRowMajor,
496  ( CBLAS_SIDE ) SIDE, ( CBLAS_UPLO ) UPLO,
497  ( CBLAS_TRANSPOSE ) TRANS_A,
498  ( CBLAS_DIAG ) DIAG,
499  tile_size_m,
500  tile_size_k,
501  alpha, TILE_A[ki][ki], tile_size_k,
502  TILE_B[mi][ki], tile_size_k );
503  }
504 
505  for ( ni = ki + 1; ni < nt; ni++ )
506  {
507  tile_size_n = ddss_tile_size( N, ni );
508  for ( mmi = 0; mmi < mt; mmi++ )
509  {
510  tile_size_mm = ddss_tile_size( M, mmi );
511 
512  #pragma oss task in( TILE_A[ki][ni] ) \
513  in( TILE_B[mmi][ki] ) \
514  inout( TILE_B[mmi][ni] ) \
515  shared( TILE_A, TILE_B ) \
516  firstprivate( ki, ni, mmi ) \
517  label( dgemm )
518  cblas_dgemm( CblasRowMajor,
519  CblasNoTrans, CblasNoTrans,
520  tile_size_mm,
521  tile_size_n,
522  tile_size_k,
523  -1.0, TILE_B[mmi][ki], tile_size_k,
524  TILE_A[ki][ni], tile_size_n,
525  alpha, TILE_B[mmi][ni], tile_size_n );
526  }
527  }
528  }
529  }
530  // --SIDE = Right & UPLO = Upper & TRANS_A = Trans--
531  else
532  {
533  for ( ki = nt - 1; ki >= 0; ki-- )
534  {
535  tile_size_k = ddss_tile_size( N, ki );
536  if ( ki == nt - 1 )
537  {
538  alpha = ALPHA;
539  }
540  else
541  {
542  alpha = 1.0;
543  }
544 
545  for ( mi = 0; mi < mt; mi++ )
546  {
547  tile_size_m = ddss_tile_size( M, mi );
548 
549  #pragma oss task in( TILE_A[ki][ki] ) \
550  inout( TILE_B[mi][ki] ) \
551  shared( TILE_A, TILE_B ) \
552  firstprivate( ki, mi ) \
553  label( dtrsm )
554  cblas_dtrsm( CblasRowMajor,
555  ( CBLAS_SIDE ) SIDE, ( CBLAS_UPLO ) UPLO,
556  ( CBLAS_TRANSPOSE ) TRANS_A,
557  ( CBLAS_DIAG ) DIAG,
558  tile_size_m,
559  tile_size_k,
560  alpha, TILE_A[ki][ki], tile_size_k,
561  TILE_B[mi][ki], tile_size_k );
562  }
563 
564  for ( ni = ki - 1; ni >= 0; ni-- )
565  {
566  tile_size_n = ddss_tile_size( N, ni );
567  for ( mmi = 0; mmi < mt; mmi++ )
568  {
569  tile_size_mm = ddss_tile_size( M, mmi );
570 
571  #pragma oss task in( TILE_A[ni][ki] ) \
572  in( TILE_B[mmi][ki] ) \
573  inout( TILE_B[mmi][ni] ) \
574  shared( TILE_A, TILE_B ) \
575  firstprivate( ki, ni, mmi ) \
576  label( dgemm )
577  cblas_dgemm( CblasRowMajor,
578  CblasNoTrans, CblasTrans,
579  tile_size_mm,
580  tile_size_n,
581  tile_size_k,
582  -1.0, TILE_B[mmi][ki], tile_size_k,
583  TILE_A[ni][ki], tile_size_k,
584  alpha, TILE_B[mmi][ni], tile_size_n );
585  }
586  }
587  }
588  }
589  }
590 
591  // --SIDE = Right & UPLO = Lower & TRANS_A = NoTrans--
592  else
593  {
594  if ( TRANS_A == NoTrans )
595  {
596  for ( ki = nt - 1; ki >= 0; ki-- )
597  {
598  tile_size_k = ddss_tile_size( N, ki );
599  if ( ki == nt - 1 )
600  {
601  alpha = ALPHA;
602  }
603  else
604  {
605  alpha = 1.0;
606  }
607 
608  for ( mi = 0; mi < mt; mi++ )
609  {
610  tile_size_m = ddss_tile_size( M, mi );
611 
612  #pragma oss task in( TILE_A[ki][ki] ) \
613  inout( TILE_B[mi][ki] ) \
614  shared( TILE_A, TILE_B ) \
615  firstprivate( ki, mi ) \
616  label( dtrsm )
617  cblas_dtrsm( CblasRowMajor,
618  ( CBLAS_SIDE ) SIDE, ( CBLAS_UPLO ) UPLO,
619  ( CBLAS_TRANSPOSE ) TRANS_A,
620  ( CBLAS_DIAG ) DIAG,
621  tile_size_m,
622  tile_size_k,
623  alpha, TILE_A[ki][ki], tile_size_k,
624  TILE_B[mi][ki], tile_size_k );
625  }
626 
627  for ( ni = ki - 1; ni >= 0; ni-- )
628  {
629  tile_size_n = ddss_tile_size( N, ni );
630  for ( mmi = 0; mmi < mt; mmi++ )
631  {
632  tile_size_mm = ddss_tile_size( M, mmi );
633 
634  #pragma oss task in( TILE_A[ki][ni] ) \
635  in( TILE_B[mmi][ki] ) \
636  inout( TILE_B[mmi][ni] ) \
637  shared( TILE_A, TILE_B ) \
638  firstprivate( ki, ni, mmi ) \
639  label( dgemm )
640  cblas_dgemm( CblasRowMajor,
641  CblasNoTrans, CblasNoTrans,
642  tile_size_mm,
643  tile_size_n,
644  tile_size_k,
645  -1.0, TILE_B[mmi][ki], tile_size_k,
646  TILE_A[ki][ni], tile_size_n,
647  alpha, TILE_B[mmi][ni], tile_size_n );
648  }
649  }
650  }
651  }
652 
653  // --SIDE = Right & UPLO = Lower & TRANS_A = Trans--
654  else
655  {
656  for ( ki = 0; ki < nt; ki++ )
657  {
658  tile_size_k = ddss_tile_size( N, ki );
659  if ( ki == 0 )
660  {
661  alpha = ALPHA;
662  }
663  else
664  {
665  alpha = 1.0;
666  }
667 
668  for ( mi = 0; mi < mt; mi++ )
669  {
670  tile_size_m = ddss_tile_size( M, mi );
671 
672  #pragma oss task in( TILE_A[ki][ki] ) \
673  inout( TILE_B[mi][ki] ) \
674  shared( TILE_A,TILE_B ) \
675  firstprivate( ki, mi ) \
676  label( dtrsm )
677  cblas_dtrsm( CblasRowMajor,
678  ( CBLAS_SIDE ) SIDE, ( CBLAS_UPLO ) UPLO,
679  ( CBLAS_TRANSPOSE ) TRANS_A,
680  ( CBLAS_DIAG ) DIAG,
681  tile_size_m,
682  tile_size_k,
683  alpha, TILE_A[ki][ki], tile_size_k,
684  TILE_B[mi][ki], tile_size_k );
685  }
686 
687  for ( ni = ki + 1; ni < nt; ni++ )
688  {
689  tile_size_n = ddss_tile_size( N, ni );
690  for ( mmi = 0; mmi < mt; mmi++ )
691  {
692  tile_size_mm = ddss_tile_size( M, mmi );
693 
694  #pragma oss task in( TILE_A[ni][ki] ) \
695  in( TILE_B[mmi][ki] ) \
696  inout( TILE_B[mmi][ni] ) \
697  shared( TILE_A, TILE_B ) \
698  firstprivate( ki, ni, mmi ) \
699  label( dgemm )
700  cblas_dgemm( CblasRowMajor,
701  CblasNoTrans, CblasTrans,
702  tile_size_mm,
703  tile_size_n,
704  tile_size_k,
705  -1.0, TILE_B[mmi][ki], tile_size_k,
706  TILE_A[ni][ki], tile_size_k,
707  alpha, TILE_B[mmi][ni], tile_size_n );
708  }
709  }
710  }
711  }
712  }
713  }
714 
715  /*********************************************
716  --From tiled data layout to flat data layout--
717  *********************************************/
718 
719  ddss_dtiled2flat( M, N, B, LDB, Bm, Bn, TILE_B );
720 
721  // --Tile matrices free--
722  free( TILE_A );
723  free( TILE_B );
724 
725  return Success;
726 
727 }
728 
void ddss_dflat2tiled(int M, int N, double *A, int LDA, int MT, int NT, double(*TILE_A)[NT][TILE_SIZE *TILE_SIZE])
void ddss_dsymflat2tiled(int M, int N, double *A, int LDA, int MT, int NT, double(*TILE_A)[NT][TILE_SIZE *TILE_SIZE], enum DDSS_UPLO UPLO)
enum LASS_RETURN kdtrsm(enum DDSS_SIDE SIDE, enum DDSS_UPLO UPLO, enum DDSS_TRANS TRANS_A, enum DDSS_DIAG DIAG, int M, int N, const double ALPHA, double *A, int LDA, double *B, int LDB)
Definition: kdtrsm.c:123
void ddss_dtiled2flat(int M, int N, double *A, int LDA, int MT, int NT, double(*TILE_A)[NT][TILE_SIZE *TILE_SIZE])
int ddss_tile_size(int M, int MT)
Definition: ddss_tile.c:52