21#include "fasp_functs.h"
47 if ( a == 1.0 )
return;
54 INT myid, mybegin, myend, nthreads;
57 nthreads = fasp_get_num_threads();
63#pragma omp parallel private(myid, mybegin, myend, i)
65 myid = omp_get_thread_num();
67 for ( i = mybegin; i < myend; ++i ) x[i] *= a;
72 for ( i = 0; i < n; ++i ) x[i] *= a;
102 INT myid, mybegin, myend, nthreads;
105 nthreads = fasp_get_num_threads();
112#pragma omp parallel private(myid, mybegin, myend, i) num_threads(nthreads)
114 myid = omp_get_thread_num();
116 for ( i = mybegin; i < myend; ++i ) y[i] += x[i];
121 for ( i = 0; i < n; ++i ) y[i] += x[i];
125 else if ( a == -1.0 ) {
128#pragma omp parallel private(myid, mybegin, myend, i) num_threads(nthreads)
130 myid = omp_get_thread_num();
132 for ( i = mybegin; i < myend; ++i ) y[i] -= x[i];
137 for ( i = 0; i < n; ++i ) y[i] -= x[i];
144#pragma omp parallel private(myid, mybegin, myend, i) num_threads(nthreads)
146 myid = omp_get_thread_num();
148 for ( i = mybegin; i < myend; ++i ) y[i] += a*x[i];
153 for ( i = 0; i < n; ++i ) y[i] += a*x[i];
357 INT myid, mybegin, myend, nthreads;
360 nthreads = fasp_get_num_threads();
366#pragma omp parallel private(myid, mybegin, myend, i) num_threads(nthreads)
368 myid = omp_get_thread_num();
370 for ( i = mybegin; i < myend; ++i ) z[i] = a*x[i] + y[i];
375 for ( i = 0; i < n; ++i ) z[i] = a*x[i] + y[i];
398 z[0] = a*x[0] + y[0];
399 z[1] = a*x[1] + y[1];
401 z[2] = a*x[2] + y[2];
402 z[3] = a*x[3] + y[3];
424 z[0] = a*x[0] + y[0];
425 z[1] = a*x[1] + y[1];
426 z[2] = a*x[2] + y[2];
428 z[3] = a*x[3] + y[3];
429 z[4] = a*x[4] + y[4];
430 z[5] = a*x[5] + y[5];
432 z[6] = a*x[6] + y[6];
433 z[7] = a*x[7] + y[7];
434 z[8] = a*x[8] + y[8];
456 z[0] = a*x[0] + y[0];
457 z[1] = a*x[1] + y[1];
458 z[2] = a*x[2] + y[2];
459 z[3] = a*x[3] + y[3];
460 z[4] = a*x[4] + y[4];
462 z[5] = a*x[5] + y[5];
463 z[6] = a*x[6] + y[6];
464 z[7] = a*x[7] + y[7];
465 z[8] = a*x[8] + y[8];
466 z[9] = a*x[9] + y[9];
468 z[10] = a*x[10] + y[10];
469 z[11] = a*x[11] + y[11];
470 z[12] = a*x[12] + y[12];
471 z[13] = a*x[13] + y[13];
472 z[14] = a*x[14] + y[14];
474 z[15] = a*x[15] + y[15];
475 z[16] = a*x[16] + y[16];
476 z[17] = a*x[17] + y[17];
477 z[18] = a*x[18] + y[18];
478 z[19] = a*x[19] + y[19];
480 z[20] = a*x[20] + y[20];
481 z[21] = a*x[21] + y[21];
482 z[22] = a*x[22] + y[22];
483 z[23] = a*x[23] + y[23];
484 z[24] = a*x[24] + y[24];
506 z[0] = a*x[0] + y[0];
507 z[1] = a*x[1] + y[1];
508 z[2] = a*x[2] + y[2];
509 z[3] = a*x[3] + y[3];
510 z[4] = a*x[4] + y[4];
511 z[5] = a*x[5] + y[5];
512 z[6] = a*x[6] + y[6];
514 z[7] = a*x[7] + y[7];
515 z[8] = a*x[8] + y[8];
516 z[9] = a*x[9] + y[9];
517 z[10] = a*x[10] + y[10];
518 z[11] = a*x[11] + y[11];
519 z[12] = a*x[12] + y[12];
520 z[13] = a*x[13] + y[13];
522 z[14] = a*x[14] + y[14];
523 z[15] = a*x[15] + y[15];
524 z[16] = a*x[16] + y[16];
525 z[17] = a*x[17] + y[17];
526 z[18] = a*x[18] + y[18];
527 z[19] = a*x[19] + y[19];
528 z[20] = a*x[20] + y[20];
530 z[21] = a*x[21] + y[21];
531 z[22] = a*x[22] + y[22];
532 z[23] = a*x[23] + y[23];
533 z[24] = a*x[24] + y[24];
534 z[25] = a*x[25] + y[25];
535 z[26] = a*x[26] + y[26];
536 z[27] = a*x[27] + y[27];
538 z[28] = a*x[28] + y[28];
539 z[29] = a*x[29] + y[29];
540 z[30] = a*x[30] + y[30];
541 z[31] = a*x[31] + y[31];
542 z[32] = a*x[32] + y[32];
543 z[33] = a*x[33] + y[33];
544 z[34] = a*x[34] + y[34];
546 z[35] = a*x[35] + y[35];
547 z[36] = a*x[36] + y[36];
548 z[37] = a*x[37] + y[37];
549 z[38] = a*x[38] + y[38];
550 z[39] = a*x[39] + y[39];
551 z[40] = a*x[40] + y[40];
552 z[41] = a*x[41] + y[41];
554 z[42] = a*x[42] + y[42];
555 z[43] = a*x[43] + y[43];
556 z[44] = a*x[44] + y[44];
557 z[45] = a*x[45] + y[45];
558 z[46] = a*x[46] + y[46];
559 z[47] = a*x[47] + y[47];
560 z[48] = a*x[48] + y[48];
590 INT myid, mybegin, myend, nthreads;
593 nthreads = fasp_get_num_threads();
599#pragma omp parallel private(myid, mybegin, myend, i) num_threads(nthreads)
601 myid = omp_get_thread_num();
603 for ( i = mybegin; i < myend; ++i ) y[i] = a*x[i] + b*y[i];
608 for ( i = 0; i < n; ++i ) y[i] = a*x[i] + b*y[i];
631 register REAL onenorm = 0.0;
635#pragma omp parallel for reduction(+:onenorm) private(i)
637 for ( i = 0; i < n; ++i ) onenorm +=
ABS(x[i]);
660 register REAL twonorm = 0.0;
664#pragma omp parallel for reduction(+:twonorm) private(i)
666 for ( i = 0; i < n; ++i ) twonorm += x[i] * x[i];
668 return sqrt(twonorm);
690 register REAL infnorm = 0.0;
694 INT myid, mybegin, myend, nthreads;
697 nthreads = fasp_get_num_threads();
703 REAL infnorm_loc = 0.0;
704#pragma omp parallel firstprivate(infnorm_loc) private(myid, mybegin, myend, i)
706 myid = omp_get_thread_num();
708 for ( i = mybegin; i < myend; ++i )
709 infnorm_loc =
MAX( infnorm_loc,
ABS(x[i]) );
711 if ( infnorm_loc > infnorm ) {
713 infnorm =
MAX( infnorm_loc, infnorm );
719 for ( i = 0; i < n; ++i ) infnorm =
MAX( infnorm,
ABS(x[i]) );
746 register REAL value = 0.0;
755#pragma omp parallel for reduction(+:value) private(i)
757 for ( i = 0; i < n; ++i ) value += x[i]*y[i];
760 for ( i = 0; i < n; ++i ) value += x[i]*y[i];
void fasp_get_start_end(const INT procid, const INT nprocs, const INT n, INT *start, INT *end)
Assign Load to each thread.
REAL fasp_blas_darray_dotprod(const INT n, const REAL *x, const REAL *y)
Inner product of two arraies x and y.
void fasp_blas_darray_axpy_nc3(const REAL a, const REAL *x, REAL *y)
y = a*x + y, length of x and y should be 3
REAL fasp_blas_darray_norminf(const INT n, const REAL *x)
Linf norm of array x.
void fasp_blas_darray_axpyz(const INT n, const REAL a, const REAL *x, const REAL *y, REAL *z)
z = a*x + y
REAL fasp_blas_darray_norm1(const INT n, const REAL *x)
L1 norm of array x.
void fasp_blas_darray_axpy_nc7(const REAL a, const REAL *x, REAL *y)
y = a*x + y, length of x and y should be 7
void fasp_blas_darray_axpy_nc5(const REAL a, const REAL *x, REAL *y)
y = a*x + y, length of x and y should be 5
void fasp_blas_darray_axpby(const INT n, const REAL a, const REAL *x, const REAL b, REAL *y)
y = a*x + b*y
void fasp_blas_darray_axpyz_nc7(const REAL a, const REAL *x, const REAL *y, REAL *z)
z = a*x + y, length of x, y and z should be 7
void fasp_blas_darray_axpy_nc2(const REAL a, const REAL *x, REAL *y)
y = a*x + y, length of x and y should be 2
void fasp_blas_darray_axpyz_nc2(const REAL a, const REAL *x, const REAL *y, REAL *z)
z = a*x + y, length of x, y and z should be 2
REAL fasp_blas_darray_norm2(const INT n, const REAL *x)
L2 norm of array x.
void fasp_blas_darray_ax(const INT n, const REAL a, REAL *x)
x = a*x
void fasp_blas_darray_axpyz_nc3(const REAL a, const REAL *x, const REAL *y, REAL *z)
z = a*x + y, length of x, y and z should be 3
void fasp_blas_darray_axpy(const INT n, const REAL a, const REAL *x, REAL *y)
y = a*x + y
void fasp_blas_darray_axpyz_nc5(const REAL a, const REAL *x, const REAL *y, REAL *z)
z = a*x + y, length of x, y and z should be 5
Main header file for the FASP project.
#define SHORT
FASP integer and floating point numbers.
#define MAX(a, b)
Definition of max, min, abs.
#define TRUE
Definition of logic type.