I am calculating the String Order parameter and I am trying to create a function like the following.

Cplx C_STRING(const System S, const ITensor n, const int ni_site, const int nj_site)
{
auto sites = S.Sites;
auto psi0 = S.psi0;
auto dn_k = ITensor(sites(1));
int nbar = 1;
for (int k = ni_site; k <= nj_site; k++)
{
auto n_k = op(sites, “N”, k);
auto nbar_k = ITensor(n_k.inds());
nbar_k.fill(nbar);
dn_k += (n_k - nbar_k);
}
auto dN = expHermitian(M_PI * dn_k, 1_i);
auto nj = SingleSiteOp(n, psi0, nj_site);
auto _nj = dN * nj; // not working
}

The form of which is

where in my calculation, it can be assumed that j=i and j=r=j. I have created that exponential operator in ITensor form, but in the last line of code, I am unable to apply that operator to MPS. How to apply that ITensor to MPS? Should I need to convert to MPO, how to do that?
In my code

MPS SingleSiteOp(const ITensor A, MPS x, const int site)
{
auto newA = A * x(site);
newA.noPrime();
x.set(site, newA);
return x;
}

Thanks for your reply. But I do not have MPO. Instead, I have the exponential of sum of many single-site operators(dN). I want to apply that dN to MPS. My question is how to do that. Should I need to convert that ITensor to MPO and then apply, or is there any more efficient way to do that?

If you have an exponential sum of operators, then that is a challenging situation. It’s possible that the sum can be compressed into an MPO, but whether there is an efficient way to do that depends very much on details of the sum and the operators.

Going back to the formula in your first post above, are you saying that this exponentially many operators corresponds to the string operator which is an exponential of a sum of local operators? In that case, yes there is an efficient MPO expression for it.

Could you write below, briefly, the mathematical expression for the operator you want to make as an MPO? Then I can probably tell you if it is possible.

My main problem is with applying e^{i\pi\sum_{i\le k<j}}\delta n_k this term on MPS. Here \delta n_k=n_k- \overline{n}, the number fluctuation from average number \overline{n}=1. I have thought in the following way

I can create all single site operators \delta n_k, but summing them up is difficult as it will throw errors like the index does not match. But if I can do that, I can use expHermitian. The next step is to apply that ITensor to the MPS. It is no longer a single site operator but a sum of many of them.

I can create an autoMPO, but how to exponentiate it because using toExpH with tau as 1 has a very poor approximation, and things diverge.

The way you should think about it is that \sum_{i \leq k < j} \delta n_k is a sum of local, single site operators, which all commute with each other. So then the exponential of this sum is a product of exponentials, namely:

\exp(\sum_{i \leq k < j} \delta n_k) = \prod_{i \leq k < j} e^{\delta n_k}

So all you really need to do is to make the individual e^{\delta n_k} operators and your operator is just a product of them. So you can apply them to your state one at a time using the code patter in this Code Formula: http://itensor.org/docs.cgi?vers=cppv3&page=classes/mps_mpo_algs

Does that plan give you what you need? As far as the constant \bar{n} shift goes, I think you can divide that up between the different operators since it is just a classical c-number constant.

Thanks for your reply. It is working now, at least with no errors, but I need to verify the result. Anyway, I have two small doubts.

I can creat \delta n_k as following

auto n_k = op(sites, “N”, k);
auto nbar_k = nbar * op(sites, “Id”, k);
auto dn_k = (n_k - nbar_k);

I can then creat the exponential operator like auto exp_dN = expHermitian(M_PI * dn_k, 1_i);. Now I need to do SingleSiteOp(exp_dN, psi_j, nk_site);. In my case the function

MPS SingleSiteOp(const ITensor A, MPS x, const int site)
{
auto newA = A * x(site);
newA.noPrime();
x.set(site, newA);
return x;
}

Here k, nk_site all refer to k^{th} site. This process needs to be done for all k. I hope this is what you have suggested to me!
2. In my previous reply I mentioned toExpH, which approximates the exponential like toExpH(ampo,tau) up to term \tau^2. It means if I put \tau as 1, it will not be a good approximation. But what if my single particle operators are commuting themselves, with which I have made ampo, then also it is a bad approximation!!