libflame  revision_anchor
Functions
bl1_trmm.c File Reference

(r)

Functions

void bl1_strmm (side1_t side, uplo1_t uplo, trans1_t trans, diag1_t diag, int m, int n, float *alpha, float *a, int a_rs, int a_cs, float *b, int b_rs, int b_cs)
 
void bl1_dtrmm (side1_t side, uplo1_t uplo, trans1_t trans, diag1_t diag, int m, int n, double *alpha, double *a, int a_rs, int a_cs, double *b, int b_rs, int b_cs)
 
void bl1_ctrmm (side1_t side, uplo1_t uplo, trans1_t trans, diag1_t diag, int m, int n, scomplex *alpha, scomplex *a, int a_rs, int a_cs, scomplex *b, int b_rs, int b_cs)
 
void bl1_ztrmm (side1_t side, uplo1_t uplo, trans1_t trans, diag1_t diag, int m, int n, dcomplex *alpha, dcomplex *a, int a_rs, int a_cs, dcomplex *b, int b_rs, int b_cs)
 
void bl1_strmm_blas (side1_t side, uplo1_t uplo, trans1_t trans, diag1_t diag, int m, int n, float *alpha, float *a, int lda, float *b, int ldb)
 
void bl1_dtrmm_blas (side1_t side, uplo1_t uplo, trans1_t trans, diag1_t diag, int m, int n, double *alpha, double *a, int lda, double *b, int ldb)
 
void bl1_ctrmm_blas (side1_t side, uplo1_t uplo, trans1_t trans, diag1_t diag, int m, int n, scomplex *alpha, scomplex *a, int lda, scomplex *b, int ldb)
 
void bl1_ztrmm_blas (side1_t side, uplo1_t uplo, trans1_t trans, diag1_t diag, int m, int n, dcomplex *alpha, dcomplex *a, int lda, dcomplex *b, int ldb)
 

Function Documentation

◆ bl1_ctrmm()

void bl1_ctrmm ( side1_t  side,
uplo1_t  uplo,
trans1_t  trans,
diag1_t  diag,
int  m,
int  n,
scomplex alpha,
scomplex a,
int  a_rs,
int  a_cs,
scomplex b,
int  b_rs,
int  b_cs 
)
220 {
221  int m_save = m;
222  int n_save = n;
223  scomplex* a_save = a;
224  scomplex* b_save = b;
225  int a_rs_save = a_rs;
226  int a_cs_save = a_cs;
227  int b_rs_save = b_rs;
228  int b_cs_save = b_cs;
229  scomplex* a_conj;
230  int dim_a;
231  int lda, inca;
232  int ldb, incb;
233  int lda_conj, inca_conj;
234  int a_was_copied;
235 
236  // Return early if possible.
237  if ( bl1_zero_dim2( m, n ) ) return;
238 
239  // If necessary, allocate, initialize, and use a temporary contiguous
240  // copy of each matrix rather than the original matrices.
241  bl1_set_dim_with_side( side, m, n, &dim_a );
242  bl1_ccreate_contigmr( uplo,
243  dim_a,
244  dim_a,
245  a_save, a_rs_save, a_cs_save,
246  &a, &a_rs, &a_cs );
247 
249  n,
250  b_save, b_rs_save, b_cs_save,
251  &b, &b_rs, &b_cs );
252 
253  // Figure out whether A was copied to contiguous memory. This is used to
254  // prevent redundant copying.
255  a_was_copied = ( a != a_save );
256 
257  // Initialize with values assuming column-major storage.
258  lda = a_cs;
259  inca = a_rs;
260  ldb = b_cs;
261  incb = b_rs;
262 
263  // Adjust the parameters based on the storage of each matrix.
264  if ( bl1_is_col_storage( b_rs, b_cs ) )
265  {
266  if ( bl1_is_col_storage( a_rs, a_cs ) )
267  {
268  // requested operation: B_c := tr( uplo( A_c ) ) * B_c
269  // effective operation: B_c := tr( uplo( A_c ) ) * B_c
270  }
271  else // if ( bl1_is_row_storage( a_rs, a_cs ) )
272  {
273  // requested operation: B_c := tr( uplo( A_r ) ) * B_c
274  // effective operation: B_c := tr( ~uplo( A_c ) )^T * B_c
275  bl1_swap_ints( lda, inca );
276 
277  bl1_toggle_uplo( uplo );
278  bl1_toggle_trans( trans );
279  }
280  }
281  else // if ( bl1_is_row_storage( b_rs, b_cs ) )
282  {
283  if ( bl1_is_col_storage( a_rs, a_cs ) )
284  {
285  // requested operation: B_r := tr( uplo( A_c ) ) * B_r
286  // effective operation: B_c := B_c * tr( uplo( A_c ) )^T
287  bl1_swap_ints( ldb, incb );
288 
289  bl1_swap_ints( m, n );
290 
291  bl1_toggle_side( side );
292  bl1_toggle_trans( trans );
293  }
294  else // if ( bl1_is_row_storage( a_rs, a_cs ) )
295  {
296  // requested operation: B_r := tr( uplo( A_r ) ) * B_r
297  // effective operation: B_c := B_c * tr( ~uplo( A_c ) )
298  bl1_swap_ints( ldb, incb );
299  bl1_swap_ints( lda, inca );
300 
301  bl1_swap_ints( m, n );
302 
303  bl1_toggle_uplo( uplo );
304  bl1_toggle_side( side );
305  }
306  }
307 
308  // Initialize with values assuming that trans is not conjnotrans.
309  a_conj = a;
310  lda_conj = lda;
311  inca_conj = inca;
312 
313  // We want to handle the conjnotrans case. The easiest way to do so is
314  // by making a conjugated copy of A.
315  if ( bl1_is_conjnotrans( trans ) && !a_was_copied )
316  {
317  int dim_a;
318 
319  bl1_set_dim_with_side( side, m, n, &dim_a );
320 
321  a_conj = bl1_callocm( dim_a, dim_a );
322  lda_conj = dim_a;
323  inca_conj = 1;
324 
325  bl1_ccopymrt( uplo,
327  dim_a,
328  dim_a,
329  a, inca, lda,
330  a_conj, inca_conj, lda_conj );
331  }
332  else if ( bl1_is_conjnotrans( trans ) && a_was_copied )
333  {
334  int dim_a;
335 
336  bl1_set_dim_with_side( side, m, n, &dim_a );
337 
338  bl1_cconjmr( uplo,
339  dim_a,
340  dim_a,
341  a_conj, inca_conj, lda_conj );
342  }
343 
344 
345  bl1_ctrmm_blas( side,
346  uplo,
347  trans,
348  diag,
349  m,
350  n,
351  alpha,
352  a_conj, lda_conj,
353  b, ldb );
354 
355  if ( bl1_is_conjnotrans( trans ) && !a_was_copied )
356  bl1_cfree( a_conj );
357 
358  // Free any temporary contiguous matrices, copying the result back to
359  // the original matrix.
360  bl1_cfree_contigm( a_save, a_rs_save, a_cs_save,
361  &a, &a_rs, &a_cs );
362 
363  bl1_cfree_saved_contigm( m_save,
364  n_save,
365  b_save, b_rs_save, b_cs_save,
366  &b, &b_rs, &b_cs );
367 }
void bl1_cconjmr(uplo1_t uplo, int m, int n, scomplex *a, int a_rs, int a_cs)
Definition: bl1_conjmr.c:23
void bl1_ccopymrt(uplo1_t uplo, trans1_t trans, int m, int n, scomplex *a, int a_rs, int a_cs, scomplex *b, int b_rs, int b_cs)
Definition: bl1_copymrt.c:223
void bl1_ctrmm_blas(side1_t side, uplo1_t uplo, trans1_t trans, diag1_t diag, int m, int n, scomplex *alpha, scomplex *a, int lda, scomplex *b, int ldb)
Definition: bl1_trmm.c:614
int bl1_is_col_storage(int rs, int cs)
Definition: bl1_is.c:90
int bl1_is_conjnotrans(trans1_t trans)
Definition: bl1_is.c:25
int bl1_zero_dim2(int m, int n)
Definition: bl1_is.c:118
scomplex * bl1_callocm(unsigned int m, unsigned int n)
Definition: bl1_allocm.c:40
void bl1_cfree_contigm(scomplex *a_save, int a_rs_save, int a_cs_save, scomplex **a, int *a_rs, int *a_cs)
Definition: bl1_free_contigm.c:45
void bl1_cfree(scomplex *p)
Definition: bl1_free.c:40
void bl1_ccreate_contigm(int m, int n, scomplex *a_save, int a_rs_save, int a_cs_save, scomplex **a, int *a_rs, int *a_cs)
Definition: bl1_create_contigm.c:81
void bl1_ccreate_contigmr(uplo1_t uplo, int m, int n, scomplex *a_save, int a_rs_save, int a_cs_save, scomplex **a, int *a_rs, int *a_cs)
Definition: bl1_create_contigmr.c:77
void bl1_cfree_saved_contigm(int m, int n, scomplex *a_save, int a_rs_save, int a_cs_save, scomplex **a, int *a_rs, int *a_cs)
Definition: bl1_free_saved_contigm.c:59
void bl1_set_dim_with_side(side1_t side, int m, int n, int *dim_new)
Definition: bl1_set_dims.c:27
@ BLIS1_CONJ_NO_TRANSPOSE
Definition: blis_type_defs.h:56
Definition: blis_type_defs.h:133

References bl1_callocm(), bl1_cconjmr(), bl1_ccopymrt(), bl1_ccreate_contigm(), bl1_ccreate_contigmr(), bl1_cfree(), bl1_cfree_contigm(), bl1_cfree_saved_contigm(), bl1_ctrmm_blas(), bl1_is_col_storage(), bl1_is_conjnotrans(), bl1_set_dim_with_side(), bl1_zero_dim2(), and BLIS1_CONJ_NO_TRANSPOSE.

Referenced by bl1_ctrmmsx(), and FLA_Trmm_external().

◆ bl1_ctrmm_blas()

void bl1_ctrmm_blas ( side1_t  side,
uplo1_t  uplo,
trans1_t  trans,
diag1_t  diag,
int  m,
int  n,
scomplex alpha,
scomplex a,
int  lda,
scomplex b,
int  ldb 
)
615 {
616 #ifdef BLIS1_ENABLE_CBLAS_INTERFACES
617  enum CBLAS_ORDER cblas_order = CblasColMajor;
618  enum CBLAS_SIDE cblas_side;
619  enum CBLAS_UPLO cblas_uplo;
620  enum CBLAS_TRANSPOSE cblas_trans;
621  enum CBLAS_DIAG cblas_diag;
622 
623  bl1_param_map_to_netlib_side( side, &cblas_side );
624  bl1_param_map_to_netlib_uplo( uplo, &cblas_uplo );
625  bl1_param_map_to_netlib_trans( trans, &cblas_trans );
626  bl1_param_map_to_netlib_diag( diag, &cblas_diag );
627 
628  cblas_ctrmm( cblas_order,
629  cblas_side,
630  cblas_uplo,
631  cblas_trans,
632  cblas_diag,
633  m,
634  n,
635  alpha,
636  a, lda,
637  b, ldb );
638 #else
639  char blas_side;
640  char blas_uplo;
641  char blas_trans;
642  char blas_diag;
643 
644  bl1_param_map_to_netlib_side( side, &blas_side );
645  bl1_param_map_to_netlib_uplo( uplo, &blas_uplo );
646  bl1_param_map_to_netlib_trans( trans, &blas_trans );
647  bl1_param_map_to_netlib_diag( diag, &blas_diag );
648 
649  F77_ctrmm( &blas_side,
650  &blas_uplo,
651  &blas_trans,
652  &blas_diag,
653  &m,
654  &n,
655  alpha,
656  a, &lda,
657  b, &ldb );
658 #endif
659 }
void F77_ctrmm(char *side, char *uplo, char *transa, char *diag, int *m, int *n, scomplex *alpha, scomplex *a, int *lda, scomplex *b, int *ldb)
CBLAS_ORDER
Definition: blis_prototypes_cblas.h:17
@ CblasColMajor
Definition: blis_prototypes_cblas.h:17
CBLAS_UPLO
Definition: blis_prototypes_cblas.h:19
CBLAS_TRANSPOSE
Definition: blis_prototypes_cblas.h:18
CBLAS_SIDE
Definition: blis_prototypes_cblas.h:21
CBLAS_DIAG
Definition: blis_prototypes_cblas.h:20
void cblas_ctrmm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side, const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag, const int M, const int N, const void *alpha, const void *A, const int lda, void *B, const int ldb)
void bl1_param_map_to_netlib_diag(diag1_t blis_diag, void *blas_diag)
Definition: bl1_param_map.c:95
void bl1_param_map_to_netlib_trans(trans1_t blis_trans, void *blas_trans)
Definition: bl1_param_map.c:15
void bl1_param_map_to_netlib_side(side1_t blis_side, void *blas_side)
Definition: bl1_param_map.c:71
void bl1_param_map_to_netlib_uplo(uplo1_t blis_uplo, void *blas_uplo)
Definition: bl1_param_map.c:47

References bl1_param_map_to_netlib_diag(), bl1_param_map_to_netlib_side(), bl1_param_map_to_netlib_trans(), bl1_param_map_to_netlib_uplo(), cblas_ctrmm(), CblasColMajor, and F77_ctrmm().

Referenced by bl1_ctrmm().

◆ bl1_dtrmm()

void bl1_dtrmm ( side1_t  side,
uplo1_t  uplo,
trans1_t  trans,
diag1_t  diag,
int  m,
int  n,
double *  alpha,
double *  a,
int  a_rs,
int  a_cs,
double *  b,
int  b_rs,
int  b_cs 
)
117 {
118  int m_save = m;
119  int n_save = n;
120  double* a_save = a;
121  double* b_save = b;
122  int a_rs_save = a_rs;
123  int a_cs_save = a_cs;
124  int b_rs_save = b_rs;
125  int b_cs_save = b_cs;
126  int dim_a;
127  int lda, inca;
128  int ldb, incb;
129 
130  // Return early if possible.
131  if ( bl1_zero_dim2( m, n ) ) return;
132 
133  // If necessary, allocate, initialize, and use a temporary contiguous
134  // copy of each matrix rather than the original matrices.
135  bl1_set_dim_with_side( side, m, n, &dim_a );
136  bl1_dcreate_contigmr( uplo,
137  dim_a,
138  dim_a,
139  a_save, a_rs_save, a_cs_save,
140  &a, &a_rs, &a_cs );
141 
143  n,
144  b_save, b_rs_save, b_cs_save,
145  &b, &b_rs, &b_cs );
146 
147  // Initialize with values assuming column-major storage.
148  lda = a_cs;
149  inca = a_rs;
150  ldb = b_cs;
151  incb = b_rs;
152 
153  // Adjust the parameters based on the storage of each matrix.
154  if ( bl1_is_col_storage( b_rs, b_cs ) )
155  {
156  if ( bl1_is_col_storage( a_rs, a_cs ) )
157  {
158  // requested operation: B_c := tr( uplo( A_c ) ) * B_c
159  // effective operation: B_c := tr( uplo( A_c ) ) * B_c
160  }
161  else // if ( bl1_is_row_storage( a_rs, a_cs ) )
162  {
163  // requested operation: B_c := tr( uplo( A_r ) ) * B_c
164  // effective operation: B_c := tr( ~uplo( A_c ) )^T * B_c
165  bl1_swap_ints( lda, inca );
166 
167  bl1_toggle_uplo( uplo );
168  bl1_toggle_trans( trans );
169  }
170  }
171  else // if ( bl1_is_row_storage( b_rs, b_cs ) )
172  {
173  if ( bl1_is_col_storage( a_rs, a_cs ) )
174  {
175  // requested operation: B_r := tr( uplo( A_c ) ) * B_r
176  // effective operation: B_c := B_c * tr( uplo( A_c ) )^T
177  bl1_swap_ints( ldb, incb );
178 
179  bl1_swap_ints( m, n );
180 
181  bl1_toggle_side( side );
182  bl1_toggle_trans( trans );
183  }
184  else // if ( bl1_is_row_storage( a_rs, a_cs ) )
185  {
186  // requested operation: B_r := tr( uplo( A_r ) ) * B_r
187  // effective operation: B_c := B_c * tr( ~uplo( A_c ) )
188  bl1_swap_ints( ldb, incb );
189  bl1_swap_ints( lda, inca );
190 
191  bl1_swap_ints( m, n );
192 
193  bl1_toggle_uplo( uplo );
194  bl1_toggle_side( side );
195  }
196  }
197 
198  bl1_dtrmm_blas( side,
199  uplo,
200  trans,
201  diag,
202  m,
203  n,
204  alpha,
205  a, lda,
206  b, ldb );
207 
208  // Free any temporary contiguous matrices, copying the result back to
209  // the original matrix.
210  bl1_dfree_contigm( a_save, a_rs_save, a_cs_save,
211  &a, &a_rs, &a_cs );
212 
213  bl1_dfree_saved_contigm( m_save,
214  n_save,
215  b_save, b_rs_save, b_cs_save,
216  &b, &b_rs, &b_cs );
217 }
void bl1_dtrmm_blas(side1_t side, uplo1_t uplo, trans1_t trans, diag1_t diag, int m, int n, double *alpha, double *a, int lda, double *b, int ldb)
Definition: bl1_trmm.c:567
void bl1_dfree_contigm(double *a_save, int a_rs_save, int a_cs_save, double **a, int *a_rs, int *a_cs)
Definition: bl1_free_contigm.c:29
void bl1_dcreate_contigmr(uplo1_t uplo, int m, int n, double *a_save, int a_rs_save, int a_cs_save, double **a, int *a_rs, int *a_cs)
Definition: bl1_create_contigmr.c:45
void bl1_dcreate_contigm(int m, int n, double *a_save, int a_rs_save, int a_cs_save, double **a, int *a_rs, int *a_cs)
Definition: bl1_create_contigm.c:47
void bl1_dfree_saved_contigm(int m, int n, double *a_save, int a_rs_save, int a_cs_save, double **a, int *a_rs, int *a_cs)
Definition: bl1_free_saved_contigm.c:36

References bl1_dcreate_contigm(), bl1_dcreate_contigmr(), bl1_dfree_contigm(), bl1_dfree_saved_contigm(), bl1_dtrmm_blas(), bl1_is_col_storage(), bl1_set_dim_with_side(), and bl1_zero_dim2().

Referenced by bl1_dtrmmsx(), and FLA_Trmm_external().

◆ bl1_dtrmm_blas()

void bl1_dtrmm_blas ( side1_t  side,
uplo1_t  uplo,
trans1_t  trans,
diag1_t  diag,
int  m,
int  n,
double *  alpha,
double *  a,
int  lda,
double *  b,
int  ldb 
)
568 {
569 #ifdef BLIS1_ENABLE_CBLAS_INTERFACES
570  enum CBLAS_ORDER cblas_order = CblasColMajor;
571  enum CBLAS_SIDE cblas_side;
572  enum CBLAS_UPLO cblas_uplo;
573  enum CBLAS_TRANSPOSE cblas_trans;
574  enum CBLAS_DIAG cblas_diag;
575 
576  bl1_param_map_to_netlib_side( side, &cblas_side );
577  bl1_param_map_to_netlib_uplo( uplo, &cblas_uplo );
578  bl1_param_map_to_netlib_trans( trans, &cblas_trans );
579  bl1_param_map_to_netlib_diag( diag, &cblas_diag );
580 
581  cblas_dtrmm( cblas_order,
582  cblas_side,
583  cblas_uplo,
584  cblas_trans,
585  cblas_diag,
586  m,
587  n,
588  *alpha,
589  a, lda,
590  b, ldb );
591 #else
592  char blas_side;
593  char blas_uplo;
594  char blas_trans;
595  char blas_diag;
596 
597  bl1_param_map_to_netlib_side( side, &blas_side );
598  bl1_param_map_to_netlib_uplo( uplo, &blas_uplo );
599  bl1_param_map_to_netlib_trans( trans, &blas_trans );
600  bl1_param_map_to_netlib_diag( diag, &blas_diag );
601 
602  F77_dtrmm( &blas_side,
603  &blas_uplo,
604  &blas_trans,
605  &blas_diag,
606  &m,
607  &n,
608  alpha,
609  a, &lda,
610  b, &ldb );
611 #endif
612 }
void F77_dtrmm(char *side, char *uplo, char *transa, char *diag, int *m, int *n, double *alpha, double *a, int *lda, double *b, int *ldb)
void cblas_dtrmm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side, const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag, const int M, const int N, const double alpha, const double *A, const int lda, double *B, const int ldb)

References bl1_param_map_to_netlib_diag(), bl1_param_map_to_netlib_side(), bl1_param_map_to_netlib_trans(), bl1_param_map_to_netlib_uplo(), cblas_dtrmm(), CblasColMajor, and F77_dtrmm().

Referenced by bl1_dtrmm().

◆ bl1_strmm()

void bl1_strmm ( side1_t  side,
uplo1_t  uplo,
trans1_t  trans,
diag1_t  diag,
int  m,
int  n,
float *  alpha,
float *  a,
int  a_rs,
int  a_cs,
float *  b,
int  b_rs,
int  b_cs 
)
14 {
15  int m_save = m;
16  int n_save = n;
17  float* a_save = a;
18  float* b_save = b;
19  int a_rs_save = a_rs;
20  int a_cs_save = a_cs;
21  int b_rs_save = b_rs;
22  int b_cs_save = b_cs;
23  int dim_a;
24  int lda, inca;
25  int ldb, incb;
26 
27  // Return early if possible.
28  if ( bl1_zero_dim2( m, n ) ) return;
29 
30  // If necessary, allocate, initialize, and use a temporary contiguous
31  // copy of each matrix rather than the original matrices.
32  bl1_set_dim_with_side( side, m, n, &dim_a );
34  dim_a,
35  dim_a,
36  a_save, a_rs_save, a_cs_save,
37  &a, &a_rs, &a_cs );
38 
40  n,
41  b_save, b_rs_save, b_cs_save,
42  &b, &b_rs, &b_cs );
43 
44  // Initialize with values assuming column-major storage.
45  lda = a_cs;
46  inca = a_rs;
47  ldb = b_cs;
48  incb = b_rs;
49 
50  // Adjust the parameters based on the storage of each matrix.
51  if ( bl1_is_col_storage( b_rs, b_cs ) )
52  {
53  if ( bl1_is_col_storage( a_rs, a_cs ) )
54  {
55  // requested operation: B_c := tr( uplo( A_c ) ) * B_c
56  // effective operation: B_c := tr( uplo( A_c ) ) * B_c
57  }
58  else // if ( bl1_is_row_storage( a_rs, a_cs ) )
59  {
60  // requested operation: B_c := tr( uplo( A_r ) ) * B_c
61  // effective operation: B_c := tr( ~uplo( A_c ) )^T * B_c
62  bl1_swap_ints( lda, inca );
63 
64  bl1_toggle_uplo( uplo );
65  bl1_toggle_trans( trans );
66  }
67  }
68  else // if ( bl1_is_row_storage( b_rs, b_cs ) )
69  {
70  if ( bl1_is_col_storage( a_rs, a_cs ) )
71  {
72  // requested operation: B_r := tr( uplo( A_c ) ) * B_r
73  // effective operation: B_c := B_c * tr( uplo( A_c ) )^T
74  bl1_swap_ints( ldb, incb );
75 
76  bl1_swap_ints( m, n );
77 
78  bl1_toggle_side( side );
79  bl1_toggle_trans( trans );
80  }
81  else // if ( bl1_is_row_storage( a_rs, a_cs ) )
82  {
83  // requested operation: B_r := tr( uplo( A_r ) ) * B_r
84  // effective operation: B_c := B_c * tr( ~uplo( A_c ) )
85  bl1_swap_ints( ldb, incb );
86  bl1_swap_ints( lda, inca );
87 
88  bl1_swap_ints( m, n );
89 
90  bl1_toggle_uplo( uplo );
91  bl1_toggle_side( side );
92  }
93  }
94 
95  bl1_strmm_blas( side,
96  uplo,
97  trans,
98  diag,
99  m,
100  n,
101  alpha,
102  a, lda,
103  b, ldb );
104 
105  // Free any temporary contiguous matrices, copying the result back to
106  // the original matrix.
107  bl1_sfree_contigm( a_save, a_rs_save, a_cs_save,
108  &a, &a_rs, &a_cs );
109 
110  bl1_sfree_saved_contigm( m_save,
111  n_save,
112  b_save, b_rs_save, b_cs_save,
113  &b, &b_rs, &b_cs );
114 }
void bl1_strmm_blas(side1_t side, uplo1_t uplo, trans1_t trans, diag1_t diag, int m, int n, float *alpha, float *a, int lda, float *b, int ldb)
Definition: bl1_trmm.c:520
void bl1_sfree_contigm(float *a_save, int a_rs_save, int a_cs_save, float **a, int *a_rs, int *a_cs)
Definition: bl1_free_contigm.c:13
void bl1_sfree_saved_contigm(int m, int n, float *a_save, int a_rs_save, int a_cs_save, float **a, int *a_rs, int *a_cs)
Definition: bl1_free_saved_contigm.c:13
void bl1_screate_contigmr(uplo1_t uplo, int m, int n, float *a_save, int a_rs_save, int a_cs_save, float **a, int *a_rs, int *a_cs)
Definition: bl1_create_contigmr.c:13
void bl1_screate_contigm(int m, int n, float *a_save, int a_rs_save, int a_cs_save, float **a, int *a_rs, int *a_cs)
Definition: bl1_create_contigm.c:13

References bl1_is_col_storage(), bl1_screate_contigm(), bl1_screate_contigmr(), bl1_set_dim_with_side(), bl1_sfree_contigm(), bl1_sfree_saved_contigm(), bl1_strmm_blas(), and bl1_zero_dim2().

Referenced by bl1_strmmsx(), and FLA_Trmm_external().

◆ bl1_strmm_blas()

void bl1_strmm_blas ( side1_t  side,
uplo1_t  uplo,
trans1_t  trans,
diag1_t  diag,
int  m,
int  n,
float *  alpha,
float *  a,
int  lda,
float *  b,
int  ldb 
)
521 {
522 #ifdef BLIS1_ENABLE_CBLAS_INTERFACES
523  enum CBLAS_ORDER cblas_order = CblasColMajor;
524  enum CBLAS_SIDE cblas_side;
525  enum CBLAS_UPLO cblas_uplo;
526  enum CBLAS_TRANSPOSE cblas_trans;
527  enum CBLAS_DIAG cblas_diag;
528 
529  bl1_param_map_to_netlib_side( side, &cblas_side );
530  bl1_param_map_to_netlib_uplo( uplo, &cblas_uplo );
531  bl1_param_map_to_netlib_trans( trans, &cblas_trans );
532  bl1_param_map_to_netlib_diag( diag, &cblas_diag );
533 
534  cblas_strmm( cblas_order,
535  cblas_side,
536  cblas_uplo,
537  cblas_trans,
538  cblas_diag,
539  m,
540  n,
541  *alpha,
542  a, lda,
543  b, ldb );
544 #else
545  char blas_side;
546  char blas_uplo;
547  char blas_trans;
548  char blas_diag;
549 
550  bl1_param_map_to_netlib_side( side, &blas_side );
551  bl1_param_map_to_netlib_uplo( uplo, &blas_uplo );
552  bl1_param_map_to_netlib_trans( trans, &blas_trans );
553  bl1_param_map_to_netlib_diag( diag, &blas_diag );
554 
555  F77_strmm( &blas_side,
556  &blas_uplo,
557  &blas_trans,
558  &blas_diag,
559  &m,
560  &n,
561  alpha,
562  a, &lda,
563  b, &ldb );
564 #endif
565 }
void F77_strmm(char *side, char *uplo, char *transa, char *diag, int *m, int *n, float *alpha, float *a, int *lda, float *b, int *ldb)
void cblas_strmm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side, const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag, const int M, const int N, const float alpha, const float *A, const int lda, float *B, const int ldb)

References bl1_param_map_to_netlib_diag(), bl1_param_map_to_netlib_side(), bl1_param_map_to_netlib_trans(), bl1_param_map_to_netlib_uplo(), cblas_strmm(), CblasColMajor, and F77_strmm().

Referenced by bl1_strmm().

◆ bl1_ztrmm()

void bl1_ztrmm ( side1_t  side,
uplo1_t  uplo,
trans1_t  trans,
diag1_t  diag,
int  m,
int  n,
dcomplex alpha,
dcomplex a,
int  a_rs,
int  a_cs,
dcomplex b,
int  b_rs,
int  b_cs 
)
370 {
371  int m_save = m;
372  int n_save = n;
373  dcomplex* a_save = a;
374  dcomplex* b_save = b;
375  int a_rs_save = a_rs;
376  int a_cs_save = a_cs;
377  int b_rs_save = b_rs;
378  int b_cs_save = b_cs;
379  dcomplex* a_conj;
380  int dim_a;
381  int lda, inca;
382  int ldb, incb;
383  int lda_conj, inca_conj;
384  int a_was_copied;
385 
386  // Return early if possible.
387  if ( bl1_zero_dim2( m, n ) ) return;
388 
389  // If necessary, allocate, initialize, and use a temporary contiguous
390  // copy of each matrix rather than the original matrices.
391  bl1_set_dim_with_side( side, m, n, &dim_a );
392  bl1_zcreate_contigmr( uplo,
393  dim_a,
394  dim_a,
395  a_save, a_rs_save, a_cs_save,
396  &a, &a_rs, &a_cs );
397 
399  n,
400  b_save, b_rs_save, b_cs_save,
401  &b, &b_rs, &b_cs );
402 
403  // Figure out whether A was copied to contiguous memory. This is used to
404  // prevent redundant copying.
405  a_was_copied = ( a != a_save );
406 
407  // Initialize with values assuming column-major storage.
408  lda = a_cs;
409  inca = a_rs;
410  ldb = b_cs;
411  incb = b_rs;
412 
413  // Adjust the parameters based on the storage of each matrix.
414  if ( bl1_is_col_storage( b_rs, b_cs ) )
415  {
416  if ( bl1_is_col_storage( a_rs, a_cs ) )
417  {
418  // requested operation: B_c := tr( uplo( A_c ) ) * B_c
419  // effective operation: B_c := tr( uplo( A_c ) ) * B_c
420  }
421  else // if ( bl1_is_row_storage( a_rs, a_cs ) )
422  {
423  // requested operation: B_c := tr( uplo( A_r ) ) * B_c
424  // effective operation: B_c := tr( ~uplo( A_c ) )^T * B_c
425  bl1_swap_ints( lda, inca );
426 
427  bl1_toggle_uplo( uplo );
428  bl1_toggle_trans( trans );
429  }
430  }
431  else // if ( bl1_is_row_storage( b_rs, b_cs ) )
432  {
433  if ( bl1_is_col_storage( a_rs, a_cs ) )
434  {
435  // requested operation: B_r := tr( uplo( A_c ) ) * B_r
436  // effective operation: B_c := B_c * tr( uplo( A_c ) )^T
437  bl1_swap_ints( ldb, incb );
438 
439  bl1_swap_ints( m, n );
440 
441  bl1_toggle_side( side );
442  bl1_toggle_trans( trans );
443  }
444  else // if ( bl1_is_row_storage( a_rs, a_cs ) )
445  {
446  // requested operation: B_r := tr( uplo( A_r ) ) * B_r
447  // effective operation: B_c := B_c * tr( ~uplo( A_c ) )
448  bl1_swap_ints( ldb, incb );
449  bl1_swap_ints( lda, inca );
450 
451  bl1_swap_ints( m, n );
452 
453  bl1_toggle_side( side );
454  bl1_toggle_uplo( uplo );
455  }
456  }
457 
458  // Initialize with values assuming that trans is not conjnotrans.
459  a_conj = a;
460  lda_conj = lda;
461  inca_conj = inca;
462 
463  // We want to handle the conjnotrans case. The easiest way to do so is
464  // by making a conjugated copy of A.
465  if ( bl1_is_conjnotrans( trans ) && !a_was_copied )
466  {
467  int dim_a;
468 
469  bl1_set_dim_with_side( side, m, n, &dim_a );
470 
471  a_conj = bl1_zallocm( dim_a, dim_a );
472  lda_conj = dim_a;
473  inca_conj = 1;
474 
475  bl1_zcopymrt( uplo,
477  dim_a,
478  dim_a,
479  a, inca, lda,
480  a_conj, inca_conj, lda_conj );
481  }
482  else if ( bl1_is_conjnotrans( trans ) && a_was_copied )
483  {
484  int dim_a;
485 
486  bl1_set_dim_with_side( side, m, n, &dim_a );
487 
488  bl1_zconjmr( uplo,
489  dim_a,
490  dim_a,
491  a_conj, inca_conj, lda_conj );
492  }
493 
494  bl1_ztrmm_blas( side,
495  uplo,
496  trans,
497  diag,
498  m,
499  n,
500  alpha,
501  a_conj, lda_conj,
502  b, ldb );
503 
504  if ( bl1_is_conjnotrans( trans ) && !a_was_copied )
505  bl1_zfree( a_conj );
506 
507  // Free any temporary contiguous matrices, copying the result back to
508  // the original matrix.
509  bl1_zfree_contigm( a_save, a_rs_save, a_cs_save,
510  &a, &a_rs, &a_cs );
511 
512  bl1_zfree_saved_contigm( m_save,
513  n_save,
514  b_save, b_rs_save, b_cs_save,
515  &b, &b_rs, &b_cs );
516 }
void bl1_zconjmr(uplo1_t uplo, int m, int n, dcomplex *a, int a_rs, int a_cs)
Definition: bl1_conjmr.c:79
void bl1_zcopymrt(uplo1_t uplo, trans1_t trans, int m, int n, dcomplex *a, int a_rs, int a_cs, dcomplex *b, int b_rs, int b_cs)
Definition: bl1_copymrt.c:328
void bl1_ztrmm_blas(side1_t side, uplo1_t uplo, trans1_t trans, diag1_t diag, int m, int n, dcomplex *alpha, dcomplex *a, int lda, dcomplex *b, int ldb)
Definition: bl1_trmm.c:661
dcomplex * bl1_zallocm(unsigned int m, unsigned int n)
Definition: bl1_allocm.c:45
void bl1_zcreate_contigm(int m, int n, dcomplex *a_save, int a_rs_save, int a_cs_save, dcomplex **a, int *a_rs, int *a_cs)
Definition: bl1_create_contigm.c:115
void bl1_zcreate_contigmr(uplo1_t uplo, int m, int n, dcomplex *a_save, int a_rs_save, int a_cs_save, dcomplex **a, int *a_rs, int *a_cs)
Definition: bl1_create_contigmr.c:109
void bl1_zfree(dcomplex *p)
Definition: bl1_free.c:45
void bl1_zfree_contigm(dcomplex *a_save, int a_rs_save, int a_cs_save, dcomplex **a, int *a_rs, int *a_cs)
Definition: bl1_free_contigm.c:61
void bl1_zfree_saved_contigm(int m, int n, dcomplex *a_save, int a_rs_save, int a_cs_save, dcomplex **a, int *a_rs, int *a_cs)
Definition: bl1_free_saved_contigm.c:82
Definition: blis_type_defs.h:138

References bl1_is_col_storage(), bl1_is_conjnotrans(), bl1_set_dim_with_side(), bl1_zallocm(), bl1_zconjmr(), bl1_zcopymrt(), bl1_zcreate_contigm(), bl1_zcreate_contigmr(), bl1_zero_dim2(), bl1_zfree(), bl1_zfree_contigm(), bl1_zfree_saved_contigm(), bl1_ztrmm_blas(), and BLIS1_CONJ_NO_TRANSPOSE.

Referenced by bl1_ztrmmsx(), and FLA_Trmm_external().

◆ bl1_ztrmm_blas()

void bl1_ztrmm_blas ( side1_t  side,
uplo1_t  uplo,
trans1_t  trans,
diag1_t  diag,
int  m,
int  n,
dcomplex alpha,
dcomplex a,
int  lda,
dcomplex b,
int  ldb 
)
662 {
663 #ifdef BLIS1_ENABLE_CBLAS_INTERFACES
664  enum CBLAS_ORDER cblas_order = CblasColMajor;
665  enum CBLAS_SIDE cblas_side;
666  enum CBLAS_UPLO cblas_uplo;
667  enum CBLAS_TRANSPOSE cblas_trans;
668  enum CBLAS_DIAG cblas_diag;
669 
670  bl1_param_map_to_netlib_side( side, &cblas_side );
671  bl1_param_map_to_netlib_uplo( uplo, &cblas_uplo );
672  bl1_param_map_to_netlib_trans( trans, &cblas_trans );
673  bl1_param_map_to_netlib_diag( diag, &cblas_diag );
674 
675  cblas_ztrmm( cblas_order,
676  cblas_side,
677  cblas_uplo,
678  cblas_trans,
679  cblas_diag,
680  m,
681  n,
682  alpha,
683  a, lda,
684  b, ldb );
685 #else
686  char blas_side;
687  char blas_uplo;
688  char blas_trans;
689  char blas_diag;
690 
691  bl1_param_map_to_netlib_side( side, &blas_side );
692  bl1_param_map_to_netlib_uplo( uplo, &blas_uplo );
693  bl1_param_map_to_netlib_trans( trans, &blas_trans );
694  bl1_param_map_to_netlib_diag( diag, &blas_diag );
695 
696  F77_ztrmm( &blas_side,
697  &blas_uplo,
698  &blas_trans,
699  &blas_diag,
700  &m,
701  &n,
702  alpha,
703  a, &lda,
704  b, &ldb );
705 #endif
706 }
void F77_ztrmm(char *side, char *uplo, char *transa, char *diag, int *m, int *n, dcomplex *alpha, dcomplex *a, int *lda, dcomplex *b, int *ldb)
void cblas_ztrmm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side, const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag, const int M, const int N, const void *alpha, const void *A, const int lda, void *B, const int ldb)

References bl1_param_map_to_netlib_diag(), bl1_param_map_to_netlib_side(), bl1_param_map_to_netlib_trans(), bl1_param_map_to_netlib_uplo(), cblas_ztrmm(), CblasColMajor, and F77_ztrmm().

Referenced by bl1_ztrmm().