Hello everyone.I am using the julia version of ITensors for my thesis project and I am currently trying to implement the operator e^{iaO} in my code where
\hat{O}=\sum_{j}X_{j}
is a sum over single-site operators and X is the pauli X.I would like to implement this operator as an MPO so that I can act on a MPS with it.My initial thought was to construct the sum over local operators as an MPO first and then take the exponential of this MPO but there doesn’t seem to be a function which calculates the exponential of MPO’s.My next thought was to write the exponential of the sum as a product of exponentials and then turn this product into a MPO.I know that for sums of local operators there is the OpSum() object but is there something similar implemented for products of local operators?
Since all of the operators X_j in your sum act on single sites, the best way to make the exponential of that sum is to just write it as a product of exponentials:
can be computed up-front and is just a particular matrix.
Then you can manually make an MPO like M = MPO(N) and set each tensor of this MPO to equal the correct matrix:
sites = siteinds(...) # put the correct arguments here for the type of sites you want
M = MPO(N)
eaX = make_exponential_operator_matrix(a)
for j=1:N
M[j] = ITensor(eaX, sites[j], sites[j]')
end
where make_exponential_operator_matrix is a function you could make which defines the operator matrix that you will put on every site.
@ntsantid - Miles has already described a nice approach, but below is my code example with OpSum().
First you can define a custom op. that defines the matrix form of the operator e^{iaX} using what’s described in the documentation Physics (SiteType) System Examples · ITensors.jl. Say you call this op eX. Then do the following (supposing your range of sites over which O acts is [a,b])
os1 = OpSum()
os1 += eX,a
for i=a+1:b
os1 *= eX,i
end
M = MPO(os1,sites) # where sites=siteinds("Qubit",N) or whatever else.
I did a similar thing where I define a custom operator for e^{-iX} but since I am dealing with single site operators I thought it would be easier to construct each term of the product as an MPO using OpSum() and then apply these MPOs to my state sequentially.Here is the code that I am using which is for a=0.2 :
using ITensors
using ITensorMPS
ITensors.op(::OpName"expix",::SiteType"Qubit")=exp(-im*0.2*[0 1; 1 0])
function apply_xerror(sites,N,psi)
phi=copy(psi)
for i in 1:N
errx=OpSum()
errx+="expix",i
phi=apply(MPO(errx,sites),phi)
end
return phi
end
Let me know if I am doing anything wrong,since I don’t have much experience with Itensor.I would also like to use a as parameter at some point,so do you know whether defining custom operators inside a function works or should it be a global definition?
I just want to mention that some of the code above (in both of your answers) seems to be conflating the idea of a sum of local operators and a product of local operators. I.e.
e^{ia X_1} + e^{ia X_2} + \ldots
versus
e^{ia X_1} \otimes e^{ia X_2} \otimes \cdots
Note that OpSum is really meant to be used for the first case, that of a sum of operators, and not primarily second case.
You can stretch this a bit with OpSum, writing things like:
Just to ensure from you, what I wrote is constructing the same thing as you did “by hand” in your example, isn’t it ? I always thought so, because sometime ago I did what I exampled in my post for long (as in too long to define by hand individually) string operators in certain models with SPT phases and the resultant numbers always agreed with other people’s results in the literature before I went ahead with it for my own calculations… I agree it is a bit of abusing OpSum()'s purpose and at first it did feel too hacky to me but then the numbers always agreed with others in the benchmarking/test cases so it seemed to do the job.
I see what you mean @adityab999 and your solution does indeed work (as long as the variable eX is defined… I guess you defined it as a particular string)? To be honest I think the *= operator for OpSum is a bit of a hidden feature and not one I had used much. Matt had been doing some work to expand the capabilities of OpSum but I think (I could be wrong!) the *= feature was experimental and not documented.
Anyway, yes if it works I agree it’s fine for now. The main thing as always is to test the output on some examples to verify that it’s working as expected.
@miles - I’d defined the custom op eX just as e^{-i\pi X} (likewise for eZ and so on) and not as a particular string… The string was then constructed with the OpSum()… I was surprised at first that this worked, given that * operator is not the same as the , operator used to construct products of operators in terms within the OpSum() but apparently operationally it works perhaps for the reasons you alluded to above.
I just wish to mention that after thus constructing the string op and then converting it to an MPO, I only needed to take its expectation values (using inner()) for my purposes… I don’t know if using this way of making a string operator and then applying it on an MPS using apply() would also work correctly (probably would).
The main thing to be careful about is that the *= will append the operator to all previous terms in the OpSum. This is consistent with your example code but just wanted to mention this behavior.