1 #include <stdio.h>
2 #include "blas_extended.h"
ssy_commit_row(enum blas_order_type order,enum blas_uplo_type uplo,int n,float * a,int lda,float * a_vec,int row)3 void ssy_commit_row(enum blas_order_type order, enum blas_uplo_type uplo,
4 		    int n, float *a, int lda, float *a_vec, int row)
5 
6 /*
7  *  Copies the given vector into the given row of symmetric matrix a.
8  */
9 {
10   int ai;
11   int i;
12   int incai1;
13   int incai2;
14   int vi;
15   int incvi = 1;
16 
17   float a_elem;
18   float *a_i = a;
19   const float *a_vec_i = a_vec;
20 
21   if ((order == blas_colmajor && uplo == blas_upper) ||
22       (order == blas_rowmajor && uplo == blas_lower)) {
23     incai1 = 1;
24     incai2 = lda;
25     ai = lda * row;
26   } else {
27     incai1 = lda;
28     incai2 = 1;
29     ai = row;
30   }
31 
32 
33 
34 
35 
36 
37   for (i = 0, vi = 0; i < row; i++, ai += incai1, vi += incvi) {
38     a_elem = a_vec_i[vi];
39     a_i[ai] = a_elem;
40   }
41 
42   for (; i < n; i++, ai += incai2, vi += incvi) {
43     a_elem = a_vec_i[vi];
44     a_i[ai] = a_elem;
45   }
46 }
dsy_commit_row(enum blas_order_type order,enum blas_uplo_type uplo,int n,double * a,int lda,double * a_vec,int row)47 void dsy_commit_row(enum blas_order_type order, enum blas_uplo_type uplo,
48 		    int n, double *a, int lda, double *a_vec, int row)
49 
50 /*
51  *  Copies the given vector into the given row of symmetric matrix a.
52  */
53 {
54   int ai;
55   int i;
56   int incai1;
57   int incai2;
58   int vi;
59   int incvi = 1;
60 
61   double a_elem;
62   double *a_i = a;
63   const double *a_vec_i = a_vec;
64 
65   if ((order == blas_colmajor && uplo == blas_upper) ||
66       (order == blas_rowmajor && uplo == blas_lower)) {
67     incai1 = 1;
68     incai2 = lda;
69     ai = lda * row;
70   } else {
71     incai1 = lda;
72     incai2 = 1;
73     ai = row;
74   }
75 
76 
77 
78 
79 
80 
81   for (i = 0, vi = 0; i < row; i++, ai += incai1, vi += incvi) {
82     a_elem = a_vec_i[vi];
83     a_i[ai] = a_elem;
84   }
85 
86   for (; i < n; i++, ai += incai2, vi += incvi) {
87     a_elem = a_vec_i[vi];
88     a_i[ai] = a_elem;
89   }
90 }
csy_commit_row(enum blas_order_type order,enum blas_uplo_type uplo,int n,void * a,int lda,void * a_vec,int row)91 void csy_commit_row(enum blas_order_type order, enum blas_uplo_type uplo,
92 		    int n, void *a, int lda, void *a_vec, int row)
93 
94 /*
95  *  Copies the given vector into the given row of symmetric matrix a.
96  */
97 {
98   int ai;
99   int i;
100   int incai1;
101   int incai2;
102   int vi;
103   int incvi = 1;
104 
105   float a_elem[2];
106   float *a_i = (float *) a;
107   const float *a_vec_i = (float *) a_vec;
108 
109   if ((order == blas_colmajor && uplo == blas_upper) ||
110       (order == blas_rowmajor && uplo == blas_lower)) {
111     incai1 = 1;
112     incai2 = lda;
113     ai = lda * row;
114   } else {
115     incai1 = lda;
116     incai2 = 1;
117     ai = row;
118   }
119 
120   incai1 *= 2;
121   incai2 *= 2;
122   ai *= 2;
123   incvi *= 2;
124 
125   for (i = 0, vi = 0; i < row; i++, ai += incai1, vi += incvi) {
126     a_elem[0] = a_vec_i[vi];
127     a_elem[1] = a_vec_i[vi + 1];
128     a_i[ai] = a_elem[0];
129     a_i[ai + 1] = a_elem[1];
130   }
131 
132   for (; i < n; i++, ai += incai2, vi += incvi) {
133     a_elem[0] = a_vec_i[vi];
134     a_elem[1] = a_vec_i[vi + 1];
135     a_i[ai] = a_elem[0];
136     a_i[ai + 1] = a_elem[1];
137   }
138 }
zsy_commit_row(enum blas_order_type order,enum blas_uplo_type uplo,int n,void * a,int lda,void * a_vec,int row)139 void zsy_commit_row(enum blas_order_type order, enum blas_uplo_type uplo,
140 		    int n, void *a, int lda, void *a_vec, int row)
141 
142 /*
143  *  Copies the given vector into the given row of symmetric matrix a.
144  */
145 {
146   int ai;
147   int i;
148   int incai1;
149   int incai2;
150   int vi;
151   int incvi = 1;
152 
153   double a_elem[2];
154   double *a_i = (double *) a;
155   const double *a_vec_i = (double *) a_vec;
156 
157   if ((order == blas_colmajor && uplo == blas_upper) ||
158       (order == blas_rowmajor && uplo == blas_lower)) {
159     incai1 = 1;
160     incai2 = lda;
161     ai = lda * row;
162   } else {
163     incai1 = lda;
164     incai2 = 1;
165     ai = row;
166   }
167 
168   incai1 *= 2;
169   incai2 *= 2;
170   ai *= 2;
171   incvi *= 2;
172 
173   for (i = 0, vi = 0; i < row; i++, ai += incai1, vi += incvi) {
174     a_elem[0] = a_vec_i[vi];
175     a_elem[1] = a_vec_i[vi + 1];
176     a_i[ai] = a_elem[0];
177     a_i[ai + 1] = a_elem[1];
178   }
179 
180   for (; i < n; i++, ai += incai2, vi += incvi) {
181     a_elem[0] = a_vec_i[vi];
182     a_elem[1] = a_vec_i[vi + 1];
183     a_i[ai] = a_elem[0];
184     a_i[ai + 1] = a_elem[1];
185   }
186 }
187 
ssy_copy_row(enum blas_order_type order,enum blas_uplo_type uplo,int n,float * a,int lda,float * a_vec,int row)188 void ssy_copy_row(enum blas_order_type order, enum blas_uplo_type uplo, int n,
189 		  float *a, int lda, float *a_vec, int row)
190 
191 /*
192  *  Copies the given row of matrix a into the supplied vector.
193  */
194 {
195   int ai;
196   int i;
197   int incai1;
198   int incai2;
199   int vi;
200   int incvi = 1;
201 
202   float a_elem;
203   const float *a_i = a;
204   float *a_vec_i = a_vec;
205 
206   if ((order == blas_colmajor && uplo == blas_upper) ||
207       (order == blas_rowmajor && uplo == blas_lower)) {
208     incai1 = 1;
209     incai2 = lda;
210     ai = lda * row;
211   } else {
212     incai1 = lda;
213     incai2 = 1;
214     ai = row;
215   }
216 
217 
218 
219 
220 
221 
222   for (i = 0, vi = 0; i < row; i++, ai += incai1, vi += incvi) {
223     a_elem = a_i[ai];
224     a_vec_i[vi] = a_elem;
225   }
226 
227   for (; i < n; i++, ai += incai2, vi += incvi) {
228     a_elem = a_i[ai];
229     a_vec_i[vi] = a_elem;
230   }
231 }
dsy_copy_row(enum blas_order_type order,enum blas_uplo_type uplo,int n,double * a,int lda,double * a_vec,int row)232 void dsy_copy_row(enum blas_order_type order, enum blas_uplo_type uplo, int n,
233 		  double *a, int lda, double *a_vec, int row)
234 
235 /*
236  *  Copies the given row of matrix a into the supplied vector.
237  */
238 {
239   int ai;
240   int i;
241   int incai1;
242   int incai2;
243   int vi;
244   int incvi = 1;
245 
246   double a_elem;
247   const double *a_i = a;
248   double *a_vec_i = a_vec;
249 
250   if ((order == blas_colmajor && uplo == blas_upper) ||
251       (order == blas_rowmajor && uplo == blas_lower)) {
252     incai1 = 1;
253     incai2 = lda;
254     ai = lda * row;
255   } else {
256     incai1 = lda;
257     incai2 = 1;
258     ai = row;
259   }
260 
261 
262 
263 
264 
265 
266   for (i = 0, vi = 0; i < row; i++, ai += incai1, vi += incvi) {
267     a_elem = a_i[ai];
268     a_vec_i[vi] = a_elem;
269   }
270 
271   for (; i < n; i++, ai += incai2, vi += incvi) {
272     a_elem = a_i[ai];
273     a_vec_i[vi] = a_elem;
274   }
275 }
csy_copy_row(enum blas_order_type order,enum blas_uplo_type uplo,int n,void * a,int lda,void * a_vec,int row)276 void csy_copy_row(enum blas_order_type order, enum blas_uplo_type uplo, int n,
277 		  void *a, int lda, void *a_vec, int row)
278 
279 /*
280  *  Copies the given row of matrix a into the supplied vector.
281  */
282 {
283   int ai;
284   int i;
285   int incai1;
286   int incai2;
287   int vi;
288   int incvi = 1;
289 
290   float a_elem[2];
291   const float *a_i = (float *) a;
292   float *a_vec_i = (float *) a_vec;
293 
294   if ((order == blas_colmajor && uplo == blas_upper) ||
295       (order == blas_rowmajor && uplo == blas_lower)) {
296     incai1 = 1;
297     incai2 = lda;
298     ai = lda * row;
299   } else {
300     incai1 = lda;
301     incai2 = 1;
302     ai = row;
303   }
304 
305   incai1 *= 2;
306   incai2 *= 2;
307   ai *= 2;
308   incvi *= 2;
309 
310   for (i = 0, vi = 0; i < row; i++, ai += incai1, vi += incvi) {
311     a_elem[0] = a_i[ai];
312     a_elem[1] = a_i[ai + 1];
313     a_vec_i[vi] = a_elem[0];
314     a_vec_i[vi + 1] = a_elem[1];
315   }
316 
317   for (; i < n; i++, ai += incai2, vi += incvi) {
318     a_elem[0] = a_i[ai];
319     a_elem[1] = a_i[ai + 1];
320     a_vec_i[vi] = a_elem[0];
321     a_vec_i[vi + 1] = a_elem[1];
322   }
323 }
zsy_copy_row(enum blas_order_type order,enum blas_uplo_type uplo,int n,void * a,int lda,void * a_vec,int row)324 void zsy_copy_row(enum blas_order_type order, enum blas_uplo_type uplo, int n,
325 		  void *a, int lda, void *a_vec, int row)
326 
327 /*
328  *  Copies the given row of matrix a into the supplied vector.
329  */
330 {
331   int ai;
332   int i;
333   int incai1;
334   int incai2;
335   int vi;
336   int incvi = 1;
337 
338   double a_elem[2];
339   const double *a_i = (double *) a;
340   double *a_vec_i = (double *) a_vec;
341 
342   if ((order == blas_colmajor && uplo == blas_upper) ||
343       (order == blas_rowmajor && uplo == blas_lower)) {
344     incai1 = 1;
345     incai2 = lda;
346     ai = lda * row;
347   } else {
348     incai1 = lda;
349     incai2 = 1;
350     ai = row;
351   }
352 
353   incai1 *= 2;
354   incai2 *= 2;
355   ai *= 2;
356   incvi *= 2;
357 
358   for (i = 0, vi = 0; i < row; i++, ai += incai1, vi += incvi) {
359     a_elem[0] = a_i[ai];
360     a_elem[1] = a_i[ai + 1];
361     a_vec_i[vi] = a_elem[0];
362     a_vec_i[vi + 1] = a_elem[1];
363   }
364 
365   for (; i < n; i++, ai += incai2, vi += incvi) {
366     a_elem[0] = a_i[ai];
367     a_elem[1] = a_i[ai + 1];
368     a_vec_i[vi] = a_elem[0];
369     a_vec_i[vi + 1] = a_elem[1];
370   }
371 }
372 
ssy_print_matrix(float * a,int n,int lda,enum blas_order_type order,enum blas_uplo_type uplo)373 void ssy_print_matrix(float *a, int n, int lda,
374 		      enum blas_order_type order, enum blas_uplo_type uplo)
375 {
376   int ai, aij;
377   int incai, incaij1, incaij2;
378   int i, j;
379 
380   float a_elem;
381   const float *a_i = a;
382 
383   if ((order == blas_colmajor && uplo == blas_upper) ||
384       (order == blas_rowmajor && uplo == blas_lower)) {
385     incai = lda;
386     incaij1 = 1;
387     incaij2 = lda;
388   } else {
389     incai = 1;
390     incaij1 = lda;
391     incaij2 = 1;
392   }
393 
394 
395 
396 
397 
398   for (i = 0, ai = 0; i < n; i++, ai += incai) {
399 
400     for (j = 0, aij = ai; j < i; j++, aij += incaij1) {
401       a_elem = a_i[aij];
402       printf("%16.8e", a_elem);
403     }
404     for (; j < n; j++, aij += incaij2) {
405       a_elem = a_i[aij];
406       printf("%16.8e", a_elem);
407     }
408     printf("\n");
409   }
410 
411 }
dsy_print_matrix(double * a,int n,int lda,enum blas_order_type order,enum blas_uplo_type uplo)412 void dsy_print_matrix(double *a, int n, int lda,
413 		      enum blas_order_type order, enum blas_uplo_type uplo)
414 {
415   int ai, aij;
416   int incai, incaij1, incaij2;
417   int i, j;
418 
419   double a_elem;
420   const double *a_i = a;
421 
422   if ((order == blas_colmajor && uplo == blas_upper) ||
423       (order == blas_rowmajor && uplo == blas_lower)) {
424     incai = lda;
425     incaij1 = 1;
426     incaij2 = lda;
427   } else {
428     incai = 1;
429     incaij1 = lda;
430     incaij2 = 1;
431   }
432 
433 
434 
435 
436 
437   for (i = 0, ai = 0; i < n; i++, ai += incai) {
438 
439     for (j = 0, aij = ai; j < i; j++, aij += incaij1) {
440       a_elem = a_i[aij];
441       printf("%24.16e", a_elem);
442     }
443     for (; j < n; j++, aij += incaij2) {
444       a_elem = a_i[aij];
445       printf("%24.16e", a_elem);
446     }
447     printf("\n");
448   }
449 
450 }
csy_print_matrix(void * a,int n,int lda,enum blas_order_type order,enum blas_uplo_type uplo)451 void csy_print_matrix(void *a, int n, int lda,
452 		      enum blas_order_type order, enum blas_uplo_type uplo)
453 {
454   int ai, aij;
455   int incai, incaij1, incaij2;
456   int i, j;
457 
458   float a_elem[2];
459   const float *a_i = (float *) a;
460 
461   if ((order == blas_colmajor && uplo == blas_upper) ||
462       (order == blas_rowmajor && uplo == blas_lower)) {
463     incai = lda;
464     incaij1 = 1;
465     incaij2 = lda;
466   } else {
467     incai = 1;
468     incaij1 = lda;
469     incaij2 = 1;
470   }
471 
472   incai *= 2;
473   incaij1 *= 2;
474   incaij2 *= 2;
475 
476   for (i = 0, ai = 0; i < n; i++, ai += incai) {
477 
478     for (j = 0, aij = ai; j < i; j++, aij += incaij1) {
479       a_elem[0] = a_i[aij];
480       a_elem[1] = a_i[aij + 1];
481       printf("(%16.8e, %16.8e)", a_elem[0], a_elem[1]);
482     }
483     for (; j < n; j++, aij += incaij2) {
484       a_elem[0] = a_i[aij];
485       a_elem[1] = a_i[aij + 1];
486       printf("(%16.8e, %16.8e)", a_elem[0], a_elem[1]);
487     }
488     printf("\n");
489   }
490 
491 }
zsy_print_matrix(void * a,int n,int lda,enum blas_order_type order,enum blas_uplo_type uplo)492 void zsy_print_matrix(void *a, int n, int lda,
493 		      enum blas_order_type order, enum blas_uplo_type uplo)
494 {
495   int ai, aij;
496   int incai, incaij1, incaij2;
497   int i, j;
498 
499   double a_elem[2];
500   const double *a_i = (double *) a;
501 
502   if ((order == blas_colmajor && uplo == blas_upper) ||
503       (order == blas_rowmajor && uplo == blas_lower)) {
504     incai = lda;
505     incaij1 = 1;
506     incaij2 = lda;
507   } else {
508     incai = 1;
509     incaij1 = lda;
510     incaij2 = 1;
511   }
512 
513   incai *= 2;
514   incaij1 *= 2;
515   incaij2 *= 2;
516 
517   for (i = 0, ai = 0; i < n; i++, ai += incai) {
518 
519     for (j = 0, aij = ai; j < i; j++, aij += incaij1) {
520       a_elem[0] = a_i[aij];
521       a_elem[1] = a_i[aij + 1];
522       printf("(%24.16e, %24.16e)", a_elem[0], a_elem[1]);
523     }
524     for (; j < n; j++, aij += incaij2) {
525       a_elem[0] = a_i[aij];
526       a_elem[1] = a_i[aij + 1];
527       printf("(%24.16e, %24.16e)", a_elem[0], a_elem[1]);
528     }
529     printf("\n");
530   }
531 
532 }
533