OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
PB_Ctzhemm.c
Go to the documentation of this file.
1/* ---------------------------------------------------------------------
2*
3* -- PBLAS auxiliary routine (version 2.0) --
4* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
5* and University of California, Berkeley.
6* April 1, 1998
7*
8* ---------------------------------------------------------------------
9*/
10/*
11* Include files
12*/
13#include "../pblas.h"
14#include "../PBpblas.h"
15#include "../PBtools.h"
16#include "../PBblacs.h"
17#include "../PBblas.h"
18
19#ifdef __STDC__
20void PB_Ctzhemm( PBTYP_T * TYPE, char * SIDE, char * UPLO, Int M, Int N,
21 Int K, Int IOFFD, char * ALPHA, char * A, Int LDA,
22 char * BC, Int LDBC, char * BR, Int LDBR, char * CC,
23 Int LDCC, char * CR, Int LDCR )
24#else
25void PB_Ctzhemm( TYPE, SIDE, UPLO, M, N, K, IOFFD, ALPHA, A, LDA, BC,
26 LDBC, BR, LDBR, CC, LDCC, CR, LDCR )
27/*
28* .. Scalar Arguments ..
29*/
30 char * SIDE, * UPLO;
31 Int IOFFD, K, LDA, LDBC, LDBR, LDCC, LDCR, M, N;
32 char * ALPHA;
33/*
34* .. Array Arguments ..
35*/
36 PBTYP_T * TYPE;
37 char * A, * BC, * BR, * CC, * CR;
38#endif
39{
40/*
41* Purpose
42* =======
43*
44* PB_Ctzhemm performs the matrix-matrix operation
45*
46* C := alpha * A * B + C,
47*
48* or
49*
50* C := alpha * B * A + C,
51*
52* where alpha is a scalar, B and C are m by k and k by n matrices and A
53* is an m by n trapezoidal symmetric or Hermitian matrix.
54*
55* Arguments
56* =========
57*
58* TYPE (local input) pointer to a PBTYP_T structure
59* On entry, TYPE is a pointer to a structure of type PBTYP_T,
60* that contains type information (See pblas.h).
61*
62* SIDE (input) pointer to CHAR
63* On entry, SIDE specifies whether op( A ) multiplies B from
64* the left or right as follows:
65*
66* SIDE = 'L' or 'l' C := alpha * A * B + C,
67*
68* SIDE = 'R' or 'r' C := alpha * B * A + C.
69*
70* UPLO (input) pointer to CHAR
71* On entry, UPLO specifies which part of the matrix A is to be
72* referenced as follows:
73*
74* UPLO = 'L' or 'l' the lower trapezoid of A is referenced,
75*
76* UPLO = 'U' or 'u' the upper trapezoid of A is referenced,
77*
78* otherwise all of the matrix A is referenced.
79*
80* M (input) INTEGER
81* On entry, M specifies the number of rows of the matrix A. M
82* must be at least zero.
83*
84* N (input) INTEGER
85* On entry, N specifies the number of columns of the matrix A.
86* N must be at least zero.
87*
88* K (input) INTEGER
89* On entry, K specifies the number of rows of the matrices BR
90* and CR and the number of columns of the matrices BC and CC. K
91* must be at least zero.
92*
93* IOFFD (input) INTEGER
94* On entry, IOFFD specifies the position of the offdiagonal de-
95* limiting the upper and lower trapezoidal part of A as follows
96* (see the notes below):
97*
98* IOFFD = 0 specifies the main diagonal A( i, i ),
99* with i = 1 ... MIN( M, N ),
100* IOFFD > 0 specifies the subdiagonal A( i+IOFFD, i ),
101* with i = 1 ... MIN( M-IOFFD, N ),
102* IOFFD < 0 specifies the superdiagonal A( i, i-IOFFD ),
103* with i = 1 ... MIN( M, N+IOFFD ).
104*
105* ALPHA (input) pointer to CHAR
106* On entry, ALPHA specifies the scalar alpha.
107*
108* A (input) pointer to CHAR
109* On entry, A is an array of dimension (LDA,N) containing the m
110* by n matrix A. Only the trapezoidal part of A determined by
111* UPLO and IOFFD is referenced.
112*
113* LDA (input) INTEGER
114* On entry, LDA specifies the leading dimension of the array A.
115* LDA must be at least max( 1, M ).
116*
117* BC (input) pointer to CHAR
118* On entry, BC is an array of dimension (LDBC,K) containing the
119* m by k matrix BC.
120*
121* LDBC (input) INTEGER
122* On entry, LDBC specifies the leading dimension of the array
123* BC. LDBC must be at least max( 1, M ).
124*
125* BR (input) pointer to CHAR
126* On entry, BR is an array of dimension (LDBR,N) containing the
127* k by n matrix BR.
128*
129* LDBR (input) INTEGER
130* On entry, LDBR specifies the leading dimension of the array
131* BR. LDBR must be at least K.
132*
133* CC (input/output) pointer to CHAR
134* On entry, CC is an array of dimension (LDCC,K) containing the
135* m by k matrix CC. On exit, CC is overwritten by the partially
136* updated matric CC.
137*
138* LDCC (input) INTEGER
139* On entry, LDCC specifies the leading dimension of the array
140* CC. LDCC must be at least max( 1, M ).
141*
142* CR (input/output) pointer to CHAR
143* On entry, CR is an array of dimension (LDCR,N) containing the
144* k by n matrix CR. On exit, CR is overwritten by the partially
145* updated matrix CR.
146*
147* LDCR (input) INTEGER
148* On entry, LDCR specifies the leading dimension of the array
149* CR. LDCR must be at least K.
150*
151* Notes
152* =====
153* N N
154* ---------------------------- -----------
155* | d | | |
156* M | d Upper | | Upper |
157* | Lower d | |d |
158* | d | M | d |
159* ---------------------------- | d |
160* | d |
161* IOFFD < 0 | Lower d |
162* | d|
163* N | |
164* ----------- -----------
165* | d Upper|
166* | d | IOFFD > 0
167* M | d |
168* | d| N
169* | Lower | ----------------------------
170* | | | Upper |
171* | | |d |
172* | | | d |
173* | | | d |
174* | | |Lower d |
175* ----------- ----------------------------
176*
177* -- Written on April 1, 1998 by
178* Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
179*
180* ---------------------------------------------------------------------
181*/
182/*
183* .. Local Scalars ..
184*/
185 char * Calph, * one, type;
186 Int i1, j1, m1, mn, n1, size;
187 cmplx Calph8;
188 cmplx16 Calph16;
189 GEMM_T gemm;
190/* ..
191* .. Executable Statements ..
192*
193*/
194 if( ( M <= 0 ) || ( N <= 0 ) ) return;
195
196 type = TYPE->type;
197
198 if( type == SCPLX )
199 { Calph = ( (char *)(Calph8 ) ); PB_Cconjg( TYPE, ALPHA, Calph ); }
200 else if( type == DCPLX )
201 { Calph = ( (char *)(Calph16) ); PB_Cconjg( TYPE, ALPHA, Calph ); }
202 else
203 { Calph = ALPHA; }
204
205 if( Mupcase( SIDE[0] ) == CLEFT )
206 {
207 if( Mupcase( UPLO[0] ) == CLOWER )
208 {
209 size = TYPE->size; one = TYPE->one; gemm = TYPE->Fgemm;
210 mn = MAX( 0, -IOFFD );
211 if( ( n1 = MIN( mn, N ) ) > 0 )
212 {
213 gemm( C2F_CHAR( NOTRAN ), C2F_CHAR( TRAN ), &M, &K, &n1, ALPHA,
214 A, &LDA, BR, &LDBR, one, CC, &LDCC );
215 gemm( C2F_CHAR( COTRAN ), C2F_CHAR( NOTRAN ), &K, &n1, &M, Calph,
216 BC, &LDBC, A, &LDA, one, CR, &LDCR );
217 }
218 n1 = M - IOFFD;
219 if( ( n1 = MIN( n1, N ) - mn ) > 0 )
220 {
221 i1 = ( j1 = mn ) + IOFFD;
222 TYPE->Fhemm( C2F_CHAR( SIDE ), C2F_CHAR( UPLO ), &n1, &K, ALPHA,
223 Mptr( A, i1, j1, LDA, size ), &LDA, Mptr( BC, i1, 0,
224 LDBC, size ), &LDBC, one, Mptr( CC, i1, 0, LDCC,
225 size ), &LDCC );
226 if( ( m1 = M - mn - n1 - IOFFD ) > 0 )
227 {
228 i1 += n1;
229 gemm( C2F_CHAR( NOTRAN ), C2F_CHAR( TRAN ), &m1, &K, &n1, ALPHA,
230 Mptr( A, i1, j1, LDA, size ), &LDA, Mptr( BR, 0, j1, LDBR,
231 size ), &LDBR, one, Mptr( CC, i1, 0, LDCC, size ), &LDCC );
232 gemm( C2F_CHAR( COTRAN ), C2F_CHAR( NOTRAN ), &K, &n1, &m1,
233 Calph, Mptr( BC, i1, 0, LDBC, size ), &LDBC, Mptr( A, i1,
234 j1, LDA, size ), &LDA, one, Mptr( CR, 0, j1, LDCR, size ),
235 &LDCR );
236 }
237 }
238 }
239 else if( Mupcase( UPLO[0] ) == CUPPER )
240 {
241 size = TYPE->size; one = TYPE->one; gemm = TYPE->Fgemm;
242 mn = MIN( M - IOFFD, N );
243 if( ( n1 = mn - MAX( 0, -IOFFD ) ) > 0 )
244 {
245 j1 = mn - n1;
246 if( ( m1 = MAX( 0, IOFFD ) ) > 0 )
247 {
248 gemm( C2F_CHAR( NOTRAN ), C2F_CHAR( TRAN ), &m1, &K, &n1, ALPHA,
249 A, &LDA, BR, &LDBR, one, CC, &LDCC );
250 gemm( C2F_CHAR( COTRAN ), C2F_CHAR( NOTRAN ), &K, &n1, &m1,
251 Calph, BC, &LDBC, A, &LDA, one, CR, &LDCR );
252 }
253 TYPE->Fhemm( C2F_CHAR( SIDE ), C2F_CHAR( UPLO ), &n1, &K, ALPHA,
254 Mptr( A, m1, j1, LDA, size ), &LDA, Mptr( BC, m1, 0,
255 LDBC, size ), &LDBC, one, Mptr( CC, m1, 0, LDCC,
256 size ), &LDCC );
257 }
258 if( ( n1 = N - MAX( 0, mn ) ) > 0 )
259 {
260 j1 = N - n1;
261 gemm( C2F_CHAR( NOTRAN ), C2F_CHAR( TRAN ), &M, &K, &n1, ALPHA,
262 Mptr( A, 0, j1, LDA, size ), &LDA, Mptr( BR, 0, j1, LDBR,
263 size ), &LDBR, one, CC, &LDCC );
264 gemm( C2F_CHAR( COTRAN ), C2F_CHAR( NOTRAN ), &K, &n1, &M, Calph,
265 BC, &LDBC, Mptr( A, 0, j1, LDA, size ), &LDA, one, Mptr( CR,
266 0, j1, LDCR, size ), &LDCR );
267 }
268 }
269 else
270 {
271 one = TYPE->one; gemm = TYPE->Fgemm;
272 gemm( C2F_CHAR( NOTRAN ), C2F_CHAR( TRAN ), &M, &K, &N, ALPHA, A,
273 &LDA, BR, &LDBR, one, CC, &LDCC );
274 gemm( C2F_CHAR( COTRAN ), C2F_CHAR( NOTRAN ), &K, &N, &M, Calph, BC,
275 &LDBC, A, &LDA, one, CR, &LDCR );
276 }
277 }
278 else
279 {
280 if( Mupcase( UPLO[0] ) == CLOWER )
281 {
282 size = TYPE->size; one = TYPE->one; gemm = TYPE->Fgemm;
283 mn = MAX( 0, -IOFFD );
284 if( ( n1 = MIN( mn, N ) ) > 0 )
285 {
286 gemm( C2F_CHAR( NOTRAN ), C2F_CHAR( COTRAN ), &M, &K, &n1, Calph,
287 A, &LDA, BR, &LDBR, one, CC, &LDCC );
288 gemm( C2F_CHAR( TRAN ), C2F_CHAR( NOTRAN ), &K, &n1, &M, ALPHA,
289 BC, &LDBC, A, &LDA, one, CR, &LDCR );
290 }
291 n1 = M - IOFFD;
292 if( ( n1 = MIN( n1, N ) - mn ) > 0 )
293 {
294 i1 = ( j1 = mn ) + IOFFD;
295 TYPE->Fhemm( C2F_CHAR( SIDE ), C2F_CHAR( UPLO ), &K, &n1, ALPHA,
296 Mptr( A, i1, j1, LDA, size ), &LDA, Mptr( BR, 0, j1,
297 LDBR, size ), &LDBR, one, Mptr( CR, 0, j1, LDCR,
298 size ), &LDCR );
299 if( ( m1 = M - mn - n1 - IOFFD ) > 0 )
300 {
301 i1 += n1;
302 gemm( C2F_CHAR( NOTRAN ), C2F_CHAR( COTRAN ), &m1, &K, &n1,
303 Calph, Mptr( A, i1, j1, LDA, size ), &LDA, Mptr( BR, 0, j1,
304 LDBR, size ), &LDBR, one, Mptr( CC, i1, 0, LDCC, size ),
305 &LDCC );
306 gemm( C2F_CHAR( TRAN ), C2F_CHAR( NOTRAN ), &K, &n1, &m1, ALPHA,
307 Mptr( BC, i1, 0, LDBC, size ), &LDBC, Mptr( A, i1, j1, LDA,
308 size ), &LDA, one, Mptr( CR, 0, j1, LDCR, size ), &LDCR );
309 }
310 }
311 }
312 else if( Mupcase( UPLO[0] ) == CUPPER )
313 {
314 size = TYPE->size; one = TYPE->one; gemm = TYPE->Fgemm;
315 mn = MIN( M - IOFFD, N );
316 if( ( n1 = mn - MAX( 0, -IOFFD ) ) > 0 )
317 {
318 j1 = mn - n1;
319 if( ( m1 = MAX( 0, IOFFD ) ) > 0 )
320 {
321 gemm( C2F_CHAR( NOTRAN ), C2F_CHAR( COTRAN ), &m1, &K, &n1,
322 Calph, A, &LDA, BR, &LDBR, one, CC, &LDCC );
323 gemm( C2F_CHAR( TRAN ), C2F_CHAR( NOTRAN ), &K, &n1, &m1,
324 ALPHA, BC, &LDBC, A, &LDA, one, CR, &LDCR );
325 }
326 TYPE->Fhemm( C2F_CHAR( SIDE ), C2F_CHAR( UPLO ), &K, &n1, ALPHA,
327 Mptr( A, m1, j1, LDA, size ), &LDA, Mptr( BR, 0, j1,
328 LDBR, size ), &LDBR, one, Mptr( CR, 0, j1, LDCR,
329 size ), &LDCR );
330 }
331 if( ( n1 = N - MAX( 0, mn ) ) > 0 )
332 {
333 j1 = N - n1;
334 gemm( C2F_CHAR( NOTRAN ), C2F_CHAR( COTRAN ), &M, &K, &n1, Calph,
335 Mptr( A, 0, j1, LDA, size ), &LDA, Mptr( BR, 0, j1, LDBR,
336 size ), &LDBR, one, CC, &LDCC );
337 gemm( C2F_CHAR( TRAN ), C2F_CHAR( NOTRAN ), &K, &n1, &M, ALPHA, BC,
338 &LDBC, Mptr( A, 0, j1, LDA, size ), &LDA, one, Mptr( CR, 0,
339 j1, LDCR, size ), &LDCR );
340 }
341 }
342 else
343 {
344 one = TYPE->one; gemm = TYPE->Fgemm;
345 gemm( C2F_CHAR( NOTRAN ), C2F_CHAR( COTRAN ), &M, &K, &N, Calph, A,
346 &LDA, BR, &LDBR, one, CC, &LDCC );
347 gemm( C2F_CHAR( TRAN ), C2F_CHAR( NOTRAN ), &K, &N, &M, ALPHA, BC,
348 &LDBC, A, &LDA, one, CR, &LDCR );
349 }
350 }
351/*
352* End of PB_Ctzhemm
353*/
354}
Int
#define Int
Definition Bconfig.h:22
cmplx
float cmplx[2]
Definition pblas.h:136
GEMM_T
F_VOID_FCT(* GEMM_T)()
Definition pblas.h:317
cmplx16
double cmplx16[2]
Definition pblas.h:137
SCPLX
#define SCPLX
Definition pblas.h:475
C2F_CHAR
#define C2F_CHAR(a)
Definition pblas.h:125
DCPLX
#define DCPLX
Definition pblas.h:476
NOTRAN
#define NOTRAN
Definition PBblas.h:44
CLEFT
#define CLEFT
Definition PBblas.h:29
TRAN
#define TRAN
Definition PBblas.h:46
CUPPER
#define CUPPER
Definition PBblas.h:26
COTRAN
#define COTRAN
Definition PBblas.h:48
CLOWER
#define CLOWER
Definition PBblas.h:25
PB_Ctzhemm
void PB_Ctzhemm()
Mptr
#define Mptr(a_, i_, j_, lda_, siz_)
Definition PBtools.h:132
Mupcase
#define Mupcase(C)
Definition PBtools.h:83
PB_Cconjg
void PB_Cconjg()
TYPE
#define TYPE
Definition clamov.c:7
BR
#define BR
Definition const.h:73
N
#define N
Definition example_user.c:10
LDA
#define LDA
Definition example_user.c:12
MIN
#define MIN(x, y)
Definition lapacke_utils.h:49
MAX
#define MAX(x, y)
Definition lapacke_utils.h:46
CR
#define CR
Definition lf_convert_c.c:62
PBTYP_T
Definition pblas.h:330