TrioCFD 1.9.9_beta
TrioCFD documentation
Loading...
Searching...
No Matches
blas1.h
1/****************************************************************************
2* Copyright (c) 2026, CEA
3* All rights reserved.
4*
5* Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
6* 1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
7* 2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
8* 3. Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.
9*
10* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
11* IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS;
12* OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
13*
14*****************************************************************************/
15
16#ifndef _BLAS1_H_
17#define _BLAS1_H_
18
19
20
21/*! @brief Forward declarations for BLAS level-1 routines.
22 *
23 */
24
25#include <arch.h>
26
27
28#ifndef _NO_BLAS1_
29extern "C"
30/*! @brief External C declarations for BLAS level-1 routines.
31 *
32 */
33{
34
35#ifdef F77_Majuscule
36 double F77DECLARE(DASUM)(const integer *n, const double* const dx, const integer *incx);
37
38
39 void F77DECLARE(DAXPY)(const integer *n, const double * const da, const double * const dx,
40 const integer *incx, double * const dy, const integer *incy);
41
42 void F77DECLARE(DCOPY)(const integer *n,
43 const double * const dx, const integer *incx,
44 double * const dy, const integer *incy);
45
46
47 double F77DECLARE(DDOT)(const integer *n, const double * const dx, const integer *incx,
48 const double * const dy, const integer *incy);
49
50 double F77DECLARE(DNRM2)(const integer *n, const double * const dx, const integer *incx);
51
52 void F77DECLARE(DROT)(const integer *n, double * const dx, const integer *incx, double * const dy,
53 const integer *incy, const double * const c, const double * const s);
54
55 void F77DECLARE(DROTG)(double * const da, double * const db, double * const c, double * const s);
56
57 void F77DECLARE(DSCAL)(const integer *n, double * const da, double * const dx, const integer *incx);
58
59 void F77DECLARE(DSWAP)(const integer *n, double * const dx, const integer *incx, double * const dy,
60 const integer *incy);
61
62 integer F77DECLARE(IDAMAX)(const integer *n, const double * const dx, const integer *incx);
63
64 float F77DECLARE(SDOT)(const integer *n, const float *cx,
65 const integer *incx, const float *cy, const integer *incy);
66
67 void F77DECLARE(SAXPY)(const integer *n, const float *da, const float *dx,
68 const integer *incx, float *dy, const integer *incy);
69
70 void F77DECLARE(SCOPY)(const integer *n, float *dx, const integer *incx,
71 float *dy, const integer *incy);
72
73 float F77DECLARE(SASUM)(const integer *n, const float *dx, const integer *incx);
74
75 float F77DECLARE(SNRM2)(const integer *n, const float *dx, const integer *incx);
76
77 void F77DECLARE(SSCAL)(const integer *n, const double * const da, float *dx,
78 const integer *incx);
79
80 integer F77DECLARE(ISAMAX)(const integer *n, const float * const dx, const integer *incx);
81
82 void F77DECLARE(SSWAP)(const integer *n, float *dx, const integer *incx,
83 float *dy, integer *incy);
84
85 void F77DECLARE(SROT)(const integer *n,
86 double * const dx,
87 const integer *incx,
88 double * const dy,
89 const integer *incy,
90 const double * const c,
91 const double * const s);
92
93 void F77DECLARE(SROTG)(double * const da,
94 double * const db,
95 double * const c,
96 double * const s);
97#else
98
99 double F77DECLARE(dasum)(const integer *n, const double* const dx, const integer *incx);
100
101
102 void F77DECLARE(daxpy)(const integer *n, const double * const da, const double * const dx,
103 const integer *incx, double * const dy, const integer *incy);
104
105 void F77DECLARE(dcopy)(const integer *n,
106 const double * const dx, const integer *incx,
107 double * const dy, const integer *incy);
108
109
110 double F77DECLARE(ddot)(const integer *n, const double * const dx, const integer *incx,
111 const double * const dy, const integer *incy);
112
113 double F77DECLARE(dnrm2)(const integer *n, const double * const dx, const integer *incx);
114
115 void F77DECLARE(drot)(const integer *n, double * const dx, const integer *incx, double * const dy,
116 const integer *incy, const double * const c, const double * const s);
117
118 void F77DECLARE(drotg)(double * const da, double * const db, double * const c, double * const s);
119
120 void F77DECLARE(dscal)(const integer *n, double * const da, double * const dx, const integer *incx);
121
122 void F77DECLARE(dswap)(const integer *n, double * const dx, const integer *incx, double * const dy,
123 const integer *incy);
124
125 integer F77DECLARE(idamax)(const integer *n, const double * const dx, const integer *incx);
126
127 float F77DECLARE(sdot)(const integer *n, const float *cx,
128 const integer *incx, const float *cy, const integer *incy);
129
130 void F77DECLARE(saxpy)(const integer *n, const float *da, const float *dx,
131 const integer *incx, float *dy, const integer *incy);
132
133 void F77DECLARE(scopy)(const integer *n, float *dx, const integer *incx,
134 float *dy, const integer *incy);
135
136 float F77DECLARE(sasum)(const integer *n, const float *dx, const integer *incx);
137
138 float F77DECLARE(snrm2)(const integer *n, const float *dx, const integer *incx);
139
140 void F77DECLARE(sscal)(const integer *n, const float *da, float *dx,
141 const integer *incx);
142
143 integer F77DECLARE(isamax)(const integer *n, const float * const dx, const integer *incx);
144
145 void F77DECLARE(sswap)(const integer *n, float *dx, const integer *incx,
146 float *dy, integer *incy);
147
148 void F77DECLARE(srot)(const integer *n,
149 double * const dx,
150 const integer *incx,
151 double * const dy,
152 const integer *incy,
153 const double * const c,
154 const double * const s);
155
156 void F77DECLARE(srotg)(double * const da,
157 double * const db,
158 double * const c,
159 double * const s);
160
161#endif
162}
163
164/*! @brief Wrapper: calls dasum (sum of absolute values of a vector).
165 *
166 */
167#ifndef F77_Majuscule
168inline double F77NAME(DASUM)(const integer *n,
169 const double* const dx,
170 const integer *incx)
171{
172 return F77NAME(dasum)(n ,dx ,incx);
173}
174
175
176
177/*! @brief Wrapper: calls daxpy (constant times a vector plus a vector).
178 *
179 */
180inline void F77NAME(DAXPY)(const integer *n,
181 const double * const da,
182 const double * const dx,
183 const integer *incx,
184 double * const dy,
185 const integer *incy)
186{
187 F77NAME(daxpy)(n ,da, dx ,incx, dy, incy);
188}
189
190
191/*! @brief Wrapper: calls dcopy (copies a vector to another vector).
192 *
193 */
194inline void F77NAME(DCOPY)(const integer *n,
195 const double * const dx,
196 const integer *incx,
197 double * const dy,
198 const integer *incy)
199{
200 F77NAME(dcopy)(n , dx , incx, dy, incy);
201}
202
203
204
205/*! @brief Wrapper: calls ddot (dot product of two vectors).
206 *
207 */
208inline double F77NAME(DDOT)(const integer *n,
209 const double * const dx,
210 const integer *incx,
211 const double * const dy,
212 const integer *incy)
213{
214 return F77NAME(ddot)(n , dx , incx, dy, incy);
215}
216
217
218/*! @brief Wrapper: calls dnrm2 (Euclidean norm of a vector).
219 *
220 */
221inline double F77NAME(DNRM2)(const integer *n,
222 const double * const dx,
223 const integer *incx)
224{
225 return F77NAME(dnrm2)(n , dx , incx);
226}
227
228
229/*! @brief Wrapper: calls drot (applies a plane rotation to two vectors).
230 *
231 */
232inline void F77NAME(DROT)(const integer *n,
233 double * const dx,
234 const integer *incx,
235 double * const dy,
236 const integer *incy,
237 const double * const c,
238 const double * const s)
239{
240 F77NAME(drot)(n , dx , incx, dy, incy, c, s);
241}
242
243
244/*! @brief Wrapper: calls drotg (constructs a Givens plane rotation).
245 *
246 */
247inline void F77NAME(DROTG)(double * const da,
248 double * const db,
249 double * const c,
250 double * const s)
251{
252 F77NAME(drotg)(da , db , c, s);
253}
254
255
256/*! @brief Wrapper: calls dscal (scales a vector by a scalar).
257 *
258 */
259inline void F77NAME(DSCAL)(const integer *n,
260 double * const da,
261 double * const dx,
262 const integer *incx)
263{
264 F77NAME(dscal)(n , da , dx, incx);
265}
266
267
268/*! @brief Wrapper: calls dswap (swaps two vectors).
269 *
270 */
271inline void F77NAME(DSWAP)(const integer *n,
272 double * const dx,
273 const integer *incx,
274 double * const dy,
275 const integer *incy)
276{
277 F77NAME(dswap)(n , dx , incx, dy, incy);
278}
279
280
281/*! @brief Wrapper: calls idamax (finds the index of the element with maximum absolute value).
282 *
283 */
284inline integer F77NAME(IDAMAX)(const integer *n,
285 const double * const dx,
286 const integer *incx)
287{
288 return F77NAME(idamax)(n , dx , incx);
289}
290
291
292/*! @brief Wrapper: calls sdot (single-precision dot product of two vectors).
293 *
294 */
295inline float F77NAME(SDOT)(const integer *n,
296 const float *dx,
297 const integer *incx,
298 const float *dy,
299 const integer *incy)
300{
301 return F77NAME(sdot)(n , dx , incx, dy, incy);
302}
303
304
305/*! @brief Wrapper: calls saxpy (single-precision constant times a vector plus a vector).
306 *
307 */
308inline void F77NAME(SAXPY)(const integer *n,
309 const float *da,
310 const float *dx,
311 const integer *incx,
312 float *dy,
313 const integer *incy)
314{
315 F77NAME(saxpy)(n , da, dx , incx, dy, incy);
316}
317
318
319/*! @brief Wrapper: calls scopy (single-precision copy of a vector).
320 *
321 */
322inline void F77NAME(SCOPY)(const integer *n,
323 float *dx,
324 const integer *incx,
325 float *dy,
326 const integer *incy)
327{
328 F77NAME(scopy)(n , dx , incx, dy, incy);
329}
330
331
332/*! @brief Wrapper: calls sasum (single-precision sum of absolute values).
333 *
334 */
335inline float F77NAME(SASUM)(const integer *n,
336 const float *dx,
337 const integer *incx)
338{
339 return F77NAME(sasum)(n , dx , incx);
340}
341
342
343/*! @brief Wrapper: calls snrm2 (single-precision Euclidean norm of a vector).
344 *
345 */
346inline float F77NAME(SNRM2)(const integer *n,
347 const float *dx,
348 const integer *incx)
349{
350 return F77NAME(snrm2)(n , dx , incx);
351}
352
353
354/*! @brief Wrapper: calls sscal (single-precision scale a vector by a constant).
355 *
356 */
357inline void F77NAME(SSCAL)(const integer *n,
358 const float *da,
359 float *dx,
360 const integer *incx)
361{
362 F77NAME(sscal)(n , da, dx , incx);
363}
364
365
366/*! @brief Wrapper: calls isamax (single-precision index of maximum absolute value element).
367 *
368 */
369inline integer F77NAME(ISAMAX)(const integer *n,
370 const float * const dx,
371 const integer *incx)
372{
373 return F77NAME(isamax)(n , dx , incx);
374}
375
376
377/*! @brief Wrapper: calls sswap (single-precision swap of two vectors).
378 *
379 */
380inline void F77NAME(SSWAP)(const integer *n,
381 float *dx,
382 const integer *incx,
383 float *dy,
384 integer *incy)
385{
386 F77NAME(sswap)(n , dx , incx, dy, incy);
387}
388
389inline void F77NAME(SROT)(const integer *n,
390 double * const dx,
391 const integer *incx,
392 double * const dy,
393 const integer *incy,
394 const double * const c,
395 const double * const s)
396{
397 F77NAME(srot)(n,dx,incx,dy,incy,c,s);
398}
399
400inline void F77NAME(SROTG)(double * const da,
401 double * const db,
402 double * const c,
403 double * const s)
404{
405 F77NAME(srotg)(da,db,c,s);
406}
407#endif
408
409#else //_NO_BLAS1_
410
411/*! @brief Prints an error message and exits (fallback stub for missing BLAS).
412 *
413 */
414static inline void ERREUR(const char* const ch)
415{
416 // Cerr << ch << finl;
417 exit();
418}
419
420/*! @brief Fallback implementation of DASUM (sum of absolute values).
421 *
422 */
423inline double F77NAME(DASUM)(const integer *n,
424 const double* const dx,
425 const integer *incx)
426{
427 double somme=0;
428 for(integer i=0; i<*n; i+=*incx)
429 somme+=dx[i];
430 return somme;
431}
432
433
434/*! @brief Fallback implementation of DAXPY (constant times a vector plus a vector).
435 *
436 */
437inline void F77NAME(DAXPY)(const integer *n,
438 const double * const da,
439 const double * const dx,
440 const integer *incx,
441 double * const dy,
442 const integer *incy)
443{
444 int j =0;
445 for(int i=0; i<*n; i+=*incx)
446 {
447 dy[i]+=*da*dx[j];
448 j+=*incy;
449 }
450}
451
452/*! @brief Fallback implementation of DCOPY (copies a vector to another vector).
453 *
454 */
455inline void F77NAME(DCOPY)(const integer *n,
456 const double * const dx,
457 const integer *incx,
458 double * const dy,
459 const integer *incy)
460{
461 int j=0;
462 for(int i=0; i<*n; i+=*incx)
463 {
464 dy[i]=dx[j];
465 j+=*incy;
466 }
467
468}
469
470
471/*! @brief Fallback implementation of DDOT (dot product of two vectors).
472 *
473 */
474inline double F77NAME(DDOT)(const integer *n,
475 const double * const dx,
476 const integer *incx,
477 const double * const dy,
478 const integer *incy)
479{
480 int j=0;
481 double psc=0;
482 for(int i=0; i<*n; i+=*incx)
483 {
484 psc+=dx[i]*dy[j];
485 j+=*incy;
486 }
487 return psc;
488}
489
490/*! @brief Fallback implementation of DNRM2 (Euclidean norm of a vector).
491 *
492 */
493inline double F77NAME(DNRM2)(const integer *n,
494 const double * const dx,
495 const integer *incx)
496{
497 integer un=1;
498 return sqrt(F77NAME(DDOT)(n, dx, &un, dx, &un));
499}
500
501/*! @brief Fallback stub for DROT (not implemented).
502 *
503 */
504inline void F77NAME(DROT)(const integer *n,
505 double * const dx,
506 const integer *incx,
507 double * const dy,
508 const integer *incy,
509 const double * const c,
510 const double * const s)
511{
512 ERREUR("DROT : pas code!");
513}
514
515/*! @brief Fallback stub for DROTG (not implemented).
516 *
517 */
518inline void F77NAME(DROTG)(double * const da,
519 double * const db,
520 double * const c,
521 double * const s)
522{
523 ERREUR("DROTG : pas code!");
524}
525
526/*! @brief Fallback implementation of DSCAL (scales a vector by a scalar).
527 *
528 */
529inline void F77NAME(DSCAL)(const integer *n,
530 double * const da,
531 double * const dx,
532 const integer *incx)
533{
534 for(int i=0; i<*n; i+= *incx)
535 dx[i]*=*da;
536}
537
538/*! @brief Fallback implementation of DSWAP (swaps two vectors).
539 *
540 */
541inline void F77NAME(DSWAP)(const integer *n,
542 double * const dx,
543 const integer *incx,
544 double * const dy,
545 const integer *incy)
546{
547 int j=0;
548 double tmp;
549 for(int i=0; i<*n; i+= *incx)
550 {
551 tmp=dx[i];
552 dx[i]=dy[j];
553 dy[j]=tmp;
554 j+=*incy;
555 }
556
557}
558
559/*! @brief Fallback implementation of IDAMAX (index of maximum absolute value element).
560 *
561 */
562inline integer F77NAME(IDAMAX)(const integer *n,
563 const double * const dx,
564 const integer *incx)
565{
566 int imax=0;
567 double themax=dx[0];
568 for(int i=1; i<*n; i+=*incx)
569 if(dx[i]>themax)
570 {
571 imax=i;
572 themax=dx[i];
573 }
574 return imax;
575}
576
577/*! @brief Fallback stub for SDOT (not implemented).
578 *
579 */
580inline float F77NAME(SDOT)(const integer *n,
581 const float *dx,
582 const integer *incx,
583 const float *dy,
584 const integer *incy)
585{
586 ERREUR("SDOT : pas code!");
587 return 0;
588}
589
590/*! @brief Fallback stub for SAXPY (not implemented).
591 *
592 */
593inline void F77NAME(SAXPY)(const integer *n,
594 const float *da,
595 const float *dx,
596 const integer *incx,
597 float *dy,
598 const integer *incy)
599{
600 ERREUR("SAXPY : pas code!");
601}
602
603/*! @brief Fallback stub for SCOPY (not implemented).
604 *
605 */
606inline void F77NAME(SCOPY)(const integer *n,
607 float *dx,
608 const integer *incx,
609 float *dy,
610 const integer *incy)
611{
612 ERREUR("SCOPY : pas code!");
613}
614
615/*! @brief Fallback stub for SASUM (not implemented).
616 *
617 */
618inline float F77NAME(SASUM)(const integer *n,
619 const float *dx,
620 const integer *incx)
621{
622 ERREUR("SASUM : pas code!");
623 return 0;
624}
625
626/*! @brief Fallback stub for SNRM2 (not implemented).
627 *
628 */
629inline float F77NAME(SNRM2)(const integer *n,
630 const float *dx,
631 const integer *incx)
632{
633 ERREUR("SNRM2 : pas code!");
634 return 0;
635}
636
637/*! @brief Fallback stub for SSCAL (not implemented).
638 *
639 */
640inline void F77NAME(SSCAL)(const integer *n,
641 const float *da,
642 float *dx,
643 const integer *incx)
644{
645 ERREUR("SSCAL : pas code!");
646}
647
648/*! @brief Fallback stub for ISAMAX (not implemented).
649 *
650 */
651inline integer F77NAME(ISAMAX)(const integer *n,
652 const float * const dx,
653 const integer *incx)
654{
655 ERREUR("ISAMAX : pas code!");
656 return 0;
657}
658
659/*! @brief Fallback stub for SSWAP (not implemented).
660 *
661 */
662inline void F77NAME(SSWAP)(const integer *n,
663 float *dx,
664 const integer *incx,
665 float *dy,
666 integer *incy)
667{
668 ERREUR("SSWAP : pas code!");
669}
670#endif //_NO_BLAS1_
671
672#endif //_BLAS1_H_
673