Site operators in denominators are accepted silently; expand powers in n-site terms

Open
  • Is something wrong?
#266 0 comments 0 reactions 1 assignee View on GitHub

Maintainers usually reply within 3 days

@Ooolab is already working on this.

Since Oct 28, 2026.

  • #318 by @marcusr2ML — open

Assessment

This issue has not been assessed yet.

Description

area:dmrg bug

Two follow-ups for site operators inside n-site terms (NSiteOperator::split, src/alps/model/include/alps/model/nsiteoperator.h), and one for the native site-term path.

1. Division by a site operator is accepted silently

Dividing by a site operator is never meaningful (Sz and n have zero eigenvalues), and no model, tutorial or example in the repository does it. Every / in a term on master divides by a number or a parameter (J#/2, n_up(i)-1/2, n(i)*(n(i)-1)/2). But the parser accepts it without an error and uses the operator in its place:

  • <SITETERM site="i">1/Sz(i)</SITETERM> (spin, local_S=1.5, fulldiag) gives the same spectrum as Sz(i).
  • <NSITETERM sites="a b">1/Sz(a)*Sz(b)</NSITETERM> and Sz(b)/Sz(a) give the same spectrum as Sz(a)*Sz(b).

A mistyped term, for example / for *, therefore runs without a message. Both paths should refuse a site operator in a denominator with an error. Division by numbers and parameters must keep working. <BONDTERM> uses the same evaluator and has not been checked yet.

2. Powers of site operators in n-site terms: expand instead of refusing

Since 5396639f8, split refuses any power of a site operator (Sz(a)^2, (Sz(a)*Sz(b))^2, (Sz(a)+Sz(b))^2, Sz(a)^p). Before that, a power was silently dropped (n(i)^2*n(j) acted as n(i)*n(j)). Site terms accept Sz(i)^2 for diagonal operators and refuse off-diagonal ones ("Cannot apply offdiagonal operator inside function argument or power"). So the same input works as a <SITETERM> and fails as an <NSITETERM>.

Proposal: in split, rewrite a positive integer power X^k as X*X*...*X before splitting, then flatten(). This is exact for diagonal and off-diagonal operators alike. Things to get right, each with a test against the written-out product:

  • Order: (A+B)^2 must become AA + AB + BA + BB, not A^2 + 2AB + B^2. Operators on one site need not commute, and fermions pick up signs. Repeating the factor and flattening keeps the order.
  • Coefficients inside the power: (K*n(a))^2 gives K^2.
  • Exponents: ^p expands only if p evaluates to a positive integer. Non-integer and negative powers are still refused. ^0 becomes 1.
  • Nothing dropped: a power anywhere in a term (inside parentheses, in a sum, with a parameter exponent) is either expanded or refused, never passed through as a single operator. That was the bug fixed in 5396639f8.

dmrg and ED both build their terms from split, so both get the change.

Dominant language
C++
Stars
32
Forks
10
Avg merge
4d 5h
Merged PRs (30d)
52

Getting set up

First steps

  1. Read the whole issue, then the project's contributing guide.
  2. Comment on the issue to say you are picking it up — it saves two people doing the same work.
  3. Fork the repository and make your change on a branch.
  4. Open a pull request that references the issue number.

More from ALPSim/ALPS

All issues in ALPSim/ALPS

Similar issues

More C++ issues

Get new issues in your inbox

A short digest of beginner-friendly GitHub issues.