Go to the documentation of this file.
23 #ifndef O2SCL_VECTOR_H
24 #define O2SCL_VECTOR_H
54 #include <gsl/gsl_vector.h>
55 #include <gsl/gsl_sys.h>
56 #include <gsl/gsl_matrix.h>
58 #include <o2scl/misc.h>
59 #include <o2scl/uniform_grid.h>
73 d=(
const double *)m->data;
76 double operator[](
size_t i)
const {
96 d=(
const double *)m->data;
101 double operator()(
size_t i,
size_t j)
const {
104 size_t size1()
const {
107 size_t size2()
const {
123 template<
class vec_t,
class vec2_t>
126 if (dest.size()<N) dest.resize(N);
131 for(i=m;i+3<N;i+=4) {
152 template<
class vec_t,
class vec2_t>
158 for(i=m;i+3<N;i+=4) {
175 template<
class mat_t,
class mat2_t>
177 size_t m=src.size1();
178 size_t n=src.size2();
179 if (dest.size1()<m || dest.size2()<n) dest.resize(m,n);
180 for(
size_t i=0;i<m;i++) {
181 for(
size_t j=0;j<n;j++) {
200 template<
class mat_t,
class mat2_t>
202 for(
size_t i=0;i<M;i++) {
203 for(
size_t j=0;j<N;j++) {
221 template<
class mat_t,
class mat2_t>
223 size_t m=src.size1();
224 size_t n=src.size2();
225 if (dest.size1()<n || dest.size2()<m) dest.resize(n,m);
226 for(
size_t i=0;i<m;i++) {
227 for(
size_t j=0;j<n;j++) {
242 template<
class mat_t,
class mat2_t>
244 for(
size_t i=0;i<m;i++) {
245 for(
size_t j=0;j<n;j++) {
261 template<
class mat_t,
class data_t>
263 size_t m=src.size1();
264 size_t n=src.size2();
268 for(
size_t i=0;i<n;i++) {
269 for(
size_t j=0;j<n;j++) {
286 template<
class mat_t,
class data_t>
289 for(
size_t i=0;i<m;i++) {
290 for(
size_t j=0;j<n;j++) {
304 size_t m=src.size1();
305 size_t n=src.size2();
307 for(
size_t i=0;i<m;i++) {
308 for(
size_t j=i+1;j<n;j++) {
309 if (src(i,j)!=0.0)
return false;
318 size_t m=src.size1();
319 size_t n=src.size2();
321 for(
size_t j=0;j<n;j++) {
322 for(
size_t i=j+1;i<m;i++) {
323 if (src(i,j)!=0.0)
return false;
333 size_t m=src.size1();
334 size_t n=src.size2();
335 for(
size_t i=0;i<m;i++) {
336 for(
size_t j=i+1;j<n;j++) {
347 size_t m=src.size1();
348 size_t n=src.size2();
349 for(
size_t j=0;j<n;j++) {
350 for(
size_t i=j+1;i<m;i++) {
363 for(
size_t i=0;i<m;i++) {
364 for(
size_t j=i+1;j<n;j++) {
365 if (src(i,j)!=0.0)
return false;
377 for(
size_t j=0;j<n;j++) {
378 for(
size_t i=j+1;i<m;i++) {
379 if (src(i,j)!=0.0)
return false;
390 for(
size_t i=0;i<m;i++) {
391 for(
size_t j=i+1;j<n;j++) {
403 for(
size_t j=0;j<n;j++) {
404 for(
size_t i=j+1;i<m;i++) {
419 template<
class vec_t,
class vec2_t,
class data_t>
428 for(i=m;i+3<N;i+=4) {
455 template<
class vec_t,
class vec2_t,
class data_t>
458 if (v2.size()<N) N=v2.size();
466 for(i=m;i+3<N;i+=4) {
489 template<
class vec_t,
class vec2_t>
491 return vector_swap<vec_t,vec2_t,double>(N,v1,v2);
505 template<
class vec_t,
class vec2_t>
507 return vector_swap<vec_t,vec2_t,double>(v1,v2);
515 template<
class vec_t,
class data_t>
530 template<
class vec_t>
532 return vector_swap<vec_t,double>(v,i,j);
541 template<
class mat_t,
class mat2_t,
class data_t>
544 for(
size_t i=0;i<M;i++) {
545 for(
size_t j=0;j<N;j++) {
560 template<
class mat_t,
class mat2_t,
class data_t>
562 return matrix_swap<mat_t,mat2_t,double>(M,N,m1,m2);
570 template<
class mat_t,
class data_t>
571 void matrix_swap(mat_t &m,
size_t i1,
size_t j1,
size_t i2,
size_t j2) {
572 data_t temp=m(i1,j1);
583 template<
class mat_t>
585 size_t i2,
size_t j2) {
586 return matrix_swap<mat_t,double>(m,i1,j1,i2,j2);
594 template<
class mat_t,
class data_t>
597 for(
size_t i=0;i<M;i++) {
611 template<
class mat_t>
613 return matrix_swap_cols<mat_t,double>(M,m,j1,j2);
622 template<
class mat_t,
class data_t>
625 for(
size_t j=0;j<N;j++) {
639 template<
class mat_t>
641 return matrix_swap_rows<mat_t,double>(N,m,i1,i2);
649 template<
class vec_t,
class data_t>
657 if (j<n && data[j] < data[j+1]) j++;
658 if (!(v < data[j]))
break;
704 template<
class vec_t,
class data_t>
717 sort_downheap<vec_t,data_t>(data,N,k);
725 sort_downheap<vec_t,data_t>(data,N,0);
733 template<
class vec_t,
class vec_
size_t>
737 const size_t pki = order[k];
742 if (j < N && data[order[j]] < data[order[j + 1]]) {
747 if (!(data[pki] < data[order[j]])) {
799 template<
class vec_t,
class vec_
size_t>
808 for (i = 0 ; i < n ; i++) {
822 sort_index_downheap<vec_t,vec_size_t>(N,data,order,k);
828 size_t tmp = order[0];
835 sort_index_downheap<vec_t,vec_size_t>(N,data,order,0);
879 template<
class vec_t,
class vec_
size_t>
894 template<
class vec_t>
896 return vector_sort<vec_t,double>(n,data);
919 template<
class vec_t,
class data_t>
922 O2SCL_ERR2(
"Subset length greater than size in ",
932 data_t xbound=data[0];
936 for(
size_t i=1;i<n;i++) {
940 }
else if (xi>=xbound) {
944 for(i1=j-1;i1>0;i1--) {
945 if (xi>smallest[i1-1])
break;
946 smallest[i1]=smallest[i1-1];
949 xbound=smallest[j-1];
970 template<
class vec_t,
class data_t>
972 size_t n=data.size();
973 if (smallest.size()<k) smallest.resize(k);
992 template<
class vec_t,
class data_t,
class vec_
size_t>
996 O2SCL_ERR2(
"Subset length greater than size in ",
997 "function vector_smallest_index().",
exc_einval);
1001 "function vector_smallest_index().",
exc_einval);
1008 data_t xbound=data[0];
1012 for(
size_t i=1;i<n;i++) {
1016 }
else if (xi>=xbound) {
1020 for(i1=j-1;i1>0;i1--) {
1021 if (xi>data[index[i1-1]])
break;
1022 index[i1]=index[i1-1];
1025 xbound=data[index[j-1]];
1032 template<
class vec_t,
class data_t,
class vec_
size_t>
1034 vec_size_t &index) {
1035 size_t n=data.size();
1036 if (index.size()<k) index.resize(k);
1037 return o2scl::vector_smallest_index<vec_t,data_t,vec_size_t>
1057 template<
class vec_t,
class data_t>
1060 O2SCL_ERR2(
"Subset length greater than size in ",
1070 data_t xbound=data[0];
1074 for(
size_t i=1;i<n;i++) {
1078 }
else if (xi<=xbound) {
1082 for(i1=j-1;i1>0;i1--) {
1083 if (xi<largest[i1-1])
break;
1084 largest[i1]=largest[i1-1];
1087 xbound=largest[j-1];
1108 template<
class vec_t,
class data_t>
1110 size_t n=data.size();
1111 if (largest.size()<k) largest.resize(k);
1120 template<
class vec_t,
class data_t>
1127 for(
size_t i=1;i<n;i++) {
1137 template<
class vec_t,
class data_t>
1140 size_t n=data.size();
1145 for(
size_t i=1;i<n;i++) {
1156 template<
class vec_t,
class data_t>
1164 for(
size_t i=1;i<n;i++) {
1175 template<
class vec_t,
class data_t>
1184 for(
size_t i=1;i<n;i++) {
1195 template<
class vec_t,
class data_t>
1202 for(
size_t i=1;i<n;i++) {
1212 template<
class vec_t,
class data_t>
1215 size_t n=data.size();
1220 for(
size_t i=1;i<n;i++) {
1231 template<
class vec_t,
class data_t>
1239 for(
size_t i=1;i<n;i++) {
1250 template<
class vec_t,
class data_t>
1259 for(
size_t i=1;i<n;i++) {
1271 template<
class vec_t,
class data_t>
1273 data_t &min, data_t &max) {
1280 for(
size_t i=1;i<n;i++) {
1294 template<
class vec_t,
class data_t>
1296 size_t &ix_min,
size_t &ix_max) {
1305 for(
size_t i=1;i<n;i++) {
1321 template<
class vec_t,
class data_t>
1323 size_t &ix_min, data_t &min,
1324 size_t &ix_max, data_t &max) {
1333 for(
size_t i=1;i<n;i++) {
1351 template<
class vec_t,
class data_t>
1353 size_t ix=vector_max_index<vec_t,data_t>(n,data);
1355 return quadratic_extremum_y<data_t>(0,1,2,data[0],data[1],data[2]);
1356 }
else if (ix==n-1) {
1357 return quadratic_extremum_y<data_t>
1358 (n-3,n-2,n-1,data[n-3],data[n-2],data[n-1]);
1360 return quadratic_extremum_y<data_t>
1361 (ix-1,ix,ix+1,data[ix-1],data[ix],data[ix+1]);
1366 template<
class vec_t,
class data_t>
1368 size_t ix=vector_max_index<vec_t,data_t>(n,y);
1369 if (ix==0 || ix==n-1)
return y[ix];
1370 return quadratic_extremum_y<data_t>(x[ix-1],x[ix],x[ix+1],
1371 y[ix-1],y[ix],y[ix+1]);
1376 template<
class vec_t,
class data_t>
1378 size_t ix=vector_max_index<vec_t,data_t>(n,y);
1379 if (ix==0 || ix==n-1)
return y[ix];
1380 return quadratic_extremum_x<data_t>(x[ix-1],x[ix],x[ix+1],
1381 y[ix-1],y[ix],y[ix+1]);
1386 template<
class vec_t,
class data_t>
1388 size_t ix=vector_min_index<vec_t,data_t>(n,data);
1390 return quadratic_extremum_y<data_t>(0,1,2,data[0],data[1],data[2]);
1391 }
else if (ix==n-1) {
1392 return quadratic_extremum_y<data_t>
1393 (n-3,n-2,n-1,data[n-3],data[n-2],data[n-1]);
1395 return quadratic_extremum_y<data_t>
1396 (ix-1,ix,ix+1,data[ix-1],data[ix],data[ix+1]);
1401 template<
class vec_t,
class data_t>
1403 size_t ix=vector_min_index<vec_t,data_t>(n,y);
1404 if (ix==0 || ix==n-1)
return y[ix];
1405 return quadratic_extremum_y<data_t>(x[ix-1],x[ix],x[ix+1],
1406 y[ix-1],y[ix],y[ix+1]);
1411 template<
class vec_t,
class data_t>
1413 size_t ix=vector_min_index<vec_t,data_t>(n,y);
1414 if (ix==0 || ix==n-1)
return y[ix];
1415 return quadratic_extremum_x<data_t>(x[ix-1],x[ix],x[ix+1],
1416 y[ix-1],y[ix],y[ix+1]);
1424 template<
class mat_t,
class data_t>
1428 std::string str=((std::string)
"Matrix with zero size (")+
1430 "matrix_max_value().";
1433 data_t max=data(0,0);
1434 for(
size_t i=0;i<m;i++) {
1435 for(
size_t j=0;j<n;j++) {
1436 if (data(i,j)>max) {
1446 template<
class mat_t,
class data_t> data_t
1448 size_t m=data.size1();
1449 size_t n=data.size2();
1451 std::string str=((std::string)
"Matrix with zero size (")+
1453 "matrix_max_value().";
1456 data_t max=data(0,0);
1457 for(
size_t i=0;i<m;i++) {
1458 for(
size_t j=0;j<n;j++) {
1459 if (data(i,j)>max) {
1469 template<
class mat_t>
double
1471 size_t m=data.size1();
1472 size_t n=data.size2();
1474 std::string str=((std::string)
"Matrix with zero size (")+
1476 "matrix_max_value_double().";
1479 double max=data(0,0);
1480 for(
size_t i=0;i<m;i++) {
1481 for(
size_t j=0;j<n;j++) {
1482 if (data(i,j)>max) {
1493 template<
class mat_t,
class data_t>
1495 size_t &i_max,
size_t &j_max, data_t &max) {
1498 std::string str=((std::string)
"Matrix with zero size (")+
1500 "matrix_max_index().";
1506 for(
size_t i=0;i<m;i++) {
1507 for(
size_t j=0;j<n;j++) {
1508 if (data(i,j)>max) {
1521 template<
class mat_t,
class data_t>
1523 size_t &i_max,
size_t &j_max, data_t &max) {
1525 size_t m=data.size1();
1526 size_t n=data.size2();
1528 std::string str=((std::string)
"Matrix with zero size (")+
1530 "matrix_max_index().";
1536 for(
size_t i=0;i<m;i++) {
1537 for(
size_t j=0;j<n;j++) {
1538 if (data(i,j)>max) {
1550 template<
class mat_t,
class data_t>
1554 std::string str=((std::string)
"Matrix with zero size (")+
1556 "matrix_min_value().";
1559 data_t min=data(0,0);
1560 for(
size_t i=0;i<m;i++) {
1561 for(
size_t j=0;j<n;j++) {
1562 if (data(i,j)<min) {
1572 template<
class mat_t,
class data_t>
1575 size_t m=data.size1();
1576 size_t n=data.size2();
1578 std::string str=((std::string)
"Matrix with zero size (")+
1580 "matrix_min_value().";
1583 data_t min=data(0,0);
1584 for(
size_t i=0;i<m;i++) {
1585 for(
size_t j=0;j<n;j++) {
1586 if (data(i,j)<min) {
1596 template<
class mat_t>
1599 size_t m=data.size1();
1600 size_t n=data.size2();
1602 std::string str=((std::string)
"Matrix with zero size (")+
1604 "matrix_min_value().";
1607 double min=data(0,0);
1608 for(
size_t i=0;i<m;i++) {
1609 for(
size_t j=0;j<n;j++) {
1610 if (data(i,j)<min) {
1621 template<
class mat_t,
class data_t>
1623 size_t &i_min,
size_t &j_min, data_t &min) {
1626 std::string str=((std::string)
"Matrix with zero size (")+
1628 "matrix_min_index().";
1634 for(
size_t i=0;i<m;i++) {
1635 for(
size_t j=0;j<n;j++) {
1636 if (data(i,j)<min) {
1649 template<
class mat_t,
class data_t>
1651 size_t &i_min,
size_t &j_min, data_t &min) {
1653 size_t m=data.size1();
1654 size_t n=data.size2();
1656 std::string str=((std::string)
"Matrix with zero size (")+
1658 "matrix_min_index().";
1664 for(
size_t i=0;i<m;i++) {
1665 for(
size_t j=0;j<n;j++) {
1666 if (data(i,j)<min) {
1678 template<
class mat_t,
class data_t>
1680 data_t &min, data_t &max) {
1687 for(
size_t i=0;i<n;i++) {
1688 for(
size_t j=0;j<m;j++) {
1689 if (data(i,j)<min) {
1691 }
else if (data(i,j)>max) {
1701 template<
class mat_t,
class data_t>
1703 data_t &min, data_t &max) {
1704 return matrix_minmax<mat_t,data_t>
1705 (data.size1(),data.size2(),data,min,max);
1711 template<
class mat_t,
class data_t>
1713 size_t &i_min,
size_t &j_min, data_t &min,
1714 size_t &i_max,
size_t &j_max, data_t &max) {
1726 for(
size_t i=0;i<n;i++) {
1727 for(
size_t j=0;j<m;j++) {
1728 if (data(i,j)<min) {
1732 }
else if (data(i,j)>max) {
1744 template<
class mat_t,
class data_t>
1748 for(
size_t i=0;i<m;i++) {
1749 for(
size_t j=0;j<n;j++) {
1758 template<
class mat_t,
class data_t>
1760 return matrix_sum<mat_t,data_t>(data.size1(),data.size2(),data);
1779 template<
class vec_t>
1782 O2SCL_ERR(
"Empty vector in function vector_lookup().",
1787 while(!std::isfinite(x[i]) && i<n-1) i++;
1793 double best=x[i], bdiff=fabs(x[i]-x0);
1795 if (std::isfinite(x[i]) && fabs(x[i]-x0)<bdiff) {
1798 bdiff=fabs(x[i]-x0);
1816 template<
class vec_t>
1826 template<
class mat_t>
1828 double x0,
size_t &i,
size_t &j) {
1830 O2SCL_ERR(
"Empty matrix in matrix_lookup().",
1834 bool found_one=
false;
1835 for(
size_t i2=0;i2<m;i2++) {
1836 for(
size_t j2=0;j2<n;j2++) {
1837 if (std::isfinite(A(i,j))) {
1838 if (found_one==
false) {
1839 dist=fabs(A(i,j)-x0);
1844 if (fabs(A(i,j)-x0)<dist) {
1845 dist=fabs(A(i,j)-x0);
1853 if (found_one==
false) {
1865 template<
class mat_t>
1867 double x0,
size_t &i,
size_t &j) {
1905 template<
class vec_t,
class data_t>
1907 size_t lo,
size_t hi) {
1909 O2SCL_ERR2(
"Low and high indexes backwards in ",
1910 "function vector_bsearch_inc().",
exc_einval);
1950 template<
class vec_t,
class data_t>
1952 size_t lo,
size_t hi) {
1954 O2SCL_ERR2(
"Low and high indexes backwards in ",
1955 "function vector_bsearch_dec().",
exc_einval);
1976 template<
class vec_t,
class data_t>
1978 size_t lo,
size_t hi) {
1979 if (x[lo]<x[hi-1]) {
1980 return vector_bsearch_inc<vec_t,data_t>(x0,x,lo,hi);
1982 return vector_bsearch_dec<vec_t,data_t>(x0,x,lo,hi);
1992 template<
class vec_t,
class data_t>
1996 if (x[lo]<x[hi-1]) {
1997 return vector_bsearch_inc<vec_t,data_t>(x0,x,lo,hi);
1999 return vector_bsearch_dec<vec_t,data_t>(x0,x,lo,hi);
2018 template<
class vec_t>
2023 if (data[0]==data[1]) {
2025 }
else if (data[0]<data[1]) {
2035 for(
size_t i=0;i<n-1 && done==
false;i++) {
2036 if (data[i]!=data[i+1]) {
2050 if (data[start]<data[start+1])
return 1;
2057 if (data[start]>data[start+1]) inc=
false;
2060 for(
size_t i=start+1;i<n-1;i++) {
2062 if (data[i]>data[i+1])
return 0;
2068 for(
size_t i=start+1;i<n-1;i++) {
2069 if (data[i]<data[i+1])
return 0;
2100 template<
class vec_t>
2108 if (data[0]==data[1]) {
2110 }
else if (data[0]>data[1]) {
2115 for(
size_t i=1;i<n-1;i++) {
2117 if (data[i]>=data[i+1])
return 0;
2123 for(
size_t i=1;i<n-1;i++) {
2124 if (data[i]<=data[i+1])
return 0;
2138 template<
class vec_t>
2151 template<
class vec_t>
2153 for(
size_t i=0;i<n;i++) {
2154 if (!std::isfinite(data[i]))
return false;
2179 template<
class vec_t,
class data_t>
2183 for(
size_t i=0;i<n;i++) {
2192 template<
class vec_t,
class rvec_t>
2195 size_t n=v_data.size();
2196 v_diffs.resize(n-1);
2197 for(
size_t i=0;i<n-1;i++) {
2198 v_diffs[i]=v_data[i+1]-v_data[i];
2209 template<
class vec_t,
class data_t> data_t
vector_sum(vec_t &data) {
2211 for(
size_t i=0;i<data.size();i++) {
2225 for(
size_t i=0;i<n;i++) {
2239 for(
size_t i=0;i<data.size();i++) {
2251 template<
class vec_t,
class data_t>
2259 }
else if (n == 1) {
2263 for (
size_t i = 0; i < n; i++) {
2264 const data_t xx = x[i];
2267 const data_t ax = fabs(xx);
2270 ssq = 1.0 + ssq * (scale / ax) * (scale / ax);
2273 ssq += (ax / scale) * (ax / scale);
2279 return scale * sqrt(ssq);
2285 template<
class vec_t,
class data_t> data_t
vector_norm(
const vec_t &x) {
2286 return vector_norm<vec_t,data_t>(x.size(),x);
2295 template<
class vec_t>
2303 }
else if (n == 1) {
2307 for (
size_t i = 0; i < n; i++) {
2308 const double xx = x[i];
2311 const double ax = fabs(xx);
2314 ssq = 1.0 + ssq * (scale / ax) * (scale / ax);
2317 ssq += (ax / scale) * (ax / scale);
2323 return scale * sqrt(ssq);
2329 return vector_norm_double<vec_t>(x.size(),x);
2337 template<
class vec_t,
class data_t>
2339 for(
size_t i=0;i<N;i++) {
2347 template<
class vec_t,
class data_t>
2355 template<
class mat_t,
class data_t>
2357 for(
size_t i=0;i<M;i++) {
2358 for(
size_t j=0;j<N;j++) {
2367 template<
class mat_t,
class data_t>
2378 template<
class vec_t,
class vec2_t>
2380 if (src.size()==0) {
2383 if (iout>=src.size()) {
2386 dest.resize(src.size()-1);
2388 for(
size_t i=0;i<src.size();i++) {
2402 template<
class vec_t,
class vec2_t>
2404 size_t iout, vec2_t &dest) {
2414 for(
size_t i=0;i<sz;i++) {
2436 template<
class vec_t,
class data_t>
2439 data_t *tmp=
new data_t[n];
2440 for(
size_t i=0;i<n;i++) {
2441 tmp[i]=data[(i+k)%n];
2443 for(
size_t i=0;i<n;i++) {
2455 template<
class vec_t,
class data_t>
2459 for(
size_t i=0;i<n/2;i++) {
2461 data[n-1-i]=data[i];
2472 template<
class vec_t,
class data_t>
2475 size_t n=data.size();
2477 for(
size_t i=0;i<n/2;i++) {
2479 data[n-1-i]=data[i];
2490 template<
class vec_t>
2494 for(
size_t i=0;i<n/2;i++) {
2496 data[n-1-i]=data[i];
2509 size_t n=data.size();
2511 for(
size_t i=0;i<n/2;i++) {
2513 data[n-1-i]=data[i];
2527 data_t operator[](
size_t &i)
const {
2592 template<
class mat_t,
class mat_row_t>
2594 return mat_row_t(M,row);
2662 mat.resize(1,n_cols);
2675 mat.resize(1,n_cols);
2694 O2SCL_ERR(
"No matrix in matrix_row_gen_ctor::operator[].",
2703 O2SCL_ERR(
"No matrix in matrix_row_gen_ctor::operator[].",
2791 return m_.size1()-1;
2848 return m_.size2()-1;
2895 const double &
operator()(
size_t row,
size_t col)
const;
2899 size_t size1()
const;
2903 size_t size2()
const;
2916 const double &
operator()(
size_t row,
size_t col)
const;
2925 size_t size1()
const;
2929 size_t size2()
const;
2941 template<
class vec1_t,
class vec2_t=std::vector<vec1_t> >
2956 std::swap(t1.
vvp,t2.
vvp);
2974 if (
vvp==0)
return 0;
2981 if (
vvp==0)
return 0;
2982 if (
vvp->size()==0)
return 0;
2983 return (*
vvp)[0].size();
2992 "matrix_view_vec_vec::operator().",
2995 if (row>=
vvp->size()) {
2997 "matrix_view_vec_vec::operator().",
3000 if (col>=(*
vvp)[row].size()) {
3002 "matrix_view_vec_vec::operator().",
3005 return (*
vvp)[row][col];
3014 "matrix_view_vec_vec::operator().",
3017 if (row>=
vvp->size()) {
3019 "matrix_view_vec_vec::operator().",
3022 if (col>=(*
vvp)[row].size()) {
3023 std::cout << row <<
" " << col <<
" "
3024 << (*vvp)[row].size() << std::endl;
3026 "matrix_view_vec_vec::operator().",
3029 return (*
vvp)[row][col];
3050 template<
class mat_t,
class mat_column_t>
3052 return mat_column_t(M,column);
3144 template<
class vec_t>
3146 bool endline=
false) {
3150 if (endline) os << std::endl;
3154 for(
size_t i=0;i<n-1;i++) os << v[i] <<
" ";
3156 if (endline) os << std::endl;
3171 template<
class vec_t>
3172 void vector_out(std::ostream &os,
const vec_t &v,
bool endline=
false) {
3179 for(
size_t i=0;i<n-1;i++) os << v[i] <<
" ";
3181 if (endline) os << std::endl;
3187 template<
class vec_t,
class data_t>
3189 g.template vector<vec_t>(v);
3194 template<
class mat_t>
3196 for(
size_t i=0;i<M;i++) {
3197 for(
size_t j=0;j<N;j++) {
3198 if (i==j) m(i,j)=1.0;
3206 template<
class mat_t>
3221 (dat_t *v,
size_t start,
size_t last) {
3231 (
const dat_t *v,
size_t start,
size_t last) {
3247 size_t start,
size_t last) {
3250 (v,boost::numeric::ublas::range(start,last));
3265 size_t start,
size_t last) {
3268 (v,boost::numeric::ublas::range(start,last));
3284 size_t start,
size_t last) {
3287 (v,boost::numeric::ublas::range(start,last));
3300 template<
class dat_t>
3307 size_t start,
size_t last) {
3311 (v,boost::numeric::ublas::range(start,last));
3324 template<
class dat_t>
3331 size_t start,
size_t last) {
3335 (v,boost::numeric::ublas::range(start,last));
3348 template<
class dat_t>
3355 size_t start,
size_t last) {
3359 (v,boost::numeric::ublas::range(start,last));
3372 template<
class dat_t>
3379 size_t start,
size_t last) {
3383 (v,boost::numeric::ublas::range(start,last));
3411 #if !O2SCL_NO_RANGE_CHECK
3413 O2SCL_ERR2(
"End before beginning in vector_range_gen::",
3414 "vector_range_gen(vec_t,size_t,size_t)",
3422 size_t last) :
v_(v2.
v_),
3424 #if !O2SCL_NO_RANGE_CHECK
3426 O2SCL_ERR2(
"End before beginning in vector_range_gen::",
3427 "vector_range_gen(vector_range_gen,size_t,size_t)",
3430 if (last>v2.
last_) {
3431 O2SCL_ERR2(
"End beyond end of previous vector in vector_range_gen::",
3432 "vector_range_gen(vector_range_gen,size_t,size_t)",
3445 #if !O2SCL_NO_RANGE_CHECK
3447 O2SCL_ERR(
"Index out of range in vector_range_gen::operator[].",
3456 #if !O2SCL_NO_RANGE_CHECK
3468 template<
class vec_t>
class const_vector_range_gen {
3486 #if !O2SCL_NO_RANGE_CHECK
3488 O2SCL_ERR2(
"End before beginning in vector_range_gen::",
3489 "vector_range_gen(vec_t,size_t,size_t)",
3497 size_t last) :
v_(v2.
v_),
3499 #if !O2SCL_NO_RANGE_CHECK
3501 O2SCL_ERR2(
"End before beginning in vector_range_gen::",
3502 "vector_range_gen(vector_range_gen,size_t,size_t)",
3505 if (last>v2.last_) {
3506 O2SCL_ERR2(
"End beyond end of previous vector in vector_range_gen::",
3507 "vector_range_gen(vector_range_gen,size_t,size_t)",
3515 size_t last) :
v_(v2.
v_),
3517 #if !O2SCL_NO_RANGE_CHECK
3519 O2SCL_ERR2(
"End before beginning in vector_range_gen::",
3520 "vector_range_gen(vector_range_gen,size_t,size_t)",
3523 if (last>v2.
last_) {
3524 O2SCL_ERR2(
"End beyond end of previous vector in vector_range_gen::",
3525 "vector_range_gen(vector_range_gen,size_t,size_t)",
3538 #if !O2SCL_NO_RANGE_CHECK
3551 template<
class data_t> vector_range_gen<std::vector<data_t> >
3559 template<
class data_t>
const const_vector_range_gen<std::vector<data_t> >
3568 template<
class data_t>
const const_vector_range_gen<std::vector<data_t> >
3577 template<
class vec_t> vector_range_gen<vec_t>
3585 template<
class vec_t>
const const_vector_range_gen<vec_t>
3587 size_t start,
size_t last) {
3594 template<
class vec_t>
const const_vector_range_gen<vec_t>
3596 size_t start,
size_t last) {
3603 template<
class vec_t>
const const_vector_range_gen<vec_t>
3605 size_t start,
size_t last) {
3620 template<
class dat_t> std::vector<dat_t>
3622 return std::vector<dat_t> (v.begin()+start,v.begin()+last);
3636 template<
class dat_t>
const std::vector<dat_t>
3639 return std::vector<dat_t> (v.begin()+start,v.begin()+last);
3646 #if defined (O2SCL_COND_FLAG) || defined (DOXYGEN)
3648 #if defined (O2SCL_ARMA) || defined (DOXYGEN)
3649 #include <armadillo>
3661 template<> arma::subview_row<double>
3662 matrix_row<arma::mat,arma::subview_row<double> >
3663 (arma::mat &M,
size_t row);
3666 template<> arma::subview_col<double>
3667 matrix_column<arma::mat,arma::subview_col<double> >
3668 (arma::mat &M,
size_t column);
3675 #if defined (O2SCL_EIGEN) || defined (DOXYGEN)
3676 #include <eigen3/Eigen/Dense>
3683 double matrix_max(
const Eigen::MatrixXd &data);
3686 double matrix_min(
const Eigen::MatrixXd &data);
3689 template<> Eigen::MatrixXd::RowXpr
3690 matrix_row<Eigen::MatrixXd,Eigen::MatrixXd::RowXpr>
3691 (Eigen::MatrixXd &M,
size_t row);
3694 template<> Eigen::MatrixXd::ColXpr
3695 matrix_column<Eigen::MatrixXd,Eigen::MatrixXd::ColXpr>
3696 (Eigen::MatrixXd &M,
size_t column);
3705 #include <o2scl/vector_special.h>
3752 template<
class T,
class F,
class A>
class matrix {
matrix_view_omit_row(mat_t &m, size_t row_omit)
Create.
const double & operator()(size_t i, size_t j) const
Return a const reference.
The default matrix type from uBlas.
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.
size_t last_
The end() iterator.
friend void swap(matrix_view_vec_vec &t1, matrix_view_vec_vec &t2)
Swap method.
matrix_row_gen_ctor(size_t n_cols=0)
Create a row object from row row of matrix m.
double & operator[](size_t i)
Return a reference to the ith column of the selected row.
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.
Matrix row object with a constructor and resize method.
mat_t & m_
A reference to the original matrix.
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.
vec_t & v_
A reference to the original vector.
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.
The default vector type from uBlas.
size_t column_
The selected column.
size_t size() const
Return the vector size.
size_t size() const
Get the size of the vector.
size_t size() const
Return the vector size.
size_t size1() const
Return the number of rows.
Construct a view of a matrix omtting a specified row.
size_t start_
The index offset.
void matrix_copy(mat_t &src, mat2_t &dest)
Simple matrix copy.
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.
@ exc_efailed
generic failure
const double & operator()(size_t row, size_t col) const
Return a reference to the element at row row and column col.
double & operator()(size_t i, size_t j)
Return a reference.
void matrix_transpose(mat_t &src, mat2_t &dest)
Simple transpose.
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
bool matrix_is_lower(mat_t &src)
Simple test that a matrix is lower triangular.
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.
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 & operator[](size_t i)
Return a reference ith element.
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.
double matrix_max(const arma::mat &data)
Armadillo version of matrix_max()
Generic object which represents a row of a const matrix.
mat_t & m_
A reference to the original matrix.
The main O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl names...
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.
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.
vector_index_vector_size(size_t n_)
Create an index vector with size n_.
A simple convenience wrapper for GSL vector objects.
data_t vector_min_value(size_t n, const vec_t &data)
Compute the minimum of the first n elements of a vector.
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.
A simple matrix view object.
mat_column_t matrix_column(mat_t &M, size_t column)
Construct a column of a matrix.
data_t vector_max_value(size_t n, const vec_t &data)
Compute the maximum of the first n elements of a vector.
Placeholder documentation of some related Boost objects.
const double & operator[](size_t i) const
Return a const reference to the ith row of the selected column.
void vector_copy(const vec_t &src, vec2_t &dest)
Simple vector copy.
const_matrix_row_gen(const mat_t &m, size_t row)
Create a row object from row row of matrix m.
const mat_t & m_
A reference to the original matrix.
double & operator()(size_t i, size_t j)
Return a reference.
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.
data_t vector_min_quad_loc(size_t n, const vec_t &x, const vec_t &y)
Location of vector minimum by quadratic fit.
double & operator()(size_t i, size_t j)
Return a reference to the ith column of the selected row.
matrix_column_gen(mat_t &m, size_t column)
Create a column object from column column of matrix m.
double matrix_min(const arma::mat &data)
Armadillo version of matrix_min()
size_t size1() const
Return the number of rows.
const_matrix_column_gen(const mat_t &m, size_t column)
Create a column object from column column of matrix m.
size_t size2() const
Return the number of columns.
const_vector_range_gen(const vec_t &v, size_t start, size_t last)
Create an object starting with index start in vector v.
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.
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.
size_t size1() const
Return the number of rows.
matrix_row_gen_ctor(mat_t &m, size_t row)
Create a row object from row row of matrix m.
vector_range_gen(vec_t &v, size_t start, size_t last)
Create an object starting with index start in vector v.
View a o2scl::table object as a matrix.
matrix_view_vec_vec(vec2_t &vv)
Create a matrix view object from the specified table and list of rows.
vector_range_gen(const vector_range_gen &v2, size_t start, size_t last)
Create an object from a previously constructed range object.
size_t row_
The selected row.
mat_t * mp
A pointer to the matrix.
size_t column_
The selected column.
void matrix_make_upper(mat_t &src)
Make a matrix upper triangular by setting the lower triangular entries to zero.
void vector_sort_index(size_t n, const vec_t &data, vec_size_t &order)
Create a permutation which sorts the first n elements of a vector (in increasing order)
Construct a view of a matrix omitting one specified column.
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.
mat_row_t matrix_row(mat_t &M, size_t row)
Construct a row of a matrix.
size_t size2() const
Return the number of columns.
data_t matrix_sum(size_t m, size_t n, const mat_t &data)
Compute the sum of matrix elements.
const double & operator[](size_t i) const
Return a const reference ith element.
@ exc_einval
invalid argument supplied by user
void vector_grid(uniform_grid< data_t > g, vec_t &v)
Fill a vector with a specified grid.
matrix_view_omit_column(mat_t &m, size_t column_omit)
Create.
const double & operator[](size_t i) const
Return a const reference to the ith column of the selected row.
Experimental vector range object.
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...
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.
double matrix_min_value_double(const mat_t &data)
Compute the minimum of a matrix.
void sort_downheap(vec_t &data, size_t n, size_t k)
Provide a downheap() function for vector_sort()
Generic object which represents a row of a matrix.
bool matrix_is_upper(mat_t &src)
Simple test that a matrix is upper triangular.
vec2_t * vvp
Pointer to the table.
data_t vector_max_quad_loc(size_t n, const vec_t &x, const vec_t &y)
Location of vector maximum by quadratic fit.
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.
const double & operator[](size_t i) const
Return a const reference to the ith row of the selected column.
size_t row_
The selected row.
void resize(size_t n_)
Resize the index vector.
double & operator()(size_t row, size_t col)
Return a reference to the element at row row and column col.
void vector_sort(size_t n, vec_t &data)
Sort a vector (in increasing order)
Generic object which represents a column of a const matrix.
size_t size2() const
Return the number of columns.
data_t vector_max_quad(size_t n, const vec_t &data)
Maximum of vector by quadratic fit.
const double & operator()(size_t row, size_t col) const
Return a reference to the element at row row and column col.
data_t operator[](size_t &i) const
Obtain the element with index i.
matrix_view_transpose(mat_t &m)
Create a row object from row row of matrix m.
const dat_t * const_vector_range(const dat_t *v, size_t start, size_t last)
Vector range function for const pointers.
mat_t & m_
A reference to the original matrix.
void matrix_make_lower(mat_t &src)
Make a matrix lower triangular by setting the upper triangular entries to zero.
void vector_reverse_double(size_t n, vec_t &data)
Reverse the first n elements in a vector of double precision numbers.
data_t matrix_min_value(size_t m, size_t n, const mat_t &data)
Compute the minimum of a matrix.
const double & operator[](size_t i) const
Return a const reference to the ith column of the selected row.
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.
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 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()
mat_t & m_
A reference to the original matrix.
bool vector_is_finite(size_t n, vec_t &data)
Test if the first n elements of a vector are finite.
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.
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
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.
const double & operator()(size_t row, size_t col) const
Return a reference to the element at row row and column col.
const mat_t & m_
A reference to the original matrix.
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.
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.
std::string itos(int x)
Convert an integer to a string.
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.
Construct a view of the transpose of a matrix.
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.
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 resize(size_t n_cols=0)
Resize.
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.
mat_t & m_
A reference to the original matrix.
A simple convenience wrapper for GSL matrix objects.
void matrix_set_all(size_t M, size_t N, mat_t &src, data_t val)
Set the first (M,N) entries in a matrix to a particular value.
Index vector with a size method.
void vector_sort_double(size_t n, vec_t &data)
Sort a vector of doubles (in increasing order)
void vector_smallest_index(size_t n, const 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.
Experimental const vector range object.
dat_t * vector_range(dat_t *v, size_t start, size_t last)
Vector range function for pointers.
size_t row_
The selected row.
size_t size2() const
Return the number of columns.
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
Generic object which represents a column of a matrix.
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 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_sum(size_t n, vec_t &data)
Compute the sum of the first n elements of a vector.
data_t vector_min_quad(size_t n, const vec_t &data)
Minimum of vector by quadratic fit.
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.
void vector_reverse(size_t n, vec_t &data)
Reverse the first n elements of a vector.
size_t size1() const
Return the number of rows.
A simple matrix view object.
const double & operator()(size_t i, size_t j) const
Return a const reference to the ith column of the selected row.
size_t start_
The index offset.
double matrix_max_value_double(const mat_t &data)
Compute the maximum of a matrix.
size_t size1() const
Return the number of rows.
void vector_diffs(const vec_t &v_data, rvec_t &v_diffs)
Create a new vector containing the differences between adjacent entries.
size_t size() const
Return size.
size_t size1() const
Return the number of rows.
void matrix_swap(size_t M, size_t N, mat_t &v1, mat2_t &v2)
Swap of the first elements in two matrices.
const vec_t & v_
A reference to the original vector.
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.
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.
size_t size2() const
Return the number of columns.
const double & operator[](size_t i) const
Return a const reference to the ith column of the selected row.
double & operator[](size_t i)
Return a reference to the ith column of the selected row.
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_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.
const double & operator()(size_t i, size_t j) const
Return a const reference.
double & operator[](size_t i)
Return a reference to the ith row of the selected column.
matrix_row_gen(mat_t &m, size_t row)
Create a row object from row row of matrix m.
std::string szttos(size_t x)
Convert a size_t to a string.
mat_t mat
A matrix to point to.
const double & operator[](size_t i) const
Return a const reference ith element.
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.
size_t size2() const
Return the number of columns.
void vector_swap(size_t N, vec_t &v1, vec2_t &v2)
Swap the first N elements of two vectors.
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.
size_t last_
The end() iterator.
Documentation generated with Doxygen. Provided under the
GNU Free Documentation License (see License Information).