contract a gate with a MPO

Hi,

I have a question about contracting gates with MPO. Suppose I have a MPO A defined on sites s. Then for each site i, A[i] has two site indexes. How to create a gate and contract it with a specific gate?

I tried

n = 4
s = siteinds("Qubit",n)
A = randomMPO(s)
g = op("X",s[1])

Now I want to act this operator on A. I find that g*A doesn’t work but apply(g,A) seems to work. But I am not sure which site index of A[1] is contracted with operator g.

Thanks

apply is a good way to apply the ITensor g to the MPO A, and it will apply g to the first MPO tensor A[1] and update this tensor for you. It is able to figure out that this is the right tensor since it has a matching index s[1].

To clarify your question: is it whether the ITensor g will be contracted with the s[1]' Index of A[1] or with the s[1] Index of A[1]? I can see why this isn’t obvious.

The answer is that it will be contracted with the primed index s[1]' of A[1] and then that tensor of the MPO will be updated. You can think of what apply is doing as being like a matrix product where each operator (g and A) are mappings from the unprimed index space to the primed index space. Then apply composes these two operators, conceptually multiplying them like matrices.

Here is a diagram of what happens when apply is called for this case:

1 Like

Hi Miles, thanks a lot for the answer. I have a further question, how to contract operator g with the index s[i] rather than s[i]'?. Basically what I am thinking is acting operator gon the two sides of A like g A g^{\dagger}.

Oh I guess the answer is simply apply(A,dag(g)) if I want to get Ag^{\dagger}.

It’s a good guess, but actually apply(::MPO,::ITensor) isn’t defined so that won’t work.

We actually have a keyword argument to apply that will do exactly what you are looking for, which is to compute gAg^\dagger. This argument is called apply_dag and you can use it like this:

B = apply(g,A; apply_dag=true)

The result will be B = gAg^\dagger.

One reason that computing Ag^\dagger isn’t as simple as apply(A,dag(g)) is that the dag(g) function, perhaps confusingly, doesn’t compute g^\dagger in the operator sense. The dag function does conjugate all of the elements of an ITensor and reverse all of the Index arrows, but it has no effect on the prime levels of the indices. So to compute the Hermitian conjugate of a single-site operator ITensor, you need to do swapprime(dag(g),0,1). Perhaps that’s something we should improve in the future, like adding a hermitian_conj(::ITensor) function or similar.

1 Like

Thanks! This is exactly what I need. All these work for double sites operators right?

Yes they should work for two-site or multi-site operators or gates (ITensors) as well.

1 Like