23 #ifndef O2SCL_VECTOR_H 24 #define O2SCL_VECTOR_H 57 #include <gsl/gsl_vector.h> 58 #include <gsl/gsl_sys.h> 59 #include <gsl/gsl_matrix.h> 61 #include <o2scl/misc.h> 62 #include <o2scl/uniform_grid.h> 76 d=(
const double *)m->data;
79 double operator[](
size_t i)
const {
99 d=(
const double *)m->data;
104 double operator()(
size_t i,
size_t j)
const {
126 template<
class vec_t,
class vec2_t>
129 if (dest.size()<N) dest.resize(N);
134 for(i=m;i+3<N;i+=4) {
155 template<
class vec_t,
class vec2_t>
161 for(i=m;i+3<N;i+=4) {
178 template<
class mat_t,
class mat2_t>
180 size_t m=src.size1();
181 size_t n=src.size2();
182 if (dest.size1()<m || dest.size2()<n) dest.resize(m,n);
183 for(
size_t i=0;i<m;i++) {
184 for(
size_t j=0;j<n;j++) {
203 template<
class mat_t,
class mat2_t>
205 for(
size_t i=0;i<M;i++) {
206 for(
size_t j=0;j<N;j++) {
224 template<
class mat_t,
class mat2_t>
226 size_t m=src.size1();
227 size_t n=src.size2();
228 if (dest.size1()<n || dest.size2()<m) dest.resize(n,m);
229 for(
size_t i=0;i<m;i++) {
230 for(
size_t j=0;j<n;j++) {
245 template<
class mat_t,
class mat2_t>
247 for(
size_t i=0;i<m;i++) {
248 for(
size_t j=0;j<n;j++) {
264 template<
class mat_t,
class data_t>
266 size_t m=src.size1();
267 size_t n=src.size2();
271 for(
size_t i=0;i<n;i++) {
272 for(
size_t j=0;j<n;j++) {
289 template<
class mat_t,
class data_t>
292 for(
size_t i=0;i<m;i++) {
293 for(
size_t j=0;j<n;j++) {
307 size_t m=src.size1();
308 size_t n=src.size2();
310 for(
size_t i=0;i<m;i++) {
311 for(
size_t j=i+1;j<n;j++) {
312 if (src(i,j)!=0.0)
return false;
321 size_t m=src.size1();
322 size_t n=src.size2();
324 for(
size_t j=0;j<n;j++) {
325 for(
size_t i=j+1;i<m;i++) {
326 if (src(i,j)!=0.0)
return false;
336 size_t m=src.size1();
337 size_t n=src.size2();
338 for(
size_t i=0;i<m;i++) {
339 for(
size_t j=i+1;j<n;j++) {
350 size_t m=src.size1();
351 size_t n=src.size2();
352 for(
size_t j=0;j<n;j++) {
353 for(
size_t i=j+1;i<m;i++) {
366 for(
size_t i=0;i<m;i++) {
367 for(
size_t j=i+1;j<n;j++) {
368 if (src(i,j)!=0.0)
return false;
380 for(
size_t j=0;j<n;j++) {
381 for(
size_t i=j+1;i<m;i++) {
382 if (src(i,j)!=0.0)
return false;
393 for(
size_t i=0;i<m;i++) {
394 for(
size_t j=i+1;j<n;j++) {
406 for(
size_t j=0;j<n;j++) {
407 for(
size_t i=j+1;i<m;i++) {
422 template<
class vec_t,
class vec2_t,
class data_t>
431 for(i=m;i+3<N;i+=4) {
458 template<
class vec_t,
class vec2_t,
class data_t>
461 if (v2.size()<N) N=v2.size();
469 for(i=m;i+3<N;i+=4) {
492 template<
class vec_t,
class vec2_t>
494 return vector_swap<vec_t,vec2_t,double>(N,v1,v2);
508 template<
class vec_t,
class vec2_t>
510 return vector_swap<vec_t,vec2_t,double>(v1,v2);
518 template<
class vec_t,
class data_t>
533 template<
class vec_t>
535 return vector_swap<vec_t,double>(v,i,j);
544 template<
class mat_t,
class mat2_t,
class data_t>
547 for(
size_t i=0;i<M;i++) {
548 for(
size_t j=0;j<N;j++) {
563 template<
class mat_t,
class mat2_t,
class data_t>
565 return matrix_swap<mat_t,mat2_t,double>(M,N,m1,m2);
573 template<
class mat_t,
class data_t>
574 void matrix_swap(mat_t &m,
size_t i1,
size_t j1,
size_t i2,
size_t j2) {
575 data_t temp=m(i1,j1);
586 template<
class mat_t>
588 size_t i2,
size_t j2) {
589 return matrix_swap<mat_t,double>(m,i1,j1,i2,j2);
597 template<
class mat_t,
class data_t>
600 for(
size_t i=0;i<M;i++) {
614 template<
class mat_t>
616 return matrix_swap_cols<mat_t,double>(M,m,j1,j2);
625 template<
class mat_t,
class data_t>
628 for(
size_t j=0;j<N;j++) {
642 template<
class mat_t>
644 return matrix_swap_rows<mat_t,double>(N,m,i1,i2);
652 template<
class vec_t,
class data_t>
660 if (j<n && data[j] < data[j+1]) j++;
661 if (!(v < data[j]))
break;
707 template<
class vec_t,
class data_t>
720 sort_downheap<vec_t,data_t>(data,N,k);
728 sort_downheap<vec_t,data_t>(data,N,0);
736 template<
class vec_t,
class vec_
size_t>
740 const size_t pki = order[k];
745 if (j < N && data[order[j]] < data[order[j + 1]]) {
750 if (!(data[pki] < data[order[j]])) {
800 template<
class vec_t,
class vec_
size_t>
809 for (i = 0 ; i < n ; i++) {
823 sort_index_downheap<vec_t,vec_size_t>(N,data,order,k);
829 size_t tmp = order[0];
836 sort_index_downheap<vec_t,vec_size_t>(N,data,order,0);
851 template<
class vec_t>
853 return vector_sort<vec_t,double>(n,data);
876 template<
class vec_t,
class data_t>
879 O2SCL_ERR2(
"Subset length greater than size in ",
889 data_t xbound=data[0];
893 for(
size_t i=1;i<n;i++) {
897 }
else if (xi>=xbound) {
901 for(i1=j-1;i1>0;i1--) {
902 if (xi>smallest[i1-1])
break;
903 smallest[i1]=smallest[i1-1];
906 xbound=smallest[j-1];
927 template<
class vec_t,
class data_t>
929 size_t n=data.size();
930 if (smallest.size()<k) smallest.resize(k);
949 template<
class vec_t,
class data_t,
class vec_
size_t>
953 O2SCL_ERR2(
"Subset length greater than size in ",
954 "function vector_smallest_index().",
exc_einval);
958 "function vector_smallest_index().",
exc_einval);
965 data_t xbound=data[0];
969 for(
size_t i=1;i<n;i++) {
973 }
else if (xi>=xbound) {
977 for(i1=j-1;i1>0;i1--) {
978 if (xi>data[index[i1-1]])
break;
979 index[i1]=index[i1-1];
982 xbound=data[index[j-1]];
989 template<
class vec_t,
class data_t,
class vec_
size_t>
992 size_t n=data.size();
993 if (index.size()<k) index.resize(k);
994 return o2scl::vector_smallest_index<vec_t,data_t,vec_size_t>
1014 template<
class vec_t,
class data_t>
1017 O2SCL_ERR2(
"Subset length greater than size in ",
1027 data_t xbound=data[0];
1031 for(
size_t i=1;i<n;i++) {
1035 }
else if (xi<=xbound) {
1039 for(i1=j-1;i1>0;i1--) {
1040 if (xi<largest[i1-1])
break;
1041 largest[i1]=largest[i1-1];
1044 xbound=largest[j-1];
1065 template<
class vec_t,
class data_t>
1067 size_t n=data.size();
1068 if (largest.size()<k) largest.resize(k);
1077 template<
class vec_t,
class data_t>
1084 for(
size_t i=1;i<n;i++) {
1094 template<
class vec_t,
class data_t>
1097 size_t n=data.size();
1102 for(
size_t i=1;i<n;i++) {
1113 template<
class vec_t,
class data_t>
1121 for(
size_t i=1;i<n;i++) {
1132 template<
class vec_t,
class data_t>
1141 for(
size_t i=1;i<n;i++) {
1152 template<
class vec_t,
class data_t>
1159 for(
size_t i=1;i<n;i++) {
1169 template<
class vec_t,
class data_t>
1172 size_t n=data.size();
1177 for(
size_t i=1;i<n;i++) {
1188 template<
class vec_t,
class data_t>
1196 for(
size_t i=1;i<n;i++) {
1207 template<
class vec_t,
class data_t>
1216 for(
size_t i=1;i<n;i++) {
1228 template<
class vec_t,
class data_t>
1230 data_t &min, data_t &max) {
1237 for(
size_t i=1;i<n;i++) {
1251 template<
class vec_t,
class data_t>
1253 size_t &ix_min,
size_t &ix_max) {
1262 for(
size_t i=1;i<n;i++) {
1278 template<
class vec_t,
class data_t>
1280 size_t &ix_min, data_t &min,
1281 size_t &ix_max, data_t &max) {
1290 for(
size_t i=1;i<n;i++) {
1308 template<
class vec_t,
class data_t>
1310 size_t ix=vector_max_index<vec_t,data_t>(n,data);
1312 return quadratic_extremum_y<data_t>(0,1,2,data[0],data[1],data[2]);
1313 }
else if (ix==n-1) {
1314 return quadratic_extremum_y<data_t>
1315 (n-3,n-2,n-1,data[n-3],data[n-2],data[n-1]);
1317 return quadratic_extremum_y<data_t>
1318 (ix-1,ix,ix+1,data[ix-1],data[ix],data[ix+1]);
1323 template<
class vec_t,
class data_t>
1325 size_t ix=vector_max_index<vec_t,data_t>(n,y);
1326 if (ix==0 || ix==n-1)
return y[ix];
1327 return quadratic_extremum_y<data_t>(x[ix-1],x[ix],x[ix+1],
1328 y[ix-1],y[ix],y[ix+1]);
1333 template<
class vec_t,
class data_t>
1335 size_t ix=vector_max_index<vec_t,data_t>(n,y);
1336 if (ix==0 || ix==n-1)
return y[ix];
1337 return quadratic_extremum_x<data_t>(x[ix-1],x[ix],x[ix+1],
1338 y[ix-1],y[ix],y[ix+1]);
1343 template<
class vec_t,
class data_t>
1345 size_t ix=vector_min_index<vec_t,data_t>(n,data);
1347 return quadratic_extremum_y<data_t>(0,1,2,data[0],data[1],data[2]);
1348 }
else if (ix==n-1) {
1349 return quadratic_extremum_y<data_t>
1350 (n-3,n-2,n-1,data[n-3],data[n-2],data[n-1]);
1352 return quadratic_extremum_y<data_t>
1353 (ix-1,ix,ix+1,data[ix-1],data[ix],data[ix+1]);
1358 template<
class vec_t,
class data_t>
1360 size_t ix=vector_min_index<vec_t,data_t>(n,y);
1361 if (ix==0 || ix==n-1)
return y[ix];
1362 return quadratic_extremum_y<data_t>(x[ix-1],x[ix],x[ix+1],
1363 y[ix-1],y[ix],y[ix+1]);
1368 template<
class vec_t,
class data_t>
1370 size_t ix=vector_min_index<vec_t,data_t>(n,y);
1371 if (ix==0 || ix==n-1)
return y[ix];
1372 return quadratic_extremum_x<data_t>(x[ix-1],x[ix],x[ix+1],
1373 y[ix-1],y[ix],y[ix+1]);
1381 template<
class mat_t,
class data_t>
1385 std::string str=((std::string)
"Matrix with zero size (")+
1387 "matrix_max_value().";
1390 data_t max=data(0,0);
1391 for(
size_t i=0;i<m;i++) {
1392 for(
size_t j=0;j<n;j++) {
1393 if (data(i,j)>max) {
1403 template<
class mat_t,
class data_t> data_t
1405 size_t m=data.size1();
1406 size_t n=data.size2();
1408 std::string str=((std::string)
"Matrix with zero size (")+
1410 "matrix_max_value().";
1413 data_t max=data(0,0);
1414 for(
size_t i=0;i<m;i++) {
1415 for(
size_t j=0;j<n;j++) {
1416 if (data(i,j)>max) {
1426 template<
class mat_t>
double 1428 size_t m=data.size1();
1429 size_t n=data.size2();
1431 std::string str=((std::string)
"Matrix with zero size (")+
1433 "matrix_max_value_double().";
1436 double max=data(0,0);
1437 for(
size_t i=0;i<m;i++) {
1438 for(
size_t j=0;j<n;j++) {
1439 if (data(i,j)>max) {
1450 template<
class mat_t,
class data_t>
1452 size_t &i_max,
size_t &j_max, data_t &max) {
1455 std::string str=((std::string)
"Matrix with zero size (")+
1457 "matrix_max_index().";
1463 for(
size_t i=0;i<m;i++) {
1464 for(
size_t j=0;j<n;j++) {
1465 if (data(i,j)>max) {
1478 template<
class mat_t,
class data_t>
1480 size_t &i_max,
size_t &j_max, data_t &max) {
1482 size_t m=data.size1();
1483 size_t n=data.size2();
1485 std::string str=((std::string)
"Matrix with zero size (")+
1487 "matrix_max_index().";
1493 for(
size_t i=0;i<m;i++) {
1494 for(
size_t j=0;j<n;j++) {
1495 if (data(i,j)>max) {
1507 template<
class mat_t,
class data_t>
1511 std::string str=((std::string)
"Matrix with zero size (")+
1513 "matrix_min_value().";
1516 data_t min=data(0,0);
1517 for(
size_t i=0;i<m;i++) {
1518 for(
size_t j=0;j<n;j++) {
1519 if (data(i,j)<min) {
1529 template<
class mat_t,
class data_t>
1532 size_t m=data.size1();
1533 size_t n=data.size2();
1535 std::string str=((std::string)
"Matrix with zero size (")+
1537 "matrix_min_value().";
1540 data_t min=data(0,0);
1541 for(
size_t i=0;i<m;i++) {
1542 for(
size_t j=0;j<n;j++) {
1543 if (data(i,j)<min) {
1553 template<
class mat_t>
1556 size_t m=data.size1();
1557 size_t n=data.size2();
1559 std::string str=((std::string)
"Matrix with zero size (")+
1561 "matrix_min_value().";
1564 double min=data(0,0);
1565 for(
size_t i=0;i<m;i++) {
1566 for(
size_t j=0;j<n;j++) {
1567 if (data(i,j)<min) {
1578 template<
class mat_t,
class data_t>
1580 size_t &i_min,
size_t &j_min, data_t &min) {
1583 std::string str=((std::string)
"Matrix with zero size (")+
1585 "matrix_min_index().";
1591 for(
size_t i=0;i<m;i++) {
1592 for(
size_t j=0;j<n;j++) {
1593 if (data(i,j)<min) {
1606 template<
class mat_t,
class data_t>
1608 size_t &i_min,
size_t &j_min, data_t &min) {
1610 size_t m=data.size1();
1611 size_t n=data.size2();
1613 std::string str=((std::string)
"Matrix with zero size (")+
1615 "matrix_min_index().";
1621 for(
size_t i=0;i<m;i++) {
1622 for(
size_t j=0;j<n;j++) {
1623 if (data(i,j)<min) {
1635 template<
class mat_t,
class data_t>
1637 data_t &min, data_t &max) {
1644 for(
size_t i=0;i<n;i++) {
1645 for(
size_t j=0;j<m;j++) {
1646 if (data(i,j)<min) {
1648 }
else if (data(i,j)>max) {
1658 template<
class mat_t,
class data_t>
1660 data_t &min, data_t &max) {
1661 return matrix_minmax<mat_t,data_t>
1662 (data.size1(),data.size2(),data,min,max);
1668 template<
class mat_t,
class data_t>
1670 size_t &i_min,
size_t &j_min, data_t &min,
1671 size_t &i_max,
size_t &j_max, data_t &max) {
1683 for(
size_t i=0;i<n;i++) {
1684 for(
size_t j=0;j<m;j++) {
1685 if (data(i,j)<min) {
1689 }
else if (data(i,j)>max) {
1701 template<
class mat_t,
class data_t>
1705 for(
size_t i=0;i<m;i++) {
1706 for(
size_t j=0;j<n;j++) {
1715 template<
class mat_t,
class data_t>
1717 return matrix_sum<mat_t,data_t>(data.size1(),data.size2(),data);
1736 template<
class vec_t>
1739 O2SCL_ERR(
"Empty vector in function vector_lookup().",
1744 while(!std::isfinite(x[i]) && i<n-1) i++;
1750 double best=x[i], bdiff=fabs(x[i]-x0);
1752 if (std::isfinite(x[i]) && fabs(x[i]-x0)<bdiff) {
1755 bdiff=fabs(x[i]-x0);
1773 template<
class vec_t>
1783 template<
class mat_t>
1785 double x0,
size_t &i,
size_t &j) {
1787 O2SCL_ERR(
"Empty matrix in matrix_lookup().",
1791 bool found_one=
false;
1792 for(
size_t i2=0;i2<m;i2++) {
1793 for(
size_t j2=0;j2<n;j2++) {
1794 if (std::isfinite(A(i,j))) {
1795 if (found_one==
false) {
1796 dist=fabs(A(i,j)-x0);
1801 if (fabs(A(i,j)-x0)<dist) {
1802 dist=fabs(A(i,j)-x0);
1810 if (found_one==
false) {
1822 template<
class mat_t>
1824 double x0,
size_t &i,
size_t &j) {
1862 template<
class vec_t,
class data_t>
1864 size_t lo,
size_t hi) {
1866 O2SCL_ERR2(
"Low and high indexes backwards in ",
1867 "function vector_bsearch_inc().",
exc_einval);
1907 template<
class vec_t,
class data_t>
1909 size_t lo,
size_t hi) {
1911 O2SCL_ERR2(
"Low and high indexes backwards in ",
1912 "function vector_bsearch_dec().",
exc_einval);
1933 template<
class vec_t,
class data_t>
1935 size_t lo,
size_t hi) {
1936 if (x[lo]<x[hi-1]) {
1937 return vector_bsearch_inc<vec_t,data_t>(x0,x,lo,hi);
1939 return vector_bsearch_dec<vec_t,data_t>(x0,x,lo,hi);
1949 template<
class vec_t,
class data_t>
1953 if (x[lo]<x[hi-1]) {
1954 return vector_bsearch_inc<vec_t,data_t>(x0,x,lo,hi);
1956 return vector_bsearch_dec<vec_t,data_t>(x0,x,lo,hi);
1975 template<
class vec_t>
1980 if (data[0]==data[1]) {
1982 }
else if (data[0]<data[1]) {
1992 for(
size_t i=0;i<n-1 && done==
false;i++) {
1993 if (data[i]!=data[i+1]) {
2007 if (data[start]<data[start+1])
return 1;
2014 if (data[start]>data[start+1]) inc=
false;
2017 for(
size_t i=start+1;i<n-1;i++) {
2019 if (data[i]>data[i+1])
return 0;
2025 for(
size_t i=start+1;i<n-1;i++) {
2026 if (data[i]<data[i+1])
return 0;
2057 template<
class vec_t>
2065 if (data[0]==data[1]) {
2067 }
else if (data[0]>data[1]) {
2072 for(
size_t i=1;i<n-1;i++) {
2074 if (data[i]>=data[i+1])
return 0;
2080 for(
size_t i=1;i<n-1;i++) {
2081 if (data[i]<=data[i+1])
return 0;
2095 template<
class vec_t>
2108 template<
class vec_t>
2110 for(
size_t i=0;i<n;i++) {
2111 if (!std::isfinite(data[i]))
return false;
2136 template<
class vec_t,
class data_t>
2140 for(
size_t i=0;i<n;i++) {
2151 template<
class vec_t,
class data_t> data_t
vector_sum(vec_t &data) {
2153 for(
size_t i=0;i<data.size();i++) {
2167 for(
size_t i=0;i<n;i++) {
2181 for(
size_t i=0;i<data.size();i++) {
2193 template<
class vec_t,
class data_t>
2201 }
else if (n == 1) {
2205 for (
size_t i = 0; i < n; i++) {
2206 const data_t xx = x[i];
2209 const data_t ax = fabs(xx);
2212 ssq = 1.0 + ssq * (scale / ax) * (scale / ax);
2215 ssq += (ax / scale) * (ax / scale);
2221 return scale * sqrt(ssq);
2227 template<
class vec_t,
class data_t> data_t
vector_norm(
const vec_t &x) {
2228 return vector_norm<vec_t,data_t>(x.size(),x);
2237 template<
class vec_t>
2245 }
else if (n == 1) {
2249 for (
size_t i = 0; i < n; i++) {
2250 const double xx = x[i];
2253 const double ax = fabs(xx);
2256 ssq = 1.0 + ssq * (scale / ax) * (scale / ax);
2259 ssq += (ax / scale) * (ax / scale);
2265 return scale * sqrt(ssq);
2271 return vector_norm_double<vec_t>(x.size(),x);
2279 template<
class vec_t,
class data_t>
2281 for(
size_t i=0;i<N;i++) {
2289 template<
class vec_t,
class data_t>
2298 template<
class vec_t,
class vec2_t>
2300 if (src.size()==0) {
2303 if (iout>=src.size()) {
2306 dest.resize(src.size()-1);
2308 for(
size_t i=0;i<src.size();i++) {
2330 template<
class vec_t,
class data_t>
2333 data_t *tmp=
new data_t[n];
2334 for(
size_t i=0;i<n;i++) {
2335 tmp[i]=data[(i+k)%n];
2337 for(
size_t i=0;i<n;i++) {
2349 template<
class vec_t,
class data_t>
2353 for(
size_t i=0;i<n/2;i++) {
2355 data[n-1-i]=data[i];
2366 template<
class vec_t,
class data_t>
2369 size_t n=data.size();
2371 for(
size_t i=0;i<n/2;i++) {
2373 data[n-1-i]=data[i];
2384 template<
class vec_t>
2388 for(
size_t i=0;i<n/2;i++) {
2390 data[n-1-i]=data[i];
2403 size_t n=data.size();
2405 for(
size_t i=0;i<n/2;i++) {
2407 data[n-1-i]=data[i];
2429 template<
class mat_t,
class mat_row_t>
2431 return mat_row_t(M,row);
2485 template<
class mat_t,
class mat_column_t>
2487 return mat_column_t(M,column);
2507 double &operator[](
size_t i) {
2508 return m_(i,column_);
2510 const double &operator[](
size_t i)
const {
2511 return m_(i,column_);
2527 template<
class vec_t>
2529 bool endline=
false) {
2534 for(
size_t i=0;i<n-1;i++) os << v[i] <<
" ";
2536 if (endline) os << std::endl;
2551 template<
class vec_t>
2552 void vector_out(std::ostream &os,
const vec_t &v,
bool endline=
false) {
2559 for(
size_t i=0;i<n-1;i++) os << v[i] <<
" ";
2561 if (endline) os << std::endl;
2567 template<
class vec_t,
class data_t>
2569 g.template vector<vec_t>(v);
2574 template<
class mat_t>
2576 for(
size_t i=0;i<M;i++) {
2577 for(
size_t j=0;j<N;j++) {
2578 if (i==j) m(i,j)=1.0;
2593 (dat_t *v,
size_t start,
size_t last) {
2603 (
const dat_t *v,
size_t start,
size_t last) {
2619 size_t start,
size_t last) {
2622 (v,boost::numeric::ublas::range(start,last));
2637 size_t start,
size_t last) {
2640 (v,boost::numeric::ublas::range(start,last));
2656 size_t start,
size_t last) {
2659 (v,boost::numeric::ublas::range(start,last));
2672 template<
class dat_t>
2679 size_t start,
size_t last) {
2683 (v,boost::numeric::ublas::range(start,last));
2696 template<
class dat_t>
2703 size_t start,
size_t last) {
2707 (v,boost::numeric::ublas::range(start,last));
2720 template<
class dat_t>
2727 size_t start,
size_t last) {
2731 (v,boost::numeric::ublas::range(start,last));
2744 template<
class dat_t>
2751 size_t start,
size_t last) {
2755 (v,boost::numeric::ublas::range(start,last));
2782 start_(start), last_(last) {
2783 #if !O2SCL_NO_RANGE_CHECK 2785 O2SCL_ERR2(
"End before beginning in vector_range_gen::",
2786 "vector_range_gen(vec_t,size_t,size_t)",
2794 size_t last) : v_(v2.v_),
2795 start_(start+v2.start_), last_(last+v2.start_) {
2796 #if !O2SCL_NO_RANGE_CHECK 2798 O2SCL_ERR2(
"End before beginning in vector_range_gen::",
2799 "vector_range_gen(vector_range_gen,size_t,size_t)",
2802 if (last>v2.
last_) {
2803 O2SCL_ERR2(
"End beyond end of previous vector in vector_range_gen::",
2804 "vector_range_gen(vector_range_gen,size_t,size_t)",
2812 return last_-start_;
2817 #if !O2SCL_NO_RANGE_CHECK 2818 if (i+start_>=last_) {
2819 O2SCL_ERR(
"Index out of range in vector_range_gen::operator[].",
2823 return v_[i+start_];
2828 #if !O2SCL_NO_RANGE_CHECK 2829 if (i+start_>=last_) {
2834 return v_[i+start_];
2857 start_(start), last_(last) {
2858 #if !O2SCL_NO_RANGE_CHECK 2860 O2SCL_ERR2(
"End before beginning in vector_range_gen::",
2861 "vector_range_gen(vec_t,size_t,size_t)",
2869 size_t last) : v_(v2.v_),
2870 start_(start+v2.start_), last_(last+v2.start_) {
2871 #if !O2SCL_NO_RANGE_CHECK 2873 O2SCL_ERR2(
"End before beginning in vector_range_gen::",
2874 "vector_range_gen(vector_range_gen,size_t,size_t)",
2877 if (last>v2.
last_) {
2878 O2SCL_ERR2(
"End beyond end of previous vector in vector_range_gen::",
2879 "vector_range_gen(vector_range_gen,size_t,size_t)",
2887 size_t last) : v_(v2.v_),
2888 start_(start+v2.start_), last_(last+v2.start_) {
2889 #if !O2SCL_NO_RANGE_CHECK 2891 O2SCL_ERR2(
"End before beginning in vector_range_gen::",
2892 "vector_range_gen(vector_range_gen,size_t,size_t)",
2895 if (last>v2.
last_) {
2896 O2SCL_ERR2(
"End beyond end of previous vector in vector_range_gen::",
2897 "vector_range_gen(vector_range_gen,size_t,size_t)",
2905 return last_-start_;
2910 #if !O2SCL_NO_RANGE_CHECK 2911 if (i+start_>=last_) {
2916 return v_[i+start_];
2959 size_t start,
size_t last) {
2968 size_t start,
size_t last) {
2977 size_t start,
size_t last) {
2992 template<
class dat_t> std::vector<dat_t>
2994 return std::vector<dat_t> (v.begin()+start,v.begin()+last);
3008 template<
class dat_t>
const std::vector<dat_t>
3011 return std::vector<dat_t> (v.begin()+start,v.begin()+last);
3018 #if defined (O2SCL_COND_FLAG) || defined (DOXYGEN) 3020 #if defined (O2SCL_ARMA) || defined (DOXYGEN) 3021 #include <armadillo> 3033 template<> arma::subview_row<double>
3034 matrix_row<arma::mat,arma::subview_row<double> >
3035 (arma::mat &M,
size_t row);
3038 template<> arma::subview_col<double>
3039 matrix_column<arma::mat,arma::subview_col<double> >
3040 (arma::mat &M,
size_t column);
3047 #if defined (O2SCL_EIGEN) || defined (DOXYGEN) 3048 #include <eigen3/Eigen/Dense> 3055 double matrix_max(
const Eigen::MatrixXd &data);
3058 double matrix_min(
const Eigen::MatrixXd &data);
3061 template<> Eigen::MatrixXd::RowXpr
3062 matrix_row<Eigen::MatrixXd,Eigen::MatrixXd::RowXpr>
3063 (Eigen::MatrixXd &M,
size_t row);
3066 template<> Eigen::MatrixXd::ColXpr
3067 matrix_column<Eigen::MatrixXd,Eigen::MatrixXd::ColXpr>
3068 (Eigen::MatrixXd &M,
size_t column);
3077 #include <o2scl/vector_special.h> 3124 template<
class T,
class F,
class A>
class matrix {
data_t vector_min_quad(size_t n, const vec_t &data)
Minimum of vector by quadratic fit.
size_t vector_lookup(size_t n, const vec_t &x, double x0)
Lookup the value x0 in the first n elements of vector x.
size_t vector_min_index(size_t n, const vec_t &data)
Compute the index which holds the minimum of the first n elements of a vector.
Experimental vector range object.
data_t vector_min_value(size_t n, const vec_t &data)
Compute the minimum of the first n elements of a vector.
data_t matrix_sum(size_t m, size_t n, const mat_t &data)
Compute the sum of matrix elements.
void vector_copy_jackknife(const vec_t &src, size_t iout, vec2_t &dest)
From a given vector, create a new vector by removing a specified element.
data_t vector_max_quad_loc(size_t n, const vec_t &x, const vec_t &y)
Location of vector maximum by quadratic fit.
data_t vector_min_quad_loc(size_t n, const vec_t &x, const vec_t &y)
Location of vector minimum by quadratic fit.
The main O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl names...
void matrix_swap_rows(size_t N, mat_t &m, size_t i1, size_t i2)
Swap the first N columns of two rows in a matrix.
void matrix_make_upper(mat_t &src)
Make a matrix upper triangular by setting the lower triangular entries to zero.
bool vector_is_finite(size_t n, vec_t &data)
Test if the first n elements of a vector are finite.
Placeholder documentation of some related Boost objects.
size_t vector_bsearch_dec(const data_t x0, const vec_t &x, size_t lo, size_t hi)
Binary search a part of an decreasing vector for x0.
void matrix_max_index(size_t m, size_t n, const mat_t &data, size_t &i_max, size_t &j_max, data_t &max)
Compute the maximum of a matrix and return the indices of the maximum element.
const_vector_range_gen(const vec_t &v, size_t start, size_t last)
Create an object starting with index start in vector v.
Generic object which represents a column of a matrix.
bool matrix_is_upper(mat_t &src)
Simple test that a matrix is upper triangular.
bool matrix_is_lower(mat_t &src)
Simple test that a matrix is lower triangular.
double & operator[](size_t i)
Return a reference ith element.
void vector_smallest_index(size_t n, vec_t &data, size_t k, vec_size_t &index)
Find the indexes of the k smallest entries among the first n entries of a vector. ...
size_t vector_bsearch_inc(const data_t x0, const vec_t &x, size_t lo, size_t hi)
Binary search a part of an increasing vector for x0.
A simple convenience wrapper for GSL vector objects.
invalid argument supplied by user
Generic object which represents a row of a matrix.
vector_range_gen< vec_t > vector_range(vector_range_gen< vec_t > &v, size_t start, size_t last)
Recursively create a o2scl::vector_range_gen object from a vector range.
void vector_reverse(size_t n, vec_t &data)
Reverse the first n elements of a vector.
size_t last_
The end() iterator.
void vector_rotate(size_t n, vec_t &data, size_t k)
"Rotate" a vector so that the kth element is now the beginning
void vector_swap_double(size_t N, vec_t &v1, vec2_t &v2)
Swap of of the first N elements of two double-precision vectors.
size_t last_
The end() iterator.
double matrix_min_value_double(const mat_t &data)
Compute the minimum of a matrix.
The default matrix type from uBlas.
void vector_min(size_t n, const vec_t &data, size_t &index, data_t &val)
Compute the minimum of the first n elements of a vector.
double vector_norm_double(size_t n, const vec_t &x)
Compute the norm of the first n entries of a vector of double precision numbers.
void matrix_make_lower(mat_t &src)
Make a matrix lower triangular by setting the upper triangular entries to zero.
double matrix_max_value_double(const mat_t &data)
Compute the maximum of a matrix.
data_t vector_sum(size_t n, vec_t &data)
Compute the sum of the first n elements of a vector.
void vector_swap(size_t N, vec_t &v1, vec2_t &v2)
Swap the first N elements of two vectors.
void vector_largest(size_t n, vec_t &data, size_t k, vec_t &largest)
Find the k largest entries of the first n elements of a vector.
std::vector< dat_t > vector_range_copy(const std::vector< dat_t > &v, size_t start, size_t last)
Vector range function template for std::vector
const_vector_range_gen(const vector_range_gen< vec_t > &v2, size_t start, size_t last)
Create an object from a previously constructed range object.
matrix_row_gen(mat_t &m, size_t row)
Create a row object from row row of matrix m.
vec_t & v_
A reference to the original vector.
void matrix_swap_double(size_t M, size_t N, mat_t &m1, mat2_t &m2)
Swap of the first elements in two double-precision matrices.
void vector_sort_index(size_t n, const vec_t &data, vec_size_t &order)
Create a permutation which sorts a vector (in increasing order)
void matrix_swap_cols(size_t M, mat_t &m, size_t j1, size_t j2)
Swap the first M rows of two columns in a matrix.
void vector_copy(const vec_t &src, vec2_t &dest)
Simple vector copy.
void matrix_copy(mat_t &src, mat2_t &dest)
Simple matrix copy.
size_t size() const
Return the vector size.
void matrix_minmax_index(size_t n, size_t m, const mat_t &data, size_t &i_min, size_t &j_min, data_t &min, size_t &i_max, size_t &j_max, data_t &max)
Compute the minimum and maximum of a matrix and return their locations.
data_t vector_max_value(size_t n, const vec_t &data)
Compute the maximum of the first n elements of a vector.
const vec_t & v_
A reference to the original vector.
dat_t * vector_range(dat_t *v, size_t start, size_t last)
Vector range function for pointers.
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
const_vector_range_gen(const const_vector_range_gen &v2, size_t start, size_t last)
Create an object from a previously constructed range object.
void vector_minmax_value(size_t n, vec_t &data, data_t &min, data_t &max)
Compute the minimum and maximum of the first n elements of a vector.
void matrix_set_identity(size_t M, size_t N, mat_t &m)
Set a matrix to unity on the diagonal and zero otherwise.
data_t vector_max_quad(size_t n, const vec_t &data)
Maximum of vector by quadratic fit.
vector_range_gen(const vector_range_gen &v2, size_t start, size_t last)
Create an object from a previously constructed range object.
mat_column_t matrix_column(mat_t &M, size_t column)
Construct a column of a matrix.
size_t start_
The index offset.
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
mat_row_t matrix_row(mat_t &M, size_t row)
Construct a row of a matrix.
void vector_out(std::ostream &os, size_t n, const vec_t &v, bool endline=false)
Output the first n elements of a vector to a stream, os.
void vector_smallest(size_t n, vec_t &data, size_t k, vec_t &smallest)
Find the k smallest entries of the first n elements of a vector.
void vector_grid(uniform_grid< data_t > g, vec_t &v)
Fill a vector with a specified grid.
void vector_minmax(size_t n, vec_t &data, size_t &ix_min, data_t &min, size_t &ix_max, data_t &max)
Compute the minimum and maximum of the first n elements of a vector.
data_t vector_norm(size_t n, const vec_t &x)
Compute the norm of the first n entries of a vector of floating-point (single or double precision) nu...
A simple convenience wrapper for GSL matrix objects.
size_t vector_max_index(size_t n, const vec_t &data)
Compute the index which holds the maximum of the first n elements of a vector.
void vector_max(size_t n, const vec_t &data, size_t &index, data_t &val)
Compute the maximum of the first n elements of a vector.
const double & operator[](size_t i) const
Return a const reference to the ith column of the selected row.
void vector_minmax_index(size_t n, vec_t &data, size_t &ix_min, size_t &ix_max)
Compute the minimum and maximum of the first n elements of a vector.
double matrix_max(const arma::mat &data)
Armadillo version of matrix_max()
const double & operator[](size_t i) const
Return a const reference ith element.
const double & operator[](size_t i) const
Return a const reference ith element.
double vector_sum_double(size_t n, vec_t &data)
Compute the sum of the first n elements of a vector of double-precision numbers.
vector_range_gen(vec_t &v, size_t start, size_t last)
Create an object starting with index start in vector v.
void matrix_minmax(size_t n, size_t m, const mat_t &data, data_t &min, data_t &max)
Compute the minimum and maximum of a matrix.
Experimental const vector range object.
size_t vector_bsearch(const data_t x0, const vec_t &x, size_t lo, size_t hi)
Binary search a part of a monotonic vector for x0.
void matrix_transpose(mat_t &src, mat2_t &dest)
Simple transpose.
void sort_index_downheap(size_t N, const vec_t &data, vec_size_t &order, size_t k)
Provide a downheap() function for vector_sort_index()
void vector_reverse_double(size_t n, vec_t &data)
Reverse the first n elements in a vector of double precision numbers.
int vector_is_monotonic(size_t n, vec_t &data)
Test if the first n elements of a vector are monotonic and increasing or decreasing.
int vector_is_strictly_monotonic(size_t n, vec_t &data)
Test if the first n elements of a vector are strictly monotonic and determine if they are increasing ...
void matrix_swap_cols_double(size_t M, mat_t &m, size_t j1, size_t j2)
Swap the first M rows of two columns in a double-precision matrix.
void matrix_min_index(size_t n, size_t m, const mat_t &data, size_t &i_min, size_t &j_min, data_t &min)
Compute the minimum of a matrix and return the indices of the minimum element.
double matrix_min(const arma::mat &data)
Armadillo version of matrix_min()
void matrix_swap(size_t M, size_t N, mat_t &v1, mat2_t &v2)
Swap of the first elements in two matrices.
void matrix_swap_rows_double(size_t N, mat_t &m, size_t i1, size_t i2)
Swap the first N columns of two rows in a double-precision matrix.
void vector_set_all(size_t N, vec_t &src, data_t &val)
Set the first N entries in a vector to a particular value.
data_t matrix_max_value(size_t m, const size_t n, const mat_t &data)
Compute the maximum of the lower-left part of a matrix.
size_t start_
The index offset.
std::string itos(int x)
Convert an integer to a string.
std::string szttos(size_t x)
Convert a size_t to a string.
const dat_t * const_vector_range(const dat_t *v, size_t start, size_t last)
Vector range function for const pointers.
void matrix_lookup(size_t m, size_t n, const mat_t &A, double x0, size_t &i, size_t &j)
Lookup an element in the first $(m,n)$ entries in a matrix.
double & operator[](size_t i)
Return a reference to the ith column of the selected row.
data_t matrix_min_value(size_t m, size_t n, const mat_t &data)
Compute the minimum of a matrix.
void vector_sort_double(size_t n, vec_t &data)
Sort a vector of doubles (in increasing order)
The default vector type from uBlas.
size_t size() const
Return the vector size.
void vector_sort(size_t n, vec_t &data)
Sort a vector (in increasing order)
void sort_downheap(vec_t &data, size_t n, size_t k)
Provide a downheap() function for vector_sort()