0 votes
asked by (190 points)


for many different reasons, i.e. 1-site TDVP algorithm, it is useful to QR decompose the MPS tensor A=QR of at position b. Since QR is not implemented in ITensor one may still use SVD, namely decompose psi.A(b) into USV and identify Q=U and SV= R. However most of the examples in ITensor are for two-site svd decomposition. I tried to implement it for one site as follow

  auto ir1 = commonIndex(psi.A(b-1),psi.A(b),Link);
  auto ir2 = findtype(psi.A(b) , Site);
  ITensor U(ir1,ir2),V,S;

However when giving


a segmentation fault is produced. The reason is that the r=3 tensor U has not its indices ordered as in psi.A(b). U has index (Link,Site,Link) while psi.A(b) has (Site, Link,Link).
Do you know how to quickly solve the problem? Or maybe there is a better way to do the decomposition. Thanks a lot!

1 Answer

0 votes
answered by (3.6k points)

Hi Jacopo,

That sounds like a good approach for doing what you want to do. In the code above, is phi a copy of psi.A(b)?

I am not sure why it would be seg faulting. The reason would not be that U and psi.A(b) have different ordering of indices, since ITensor doesn't care what the ordering of the indices are. Could you run your code in debug mode and see if it outputs any more specific errors?

In order to help debug it would be helpful to have a minimal working example of the code that is failing.

Here is a minimal example of how this approach would work (this runs correctly for me):

auto N = 10;
auto sites = SpinHalf(N);
auto psi = MPS(sites);
auto b = 4;
// Set the center site gauge with ITensor
// Move the center gauge manually
auto lb = linkInd(psi,b);
ITensor U,S,V(lb);
psi.Aref(b+1) *= S*V;

linkInd(psi,b) is a shorthand for commonIndex(psi.A(b),psi.A(b+1),Link). Also note that if you put indices only on the V tensor, the svd function will use those to determine how to reshape the input tensor.

Note that a more efficient alternative to using svd, if you don't need extremely high accuracy of your singular values, is to use denmatDecomp ( http://itensor.org/docs.cgi?page=classes/decomp ), which directly computes T = U*B where B = S*V.


commented by (190 points)
Dear Matt,
thanks a lot, the code you wrote works indeed very well. I guess the error I was getting was because I was imposing the index lb on the matrix U instead of the matrix V (and yes, phi was a copy of psi.A(b) ). Thanks a lot for your time and apologies for the numerous questions!
commented by (3.6k points)
Glad to hear it works!

It would also be fine specifying the indices of `U`, I was just showing a slightly simpler way since there is only one index that goes onto `V`. Perhaps there was another small bug in your code that rewriting managed to fix.

Welcome to ITensor Support Q&A, where you can ask questions and receive answers from other members of the community.

Formatting Tips:
  • To format code, indent by four spaces
  • To format inline LaTeX, surround it by @@ on both sides
  • To format LaTeX on its own line, surround it by $$ above and below
  • For LaTeX, it may be necessary to backslash-escape underscore characters to obtain proper formatting. So for example writing \sum\_i to represent a sum over i.
If you cannot register due to firewall issues (e.g. you cannot see the capcha box) please email Miles Stoudenmire to ask for an account.

To report ITensor bugs, please use the issue tracker.