23 #ifndef O2SCL_TENSOR_GRID_H 24 #define O2SCL_TENSOR_GRID_H 36 #include <gsl/gsl_matrix.h> 37 #include <gsl/gsl_ieee_utils.h> 39 #include <o2scl/err_hnd.h> 40 #include <o2scl/interp.h> 41 #include <o2scl/tensor.h> 42 #include <o2scl/table3d.h> 52 template<
class vec_t,
class vec_
size_t>
55 template<
class vec_t,
class vec_
size_t>
60 #ifndef DOXYGEN_NO_O2NS 64 typedef boost::numeric::ublas::range ub_range;
148 template<
class vec_t=std::vector<
double>,
149 class vec_
size_t=std::vector<
size_t> >
class tensor_grid :
150 public tensor<double,vec_t,vec_size_t> {
154 #ifndef DOXYGEN_INTERNAL 158 #ifdef O2SCL_NEVER_DEFINED 189 template<
class size_vec_t>
191 tensor<double,vec_t,vec_size_t>(rank,dim) {
195 for(
size_t i=0;i<this->rk;i++) {
197 O2SCL_ERR((((std::string)
"Requested zero size with non-zero ")+
198 "rank for index "+
szttos(i)+
" in tensor_grid::"+
199 "tensor_grid(size_t,size_vec_t)").c_str(),
209 tensor<double,vec_t,vec_size_t>() {
213 for(
size_t j=0;j<this->rk;j++) {
214 this->size.push_back(ugs[j].get_npoints());
215 tot*=ugs[j].get_npoints();
217 this->data.resize(tot);
235 if (this->rk>0 && grid_set) {
237 for(
size_t i=0;i<this->rk;i++) {
241 if (tot2!=grid.size()) {
242 O2SCL_ERR2(
"Value grid_set is true but grid vector size ",
243 "is wrong in tensor_grid::is_valid().",
248 if (!grid_set && grid.size()>0) {
249 O2SCL_ERR2(
"Value grid_set is false but grid vector size ",
250 "is not zero in tensor_grid::is_valid().",
291 template<
class vec2_t>
292 void set_val(
const vec2_t &grdp,
double val) {
295 vec_size_t index(this->rk);
296 for(
size_t i=0;i<this->rk;i++) {
297 index[i]=lookup_grid(i,grdp[i]);
302 for(
size_t i=1;i<this->rk;i++) {
320 template<
class vec2_t,
class vec3_t>
321 void set_val(
const vec2_t &grdp,
double val, vec3_t &closest) {
324 vec_size_t index(this->rk);
325 for(
size_t i=0;i<this->rk;i++) {
326 index[i]=lookup_grid_val(i,grdp[i],closest[i]);
331 for(
size_t i=1;i<this->rk;i++) {
346 template<
class vec2_t>
double get_val(
const vec2_t &gridp) {
349 vec_size_t index(this->rk);
350 for(
size_t i=0;i<this->rk;i++) {
351 index[i]=lookup_grid(i,gridp[i]);
356 for(
size_t i=1;i<this->rk;i++) {
362 return this->data[ix];
371 template<
class vec2_t,
class vec3_t>
372 double get_val(
const vec2_t &gridp, vec3_t &closest) {
375 vec_size_t index(this->rk);
376 for(
size_t i=0;i<this->rk;i++) {
377 index[i]=lookup_grid_val(i,gridp[i],closest[i]);
382 for(
size_t i=1;i<this->rk;i++) {
388 return this->data[ix];
409 template<
class size_vec2_t>
410 void resize(
size_t rank,
const size_vec2_t &dim) {
413 for(
size_t i=0;i<rank;i++) {
415 O2SCL_ERR((((std::string)
"Requested zero size with non-zero ")+
416 "rank for index "+
szttos(i)+
" in tensor_grid::"+
423 this->size.resize(this->rk);
430 this->data.resize(0);
434 for(
size_t i=0;i<this->rk;i++) {
435 this->size[i]=dim[i];
438 this->data.resize(tot);
471 template<
class vec2_t>
474 O2SCL_ERR2(
"Tried to set grid for empty tensor in ",
475 "tensor_grid::set_grid_packed().",
exc_einval);
478 for(
size_t i=0;i<this->rk;i++) ngrid+=this->size[i];
480 for(
size_t i=0;i<ngrid;i++) {
489 template<
class vec_vec_t>
492 O2SCL_ERR2(
"Tried to set grid for empty tensor in ",
496 for(
size_t i=0;i<this->rk;i++) ngrid+=this->size[i];
499 for(
size_t i=0;i<this->rk;i++) {
500 for(
size_t j=0;j<this->size[i];j++) {
501 grid[k]=grid_vecs[i][j];
513 for(
size_t i=0;i<this->rk;i++) ngrid+=this->size[i];
516 for(
size_t i=0;i<this->rk;i++) {
517 for(
size_t j=0;j<this->size[i];j++) {
528 template<
class vec2_t>
530 if (grid_set==
false) {
535 O2SCL_ERR2(
"Tried to set grid for empty tensor in ",
539 for(
size_t i=0;i<this->rk;i++) {
540 for(
size_t j=0;j<this->size[i];j++) {
552 template<
class vec2_t>
554 if (grid_set==
false) {
556 "tensor_grid::set_grid_i_func().",
exc_einval);
559 O2SCL_ERR2(
"Tried to set grid for empty tensor in ",
560 "tensor_grid::set_grid_i_func().",
exc_einval);
564 std::map<std::string,double> vars;
565 calc.
compile(func.c_str(),&vars);
568 for(
size_t i=0;i<this->rk;i++) {
569 for(
size_t j=0;j<this->size[i];j++) {
571 vars[
"i"]=((double)j);
572 grid[k]=calc.
eval(&vars);
587 O2SCL_ERR2(
"Tried to set grid for empty tensor in ",
591 for(
size_t i=0;i<this->rk;i++) ngrid+=ugs[i].get_npoints();
594 for(
size_t i=0;i<this->rk;i++) {
595 for(
size_t j=0;j<ugs[i].get_npoints();j++) {
609 template<
class rvec_t>
void copy_grid(
size_t i, rvec_t &v) {
610 v.resize(this->size[i]);
612 for(
size_t k=0;k<i;k++) istart+=this->size[k];
613 for(
size_t k=0;k<this->size[i];k++) {
622 O2SCL_ERR(
"Grid not set in tensor_grid::get_grid().",
627 " greater than or equal to rank, "+
szttos(this->rk)+
628 ", in tensor_grid::get_grid().").c_str(),
632 for(
size_t k=0;k<i;k++) istart+=this->size[k];
633 return grid[istart+j];
639 O2SCL_ERR(
"Grid not set in tensor_grid::get_grid().",
644 " greater than or equal to rank, "+
szttos(this->rk)+
645 ", in tensor_grid::get_grid().").c_str(),
649 for(
size_t k=0;k<i;k++) istart+=this->size[k];
663 " greater than or equal to rank, "+
szttos(this->rk)+
664 ", in tensor_grid::lookup_grid_val().").c_str(),
667 if (grid_set==
false) {
668 O2SCL_ERR(
"Grid not set in tensor_grid::lookup_grid_val().",
677 for(
size_t j=0;j<i;j++) istart+=this->size[j];
679 double min=fabs(grid[istart]-temp);
681 for(
size_t j=istart;j<istart+this->size[i];j++) {
682 if (fabs(grid[j]-temp)<min) {
684 min=fabs(grid[j]-temp);
694 return lookup_grid_val(i,val,val2);
707 template<
class vec2_t,
class size_vec2_t>
709 for(
size_t k=0;k<this->rk;k++) {
710 indices[k]=lookup_grid(k,vals[k]);
725 O2SCL_ERR(
"Grid not set in tensor_grid::lookup_grid_packed_val().",
731 ", in tensor_grid::lookup_grid_packed_val().").c_str(),
735 for(
size_t j=0;j<i;j++) istart+=this->size[j];
737 double min=fabs(grid[istart]-val);
739 for(
size_t j=istart;j<istart+this->size[i];j++) {
740 if (fabs(grid[j]-val)<min) {
742 min=fabs(grid[j]-val);
754 return lookup_grid_packed_val(i,val,val2);
763 template<
class size_vec2_t,
class vec2_t>
766 if (this->rk<1+ifix.size()) {
768 "tensor_grid::copy_slice_interp().",
771 if (ifix.size()!=vals.size()) {
772 O2SCL_ERR2(
"Mismatch between indices and values in ",
773 "tensor_grid::copy_slice_interp().",
778 size_t rank_new=this->rk-ifix.size();
781 std::vector<size_t> sz_new;
782 std::vector<std::vector<double> > grid_new;
783 for(
size_t i=0;i<this->rk;i++) {
785 for(
size_t j=0;j<ifix.size();j++) {
786 if (ifix[j]==i) found=
true;
789 sz_new.push_back(this->get_size(i));
790 std::vector<double> grid_temp;
791 for(
size_t j=0;j<this->get_size(i);j++) {
792 grid_temp.push_back(this->get_grid(i,j));
794 grid_new.push_back(grid_temp);
803 std::vector<size_t> ix_new(rank_new);
804 std::vector<double> point_old(this->rk);
813 for(
size_t j=0;j<this->rk;j++) {
815 for(
size_t k=0;k<ifix.size();k++) {
816 if (ifix[k]==j) ix_found=k;
819 point_old[j]=this->get_grid(j,ix_new[j]);
821 point_old[j]=vals[ix_found];
839 void convert_table3d_sum
840 (
size_t ix_x,
size_t ix_y,
table3d &tab, std::string x_name=
"x",
841 std::string y_name=
"y", std::string slice_name=
"z") {
847 if (nx==0 && ny==0) {
849 if (x_name.length()==0) x_name=
"x";
850 if (y_name.length()==0) y_name=
"y";
854 std::vector<double> grid_x, grid_y;
855 copy_grid(ix_x,grid_x);
856 copy_grid(ix_y,grid_y);
857 tab.
set_xy(
"x",grid_x.size(),grid_x,
858 "y",grid_y.size(),grid_y);
864 if (nx!=this->size[ix_x] || ny!=this->size[ix_y]) {
866 "tensor_grid::convert_table3d_sum().",
exc_einval);
871 std::vector<size_t> ix;
872 for(
size_t i=0;i<this->total_size();i++) {
873 this->unpack_index(i,ix);
874 tab.
set(ix[ix_x],ix[ix_y],slice_name,
875 tab.
get(ix[ix_x],ix[ix_y],slice_name)+
906 template<
class size_vec2_t>
908 table3d &tab, std::string slice_name=
"z") {
910 if (ix_x>=this->rk || ix_y>=this->rk || ix_x==ix_y) {
911 O2SCL_ERR2(
"Either indices greater than rank or x and y ind",
912 "ices equal in tensor_grid::copy_table3d_align().",
921 if (nx!=this->size[ix_x] || ny!=this->size[ix_y]) {
923 "tensor_grid::copy_table3d_align().",
exc_einval);
931 for(
size_t i=0;i<nx;i++) {
932 for(
size_t j=0;j<ny;j++) {
935 double val=this->
get(index);
936 tab.
set(i,j,slice_name,val);
946 template<
class size_vec2_t>
947 void copy_table3d_align_setxy
948 (
size_t ix_x,
size_t ix_y, size_vec2_t &index,
949 table3d &tab, std::string x_name=
"x", std::string y_name=
"y",
950 std::string slice_name=
"z") {
956 if (nx==0 && ny==0) {
958 if (x_name.length()==0) x_name=
"x";
959 if (y_name.length()==0) y_name=
"y";
963 std::vector<double> grid_x, grid_y;
964 copy_grid(ix_x,grid_x);
965 copy_grid(ix_y,grid_y);
966 tab.
set_xy(
"x",grid_x.size(),grid_x,
967 "y",grid_y.size(),grid_y);
972 copy_table3d_align(ix_x,ix_y,index,tab,slice_name);
997 template<
class size_vec2_t>
999 table3d &tab, std::string slice_name=
"z") {
1001 if (ix_x>=this->rk || ix_y>=this->rk || ix_x==ix_y) {
1002 O2SCL_ERR2(
"Either indices greater than rank or x and y ",
1003 "indices equal in tensor_grid::copy_table3d_interp().",
1011 if (nx==0 && ny==0) {
1013 return copy_table3d_align(ix_x,ix_y,index,tab,slice_name);
1017 std::vector<double> vals(this->rk);
1018 for(
size_t i=0;i<this->rk;i++) {
1019 if (i!=ix_x && i!=ix_y) vals[i]=this->get_grid(i,index[i]);
1027 for(
size_t i=0;i<nx;i++) {
1028 for(
size_t j=0;j<ny;j++) {
1040 template<
class vec2_t>
1043 std::string slice_name=
"z",
1046 if (ix_x>=this->rk || ix_y>=this->rk || ix_x==ix_y) {
1047 O2SCL_ERR2(
"Either indices greater than rank or x and y ",
1048 "indices equal in tensor_grid::copy_table3d_interp().",
1051 if (values.size()!=this->rk) {
1052 O2SCL_ERR2(
"Values array not equal to rank ",
1053 "in tensor_grid::copy_table3d_interp_values().",
1071 for(
size_t i=0;i<nx;i++) {
1072 for(
size_t j=0;j<ny;j++) {
1077 std::cout <<
"At location values: ";
1078 for(
size_t k=0;k<values.size();k++) {
1079 std::cout << values[k] <<
" ";
1081 std::cout <<
"Interpolated to get: " 1082 << i <<
" " << j <<
" " << slice_name <<
" " 1098 template<
class vec2_t>
1099 void copy_table3d_interp_values_setxy
1100 (
size_t ix_x,
size_t ix_y, vec2_t &values,
table3d &tab,
1101 std::string x_name=
"x", std::string y_name=
"y",
1102 std::string slice_name=
"z") {
1108 if (nx==0 && ny==0) {
1110 if (x_name.length()==0) x_name=
"x";
1111 if (y_name.length()==0) y_name=
"y";
1115 std::vector<double> grid_x, grid_y;
1116 copy_grid(ix_x,grid_x);
1117 copy_grid(ix_y,grid_y);
1118 tab.
set_xy(
"x",grid_x.size(),grid_x,
1119 "y",grid_y.size(),grid_y);
1124 copy_table3d_interp_values(ix_x,ix_y,values,tab,slice_name);
1174 template<
class range_t=ub_range,
1175 class data_range_t=ubvector_range,
1176 class index_range_t=ubvector_size_t_range>
1183 interp_t si(this->size[0],grid,this->data,itype);
1184 return si.eval(vals[0]);
1190 for(
size_t i=1;i<this->rk;i++) ss*=this->size[i];
1193 std::vector<vec_t> yvec(ss);
1194 std::vector<interp_t *> si(ss);
1195 for(
size_t i=0;i<ss;i++) {
1196 yvec[i].resize(this->size[0]);
1201 index_range_t size_new(this->size,ub_range(1,this->rk));
1202 tdat.
resize(this->rk-1,size_new);
1205 data_range_t grid_new(grid,ub_range(this->size[0],grid.size()));
1209 vec_size_t co(this->rk);
1210 for(
size_t i=0;i<this->rk;i++) co[i]=0;
1215 while(done==
false) {
1218 for(
size_t i=0;i<this->size[0];i++) {
1220 yvec[cnt][i]=this->
get(co);
1223 si[cnt]=
new interp_t(this->size[0],grid,yvec[cnt],itype);
1225 index_range_t co2(co,ub_range(1,this->rk));
1226 tdat.
set(co2,si[cnt]->eval(vals[0]));
1232 for(
int j=((
int)this->rk)-1;j>0;j--) {
1233 if (co[j]>=this->size[j]) {
1240 if (cnt==ss) done=
true;
1248 for(
size_t i=0;i<ss;i++) {
1278 std::vector<size_t> loc(this->rk);
1279 std::vector<double> gnew;
1280 for(
size_t i=0;i<this->rk;i++) {
1281 std::vector<double> grid_unpacked(this->size[i]);
1282 for(
size_t j=0;j<this->size[i];j++) {
1283 grid_unpacked[j]=grid[j+rgs];
1286 loc[i]=sv.
find(v[i]);
1287 gnew.push_back(grid_unpacked[loc[i]]);
1288 gnew.push_back(grid_unpacked[loc[i]+1]);
1294 std::vector<size_t> snew(this->rk);
1295 for(
size_t i=0;i<this->rk;i++) {
1303 std::vector<size_t> index_new(this->rk), index_old(this->rk);
1305 for(
size_t j=0;j<this->rk;j++) index_old[j]=index_new[j]+loc[j];
1306 tnew.
set(index_new,this->
get(index_old));
1324 template<
class vec2_
size_t>
1328 return this->data[0]+(this->data[1]-this->data[0])/
1329 (grid[1]-grid[0])*(v[0]-grid[0]);
1332 size_t last=this->rk-1;
1333 double frac=(v[last]-get_grid(last,0))/
1334 (get_grid(last,1)-get_grid(last,0));
1343 std::vector<size_t> index(this->rk);
1345 index[this->rk-1]=0;
1346 double val_lo=this->
get(index);
1347 index[this->rk-1]=1;
1348 double val_hi=this->
get(index);
1349 tnew.
set(index,val_lo+frac*(val_hi-val_lo));
1367 template<
class vec2_
size_t,
class vec2_t>
1372 std::vector<size_t> loc(this->rk);
1373 std::vector<double> gnew;
1374 for(
size_t i=0;i<this->size[0];i++) {
1375 gnew.push_back(grid[i]);
1379 for(
size_t i=1;i<this->rk;i++) {
1380 std::vector<double> grid_unpacked(this->size[i]);
1381 for(
size_t j=0;j<this->size[i];j++) {
1382 grid_unpacked[j]=grid[j+rgs];
1385 loc[i]=sv.
find(v[i]);
1386 gnew.push_back(grid_unpacked[loc[i]]);
1387 gnew.push_back(grid_unpacked[loc[i]+1]);
1393 std::vector<size_t> snew(this->rk);
1394 snew[0]=this->size[0];
1395 for(
size_t i=1;i<this->rk;i++) {
1403 std::vector<size_t> index_new(this->rk), index_old(this->rk);
1405 for(
size_t j=0;j<this->rk;j++) {
1406 index_old[j]=index_new[j]+loc[j];
1408 tnew.
set(index_new,this->
get(index_old));
1428 template<
class vec2_
size_t,
class vec2_t>
1432 size_t n=this->size[0];
1434 vec_size_t ix0(2), ix1(2);
1437 for(
size_t i=0;i<n;i++) {
1440 res[i]=this->
get(ix0)+(this->
get(ix1)-this->
get(ix0))/
1441 (grid[n+1]-grid[n])*(v[1]-grid[n]);
1446 size_t last=this->rk-1;
1447 double frac=(v[last]-get_grid(last,0))/
1448 (get_grid(last,1)-get_grid(last,0));
1457 std::vector<size_t> index(this->rk);
1459 index[this->rk-1]=0;
1460 double val_lo=this->
get(index);
1461 index[this->rk-1]=1;
1462 double val_hi=this->
get(index);
1463 tnew.
set(index,val_lo+frac*(val_hi-val_lo));
1482 template<
class vec2_
size_t,
class vec2_t>
1485 size_t n=this->size[ifree];
1489 std::vector<size_t> map;
1490 map.push_back(ifree);
1491 for(
size_t i=0;i<this->rk;i++) {
1499 vec_size_t loc(this->rk);
1501 for(
size_t i=0;i<this->rk;i++) {
1502 vec_t grid_unpacked(this->size[i]);
1503 for(
size_t j=0;j<this->size[i];j++) {
1504 grid_unpacked[j]=grid[j+rgs];
1508 loc[i]=sv.
find(v[i]);
1514 std::vector<double> gnew, vnew;
1515 for(
size_t new_ix=0;new_ix<this->rk;new_ix++) {
1516 for(
size_t old_ix=0;old_ix<this->rk;old_ix++) {
1517 if (map[new_ix]==old_ix) {
1518 vnew.push_back(v[old_ix]);
1519 if (old_ix==ifree) {
1520 for(
size_t j=0;j<this->size[old_ix];j++) {
1521 gnew.push_back(this->get_grid(old_ix,j));
1524 gnew.push_back(this->get_grid(old_ix,loc[old_ix]));
1525 gnew.push_back(this->get_grid(old_ix,loc[old_ix]+1));
1535 std::vector<size_t> snew;
1537 for(
size_t i=0;i<this->rk;i++) {
1549 std::vector<size_t> index_new(this->rk), index_old(this->rk);
1551 for(
size_t j=0;j<this->rk;j++) {
1552 index_old[map[j]]=index_new[j]+loc[map[j]];
1554 tnew.
set(index_new,this->
get(index_old));
1564 template<
class vecf_t,
class vecf_
size_t>
friend void o2scl_hdf::hdf_output
1578 template<
class vec_t=std::vector<
double>,
1592 this->data.resize(sz);
1593 this->grid_set=
false;
1596 #ifdef O2SCL_NEVER_DEFINED 1604 double &
get(
size_t ix1) {
1610 const double &
get(
size_t ix1)
const {
1616 void set(
size_t ix1,
double val) {
1622 template<
class range_t=ub_range,
class data_range_t=ubvector_range,
1623 class index_range_t=ubvector_size_t_range>
1626 <range_t,data_range_t,index_range_t>(&x);
1639 template<
class vec_t=std::vector<
double>,
1655 this->data.resize(tot);
1656 this->grid_set=
false;
1659 #ifdef O2SCL_NEVER_DEFINED 1667 double &
get(
size_t ix1,
size_t ix2) {
1668 size_t sz[2]={ix1,ix2};
1673 const double &
get(
size_t ix1,
size_t ix2)
const {
1674 size_t sz[2]={ix1,ix2};
1679 void set(
size_t ix1,
size_t ix2,
double val) {
1680 size_t sz[2]={ix1,ix2};
1686 template<
class range_t=ub_range,
class data_range_t=ubvector_range,
1687 class index_range_t=ubvector_size_t_range>
1689 double arr[2]={x,y};
1691 <range_t,data_range_t,index_range_t>(arr);
1696 double arr[2]={x,y};
1704 template<
class vec_t=std::vector<
double>,
1721 size_t tot=sz*sz2*sz3;
1722 this->data.resize(tot);
1723 this->grid_set=
false;
1726 #ifdef O2SCL_NEVER_DEFINED 1734 double &
get(
size_t ix1,
size_t ix2,
size_t ix3) {
1735 size_t sz[3]={ix1,ix2,ix3};
1740 const double &
get(
size_t ix1,
size_t ix2,
size_t ix3)
const {
1741 size_t sz[3]={ix1,ix2,ix3};
1746 void set(
size_t ix1,
size_t ix2,
size_t ix3,
double val) {
1747 size_t sz[3]={ix1,ix2, ix3};
1753 template<
class range_t=ub_range,
class data_range_t=ubvector_range,
1754 class index_range_t=ubvector_size_t_range>
1755 double interp(
double x,
double y,
double z) {
1756 double arr[3]={x,y,z};
1758 <range_t,data_range_t,index_range_t>(arr);
1763 double arr[3]={x,y,z};
1771 template<
class vec_t=std::vector<
double>,
1789 size_t tot=sz*sz2*sz3*sz4;
1790 this->data.resize(tot);
1791 this->grid_set=
false;
1794 #ifdef O2SCL_NEVER_DEFINED 1802 double &
get(
size_t ix1,
size_t ix2,
size_t ix3,
size_t ix4) {
1803 size_t sz[4]={ix1,ix2,ix3,ix4};
1808 const double &
get(
size_t ix1,
size_t ix2,
size_t ix3,
1810 size_t sz[4]={ix1,ix2,ix3,ix4};
1815 void set(
size_t ix1,
size_t ix2,
size_t ix3,
size_t ix4,
1817 size_t sz[4]={ix1,ix2,ix3,ix4};
1823 template<
class range_t=ub_range,
class data_range_t=ubvector_range,
1824 class index_range_t=ubvector_size_t_range>
1825 double interp(
double x,
double y,
double z,
double a) {
1826 double arr[4]={x,y,z,a};
1828 <range_t,data_range_t,index_range_t>(arr);
1833 double arr[4]={x,y,z,a};
1839 #ifndef DOXYGEN_NO_O2NS bool is_grid_set() const
Return true if the grid has been set.
double get_grid_x(size_t ix) const
Get x grid point at index ix.
void unpack_index(size_t ix, size_vec_t &index)
Unpack the single vector index into indices.
Tensor class with arbitrary dimensions with a grid.
tensor_grid()
Create an empty tensor with zero rank.
void interp_linear_power_two_vec0(vec2_size_t &v, vec2_t &res)
Perform linear interpolation assuming that the last n-1 indices can take only two values...
The main O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl names...
Rank 1 tensor with a grid.
size_t lookup_grid_packed(size_t i, double val)
Lookup internal packed grid index for point closest to val.
double interpolate(double *vals)
Interpolate values vals into the tensor, returning the result.
virtual ~tensor_grid()
Destructor.
void compile(const char *expr, std::map< std::string, double > *vars=0, bool debug=false, std::map< std::string, int > opPrec=opPrecedence)
Compile expression expr using variables specified in vars.
size_t lookup_grid_val(size_t i, const double &val, double &val2)
Lookup index for grid closest to val, returning the grid point.
size_t lookup_grid_packed_val(size_t i, double val, double &val2)
Lookup internal packed grid index for point closest to val and store closest value in val2...
size_t find(const double x0)
Search an increasing or decreasing vector for the interval containing x0
tensor_grid(size_t rank, const size_vec_t &dim)
Create a tensor of rank rank with sizes given in dim.
sanity check failed - shouldn't happen
double & get(size_t ix, size_t iy, std::string name)
Get element in slice name at location ix,iy
invalid argument supplied by user
void clear()
Clear the tensor of all data and free allocated memory.
void set_grid(std::vector< uniform_grid< double > > &ugs)
Set grid from a vector of uniform grid objects.
Rank 2 tensor with a grid.
double interp_linear(vec2_size_t &v)
Perform a linear interpolation of v into the function implied by the tensor and grid.
void is_valid() const
Check that the o2scl::tensor_grid object is valid.
void copy_table3d_align(size_t ix_x, size_t ix_y, size_vec2_t &index, table3d &tab, std::string slice_name="z")
Create a slice in a o2scl::table3d object with an aligned grid.
void resize(size_t rank, const size_vec2_t &dim)
Resize the tensor to rank rank with sizes given in dim.
void set_interp_type(size_t interp_type)
Set interpolation type for interpolate()
bool is_slice(std::string name, size_t &ix) const
Return true if slice is already present.
bool is_xy_set() const
True if the grid has been set.
void set_slice_all(std::string name, double val)
Set all of the values in slice name to val.
Rank 3 tensor with a grid.
void default_grid()
Use a default grid which just uses the index.
void copy_table3d_interp(size_t ix_x, size_t ix_y, size_vec2_t &index, table3d &tab, std::string slice_name="z")
Copy to a slice in a table3d object using interpolation.
size_t lookup_grid(size_t i, double val)
Lookup index for grid closest to val.
tensor_grid3(size_t sz, size_t sz2, size_t sz3)
Create a rank 3 tensor of size (sz,sz2,sz3)
size_t total_size() const
Returns the size of the tensor (the product of the sizes over every index)
double get_grid_y(size_t iy) const
Get y grid point at index iy.
double interp(double x)
Interpolate x and return the results.
vec_t & get_data()
Return a reference to the data (for HDF I/O)
tensor_grid(std::vector< uniform_grid< double > > &ugs)
Create a tensor with a grid defined by a set of o2scl::uniform_grid objects.
void get_size(size_t &nx, size_t &ny) const
Get the size of the slices.
double interp(double x, double y, double z)
Interpolate (x,y,z) and return the results.
double get_grid(size_t i, size_t j) const
Lookup jth value on the ith grid.
double interp(double x, double y, double z, double a)
Interpolate (x,y,z,a) and return the results.
double & get(const size_vec_t &index)
Get the element indexed by index.
double interp_linear(double x, double y)
Interpolate (x,y) and return the results.
dat_t * vector_range(dat_t *v, size_t start, size_t last)
Vector range function for pointers.
tensor_grid2(size_t sz, size_t sz2)
Create a rank 2 tensor of size (sz,sz2)
double interp(double x, double y)
Interpolate (x,y) and return the results.
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
The O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl namespace ...
tensor_grid4()
Create an empty tensor.
double interp_linear_power_two(vec2_size_t &v)
Perform linear interpolation assuming that all indices can take only two values.
Rank 4 tensor with a grid.
void is_valid() const
Check that the o2scl::tensor object is valid.
double eval(std::map< std::string, double > *vars=0)
Evalate the previously compiled expression using variables specified in vars.
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
void interp_linear_vec(vec2_size_t &v, size_t ifree, vec2_t &res)
Perform a linear interpolation of v into the tensor leaving one index free resulting in a vector...
void set_grid_packed(const vec2_t &grid_vec)
Set the grid.
void interp_linear_vec0(vec2_size_t &v, vec2_t &res)
Perform a linear interpolation of v[1] to v[n-1] resulting in a vector.
double get_val(const vec2_t &gridp, vec3_t &closest)
Get the element closest to grid point gridp, store grid values in closest and return value...
void set_grid(const vec_vec_t &grid_vecs)
Set grid from a vector of vectors of grid points.
void copy_grid(size_t i, rvec_t &v)
Copy grid for index i to vector v.
void set_grid(size_t i, size_t j, double val)
Set the jth value on the ith grid.
void copy_table3d_interp_values(size_t ix_x, size_t ix_y, vec2_t &values, table3d &tab, std::string slice_name="z", int verbose=0)
Copy to a slice in a table3d object using interpolation.
Evaluate a mathematical expression in a string.
double interp_linear(double x, double y, double z, double a)
Interpolate (x,y,z,a) and return the results.
void clear()
Clear the tensor of all data and free allocated memory.
double get_val(const vec2_t &gridp)
Get the element closest to grid point gridp.
void set_val(const vec2_t &grdp, double val, vec3_t &closest)
Set the element closest to grid point grdp to value val.
tensor_grid copy_slice_interp(size_vec2_t &ifix, vec2_t &vals)
Copy an abitrary slice by fixing 1 or more indices and use interpolation to return a new tensor_grid ...
void set(size_t ix, size_t iy, std::string name, double val)
Set element in slice name at location ix,iy to value val.
void set_grid_i_vec(size_t ix, const vec2_t &grid_vec)
Set grid for one index from a vector.
Searching class for monotonic data with caching.
size_t itype
Interpolation type.
A data structure containing one or more slices of two-dimensional data points defined on a grid...
Store data in an O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$sc...
void new_slice(std::string name)
Add a new slice.
bool grid_set
If true, the grid has been set by the user.
tensor_grid2()
Create an empty tensor.
void set(const size_vec_t &index, double val)
Set the element indexed by index to value val.
vec_t grid
A rank-sized set of arrays for the grid points.
tensor_grid4(size_t sz, size_t sz2, size_t sz3, size_t sz4)
Create a rank 4 tensor of size (sz,sz2,sz3,sz4)
void set_grid_i_func(size_t ix, std::string func)
Set grid for one index from a function.
void lookup_grid_vec(const vec2_t &vals, size_vec2_t &indices) const
Lookup indices for grid closest point to vals.
vec_size_t size
A rank-sized vector of the sizes of each dimension.
void set_val(const vec2_t &grdp, double val)
Set the element closest to grid point grdp to value val.
void set_xy(std::string x_name, size_t nx, const vec_t &x, std::string y_name, size_t ny, const vec2_t &y)
Initialize the x-y grid.
tensor_grid3()
Create an empty tensor.
Linear interpolation (GSL)
std::string szttos(size_t x)
Convert a size_t to a string.
void hdf_input(hdf_file &hf, o2scl::table< vec_t > &t, std::string name)
Input a o2scl::table object from a hdf_file.
tensor_grid1()
Create an empty tensor.
double interp_linear(double x)
Interpolate x and return the results.
bool is_size_set() const
True if the size of the table has been set.
tensor_grid1(size_t sz)
Create a rank 2 tensor of size (sz,sz2,sz3)
Tensor class with arbitrary dimensions.
double interp_linear(double x, double y, double z)
Interpolate (x,y,z) and return the results.