hp2FEM
0.1
|
00001 00010 #ifndef _GK_MKBLAS_H_ 00011 #define _GK_MKBLAS_H_ 00012 00013 00014 #define GK_MKBLAS(PRFX, TYPE, OUTTYPE) \ 00015 /*************************************************************************/\ 00016 \ 00017 /*************************************************************************/\ 00018 TYPE *PRFX ## incset(size_t n, TYPE baseval, TYPE *x)\ 00019 {\ 00020 size_t i;\ 00021 \ 00022 for (i=0; i<n; i++)\ 00023 x[i] = baseval+i;\ 00024 \ 00025 return x;\ 00026 }\ 00027 \ 00028 /*************************************************************************/\ 00029 \ 00030 /*************************************************************************/\ 00031 TYPE PRFX ## max(size_t n, TYPE *x)\ 00032 {\ 00033 size_t i, max=0; \ 00034 \ 00035 if (n <= 0) return (TYPE) 0;\ 00036 \ 00037 for (i=1; i<n; i++)\ 00038 max = (x[i] > x[max] ? i : max);\ 00039 \ 00040 return x[max];\ 00041 }\ 00042 \ 00043 \ 00044 /*************************************************************************/\ 00045 \ 00046 /*************************************************************************/\ 00047 TYPE PRFX ## min(size_t n, TYPE *x)\ 00048 {\ 00049 size_t i, min=0;\ 00050 \ 00051 if (n <= 0) return (TYPE) 0;\ 00052 \ 00053 for (i=1; i<n; i++)\ 00054 min = (x[i] < x[min] ? i : min);\ 00055 \ 00056 return x[min];\ 00057 }\ 00058 \ 00059 \ 00060 /*************************************************************************/\ 00061 \ 00062 /*************************************************************************/\ 00063 size_t PRFX ## argmax(size_t n, TYPE *x)\ 00064 {\ 00065 size_t i, max=0;\ 00066 \ 00067 for (i=1; i<n; i++)\ 00068 max = (x[i] > x[max] ? i : max);\ 00069 \ 00070 return max;\ 00071 }\ 00072 \ 00073 \ 00074 /*************************************************************************/\ 00075 \ 00076 /*************************************************************************/\ 00077 size_t PRFX ## argmin(size_t n, TYPE *x)\ 00078 {\ 00079 size_t i, min=0;\ 00080 \ 00081 for (i=1; i<n; i++)\ 00082 min = (x[i] < x[min] ? i : min);\ 00083 \ 00084 return min;\ 00085 }\ 00086 \ 00087 \ 00088 /*************************************************************************/\ 00089 \ 00090 /*************************************************************************/\ 00091 size_t PRFX ## argmax_n(size_t n, TYPE *x, size_t k)\ 00092 {\ 00093 size_t i, max_n;\ 00094 PRFX ## kv_t *cand;\ 00095 \ 00096 cand = PRFX ## kvmalloc(n, "GK_ARGMAX_N: cand");\ 00097 \ 00098 for (i=0; i<n; i++) {\ 00099 cand[i].val = i;\ 00100 cand[i].key = x[i];\ 00101 }\ 00102 PRFX ## kvsortd(n, cand);\ 00103 \ 00104 max_n = cand[k-1].val;\ 00105 \ 00106 gk_free((void *)&cand, LTERM);\ 00107 \ 00108 return max_n;\ 00109 }\ 00110 \ 00111 \ 00112 /*************************************************************************/\ 00113 \ 00114 /**************************************************************************/\ 00115 OUTTYPE PRFX ## sum(size_t n, TYPE *x, size_t incx)\ 00116 {\ 00117 size_t i;\ 00118 OUTTYPE sum = 0;\ 00119 \ 00120 for (i=0; i<n; i++, x+=incx)\ 00121 sum += (*x);\ 00122 \ 00123 return sum;\ 00124 }\ 00125 \ 00126 \ 00127 /*************************************************************************/\ 00128 \ 00129 /**************************************************************************/\ 00130 TYPE *PRFX ## scale(size_t n, TYPE alpha, TYPE *x, size_t incx)\ 00131 {\ 00132 size_t i;\ 00133 \ 00134 for (i=0; i<n; i++, x+=incx)\ 00135 (*x) *= alpha;\ 00136 \ 00137 return x;\ 00138 }\ 00139 \ 00140 \ 00141 /*************************************************************************/\ 00142 \ 00143 /**************************************************************************/\ 00144 OUTTYPE PRFX ## norm2(size_t n, TYPE *x, size_t incx)\ 00145 {\ 00146 size_t i;\ 00147 OUTTYPE partial = 0;\ 00148 \ 00149 for (i=0; i<n; i++, x+=incx)\ 00150 partial += (*x) * (*x);\ 00151 \ 00152 return (partial > 0 ? (OUTTYPE)sqrt((double)partial) : (OUTTYPE)0);\ 00153 }\ 00154 \ 00155 \ 00156 /*************************************************************************/\ 00157 \ 00158 /**************************************************************************/\ 00159 OUTTYPE PRFX ## dot(size_t n, TYPE *x, size_t incx, TYPE *y, size_t incy)\ 00160 {\ 00161 size_t i;\ 00162 OUTTYPE partial = 0.0;\ 00163 \ 00164 for (i=0; i<n; i++, x+=incx, y+=incy)\ 00165 partial += (*x) * (*y);\ 00166 \ 00167 return partial;\ 00168 }\ 00169 \ 00170 \ 00171 /*************************************************************************/\ 00172 \ 00173 /**************************************************************************/\ 00174 TYPE *PRFX ## axpy(size_t n, TYPE alpha, TYPE *x, size_t incx, TYPE *y, size_t incy)\ 00175 {\ 00176 size_t i;\ 00177 TYPE *y_in = y;\ 00178 \ 00179 for (i=0; i<n; i++, x+=incx, y+=incy)\ 00180 *y += alpha*(*x);\ 00181 \ 00182 return y_in;\ 00183 }\ 00184 00185 00186 00187 #define GK_MKBLAS_PROTO(PRFX, TYPE, OUTTYPE) \ 00188 TYPE *PRFX ## incset(size_t n, TYPE baseval, TYPE *x);\ 00189 TYPE PRFX ## max(size_t n, TYPE *x);\ 00190 TYPE PRFX ## min(size_t n, TYPE *x);\ 00191 size_t PRFX ## argmax(size_t n, TYPE *x);\ 00192 size_t PRFX ## argmin(size_t n, TYPE *x);\ 00193 size_t PRFX ## argmax_n(size_t n, TYPE *x, size_t k);\ 00194 OUTTYPE PRFX ## sum(size_t n, TYPE *x, size_t incx);\ 00195 TYPE *PRFX ## scale(size_t n, TYPE alpha, TYPE *x, size_t incx);\ 00196 OUTTYPE PRFX ## norm2(size_t n, TYPE *x, size_t incx);\ 00197 OUTTYPE PRFX ## dot(size_t n, TYPE *x, size_t incx, TYPE *y, size_t incy);\ 00198 TYPE *PRFX ## axpy(size_t n, TYPE alpha, TYPE *x, size_t incx, TYPE *y, size_t incy);\ 00199 00200 00201 #endif