Is there a way to get\set slices of a tensor?

+1 vote

For example, here's what I would hope to achieve with the following code:

Index i = Index("i", 2);
Index j = Index("j", 2);
Index k = Index("k", 2);
ITensor A = ITensor(i, j, k);
ITensor B = ITensor(i, j);
A.fill(0.0);
B.fill(1.0);
A.set(k(2), B); // the second k-slice of A is now filled with 1.0's from B


This should result in a tensor A which has two slices, with 0.0 and 1.0 for all i and j, for different k.

Looking at the docs and source, it seems that we can only set\get individual elements of a tensor, which requires a complete list of IndexVal to be passed. (Thankfully, this can now be done with a vector of IndexVals!)

The functionality here would be the same, except we're passing a partial list of IndexVal and setting not a Real\Cplx, but a Tensor.

Does this feature exist already, or is there a way to achieve it? If not, I think we should definitely include this convenient feature, as it is present and easy to use everywhere else.

selected by

The most straightforward way I can think to do it right now is the following:

auto i = Index("i", 2);
auto j = Index("j", 2);
auto k = Index("k", 2);
auto A = ITensor(i, j, k);
auto B = ITensor(i, j);
A.fill(0.0);
B.fill(1.0);
for(int ii = 1; ii <= i.m(); ii++)
for(int jj = 1; jj <= j.m(); jj++)
A.set(i(ii),j(jj),k(2), B.real(i(ii),j(jj)));


We would definitely like to support a nicer interface for this kind of slicing. However, it requires some thinking to figure out the right interface for slicing into a QN conserving/block sparse ITensor. It may be a little while before we get to it, but it is on our todo list!

commented by (240 points)
Hmm... Yeah, but that's exactly what I'm trying to avoid here: copying one by one. It will be cumbersome to use, slow to execute, and complicated to implement for non-fixed ranks. I can whip a function to do it all automatically, of course, so I'll do that for now.

It's unfortunate that this is not on the horizon yet, but it's understandable.

If I had an easy way to iterate over the non-zero elements (like PrintData does), this would be a breeze! Oh well, I'll see what I can do then.

Thanks!
commented by (10.7k points)
Looping like that shouldn't be too slow (there is no "copying", just setting elements). But yes, it has the downside that it is harder to make generic for ITensors with an arbitrary number of indices. We definitely agree we should have a better operation for doing this, and it is something that we have discussed as an important feature to add.

One thing I can point out that could make this operation slightly more convenient is that you can also index into an ITensor with integers, not just IndexVals, as long as you have "ordered" the ITensor data correctly:

A = order(A.i,j,k);    // Make sure the data of A is ordered according to Indices i,j,k
B = order(B,i,j);       // Make sure the data of B is ordered according to Indices i,j
for(int ii = 1; ii <= i.m(); ii++)
for(int jj = 1; jj <= j.m(); jj++)
A.set(ii,jj,2, B.real(ii,jj));

One of the main design choices is how to deal with more general types of slicing, for example if one of the dimensions you were slicing into was a range of values, not just a fixed value. For example, we may have to introduce an "IndexRange" object that generalizes an IndexVal, for example something like i(2,4) would represent the Index i over the range of values 2,3,4. Then, the question is, with these types of objects is there a nice way to deal with block-sparse tensors, where the indexing is not as straightforward (i.e. how do you in general make sure two IQTensor slices are compatible with each other).

A more automated way to iterate over all of the non-zero IndexVals of an ITensor/IQTensor is also good suggestion, but again requires some thinking for how to do it in general for sparse tensors. In the longer term, we are working on porting ITensor to Julia, which has a lot of sophisticated Array indexing which we hope to leverage to make some of these operations easier to implement for ITensor.
commented by (10.7k points)
For other people referring to this question: you can now use the "iterInds" function discussed in this answer to help with this operation: http://www.itensor.org/support/1361/iterating-over-non-zero-elements-of-a-sparse-itensor