00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043 #include "blas_dh.h"
00044
00045 #undef __FUNC__
00046 #define __FUNC__ "matvec_euclid_seq"
00047 void
00048 matvec_euclid_seq (int n, int *rp, int *cval, double *aval, double *x,
00049 double *y)
00050 {
00051 START_FUNC_DH int i, j;
00052 int from, to, col;
00053 double sum;
00054
00055 if (np_dh > 1)
00056 SET_V_ERROR ("only for sequential case!\n");
00057
00058 {
00059 for (i = 0; i < n; ++i)
00060 {
00061 sum = 0.0;
00062 from = rp[i];
00063 to = rp[i + 1];
00064 for (j = from; j < to; ++j)
00065 {
00066 col = cval[j];
00067 sum += (aval[j] * x[col]);
00068 }
00069 y[i] = sum;
00070 }
00071 }
00072 END_FUNC_DH}
00073
00074 #undef __FUNC__
00075 #define __FUNC__ "Axpy"
00076 void
00077 Axpy (int n, double alpha, double *x, double *y)
00078 {
00079 START_FUNC_DH int i;
00080
00081 for (i = 0; i < n; ++i)
00082 {
00083 y[i] = alpha * x[i] + y[i];
00084 }
00085 END_FUNC_DH}
00086
00087
00088 #undef __FUNC__
00089 #define __FUNC__ "CopyVec"
00090 void
00091 CopyVec (int n, double *xIN, double *yOUT)
00092 {
00093 START_FUNC_DH int i;
00094
00095 for (i = 0; i < n; ++i)
00096 {
00097 yOUT[i] = xIN[i];
00098 }
00099 END_FUNC_DH}
00100
00101
00102 #undef __FUNC__
00103 #define __FUNC__ "ScaleVec"
00104 void
00105 ScaleVec (int n, double alpha, double *x)
00106 {
00107 START_FUNC_DH int i;
00108
00109 for (i = 0; i < n; ++i)
00110 {
00111 x[i] *= alpha;
00112 }
00113 END_FUNC_DH}
00114
00115 #undef __FUNC__
00116 #define __FUNC__ "InnerProd"
00117 double
00118 InnerProd (int n, double *x, double *y)
00119 {
00120 START_FUNC_DH double result, local_result = 0.0;
00121
00122 int i;
00123
00124 for (i = 0; i < n; ++i)
00125 {
00126 local_result += x[i] * y[i];
00127 }
00128
00129 if (np_dh > 1)
00130 {
00131 MPI_Allreduce (&local_result, &result, 1, MPI_DOUBLE, MPI_SUM, comm_dh);
00132 }
00133 else
00134 {
00135 result = local_result;
00136 }
00137
00138 END_FUNC_VAL (result)}
00139
00140 #undef __FUNC__
00141 #define __FUNC__ "Norm2"
00142 double
00143 Norm2 (int n, double *x)
00144 {
00145 START_FUNC_DH double result, local_result = 0.0;
00146 int i;
00147
00148 for (i = 0; i < n; ++i)
00149 {
00150 local_result += (x[i] * x[i]);
00151 }
00152
00153 if (np_dh > 1)
00154 {
00155 MPI_Allreduce (&local_result, &result, 1, MPI_DOUBLE, MPI_SUM, comm_dh);
00156 }
00157 else
00158 {
00159 result = local_result;
00160 }
00161 result = sqrt (result);
00162 END_FUNC_VAL (result)}