google / TensorNetwork

A library for easy and efficient manipulation of tensor networks.
Apache License 2.0
1.81k stars 356 forks source link

zeros_like, ones_like, empty_like for BlockSparseTensor and symmetric backend #685

Open mganahl opened 4 years ago

mganahl commented 4 years ago

This is useful for initializing new BlockSparseTensor objects from existing ones with compatible shapes, charges and flows

mganahl commented 4 years ago

@gevenbly, as you pointed out earlier, this could potentially free us of the Index class. One case I'd like to cover though is where one needs to initialize a BlockSparseTensor from the shapes of two or more other BlockSparseTensors. This happens for example in sparse solver initialization for DMRG if one has no access to the current MPS tensor, or when initializing eigs for finding reduced density matrices of infinite systems. For the latter we could extend e.g ones_like to

res = ones_like([A,B,C],axes = [(0,),(2,4),(3,1)]) 
#res.shape has shape (A[0], B[2], B[4], C[3], C[1])
#but is not contractable with A,B,C due to mismatching flow signature
#calling e.g. res.flip_flows() would flip all flows and make it contractable

or

res = tn.contractable_ones(tensors = [A,B],axes=[(0,),(2,1)])
#res.shape=(A.shape[0], B.shape[2], B.shape[1]); 
#the flow-signature of res is such that it is contractable with A, B. 
tn.ncon([res,A,B],[[1,2,3],[1,-1],[-2,3,2]])
gevenbly commented 4 years ago

While we are here, I would also suggest it could be very useful to include a rand_like or random_like function (even though I don't believe this exists for numpy?). I am not against the extensions to one_like that you suggest above, although I do not think that it is necessary. To initialize a blocksparse tensor, one is required to provide the charge and flow information; if one wants to build a new tensor matching the dimensions of some existing tensors A and B the obvious way to do so is to use the A.charges and A.flows. So the usual way of doing things is numpy is:

# numpy arrays only require a shape to be initialized
A = np.random.rand(A_shape)
B = np.random.rand(B_shape)
# can initialize a random tensor from shape of A and B
C = np.random.rand(A.shape[0],B.shape[2])

then the following code seems entirely consistent to me:

# blocksparse arrays require charge and flow information to be initialized
A = BS.random(charges=A_charges, flows=A_flows)
B = BS.random(charges=B_charges, flows=B_flows)
# can initialize a random tensor from charges and flows of A and B
C = BS.random(charges=[A.charges[0],B.charges[2]], flows=[A.flows[0],B.flows[2]])

Here one could use np.invert on the flows or .conj() the final tensor if it is desired that the flows should be reversed. Why does there need to be a special syntax to build tensors to match the index dimensions of other tensors? The only advantage I see is that the syntax may be more compact (which comes at the expense of introducing some new behavior which is not present in numpy).

mganahl commented 4 years ago

yes, rand_like or random_like is a good addition.

The aim is to provide an interface for initializing tensors that is the same for all backends, i.e. to avoid having to explicitly access charge or flow information of BlockSparseTensor. In some situations, one could want to initialize new tensors deep within a tensor-network algorithm, not only at initialization time. This API is aimed at covering such cases.

gevenbly commented 4 years ago

Sure, but it is already the case that blocksparse backend necessarily behaves differently to other backends when it comes to tensor initialization (such that there could inevitably be many places deep within code where tensor initializations need to be changed). The extension to ones_like you suggest could mitigate the need for some of the changes between code written for the numpy and blocksparse backends (provided the user is willing to write the numpy part in a non-standard way). So I think the ones_like you suggest could be a very useful addition to the toolbox and there is no any harm in having it (my argument above was simply that it is not strictly necessary, since there already exists a straight-forward way of achieving the same goal).

On a related note, I think it would also be useful to have a tensor_expand type function. The idea is that this function expands the dimensions of an existing tensor (padding the tensor with zeros) in order to match the indices of some other existing tensors or indices, rather than creating a new tensor from scratch. So it could look something like:

C_new = tn.tensor_expand(C, tensors = [A,B],axes=[(0,),(2,1)])

This was something that I used extensively in my previous libraries (both for standard and for blocksparse algorithms). For instance, in a TEBD algorithm, I would only define the left/right density matrices rho_L and rho_R (given from contracting the infinite MPS from the left/right boundary) at the start of my code, but then just expand them to match the dimensions of the MPS tensors as necessary. From the computational perspective, this was more efficient as you always reuse the previous boundary contraction information, even if the dimensions have changed from the previous iteration. (My own philosophy is that your algorithm is probably inefficient if you need to define random tensors during each iteration. You should always be reusing the information from the previous iteration, unless you need a truly random process).