LAPACK  3.11.0
LAPACK: Linear Algebra PACKage
c_sblas2.c
1 /*
2  * Written by D.P. Manley, Digital Equipment Corporation.
3  * Prefixed "C_" to BLAS routines and their declarations.
4  *
5  * Modified by T. H. Do, 1/23/98, SGI/CRAY Research.
6  */
7 #include <stdlib.h>
8 #include "cblas.h"
9 #include "cblas_test.h"
10 
11 void F77_sgemv(CBLAS_INT *layout, char *transp, CBLAS_INT *m, CBLAS_INT *n, float *alpha,
12  float *a, CBLAS_INT *lda, float *x, CBLAS_INT *incx, float *beta,
13  float *y, CBLAS_INT *incy ) {
14 
15  float *A;
16  CBLAS_INT i,j,LDA;
17  CBLAS_TRANSPOSE trans;
18 
19  get_transpose_type(transp, &trans);
20  if (*layout == TEST_ROW_MJR) {
21  LDA = *n+1;
22  A = ( float* )malloc( (*m)*LDA*sizeof( float ) );
23  for( i=0; i<*m; i++ )
24  for( j=0; j<*n; j++ )
25  A[ LDA*i+j ]=a[ (*lda)*j+i ];
26  cblas_sgemv( CblasRowMajor, trans,
27  *m, *n, *alpha, A, LDA, x, *incx, *beta, y, *incy );
28  free(A);
29  }
30  else if (*layout == TEST_COL_MJR)
31  cblas_sgemv( CblasColMajor, trans,
32  *m, *n, *alpha, a, *lda, x, *incx, *beta, y, *incy );
33  else
34  cblas_sgemv( UNDEFINED, trans,
35  *m, *n, *alpha, a, *lda, x, *incx, *beta, y, *incy );
36 }
37 
38 void F77_sger(CBLAS_INT *layout, CBLAS_INT *m, CBLAS_INT *n, float *alpha, float *x, CBLAS_INT *incx,
39  float *y, CBLAS_INT *incy, float *a, CBLAS_INT *lda ) {
40 
41  float *A;
42  CBLAS_INT i,j,LDA;
43 
44  if (*layout == TEST_ROW_MJR) {
45  LDA = *n+1;
46  A = ( float* )malloc( (*m)*LDA*sizeof( float ) );
47 
48  for( i=0; i<*m; i++ ) {
49  for( j=0; j<*n; j++ )
50  A[ LDA*i+j ]=a[ (*lda)*j+i ];
51  }
52 
53  cblas_sger(CblasRowMajor, *m, *n, *alpha, x, *incx, y, *incy, A, LDA );
54  for( i=0; i<*m; i++ )
55  for( j=0; j<*n; j++ )
56  a[ (*lda)*j+i ]=A[ LDA*i+j ];
57  free(A);
58  }
59  else
60  cblas_sger( CblasColMajor, *m, *n, *alpha, x, *incx, y, *incy, a, *lda );
61 }
62 
63 void F77_strmv(CBLAS_INT *layout, char *uplow, char *transp, char *diagn,
64  CBLAS_INT *n, float *a, CBLAS_INT *lda, float *x, CBLAS_INT *incx) {
65  float *A;
66  CBLAS_INT i,j,LDA;
67  CBLAS_TRANSPOSE trans;
68  CBLAS_UPLO uplo;
69  CBLAS_DIAG diag;
70 
71  get_transpose_type(transp,&trans);
72  get_uplo_type(uplow,&uplo);
73  get_diag_type(diagn,&diag);
74 
75  if (*layout == TEST_ROW_MJR) {
76  LDA = *n+1;
77  A = ( float* )malloc( (*n)*LDA*sizeof( float ) );
78  for( i=0; i<*n; i++ )
79  for( j=0; j<*n; j++ )
80  A[ LDA*i+j ]=a[ (*lda)*j+i ];
81  cblas_strmv(CblasRowMajor, uplo, trans, diag, *n, A, LDA, x, *incx);
82  free(A);
83  }
84  else if (*layout == TEST_COL_MJR)
85  cblas_strmv(CblasColMajor, uplo, trans, diag, *n, a, *lda, x, *incx);
86  else {
87  cblas_strmv(UNDEFINED, uplo, trans, diag, *n, a, *lda, x, *incx);
88  }
89 }
90 
91 void F77_strsv(CBLAS_INT *layout, char *uplow, char *transp, char *diagn,
92  CBLAS_INT *n, float *a, CBLAS_INT *lda, float *x, CBLAS_INT *incx ) {
93  float *A;
94  CBLAS_INT i,j,LDA;
95  CBLAS_TRANSPOSE trans;
96  CBLAS_UPLO uplo;
97  CBLAS_DIAG diag;
98 
99  get_transpose_type(transp,&trans);
100  get_uplo_type(uplow,&uplo);
101  get_diag_type(diagn,&diag);
102 
103  if (*layout == TEST_ROW_MJR) {
104  LDA = *n+1;
105  A = ( float* )malloc( (*n)*LDA*sizeof( float ) );
106  for( i=0; i<*n; i++ )
107  for( j=0; j<*n; j++ )
108  A[ LDA*i+j ]=a[ (*lda)*j+i ];
109  cblas_strsv(CblasRowMajor, uplo, trans, diag, *n, A, LDA, x, *incx );
110  free(A);
111  }
112  else
113  cblas_strsv(CblasColMajor, uplo, trans, diag, *n, a, *lda, x, *incx );
114 }
115 void F77_ssymv(CBLAS_INT *layout, char *uplow, CBLAS_INT *n, float *alpha, float *a,
116  CBLAS_INT *lda, float *x, CBLAS_INT *incx, float *beta, float *y,
117  CBLAS_INT *incy) {
118  float *A;
119  CBLAS_INT i,j,LDA;
120  CBLAS_UPLO uplo;
121 
122  get_uplo_type(uplow,&uplo);
123 
124  if (*layout == TEST_ROW_MJR) {
125  LDA = *n+1;
126  A = ( float* )malloc( (*n)*LDA*sizeof( float ) );
127  for( i=0; i<*n; i++ )
128  for( j=0; j<*n; j++ )
129  A[ LDA*i+j ]=a[ (*lda)*j+i ];
130  cblas_ssymv(CblasRowMajor, uplo, *n, *alpha, A, LDA, x, *incx,
131  *beta, y, *incy );
132  free(A);
133  }
134  else
135  cblas_ssymv(CblasColMajor, uplo, *n, *alpha, a, *lda, x, *incx,
136  *beta, y, *incy );
137 }
138 
139 void F77_ssyr(CBLAS_INT *layout, char *uplow, CBLAS_INT *n, float *alpha, float *x,
140  CBLAS_INT *incx, float *a, CBLAS_INT *lda) {
141  float *A;
142  CBLAS_INT i,j,LDA;
143  CBLAS_UPLO uplo;
144 
145  get_uplo_type(uplow,&uplo);
146 
147  if (*layout == TEST_ROW_MJR) {
148  LDA = *n+1;
149  A = ( float* )malloc( (*n)*LDA*sizeof( float ) );
150  for( i=0; i<*n; i++ )
151  for( j=0; j<*n; j++ )
152  A[ LDA*i+j ]=a[ (*lda)*j+i ];
153  cblas_ssyr(CblasRowMajor, uplo, *n, *alpha, x, *incx, A, LDA);
154  for( i=0; i<*n; i++ )
155  for( j=0; j<*n; j++ )
156  a[ (*lda)*j+i ]=A[ LDA*i+j ];
157  free(A);
158  }
159  else
160  cblas_ssyr(CblasColMajor, uplo, *n, *alpha, x, *incx, a, *lda);
161 }
162 
163 void F77_ssyr2(CBLAS_INT *layout, char *uplow, CBLAS_INT *n, float *alpha, float *x,
164  CBLAS_INT *incx, float *y, CBLAS_INT *incy, float *a, CBLAS_INT *lda) {
165  float *A;
166  CBLAS_INT i,j,LDA;
167  CBLAS_UPLO uplo;
168 
169  get_uplo_type(uplow,&uplo);
170 
171  if (*layout == TEST_ROW_MJR) {
172  LDA = *n+1;
173  A = ( float* )malloc( (*n)*LDA*sizeof( float ) );
174  for( i=0; i<*n; i++ )
175  for( j=0; j<*n; j++ )
176  A[ LDA*i+j ]=a[ (*lda)*j+i ];
177  cblas_ssyr2(CblasRowMajor, uplo, *n, *alpha, x, *incx, y, *incy, A, LDA);
178  for( i=0; i<*n; i++ )
179  for( j=0; j<*n; j++ )
180  a[ (*lda)*j+i ]=A[ LDA*i+j ];
181  free(A);
182  }
183  else
184  cblas_ssyr2(CblasColMajor, uplo, *n, *alpha, x, *incx, y, *incy, a, *lda);
185 }
186 
187 void F77_sgbmv(CBLAS_INT *layout, char *transp, CBLAS_INT *m, CBLAS_INT *n, CBLAS_INT *kl, CBLAS_INT *ku,
188  float *alpha, float *a, CBLAS_INT *lda, float *x, CBLAS_INT *incx,
189  float *beta, float *y, CBLAS_INT *incy ) {
190 
191  float *A;
192  CBLAS_INT i,irow,j,jcol,LDA;
193  CBLAS_TRANSPOSE trans;
194 
195  get_transpose_type(transp, &trans);
196 
197  if (*layout == TEST_ROW_MJR) {
198  LDA = *ku+*kl+2;
199  A = ( float* )malloc( (*n+*kl)*LDA*sizeof( float ) );
200  for( i=0; i<*ku; i++ ){
201  irow=*ku+*kl-i;
202  jcol=(*ku)-i;
203  for( j=jcol; j<*n; j++ )
204  A[ LDA*(j-jcol)+irow ]=a[ (*lda)*j+i ];
205  }
206  i=*ku;
207  irow=*ku+*kl-i;
208  for( j=0; j<*n; j++ )
209  A[ LDA*j+irow ]=a[ (*lda)*j+i ];
210  for( i=*ku+1; i<*ku+*kl+1; i++ ){
211  irow=*ku+*kl-i;
212  jcol=i-(*ku);
213  for( j=jcol; j<(*n+*kl); j++ )
214  A[ LDA*j+irow ]=a[ (*lda)*(j-jcol)+i ];
215  }
216  cblas_sgbmv( CblasRowMajor, trans, *m, *n, *kl, *ku, *alpha,
217  A, LDA, x, *incx, *beta, y, *incy );
218  free(A);
219  }
220  else
221  cblas_sgbmv( CblasColMajor, trans, *m, *n, *kl, *ku, *alpha,
222  a, *lda, x, *incx, *beta, y, *incy );
223 }
224 
225 void F77_stbmv(CBLAS_INT *layout, char *uplow, char *transp, char *diagn,
226  CBLAS_INT *n, CBLAS_INT *k, float *a, CBLAS_INT *lda, float *x, CBLAS_INT *incx) {
227  float *A;
228  CBLAS_INT irow, jcol, i, j, LDA;
229  CBLAS_TRANSPOSE trans;
230  CBLAS_UPLO uplo;
231  CBLAS_DIAG diag;
232 
233  get_transpose_type(transp,&trans);
234  get_uplo_type(uplow,&uplo);
235  get_diag_type(diagn,&diag);
236 
237  if (*layout == TEST_ROW_MJR) {
238  LDA = *k+1;
239  A = ( float* )malloc( (*n+*k)*LDA*sizeof( float ) );
240  if (uplo == CblasUpper) {
241  for( i=0; i<*k; i++ ){
242  irow=*k-i;
243  jcol=(*k)-i;
244  for( j=jcol; j<*n; j++ )
245  A[ LDA*(j-jcol)+irow ]=a[ (*lda)*j+i ];
246  }
247  i=*k;
248  irow=*k-i;
249  for( j=0; j<*n; j++ )
250  A[ LDA*j+irow ]=a[ (*lda)*j+i ];
251  }
252  else {
253  i=0;
254  irow=*k-i;
255  for( j=0; j<*n; j++ )
256  A[ LDA*j+irow ]=a[ (*lda)*j+i ];
257  for( i=1; i<*k+1; i++ ){
258  irow=*k-i;
259  jcol=i;
260  for( j=jcol; j<(*n+*k); j++ )
261  A[ LDA*j+irow ]=a[ (*lda)*(j-jcol)+i ];
262  }
263  }
264  cblas_stbmv(CblasRowMajor, uplo, trans, diag, *n, *k, A, LDA, x, *incx);
265  free(A);
266  }
267  else
268  cblas_stbmv(CblasColMajor, uplo, trans, diag, *n, *k, a, *lda, x, *incx);
269 }
270 
271 void F77_stbsv(CBLAS_INT *layout, char *uplow, char *transp, char *diagn,
272  CBLAS_INT *n, CBLAS_INT *k, float *a, CBLAS_INT *lda, float *x, CBLAS_INT *incx) {
273  float *A;
274  CBLAS_INT irow, jcol, i, j, LDA;
275  CBLAS_TRANSPOSE trans;
276  CBLAS_UPLO uplo;
277  CBLAS_DIAG diag;
278 
279  get_transpose_type(transp,&trans);
280  get_uplo_type(uplow,&uplo);
281  get_diag_type(diagn,&diag);
282 
283  if (*layout == TEST_ROW_MJR) {
284  LDA = *k+1;
285  A = ( float* )malloc( (*n+*k)*LDA*sizeof( float ) );
286  if (uplo == CblasUpper) {
287  for( i=0; i<*k; i++ ){
288  irow=*k-i;
289  jcol=(*k)-i;
290  for( j=jcol; j<*n; j++ )
291  A[ LDA*(j-jcol)+irow ]=a[ (*lda)*j+i ];
292  }
293  i=*k;
294  irow=*k-i;
295  for( j=0; j<*n; j++ )
296  A[ LDA*j+irow ]=a[ (*lda)*j+i ];
297  }
298  else {
299  i=0;
300  irow=*k-i;
301  for( j=0; j<*n; j++ )
302  A[ LDA*j+irow ]=a[ (*lda)*j+i ];
303  for( i=1; i<*k+1; i++ ){
304  irow=*k-i;
305  jcol=i;
306  for( j=jcol; j<(*n+*k); j++ )
307  A[ LDA*j+irow ]=a[ (*lda)*(j-jcol)+i ];
308  }
309  }
310  cblas_stbsv(CblasRowMajor, uplo, trans, diag, *n, *k, A, LDA, x, *incx);
311  free(A);
312  }
313  else
314  cblas_stbsv(CblasColMajor, uplo, trans, diag, *n, *k, a, *lda, x, *incx);
315 }
316 
317 void F77_ssbmv(CBLAS_INT *layout, char *uplow, CBLAS_INT *n, CBLAS_INT *k, float *alpha,
318  float *a, CBLAS_INT *lda, float *x, CBLAS_INT *incx, float *beta,
319  float *y, CBLAS_INT *incy) {
320  float *A;
321  CBLAS_INT i,j,irow,jcol,LDA;
322  CBLAS_UPLO uplo;
323 
324  get_uplo_type(uplow,&uplo);
325 
326  if (*layout == TEST_ROW_MJR) {
327  LDA = *k+1;
328  A = ( float* )malloc( (*n+*k)*LDA*sizeof( float ) );
329  if (uplo == CblasUpper) {
330  for( i=0; i<*k; i++ ){
331  irow=*k-i;
332  jcol=(*k)-i;
333  for( j=jcol; j<*n; j++ )
334  A[ LDA*(j-jcol)+irow ]=a[ (*lda)*j+i ];
335  }
336  i=*k;
337  irow=*k-i;
338  for( j=0; j<*n; j++ )
339  A[ LDA*j+irow ]=a[ (*lda)*j+i ];
340  }
341  else {
342  i=0;
343  irow=*k-i;
344  for( j=0; j<*n; j++ )
345  A[ LDA*j+irow ]=a[ (*lda)*j+i ];
346  for( i=1; i<*k+1; i++ ){
347  irow=*k-i;
348  jcol=i;
349  for( j=jcol; j<(*n+*k); j++ )
350  A[ LDA*j+irow ]=a[ (*lda)*(j-jcol)+i ];
351  }
352  }
353  cblas_ssbmv(CblasRowMajor, uplo, *n, *k, *alpha, A, LDA, x, *incx,
354  *beta, y, *incy );
355  free(A);
356  }
357  else
358  cblas_ssbmv(CblasColMajor, uplo, *n, *k, *alpha, a, *lda, x, *incx,
359  *beta, y, *incy );
360 }
361 
362 void F77_sspmv(CBLAS_INT *layout, char *uplow, CBLAS_INT *n, float *alpha, float *ap,
363  float *x, CBLAS_INT *incx, float *beta, float *y, CBLAS_INT *incy) {
364  float *A,*AP;
365  CBLAS_INT i,j,k,LDA;
366  CBLAS_UPLO uplo;
367 
368  get_uplo_type(uplow,&uplo);
369 
370  if (*layout == TEST_ROW_MJR) {
371  LDA = *n;
372  A = ( float* )malloc( LDA*LDA*sizeof( float ) );
373  AP = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
374  if (uplo == CblasUpper) {
375  for( j=0, k=0; j<*n; j++ )
376  for( i=0; i<j+1; i++, k++ )
377  A[ LDA*i+j ]=ap[ k ];
378  for( i=0, k=0; i<*n; i++ )
379  for( j=i; j<*n; j++, k++ )
380  AP[ k ]=A[ LDA*i+j ];
381  }
382  else {
383  for( j=0, k=0; j<*n; j++ )
384  for( i=j; i<*n; i++, k++ )
385  A[ LDA*i+j ]=ap[ k ];
386  for( i=0, k=0; i<*n; i++ )
387  for( j=0; j<i+1; j++, k++ )
388  AP[ k ]=A[ LDA*i+j ];
389  }
390  cblas_sspmv( CblasRowMajor, uplo, *n, *alpha, AP, x, *incx, *beta, y,
391  *incy );
392  free(A); free(AP);
393  }
394  else
395  cblas_sspmv( CblasColMajor, uplo, *n, *alpha, ap, x, *incx, *beta, y,
396  *incy );
397 }
398 
399 void F77_stpmv(CBLAS_INT *layout, char *uplow, char *transp, char *diagn,
400  CBLAS_INT *n, float *ap, float *x, CBLAS_INT *incx) {
401  float *A, *AP;
402  CBLAS_INT i, j, k, LDA;
403  CBLAS_TRANSPOSE trans;
404  CBLAS_UPLO uplo;
405  CBLAS_DIAG diag;
406 
407  get_transpose_type(transp,&trans);
408  get_uplo_type(uplow,&uplo);
409  get_diag_type(diagn,&diag);
410 
411  if (*layout == TEST_ROW_MJR) {
412  LDA = *n;
413  A = ( float* )malloc( LDA*LDA*sizeof( float ) );
414  AP = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
415  if (uplo == CblasUpper) {
416  for( j=0, k=0; j<*n; j++ )
417  for( i=0; i<j+1; i++, k++ )
418  A[ LDA*i+j ]=ap[ k ];
419  for( i=0, k=0; i<*n; i++ )
420  for( j=i; j<*n; j++, k++ )
421  AP[ k ]=A[ LDA*i+j ];
422  }
423  else {
424  for( j=0, k=0; j<*n; j++ )
425  for( i=j; i<*n; i++, k++ )
426  A[ LDA*i+j ]=ap[ k ];
427  for( i=0, k=0; i<*n; i++ )
428  for( j=0; j<i+1; j++, k++ )
429  AP[ k ]=A[ LDA*i+j ];
430  }
431  cblas_stpmv( CblasRowMajor, uplo, trans, diag, *n, AP, x, *incx );
432  free(A); free(AP);
433  }
434  else
435  cblas_stpmv( CblasColMajor, uplo, trans, diag, *n, ap, x, *incx );
436 }
437 
438 void F77_stpsv(CBLAS_INT *layout, char *uplow, char *transp, char *diagn,
439  CBLAS_INT *n, float *ap, float *x, CBLAS_INT *incx) {
440  float *A, *AP;
441  CBLAS_INT i, j, k, LDA;
442  CBLAS_TRANSPOSE trans;
443  CBLAS_UPLO uplo;
444  CBLAS_DIAG diag;
445 
446  get_transpose_type(transp,&trans);
447  get_uplo_type(uplow,&uplo);
448  get_diag_type(diagn,&diag);
449 
450  if (*layout == TEST_ROW_MJR) {
451  LDA = *n;
452  A = ( float* )malloc( LDA*LDA*sizeof( float ) );
453  AP = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
454  if (uplo == CblasUpper) {
455  for( j=0, k=0; j<*n; j++ )
456  for( i=0; i<j+1; i++, k++ )
457  A[ LDA*i+j ]=ap[ k ];
458  for( i=0, k=0; i<*n; i++ )
459  for( j=i; j<*n; j++, k++ )
460  AP[ k ]=A[ LDA*i+j ];
461 
462  }
463  else {
464  for( j=0, k=0; j<*n; j++ )
465  for( i=j; i<*n; i++, k++ )
466  A[ LDA*i+j ]=ap[ k ];
467  for( i=0, k=0; i<*n; i++ )
468  for( j=0; j<i+1; j++, k++ )
469  AP[ k ]=A[ LDA*i+j ];
470  }
471  cblas_stpsv( CblasRowMajor, uplo, trans, diag, *n, AP, x, *incx );
472  free(A); free(AP);
473  }
474  else
475  cblas_stpsv( CblasColMajor, uplo, trans, diag, *n, ap, x, *incx );
476 }
477 
478 void F77_sspr(CBLAS_INT *layout, char *uplow, CBLAS_INT *n, float *alpha, float *x,
479  CBLAS_INT *incx, float *ap ){
480  float *A, *AP;
481  CBLAS_INT i,j,k,LDA;
482  CBLAS_UPLO uplo;
483 
484  get_uplo_type(uplow,&uplo);
485 
486  if (*layout == TEST_ROW_MJR) {
487  LDA = *n;
488  A = ( float* )malloc( LDA*LDA*sizeof( float ) );
489  AP = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
490  if (uplo == CblasUpper) {
491  for( j=0, k=0; j<*n; j++ )
492  for( i=0; i<j+1; i++, k++ )
493  A[ LDA*i+j ]=ap[ k ];
494  for( i=0, k=0; i<*n; i++ )
495  for( j=i; j<*n; j++, k++ )
496  AP[ k ]=A[ LDA*i+j ];
497  }
498  else {
499  for( j=0, k=0; j<*n; j++ )
500  for( i=j; i<*n; i++, k++ )
501  A[ LDA*i+j ]=ap[ k ];
502  for( i=0, k=0; i<*n; i++ )
503  for( j=0; j<i+1; j++, k++ )
504  AP[ k ]=A[ LDA*i+j ];
505  }
506  cblas_sspr( CblasRowMajor, uplo, *n, *alpha, x, *incx, AP );
507  if (uplo == CblasUpper) {
508  for( i=0, k=0; i<*n; i++ )
509  for( j=i; j<*n; j++, k++ )
510  A[ LDA*i+j ]=AP[ k ];
511  for( j=0, k=0; j<*n; j++ )
512  for( i=0; i<j+1; i++, k++ )
513  ap[ k ]=A[ LDA*i+j ];
514  }
515  else {
516  for( i=0, k=0; i<*n; i++ )
517  for( j=0; j<i+1; j++, k++ )
518  A[ LDA*i+j ]=AP[ k ];
519  for( j=0, k=0; j<*n; j++ )
520  for( i=j; i<*n; i++, k++ )
521  ap[ k ]=A[ LDA*i+j ];
522  }
523  free(A); free(AP);
524  }
525  else
526  cblas_sspr( CblasColMajor, uplo, *n, *alpha, x, *incx, ap );
527 }
528 
529 void F77_sspr2(CBLAS_INT *layout, char *uplow, CBLAS_INT *n, float *alpha, float *x,
530  CBLAS_INT *incx, float *y, CBLAS_INT *incy, float *ap ){
531  float *A, *AP;
532  CBLAS_INT i,j,k,LDA;
533  CBLAS_UPLO uplo;
534 
535  get_uplo_type(uplow,&uplo);
536 
537  if (*layout == TEST_ROW_MJR) {
538  LDA = *n;
539  A = ( float* )malloc( LDA*LDA*sizeof( float ) );
540  AP = ( float* )malloc( (((LDA+1)*LDA)/2)*sizeof( float ) );
541  if (uplo == CblasUpper) {
542  for( j=0, k=0; j<*n; j++ )
543  for( i=0; i<j+1; i++, k++ )
544  A[ LDA*i+j ]=ap[ k ];
545  for( i=0, k=0; i<*n; i++ )
546  for( j=i; j<*n; j++, k++ )
547  AP[ k ]=A[ LDA*i+j ];
548  }
549  else {
550  for( j=0, k=0; j<*n; j++ )
551  for( i=j; i<*n; i++, k++ )
552  A[ LDA*i+j ]=ap[ k ];
553  for( i=0, k=0; i<*n; i++ )
554  for( j=0; j<i+1; j++, k++ )
555  AP[ k ]=A[ LDA*i+j ];
556  }
557  cblas_sspr2( CblasRowMajor, uplo, *n, *alpha, x, *incx, y, *incy, AP );
558  if (uplo == CblasUpper) {
559  for( i=0, k=0; i<*n; i++ )
560  for( j=i; j<*n; j++, k++ )
561  A[ LDA*i+j ]=AP[ k ];
562  for( j=0, k=0; j<*n; j++ )
563  for( i=0; i<j+1; i++, k++ )
564  ap[ k ]=A[ LDA*i+j ];
565  }
566  else {
567  for( i=0, k=0; i<*n; i++ )
568  for( j=0; j<i+1; j++, k++ )
569  A[ LDA*i+j ]=AP[ k ];
570  for( j=0, k=0; j<*n; j++ )
571  for( i=j; i<*n; i++, k++ )
572  ap[ k ]=A[ LDA*i+j ];
573  }
574  free(A);
575  free(AP);
576  }
577  else
578  cblas_sspr2( CblasColMajor, uplo, *n, *alpha, x, *incx, y, *incy, ap );
579 }