66 o2cblas_ConjTrans=113};
86 for(
size_t i=0;i<m;i++) {
87 for(
size_t j=0;j<n;j++) {
88 if (!std::isfinite(O2SCL_IX2(data,i,j)))
return false;
111 template<
class vec_t>
double dasum(
const size_t N,
const vec_t &X) {
113 for(
size_t i=0;i<N;i++) {
124 template<
class vec_t,
class vec2_t>
125 void daxpy(
const double alpha,
const size_t N,
const vec_t &X,
134 const size_t m=N % 4;
137 O2SCL_IX(Y,i)+=alpha*O2SCL_IX(X,i);
140 for (i=m;i+3<N;i+=4) {
141 O2SCL_IX(Y,i)+=alpha*O2SCL_IX(X,i);
142 O2SCL_IX(Y,i+1)+=alpha*O2SCL_IX(X,i+1);
143 O2SCL_IX(Y,i+2)+=alpha*O2SCL_IX(X,i+2);
144 O2SCL_IX(Y,i+3)+=alpha*O2SCL_IX(X,i+3);
149 template<
class vec_t,
class vec2_t>
150 double ddot(
const size_t N,
const vec_t &X,
const vec2_t &Y) {
155 const size_t m=N % 4;
158 r+=O2SCL_IX(X,i)*O2SCL_IX(Y,i);
161 for (i=m;i+3<N;i+=4) {
162 r+=O2SCL_IX(X,i)*O2SCL_IX(Y,i);
163 r+=O2SCL_IX(X,i+1)*O2SCL_IX(Y,i+1);
164 r+=O2SCL_IX(X,i+2)*O2SCL_IX(Y,i+2);
165 r+=O2SCL_IX(X,i+3)*O2SCL_IX(Y,i+3);
183 template<
class vec_t>
double dnrm2(
const size_t N,
const vec_t &X) {
192 return fabs(O2SCL_IX(X,0));
196 const double x=O2SCL_IX(X,i);
199 const double ax=fabs(x);
202 ssq=1.0+ssq*(scale/ax)*(scale/ax);
205 ssq+=(ax/scale)*(ax/scale);
211 return scale*sqrt(ssq);
216 template<
class vec_t>
217 void dscal(
const double alpha,
const size_t N, vec_t &X) {
220 const size_t m=N % 4;
223 O2SCL_IX(X,i)*=alpha;
226 for (i=m;i+3<N;i+=4) {
227 O2SCL_IX(X,i)*=alpha;
228 O2SCL_IX(X,i+1)*=alpha;
229 O2SCL_IX(X,i+2)*=alpha;
230 O2SCL_IX(X,i+3)*=alpha;
244 template<
class mat_t,
class vec_t,
class vec2_t>
247 const size_t N,
const double alpha,
const mat_t &A,
248 const vec_t &X,
const double beta, vec2_t &Y) {
254 const int Trans=(TransA != o2cblas_ConjTrans) ? TransA : o2cblas_Trans;
256 if (M == 0 || N == 0) {
260 if (alpha == 0.0 && beta == 1.0) {
264 if (Trans == o2cblas_NoTrans) {
275 for (i=0;i<lenY;i++) {
279 }
else if (beta != 1.0) {
281 for (i=0;i<lenY;i++) {
282 O2SCL_IX(Y,iy) *= beta;
291 if ((order == o2cblas_RowMajor && Trans == o2cblas_NoTrans) ||
292 (order == o2cblas_ColMajor && Trans == o2cblas_Trans)) {
296 for (i=0;i<lenY;i++) {
299 for (j=0;j<lenX;j++) {
300 temp+=O2SCL_IX(X,ix)*O2SCL_IX2(A,i,j);
303 O2SCL_IX(Y,iy)+=alpha*temp;
307 }
else if ((order == o2cblas_RowMajor && Trans == o2cblas_Trans) ||
308 (order == o2cblas_ColMajor && Trans == o2cblas_NoTrans)) {
312 for (j=0;j<lenX;j++) {
313 const double temp=alpha*O2SCL_IX(X,ix);
316 for (i=0;i<lenY;i++) {
317 O2SCL_IX(Y,iy)+=temp*O2SCL_IX2(A,j,i);
334 template<
class mat_t,
class vec_t>
339 const size_t M,
const size_t N,
const mat_t &A, vec_t &X) {
341 const int nonunit=(Diag == o2cblas_NonUnit);
344 const int Trans=(TransA != o2cblas_ConjTrans) ? TransA : o2cblas_Trans;
352 if ((order == o2cblas_RowMajor && Trans == o2cblas_NoTrans &&
353 Uplo == o2cblas_Upper) ||
354 (order == o2cblas_ColMajor && Trans == o2cblas_Trans &&
355 Uplo == o2cblas_Lower)) {
363 O2SCL_IX(X,ix)=O2SCL_IX(X,ix)/O2SCL_IX2(A,N-1,N-1);
367 for (i=(
int)(N-1);i>0 && i--;) {
368 double tmp=O2SCL_IX(X,ix);
370 for (j=i+1;j<((int)N);j++) {
371 const double Aij=O2SCL_IX2(A,i,j);
372 tmp-=Aij*O2SCL_IX(X,jx);
376 O2SCL_IX(X,ix)=tmp/O2SCL_IX2(A,i,i);
383 }
else if ((order == o2cblas_RowMajor && Trans == o2cblas_NoTrans &&
384 Uplo == o2cblas_Lower) ||
385 (order == o2cblas_ColMajor && Trans == o2cblas_Trans &&
386 Uplo == o2cblas_Upper)) {
391 O2SCL_IX(X,ix)=O2SCL_IX(X,ix)/O2SCL_IX2(A,0,0);
394 for (i=1;i<((int)N);i++) {
395 double tmp=O2SCL_IX(X,ix);
398 const double Aij=O2SCL_IX2(A,i,j);
399 tmp-=Aij*O2SCL_IX(X,jx);
403 O2SCL_IX(X,ix)=tmp/O2SCL_IX2(A,i,i);
410 }
else if ((order == o2cblas_RowMajor && Trans == o2cblas_Trans &&
411 Uplo == o2cblas_Upper) ||
412 (order == o2cblas_ColMajor && Trans == o2cblas_NoTrans &&
413 Uplo == o2cblas_Lower)) {
420 O2SCL_IX(X,ix)=O2SCL_IX(X,ix)/O2SCL_IX2(A,0,0);
423 for (i=1;i<((int)N);i++) {
424 double tmp=O2SCL_IX(X,ix);
427 const double Aji=O2SCL_IX2(A,j,i);
428 tmp-=Aji*O2SCL_IX(X,jx);
432 O2SCL_IX(X,ix)=tmp/O2SCL_IX2(A,i,i);
439 }
else if ((order == o2cblas_RowMajor && Trans == o2cblas_Trans &&
440 Uplo == o2cblas_Lower) ||
441 (order == o2cblas_ColMajor && Trans == o2cblas_NoTrans &&
442 Uplo == o2cblas_Upper)) {
449 O2SCL_IX(X,ix)=O2SCL_IX(X,ix)/O2SCL_IX2(A,N-1,N-1);
452 for (i=(
int)(N-1);i>0 && i--;) {
453 double tmp=O2SCL_IX(X,ix);
455 for (j=i+1;j<((int)N);j++) {
456 const double Aji=O2SCL_IX2(A,j,i);
457 tmp-=Aji*O2SCL_IX(X,jx);
461 O2SCL_IX(X,ix)=tmp/O2SCL_IX2(A,i,i);
476 template<
class mat_t,
class vec_t>
481 const mat_t &A, vec_t &x) {
485 const int nonunit=(Diag == o2cblas_NonUnit);
486 const int Trans=(TransA != o2cblas_ConjTrans) ? TransA : o2cblas_Trans;
488 if ((Order == o2cblas_RowMajor && Trans == o2cblas_NoTrans &&
489 Uplo == o2cblas_Upper) ||
490 (Order == o2cblas_ColMajor && Trans == o2cblas_Trans &&
491 Uplo == o2cblas_Lower)) {
497 const size_t j_min=i+1;
498 const size_t j_max=N;
500 for (j=j_min;j<j_max;j++) {
501 temp+=O2SCL_IX(x,jx)*O2SCL_IX2(A,i,j);
505 O2SCL_IX(x,i)=temp+O2SCL_IX(x,i)*O2SCL_IX2(A,i,i);
511 }
else if ((Order == o2cblas_RowMajor && Trans == o2cblas_NoTrans &&
512 Uplo == o2cblas_Lower) ||
513 (Order == o2cblas_ColMajor && Trans == o2cblas_Trans &&
514 Uplo == o2cblas_Upper)) {
519 for (i=N;i>0 && i--;) {
521 const size_t j_min=0;
522 const size_t j_max=i;
524 for (j=j_min;j<j_max;j++) {
525 temp+=O2SCL_IX(x,jx)*O2SCL_IX2(A,i,j);
529 O2SCL_IX(x,ix)=temp+O2SCL_IX(x,ix)*O2SCL_IX2(A,i,i);
531 O2SCL_IX(x,ix)+=temp;
536 }
else if ((Order == o2cblas_RowMajor && Trans == o2cblas_Trans &&
537 Uplo == o2cblas_Upper) ||
538 (Order == o2cblas_ColMajor && Trans == o2cblas_NoTrans &&
539 Uplo == o2cblas_Lower)) {
543 for (i=N;i>0 && i--;) {
545 const size_t j_min=0;
546 const size_t j_max=i;
548 for (j=j_min;j<j_max;j++) {
549 temp+=O2SCL_IX(x,jx)*O2SCL_IX2(A,j,i);
553 O2SCL_IX(x,ix)=temp+O2SCL_IX(x,ix)*O2SCL_IX2(A,i,i);
555 O2SCL_IX(x,ix)+=temp;
560 }
else if ((Order == o2cblas_RowMajor && Trans == o2cblas_Trans &&
561 Uplo == o2cblas_Lower) ||
562 (Order == o2cblas_ColMajor && Trans == o2cblas_NoTrans &&
563 Uplo == o2cblas_Upper)) {
567 const size_t j_min=i+1;
568 const size_t j_max=N;
570 for (j=j_min;j<j_max;j++) {
571 temp+=O2SCL_IX(x,jx)*O2SCL_IX2(A,j,i);
575 O2SCL_IX(x,i)=temp+O2SCL_IX(x,i)*O2SCL_IX2(A,i,i);
582 O2SCL_ERR(
"Unrecognized operation in dtrmv().",
608 template<
class mat_t>
612 const size_t N,
const size_t K,
const double alpha,
613 const mat_t &A,
const mat_t &B,
const double beta, mat_t &C) {
619 if (alpha == 0.0 && beta == 1.0) {
631 if (Order == o2cblas_RowMajor) {
640 O2SCL_IX2(C,i,j)=0.0;
643 }
else if (beta != 1.0) {
646 O2SCL_IX2(C,i,j)*=beta;
655 TransF=(TransA == o2cblas_ConjTrans) ? o2cblas_Trans : TransA;
656 TransG=(TransB == o2cblas_ConjTrans) ? o2cblas_Trans : TransB;
658 if (TransF == o2cblas_NoTrans && TransG == o2cblas_NoTrans) {
664 const double temp=alpha*O2SCL_IX2(A,i,k);
667 O2SCL_IX2(C,i,j)+=temp*O2SCL_IX2(B,k,j);
673 }
else if (TransF == o2cblas_NoTrans && TransG == o2cblas_Trans) {
681 temp+=O2SCL_IX2(A,i,k)*O2SCL_IX2(B,j,k);
683 O2SCL_IX2(C,i,j)+=alpha*temp;
687 }
else if (TransF == o2cblas_Trans && TransG == o2cblas_NoTrans) {
691 const double temp=alpha*O2SCL_IX2(A,k,i);
694 O2SCL_IX2(C,i,j)+=temp*O2SCL_IX2(B,k,j);
700 }
else if (TransF == o2cblas_Trans && TransG == o2cblas_Trans) {
706 temp+=O2SCL_IX2(A,k,i)*O2SCL_IX2(B,j,k);
708 O2SCL_IX2(C,i,j)+=alpha*temp;
727 O2SCL_IX2(C,i,j)=0.0;
730 }
else if (beta != 1.0) {
733 O2SCL_IX2(C,i,j)*=beta;
742 TransF=(TransB == o2cblas_ConjTrans) ? o2cblas_Trans : TransB;
743 TransG=(TransA == o2cblas_ConjTrans) ? o2cblas_Trans : TransA;
745 if (TransF == o2cblas_NoTrans && TransG == o2cblas_NoTrans) {
751 const double temp=alpha*O2SCL_IX2(B,i,k);
754 O2SCL_IX2(C,i,j)+=temp*O2SCL_IX2(A,k,j);
760 }
else if (TransF == o2cblas_NoTrans && TransG == o2cblas_Trans) {
768 temp+=O2SCL_IX2(B,i,k)*O2SCL_IX2(A,j,k);
770 O2SCL_IX2(C,i,j)+=alpha*temp;
774 }
else if (TransF == o2cblas_Trans && TransG == o2cblas_NoTrans) {
778 const double temp=alpha*O2SCL_IX2(B,k,i);
781 O2SCL_IX2(C,i,j)+=temp*O2SCL_IX2(A,k,j);
787 }
else if (TransF == o2cblas_Trans && TransG == o2cblas_Trans) {
793 temp+=O2SCL_IX2(B,k,i)*O2SCL_IX2(A,j,k);
795 O2SCL_IX2(C,i,j)+=alpha*temp;
823 template<
class vec_t,
class vec2_t>
825 vec2_t &Y,
const size_t ie) {
829 if (alpha == 0.0)
return;
830 #if O2SCL_NO_RANGE_CHECK 837 const size_t m=(N-ie) % 4;
839 for (i=ie;i<ie+m;i++) {
840 O2SCL_IX(Y,i)+=alpha*O2SCL_IX(X,i);
844 O2SCL_IX(Y,i)+=alpha*O2SCL_IX(X,i);
845 O2SCL_IX(Y,i+1)+=alpha*O2SCL_IX(X,i+1);
846 O2SCL_IX(Y,i+2)+=alpha*O2SCL_IX(X,i+2);
847 O2SCL_IX(Y,i+3)+=alpha*O2SCL_IX(X,i+3);
860 template<
class vec_t,
class vec2_t>
861 double ddot_subvec(
const size_t N,
const vec_t &X,
const vec2_t &Y,
866 #if O2SCL_NO_RANGE_CHECK 873 const size_t m=(N-ie) % 4;
875 for (i=ie;i<ie+m;i++) {
876 r+=O2SCL_IX(X,i)*O2SCL_IX(Y,i);
880 r+=O2SCL_IX(X,i)*O2SCL_IX(Y,i);
881 r+=O2SCL_IX(X,i+1)*O2SCL_IX(Y,i+1);
882 r+=O2SCL_IX(X,i+2)*O2SCL_IX(Y,i+2);
883 r+=O2SCL_IX(X,i+3)*O2SCL_IX(Y,i+3);
901 template<
class vec_t>
907 #if O2SCL_NO_RANGE_CHECK 915 return fabs(O2SCL_IX(X,ie));
918 for (
size_t i=ie;i<N;i++) {
919 const double x=O2SCL_IX(X,i);
922 const double ax=fabs(x);
925 ssq=1.0+ssq*(scale/ax)*(scale/ax);
928 ssq+=(ax/scale)*(ax/scale);
934 return scale*sqrt(ssq);
946 template<
class vec_t>
950 #if O2SCL_NO_RANGE_CHECK 958 const size_t m=(N-ie) % 4;
960 for (i=ie;i<ie+m;i++) {
961 O2SCL_IX(X,i)*=alpha;
965 O2SCL_IX(X,i)*=alpha;
966 O2SCL_IX(X,i+1)*=alpha;
967 O2SCL_IX(X,i+2)*=alpha;
968 O2SCL_IX(X,i+3)*=alpha;
985 template<
class mat_t,
class vec_t>
987 const size_t ir,
const size_t ic, vec_t &y) {
989 #if O2SCL_NO_RANGE_CHECK 1001 const size_t m=(M-ir) % 4;
1003 for (i=ir;i<m+ir;i++) {
1004 O2SCL_IX(y,i)+=alpha*O2SCL_IX2(X,i,ic);
1008 O2SCL_IX(y,i)+=alpha*O2SCL_IX2(X,i,ic);
1009 O2SCL_IX(y,i+1)+=alpha*O2SCL_IX2(X,i+1,ic);
1010 O2SCL_IX(y,i+2)+=alpha*O2SCL_IX2(X,i+2,ic);
1011 O2SCL_IX(y,i+3)+=alpha*O2SCL_IX2(X,i+3,ic);
1027 template<
class mat_t,
class vec_t>
1029 const size_t ic,
const vec_t &y) {
1030 #if O2SCL_NO_RANGE_CHECK 1039 const size_t m=(M-ir) % 4;
1041 for (i=ir;i<m+ir;i++) {
1042 r+=O2SCL_IX2(X,i,ic)*O2SCL_IX(y,i);
1046 r+=O2SCL_IX2(X,i,ic)*O2SCL_IX(y,i);
1047 r+=O2SCL_IX2(X,i+1,ic)*O2SCL_IX(y,i+1);
1048 r+=O2SCL_IX2(X,i+2,ic)*O2SCL_IX(y,i+2);
1049 r+=O2SCL_IX2(X,i+3,ic)*O2SCL_IX(y,i+3);
1070 template<
class mat_t>
1078 #if O2SCL_NO_RANGE_CHECK 1087 return fabs(O2SCL_IX2(A,ir,ic));
1090 for (i=ir;i<M;i++) {
1091 const double x=O2SCL_IX2(A,i,ic);
1094 const double ax=fabs(x);
1097 ssq=1.0+ssq*(scale/ax)*(scale/ax);
1100 ssq+=(ax/scale)*(ax/scale);
1106 return scale*sqrt(ssq);
1119 template<
class mat_t>
1121 const size_t M,
const double alpha) {
1123 #if O2SCL_NO_RANGE_CHECK 1131 const size_t m=(M-ir) % 4;
1133 for (i=ir;i<m+ir;i++) {
1134 O2SCL_IX2(A,i,ic)*=alpha;
1138 O2SCL_IX2(A,i,ic)*=alpha;
1139 O2SCL_IX2(A,i+1,ic)*=alpha;
1140 O2SCL_IX2(A,i+2,ic)*=alpha;
1141 O2SCL_IX2(A,i+3,ic)*=alpha;
1157 template<
class mat_t>
1161 #if O2SCL_NO_RANGE_CHECK 1170 const size_t m=(M-ir) % 4;
1172 for (i=ir;i<m+ir;i++) {
1173 r+=fabs(O2SCL_IX2(A,i,ic));
1177 r+=fabs(O2SCL_IX2(A,i,ic));
1178 r+=fabs(O2SCL_IX2(A,i+1,ic));
1179 r+=fabs(O2SCL_IX2(A,i+2,ic));
1180 r+=fabs(O2SCL_IX2(A,i+3,ic));
1203 template<
class mat_t,
class vec_t>
1205 const size_t ir,
const size_t ic, vec_t &Y) {
1207 #if O2SCL_NO_RANGE_CHECK 1219 const size_t m=(N-ic) % 4;
1221 for (i=ic;i<m+ic;i++) {
1222 O2SCL_IX(Y,i)+=alpha*O2SCL_IX2(X,ir,i);
1226 O2SCL_IX(Y,i)+=alpha*O2SCL_IX2(X,ir,i);
1227 O2SCL_IX(Y,i+1)+=alpha*O2SCL_IX2(X,ir,i+1);
1228 O2SCL_IX(Y,i+2)+=alpha*O2SCL_IX2(X,ir,i+2);
1229 O2SCL_IX(Y,i+3)+=alpha*O2SCL_IX2(X,ir,i+3);
1249 template<
class mat_t,
class vec_t>
1251 const size_t ic,
const vec_t &Y) {
1253 #if O2SCL_NO_RANGE_CHECK 1262 const size_t m=(N-ic) % 4;
1264 for (i=ic;i<m+ic;i++) {
1265 r+=O2SCL_IX2(X,ir,i)*O2SCL_IX(Y,i);
1269 r+=O2SCL_IX2(X,ir,i)*O2SCL_IX(Y,i);
1270 r+=O2SCL_IX2(X,ir,i+1)*O2SCL_IX(Y,i+1);
1271 r+=O2SCL_IX2(X,ir,i+2)*O2SCL_IX(Y,i+2);
1272 r+=O2SCL_IX2(X,ir,i+3)*O2SCL_IX(Y,i+3);
1288 template<
class mat_t>
1297 return fabs(O2SCL_IX2(M,ir,ic));
1300 for (i=ic;i<N;i++) {
1301 const double x=O2SCL_IX2(M,ir,i);
1304 const double ax=fabs(x);
1307 ssq=1.0+ssq*(scale/ax)*(scale/ax);
1310 ssq+=(ax/scale)*(ax/scale);
1316 return scale*sqrt(ssq);
1331 template<
class mat_t>
1333 const size_t N,
const double alpha) {
1335 #if O2SCL_NO_RANGE_CHECK 1343 const size_t m=(N-ic) % 4;
1345 for (i=ic;i<m+ic;i++) {
1346 O2SCL_IX2(A,ir,i)*=alpha;
1350 O2SCL_IX2(A,ir,i)*=alpha;
1351 O2SCL_IX2(A,ir,i+1)*=alpha;
1352 O2SCL_IX2(A,ir,i+2)*=alpha;
1353 O2SCL_IX2(A,ir,i+3)*=alpha;
void daxpy_subvec(const double alpha, const size_t N, const vec_t &X, vec2_t &Y, const size_t ie)
Compute beginning with index ie and ending with index N-1.
double dasum(const size_t N, const vec_t &X)
Compute the absolute sum of vector elements.
void dtrsv(const enum o2cblas_order order, const enum o2cblas_uplo Uplo, const enum o2cblas_transpose TransA, const enum o2cblas_diag Diag, const size_t M, const size_t N, const mat_t &A, vec_t &X)
Compute .
void daxpy_subcol(const double alpha, const size_t M, const mat_t &X, const size_t ir, const size_t ic, vec_t &y)
Compute for a subcolumn of a matrix.
void dgemm(const enum o2cblas_order Order, const enum o2cblas_transpose TransA, const enum o2cblas_transpose TransB, const size_t M, const size_t N, const size_t K, const double alpha, const mat_t &A, const mat_t &B, const double beta, mat_t &C)
Compute .
bool matrix_is_finite(size_t m, size_t n, mat_t &data)
Test if the first n elements of a matrix are finite.
invalid argument supplied by user
o2cblas_uplo
Upper- or lower-triangular.
double ddot_subrow(const size_t N, const mat_t &X, const size_t ir, const size_t ic, const vec_t &Y)
Compute for a subrow of a matrix.
o2cblas_diag
Unit or generic diagonal.
double dnrm2(const size_t N, const vec_t &X)
Compute the norm of the vector X.
void dscal_subrow(mat_t &A, const size_t ir, const size_t ic, const size_t N, const double alpha)
Compute for a subrow of a matrix.
void dgemv(const enum o2cblas_order order, const enum o2cblas_transpose TransA, const size_t M, const size_t N, const double alpha, const mat_t &A, const vec_t &X, const double beta, vec2_t &Y)
Compute .
Namespace for O2scl CBLAS function templates.
void dtrmv(const enum o2cblas_order Order, const enum o2cblas_uplo Uplo, const enum o2cblas_transpose TransA, const enum o2cblas_diag Diag, const size_t N, const mat_t &A, vec_t &x)
Compute for the triangular matrix A.
o2cblas_transpose
Transpose operations.
double dasum_subcol(mat_t &A, const size_t ir, const size_t ic, const size_t M)
Compute for a subcolumn of a matrix.
double ddot(const size_t N, const vec_t &X, const vec2_t &Y)
Compute .
double dnrm2_subcol(const mat_t &A, const size_t ir, const size_t ic, const size_t M)
Compute the norm of a subcolumn of a matrix.
double dnrm2_subvec(const size_t N, const vec_t &X, const size_t ie)
Compute the norm of the vector X beginning with index ie and ending with index N-1.
double dnrm2_subrow(const mat_t &M, const size_t ir, const size_t ic, const size_t N)
Compute the norm of a subrow of a matrix.
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
double ddot_subvec(const size_t N, const vec_t &X, const vec2_t &Y, const size_t ie)
Compute beginning with index ie and ending with index N-1.
void dscal(const double alpha, const size_t N, vec_t &X)
Compute .
o2cblas_side
Left or right sided operation.
void daxpy(const double alpha, const size_t N, const vec_t &X, vec2_t &Y)
Compute .
double ddot_subcol(const size_t M, const mat_t &X, const size_t ir, const size_t ic, const vec_t &y)
Compute for a subcolumn of a matrix.
void daxpy_subrow(const double alpha, const size_t N, const mat_t &X, const size_t ir, const size_t ic, vec_t &Y)
Compute for a subrow of a matrix.
o2cblas_order
Matrix order, either column-major or row-major.
void dscal_subvec(const double alpha, const size_t N, vec_t &X, const size_t ie)
Compute beginning with index ie and ending with index N-1.
void dscal_subcol(mat_t &A, const size_t ir, const size_t ic, const size_t M, const double alpha)
Compute for a subcolumn of a matrix.