Feature request: index size scaling

[lkdvos]

In the context of specific tensor network algorithms, often the sizes of some of the tensors are not fixed, and typically one of the sizes is used as a parameter to tune the accuracy/expressiveness of some of the tensor network ansatzes. (For example, the bonddimension of a matrix product state).

It would be cool to be able to specify a contraction cost that is a symbolic polynomial in some variable, and to determine the optimal order in the limit where this specific size dominates. This effectively means that the metric of optimizing is to select the order that minimizes the exponent of that symbolic variable, and then selects the option with the smallest prefactor.

As a concrete example of what I mean, this exists in TensorOperations.jl, see for example the section on @tensoropt.

[mofeing]

mmm this could be done by refactoring the EinExpr.size Dict-field to allow symbolic expressions from Symbolics.

struct EinExpr
    ...
    size::Dict{Symbol, Union{Int, Symbolics.Num}}

It will also need some specialized methods of resource counters like flops.

Questions

1. @lkdvos What should be the behavior of comparing 2 symbolic polynomials? Or comparing 1 symbolic polynomial against 1 BigInt?
2. @jofrevalles This will probably incur into type-instability (less probable) and runtime dispatch (more probable). Can it be avoided if we are not using symbolic polynomials? If not, is the performance overhead significant?

[lkdvos]

In principle, the behavior I had in mind is to only allow a single variable, for example x, and then formally take the limit of x being "large". This means that:

x > 1
2x > x
x^2 > x
...

[mofeing]

That will require the computation of limits of symbolic expressions. I'm not so familiar with Symbolics, but lets see if it can manage computing

$$ \lim_{x\to\infty} x $$

$$ \lim_{x\to\infty} \frac{x}{x^2} $$

[lkdvos]

You could also just define a very lightweight polynomial type manually. It does not need to be able to do that much, addition, multiplication, scalar multiplication and comparison should really be sufficient, and this does not add a dependency on symbolics

[mofeing]

A lighter (or custom) alternative to Symbolics can be used, but I don't want to use a naive Expr. I will discuss alternatives next week with @bsc-quantic/software

[mofeing]

So I've been trying to see if Symbolics is up to the job.

using Symbolics
@variables χ

a = χ^2
b = χ^3

This doesn't work because it returns false always.

substitute(a < b, χ => Inf) # = false
substitute(a > b, χ => Inf) # = false

This seems to work with polynomials...

substitute(a / b, χ => Inf) # = 0.0
substitute(b / a, χ => Inf) # = Inf

...but not with asymptotically bigger functions like exp.

c = exp(χ)

substitute(a / c, χ => Inf) # = NaN
substitute(c / a, χ => Inf) # = NaN

Symbolics doesn't implement symbolic limit yet, but I found out that SciPy uses the Gruntz algorithm. Furthermore, it does have a simplified implementation for when $\lim_{\chi\rightarrow\infty}$ which is what we are looking for.

@lkdvos If you guarantee that all you need is polynomials, I can support it and implement it now. But if you need something more complex, then we need to implement this "Gruntz algorithm".

[mofeing]

ey @lkdvos I want to retake on this in order to have it for v0.6 (which I plan to release in the general registry)

sth I've thinking about is that maybe I don't really need to refactorize many things as long as you ONLY use polynomials and you represent these polynomials as "vectors" which are easily comparable

[lkdvos]

I think that should be fine, I am used to this implementation: https://github.com/Jutho/TensorOperations.jl/blob/master/src/indexnotation/poly.jl In particular, there is not really a need a priori for any functions except the basic * and +, as this is typically all that is needed for the basic costs of contractions.

[mofeing]

Okay, going to try it and will inform you again.