Exploiting Dynamic Sparsity in Einsum
Authors: Christoph Staudt
Mark Blacher
Tim Hoffmann
Kaspar Kasche
Olaf Beyersdorff
Joachim Giesen
Venue: NeurIPS 2025 poster ()
Keywords: einsum
tensor network
dynamic sparsity
tensor contraction
Primary area: infrastructure
TL;DR: Evaluating real-world einsum expressions can lead to hyper-sparse intermediate tensors that are hard to predict, we exploit this dynamic sparsity in a hybrid algorithm.

Abstract: Einsum expressions specify an output tensor in terms of several input tensors. They offer a simple yet expressive abstraction for many computational tasks in artificial intelligence and beyond. However, evaluating einsum expressions poses hard algorithmic problems that depend on the representation of the tensors. Two popular representations are multidimensional arrays and coordinate lists. The latter is a more compact representation for sparse tensors, that is, tensors where a significant proportion of the entries are zero. So far, however, most of the popular einsum implementations use the multidimensional array representation for tensors. Here, we show on a non-trivial example that, when evaluating einsum expressions, coordinate lists can be exponentially more efficient than multidimensional arrays. In practice, however, coordinate lists can also be significantly less efficient than multidimensional arrays, but it is hard to decide from the input tensors whether this will be the case. Sparsity evolves dynamically in intermediate tensors during the evaluation of an einsum expression. Therefore, we introduce a hybrid solution where the representation is switched on the fly from multidimensional arrays to coordinate lists depending on the sparsity of the remaining tensors. In our experiments on established benchmark einsum expressions, the hybrid solution is consistently competitive with or outperforms the better of the two static representations.

Exploiting Dynamic Sparsity in Einsum
Christoph Staudt", Mark Blacher, Tim Hoffmann, Kaspar Kasche, Olaf Beyersdorff, Joachim Giesen Friedrich Schiller University Jena
Abstract
Einsum expressions specify an output tensor in terms of several input tensors. They offer simple yet expressive abstraction for many computational tasks in artificial intelligence and beyond. However, evaluating einsum expressions poses hard algorithmic problems that depend on the representation of the tensors. Two popular representations are multidimensional arrays and coordinate lists. The latter is more compact representation for sparse tensors, that is, tensors where significant proportion of the entries are zero. So far, however, most of the popular einsum implementations use the multidimensional array representation for tensors. Here, we show on non-trivial example that, when evaluating einsum expressions, coordinate lists can be exponentially more efficient than multidimensional arrays. In practice, however, coordinate lists can also be significantly less efficient than multidimensional arrays. but hard decide from the input tensors whether this will be the case. Sparsity evolves dynamically in intermediate tensors during the evaluation of an einsum expression. Therefore, we introduce a hybrid solution where the representation is switched on the fly from multidimensional arrays to co- ordinate lists depending on the sparsity of the remaining tensors. our experiments established benchmark einsum expressions, the hybrid solution is consistently competitive with or outperforms the better of the two static representations.
Introduction
Einsum was introduced in NumPy [9 in 2011 and has since become the quasi- standard for specifying pure tensor expressions, meaning expressions that involve no non-linear operations tensor entries. is now an integral part of machine learning frameworks such as TensorFlow [.. PyTorch 255, and Nvidia's cuQuantum [5] library for quantum information science. So far, these frameworks support only dense tensor formats einsum. While dense formats are widely used, sparse tensor formats such coordinate lists offer potential advantages when many tensor entries are zero. They reduce the number of floating-poit poinn operations (flops) compared to dense formats and, as we demonstrate yield exponential improvements in specific cases. However, sparse formats
in this paper, may also have drawbacks. They may require more memory when tensors are not particularly sparse, and while in theory their space usage grows at most linearly, in practice the memory overhead can be considerable. Moreover, they often use a less favorable data layout for computation, which can lead significantly slower execution. Thus, choosing the appropriate tensor format is critical for efficient execution of einsum expressions.
Sparse tensor formats in einsum supported in the Python package sparse 226 and also in SQL (8]. These implementations, however, rely on static format choices. In our experiments, we show that sparsity can evolve dynamically in intermediate tensors during the evaluation of einsum expression, especially for tensor expressions that involve a large number of input tensors. Such expressions
"Corresponiing author: christoph. staudteuni-jent -jena.de
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

are not yet a primary focus of established frameworks, but are becoming increasingly important in model classes such as probabílistic circuits [T]
23 32 and probabilistic neural eircuits [38, which often rely on einsum as a computational backend, as well as in neuro-symbolic models that combine 12 Moreover, the simulation and validation of quantum 4 35 29 also benefit directly from dynamically sparse einsum. However, since sparsity can evolve dynamically, it is difficult to determine from the expression and the input tensors alone which tensor format is best for given einsum expression. Therefore, to efficiently exploit dynamic sparsity in einsum, we introduce a hybrid solution, where the representation is switched the fly from a dense to sparse format once the sparsity of the remaining tensors falls below a predefined threshold. In our experiments we show that the overhead for switching from one format to the other is relatively low. As result, the hybrid solution is consistently competitive with or outperforms the better of the two static representations.
deep learning with logical constraints T0
machine learning algorithms 17,
Our approach to einsum aligns well with a current trend in the machine learning community that is shifting from traditional static sparse modeling to dynamic sparsity [7 21 24
31
3
34 37], which
is broadly defined any kind of sparse computation or memory compression where the sparsity pattern is input-dependent or adaptive [23].
2 Einsum expressions and tensor formats
In this section, we provide the necessary background on einsum expressions and their evalua- tion, as well as dense and sparse formats for storing tensors. Einsum expressions are of the form einsum .I,,...; In - O; ......: Tn), where the T. are input tensors with index strings I,, and O is the index string of the output tensor. The index strings consist of symbols that represent tensor axes, sometimes also called dimensions or modes. The size of the index string the order of the tensor T.. Notably, input tensors can themselves be einsum expressions.
2.1 Evaluating einsum expressions
All indices that appear in at least one input tensor but not in the output tensor are summed over. Indices that occur in two or more input tensors signify a combination, and if such in- dices are also summed over, they define contraction. The simplest combination is the elemen- twise product of two vectors u and v, written as cinsum (i, + i; u, v). The simplest contrac- tion is the inner product u. v = einsum (i, +; u, v), where the output index string is empty because the result is a scalar. More complex expressions involve multiple contraction indices, such as the matrixxvxaxxxxoooor product einsum (ij.jk, k - i; A, B, v), which results in a vector. Contracting over index j first, as in einsum (ik, k + é; einsum (ij, + ik; A, B), u). requires 2- 1- [5|- |k] + 2- lil [k] floating point operations (flops), whereas contracting over & first, as in einsum (ūj,j + i; A, cinsum (jk, k + j; B, v)), requires only 2 |5 |k| + 2 |i| [k]] flops. Thus, the second contraction order is more efficient.
Computing an optimal contraction order that minimizes the number of flops is, in general, an hP-hard problem 11] Therefore, heuristics such as greedy strategies [ or graph partitioning methods [8] 30) are used in practical solvers for computing tensor contraction orders. The number of flops depends not only on the contraction order but also on the choice of data structure for representing tensors. Sparse tensor formats can significantly reduce the number of flops compared to dense formats when the tensors involved in the contractions are sparse. Note, however, that in einsum expressions involving a large number of tensors, it is not known in advance whether intermediate tensors will become sparse.
2.2 Dense and sparse tensor formats
In machine learning contexts, tensors are often equated with multidimensional arrays. Here, we take more abstract view and consider tensors as multivariate functions of the form T A - R. where A = A₁ X .... x A,, with finite sets A,, which are called tensor axes. The natural data structure for storing tensor T is multidimensional array, where T(a.,..., a,.] = T(a,, ... a,)) Alternatively, tensors can be stored as lists of coordinate entries of the form (aI,..., (n, T(a,,.... :.,)....)))) (n))). Lists have the advantage that sparse tensors, meaning tensors in which significant proportion of function values are zero, can be stored more compactly by only storing the coordinates where the function
2

---

values are nonzero. For dense tensors, however, the list format is less efficient because of the overhead incurred by explicitly storing the coordinates.
3 Exponential separation
In this section, we show that there exist cinsum expressions that can be evaluated exponentially more efficiently when using coordinate lists instead of multidimensional arrays. These expressions are represented as tensor hypernetworss, which form the underlying algebraic structure of einsum expressions. The tensor hypernetworks used to demonstrate the exponential separation between data structures stem from model counting problems. The model counting problem #SAT is of independent interest in artificial intelligence, because variety of real-world questions can be encoded very naturally and succinctly in #SAT, including problems in probabilistic reasoning [2 20 risk analysis [14 36), and explainable artificial intelligence 3 27
3.1 Tensor hypernetworks and einsum
Tensor hypernetworks define output tensor in terms of input tensors. Formally, let A (Ai...., Am be a set of axes. A tensor hypernetwork T is given by a finite set (T T,} of input tensors such that the axes of each imput tensor T; are subset A, S A of the axes set. The subset A; is encoded in the index string I; of the corresponding einsum expression. For instance, the tensor hypernetwork corresponding to the einsum expression einsum (ij. jk, k i; A, B, v) from the previous section has three axes. The first axis, i, is used only by the matrix A. The second axis, j. is shared by the matrices A and B. The third axis, k, is shared by the matrix B and the vector v. Therefore, the set of axes for the matrix A contains the first and second axes, indexed by ij. The set of axes of the matrix B contains the second and third axes, indexed by jk; and the of axes of the vector v contains the third axis, indexed by k.
The output tensor is given by combining input tensors that share common axis and by summing over all axes in A that are not part of the output tensor. Here, in order to keep the exposition simple, assume that the output tensor has no axes, that is, it is scalar, and thus we The full combination of the tensor hypernetwork T is a tensor T with axes A, such that, for a E A, T(a) II. T, (a).), where a is the projection of a onto the axes in A; The full summation of T, that is, the scalar output tensor, is the sum of all function values of T and can thus be written as
over all axes in A.
ΣT(a) - Σ -E-III I(a.). AEΛ 1=1
REA
Algorithmically, of course, it is not a good idea to first compute the full combination and then the full contraction, because the full combination results in large intermediate tensor. In practice, the tensor hypernetwork is contracted by contracting over the axes one after the other, following an optimized contraction order. To contract over an axis, one collects all tensors that share the axis, computes the full combination of these tensors, and sums it over the given axis. Finally, the collected tensors are replaced by the sum over the full combination, which is single tensor but not necessarily a scalar.
The combinatorial structure of a tensor hypernetwork can be encoded in a hypergraph: the nodes given by the input tensors in T and for each axis A; there is hyperedge that contains all input tensors T3 with A, E Aj. Contracting over an axis in a tensor hypernetwork corresponds to the contraction of the corresponding hyperedge in the hypergraph. The contraction of a hyperedge is illustrated by example in Figure|
3.2 Model counting by tensor contractions
Problems on propositional formulas naturally fit into the framework of tensor hypernetworks 6 16 A propositional formula is multivariate function of the form F {0, 1), - {0. 1} and thus tensor with axes A; = {0, 1} for i e [n].
Given propositional formula F, the satisfiability problem SAT asks if there exists E 10,1) such that F(a) = Such satisfying assignment is called model for F. The model counting problem #SAT asks to count all models, that is, all satisfying assignments.
3

---

Figure I: Left: Tensor hypernetwork for the grid formula GRIDy (cf. Section 3.3
Here, large black bullets represent input tensors (matrices), and small gray bullets represent hyperedges. The tensor hypernetwork has one hyperedge of size four (at the center), and four hyperedges of sizes three and two, respectively. Right: The tensor hypernetwork after contracting the hyperedge of size four. Four second-order tensors have been contracted into one fourth-order tensor.
Both problems, satisfiability and model counting, can be served by tensor hypernetwork contractions. This is possible because any propositional formula has representation in conjunctive normal form (CNF), that is, as conjunction of clauses. A clause is disjunction of literals, and literal is either variable x; its negation T.
Let F be a propositional formula in conjunctive normal form, clauses( be its clauses, and vars( F) be the set of variables in F. For clause C' E clauses(F), let vars(C) S vars(F) be the set of variables in C. Since any clause C' is itself a propositional formula, is also a tensor Tc. The axes of Te correspond to its variables vCrs(C). The tensor hypernetwork F that results from the CNF representation of the formula F' is given as F {Tc C' clauses(F)). The full combination of F corresponds to the formula F, and the full summation
Σ
II Te(а)),
aE{0,1)" CEclauses(F)
where die is the projection of onto the axes of Te that correspond to vars(C), is the propositional model count.
3.3
Exponential separation
We separate the two data structures, multidimensional arrays and coordinate lists, on family of grid formulas. Let [n] = [I,...,.. n), the formula GRID, has variables Tit for i, E [ri] and clauses
-iy 24+1 and Fi.j -# for i [n 1], E [n]. Bi, V Ris+1 and Ti,j V -iij+-. for i [m], j E [n 11.
That is, GRID, has 4n(n 1) clauses and n² variables. The clauses in the first line encode the condition Fi = 1e+1, and the two clauses the second line encode vi,j = di That is, by construction, satisfy GRID, all variables have take the same value. Therefore, GRID, has exactly two models.
If the clauses encoded in multidimensional arrays, then these multidimensional arrays given by the following matrices T (11) and V V- = =(d1). where, in both matrices, the columns labeled by Zz 0, and the are labeled by I4+1i 0, 1, starting in the top left corner. The encoding of the clauses -iy V ri,t+1 and Fi5 V -1. into matrices follows analogously.
Contractions are over the axes, that is, the variables. The contractions over Vijj and Ti+1 respec- tively, involve both tensors Ta,,, and Ta.,vV-urr which are contracted by an elementwise multiplication. Therefore, the tensors T-z., V2 and T2,,y-uu V-R++1., can be replaced with their elementwise product in the tensor hypernetwork encoding, which reduces the number of tensors to 2n(n - 1) and is standard practice in contraction algorithms. the multidimensional array encoding, all the 2n 1) tensors (2 x 2) -matrices with the diagonal and ( on the off-diagonal. Their coordinate list encoding is given as ((0,0, 1), (1,
Lemma Contracting the tensor hypernetwork encoding of GRID,, n 2 needs 4n² 4n 1 flops, when using the coordinate list encoding.
Proof. As we have explained before, the tensor hypernetwork contracted by contracting over the n² axes one after the other. Since all axes for the GRID,, tensor hypernetwork have size two, the

---

arguments of all tensors, including the intermediate tensors, can be encoded by bitstrings. The lengths of the bitstrings are the orders of the corresponding tensors. For the proof, we can assume that the bits are ordered by the contraction order of the corresponding axes. Therefore, the first bit always corresponds to the contraction axis.
We show by induction that any algorithm for contracting the tensor hypernetwork encoding of GRID,, keeps the following invariant: All tensors that arise during the contraction, except for the last tensor, have exactly two non-zero elements, namely, at the all-zeroes and the all-ones bitstrings, where they evaluate to The last tensor gives the value of the tensor hypernetwork, which is a scalar.
Base case. By construction, the invariant holds true for the input tensors that encode the clauses as shown above.
Inductive step. If the invariant holds true before contraction, then it also must hold true after the contraction, because, as we show below, computing the full combination for the given axis and the subsequent summation over this axis both keep the invariant.
Full combination: Assume that e tensors need to be contracted. The value of the full combination at every bitstring is given by the product of the values of the tensors that need to be contracted the projection of the bitstring onto the axes of these tensors. In particular, we only multiply elements whose first bit is equal. By our induction hypothesis all to be contracted tensors only take non-zero values at the all-zeroes and the all-ones bitstrings. Therefore, only two products are non-zero, namely the products at the all-ones and the all-zeroes bitstrings.
Summation; We need to sum all entries that differ on the first bit, but are equal on all other bits. Since the first bit is always equal the other bits, no summation occurs, until reach the last contraction, where only one bit is remaining. The summation for the last contraction always boils down to the This results in the correct model count 2.
addition
Required flops. As discussed above, the full combination of e tensors needs 2(c - 1) flops. Overall need to combine 2n(n- 1) tensors leading to 2(2n(n 1) 0011 4n(n 1) 2 flops. As only the final summation needs a flop, the total number of flops is 4n 4n 1.
The proof shows that, when intermediate tensors become large in terms of the number of axes, they become hyper-sparse. Next, show that the intermediate tensors, actually become large in terms of the number of axes. The exponential separation then follows, because all entries must be computed in the dense tensor encoding.
Lemma 2. Contracting the tensor hypernetwork encoding of GRID,, needs at least 2D flops, when using multidimensional arrays.
Proof. Mostly follows from Proposition 4.2 in Markov and Shi 222 A full proof is included in the appendix.
We summarize Lemmas 2 in the following proposition.
Proposition 1. Contracting the tensor hypernetwork encoding of GRID, needs exponentially more flops, when using multidimensional arrays instead of coordinate list encodings.
Here, the data structures can be separated exponentially, because the intermediate tensors become hyper-sparse. In general, it is hard to tell from the input tensors if intermediate tensors will become sparse enough to warrant the sparse representation, which can also incur significant overhead over the dense encoding. Even in the grid formulas the input tensors have an average sparsity of 50%, that is, only half of the entries are zero, which, in general, is not considered sparse.
4 Experiments
We first verify Proposition[[eemprrcc
by simulating the model counting instance GRID,, Since this is a synthetic result, we also analyze dynamic sparsity on recent benchmark of 168 einsum expressions [9] from areas such as probabilistic models, weighted model counting, tensorized language models, and quantum computing. Our findings show that hyper-sparse instances oceur in practice but often become evident only during execution. Sparse formats are efficient in such
5

---

cases but introduce overhead on dense instances. To address this, we propose hybrid approach that performs well across the benchmark.
All experiments were run five times and we report median runtimes. Experiments were conducted 64-core system (4 x Intel Xeon Gold 6130, 2.1 GHz) with 1.5 TB RAM. Code is available https://gtt thub. cwm/tii- group wonmtiiggroup//////oooooooooopet dense baseline, as it has been shown perform best on the einsum benchmark. Existing sparse einsum libraries cannot handle all benchmark cases due to index limitations, especially in expressions with many input tensors and large intermediate tensors, we use our own prototype implementation. We exclude ten dense instances from the einsum benchmark due to resource constraints, leaving 158 instances. Optimized contraction orders for all instances are provided by the benchmark.
saarsspatttttettiisreeee We use PyTorch [25] for the
4.1 Exponential speedup on GRID,
To confirm Proposition empirically, we executed the einsum expressions for model counting tensor hypernetworss GRID, for [2, 80] using both dense and sparse implementations. As shown in Figure 2] the exponential separation between the two tensor formats clearly visible. Moreover, the figure illustrates that, predicted by theory, the runtime for the dense tensor format grows exponentially in n, while the runtime for the sparse format grows only quadratically.
dense
Ronning Time
1000
spurse
Running: Time
0.6
0.4
500
0.0
20
40
60
20
40
60
80
Grid Site Length
Grid Site Length
Figure 2: Running times for dense and sparse tensor formats when evaluating einsum expressions for GRID, formulas with n E [2,80), 80].
4.2 Dynamic sparsity in the benchmark
The results on the synthetic GRID,, expressions show that the choice of tensor format can significantly impact the time needed to evaluate einsum expressions if tensors become sparse during the evaluation. We refer to the latter phenomenon as dynamic sparsity. To assess if dynamic sparsity also arises in practice, we analyzed the sparsity evolution, in particular the average density, during the evaluation of einsum expressions from the benchmark dataset.
Average density. In order to benefit from sparse tensor formats, the fraction of non-zero entries, and thus saved flops, needs to be small across the whole contraction sequence. Many sparse, but small contractions, will not lead to a significant speedup, if the running time is dominated by a few large and dense contractions. Therefore, we introduce avg, which denotes the averuge density of all tensors that arise in the contraction sequence up to the 7-th contraction. Here, the average density is the sum of all non-zero entries in all tensors up to the 7-th contraction, divided by the total number of elements across those tensors. The average density of the input tensors is avg.- We refer to avgTmax at the end of the contraction sequence simply as the average density of the expression.
Density evolution. We analyzed the average density of the benchmark expressions as a function o. In 104 out of the 158 considered benchmark instances, the difference between the maximum and minimum average density exceeds 0.2, indicating significant changes in sparsity over the course of their contraction sequences. Figure shoows the evolution of average density for these instances. The instances are grouped by their final average density into three main categories. For each group, the bold red line represents the median average density evolution. A larger version of this figure available in the appendix.

---

Group PS... (
Group 2: or (0.0.8)
1.0
01.0
0.9
Completion
Completion
Croup MiM... (0.0..0.) 6.8]
Group 3h: Mixal Dersit
0.5
114
0.0
0.5 Completion
0.0
0.2
0.9
Figure 3: Evolution of average tensor density across instances, grouped by their final values. Within each group, bold red line indicates the median density. The last group comprises instances transitioning from sparse (< 5%) to dense (>> 10%). An exponential scaling at both ends is used along the x-axis to show the fraction of completed contractions T/Tmax- A green dashed line marks the 5% density threshold.
The first group contains 42 instances whose initial average densities avg, vary widely, but become nearly dense toward the end of the contraction sequence. Several instances begin relatively sparse, but then transition to dense tensors as zeros are eliminated through summations.
The second group contains 28 instances that become very sparse. Five of these instances have fewer than 10 non-zero zero entries. Interestingly, none these instances starts out sparse. Their average sparsity arises only during large, flop-heavy contractions toward the end of the contraction sequence.
The third group contains 34 instances that gradually decrease in average density but maintain an average density above 5% and are split into two subgroups, Group 3a and Group 3b. Group 3b contains seven instances that exhibit more complex behavior. These instances first become sparse during the contraction sequence, but then their average density increases again. For instances that exhibit this mixed-density pattern, a sparse backend would not be optimal.
4.3 Sparse vs. dense on the benchmark dataset
The density evolution shows that there is potential for a sparse implementation to be faster than dense implementation on at least 28 instances. We confirm this by executing all instances using both the sparse and dense tensor formats. Figure 4 shows the relative speedup, (la/ta), where la is the running time for the dense format and is the running time for the sparse format. The red
7

---

line represents a speedup of one. Thus. points above the red line correspond to instances where the sparse format faster than the dense format. This appears to be the case for the majority of instances with average density below 5%. Afterwards, the dense format is mostly faster. In both cases, the speedup can reach up to 1,000x. Overall, the median slowdown of sparse vs. dense is around 17x on all instances and around 30x dense instances. Thus, using the sparse format for all benchmark instances is not viable in practice.
10a
10a
Relative Speedup
10'
100
10
other
int
10-2
small
10-3
(0.0, 0.01] (0.01, 0.05] (0.05, 0.8] (0.8, 1.0]
(0.0, 0.01] (0.01. 0.05]
Average Density
Average Density (double)
Figure 4: Relative speedup of the sparse over the dense tensor format on the benchmark dataset. The red line marks equal performance. The orange dots represent instances with an integer data type, and the green dots represent instances with fewer than 10th flops. The right plot shows the relative performance on the sparse integer instances after conversion to the double data type.
Outliers. There are some outliers where sparse is faster even though the average density is high, primarily for two reasons. First, some of the instances, marked in orange, have very few flops. Therefore, the choice of tensor format has little impact. Our sparse implementation is written purely C++, and thus has lower overhead during einsum expression evaluations, which leads slightly better running times on these instances. Second, as shown in the benchmark paper [9] Py Torch can experience significant slowdowns on integer operations, whereas our sparse implementation does not have this problem, which results in a speedup for these instances. To validate the speedup for instances with densities below 5%, we converted all input tensors from integers doubles and repeated the experiment. The speedup became smaller but persisted for sparse instances. Moreover, there are a few sparse outliers where PyTorch performs more efficiently, although their average density is below 5%.
4.4 Algorithmically exploiting dynamic sparsity
These observations show that sparsity in einsum is inherently dynamie: tensors that sparse can create dense intermediates, tensors with few nonzeros can produce hyper-sparse intermediates later. There are no sparsity estimators tailored to this setting, and while existing cardinality estimators from the database literature are likely adaptable they are often unreliable and expensive on large expressions. Rather than relying on a costly upfront estimate, our hybrid approach detects sparsity during execution and adapts the representation on the fly. Below we first describe the hybrid algorithm in detail and then evaluate effectiveness on the benchmark dataset.
Hybrid algorithm. The hybrid algorithm evaluates contractions in the dense format while the average density of the remaining (not-yet-contracted) tensors stays above chosen threshold. Once the average density falls below that threshold, the algorithm switches to the sparse format. For example, Figure T there are twelve remaining tensors on the left. After the first contraction the four center tensors are combined into single intermediate tensor, leaving nine remaining tensors

---

2x
2x
100%
2.x
2x
1x
2x
2X
100x
10'
Slowdoww
10
Slow Down/ Spood Up
10'
10'
Slow Down Speed
Specelup
Slowdown Speedup
10P
10a
10
10a
101
10!
10!
101
10
109
10"
10P
2
23
113
26
87
10
Dense Instances
Sparse Instances
Dense Instances
Spurse Instances
(a) hybrid sparse
(b) hybrid dense
Figure 5: Slowdown and speedup of hybrid algorithm on dense and sparse instances (av
0.05).
the right. Because all other tensors are unchanged by this contraction, updating the average density requires only computing the density of the new intermediate tensor and updating the running totals.
For the experiments reported here we used a 5% threshold chosen empirically from the benchmark. As indicated by Figure 4| this threshold appears to be the point where a speedup is most likely. In general, the threshold depends mostly the sparse implementation and should be chosen rather low, since this will only slightly delay the switch on instances that become hyper-sparse.
To minimize the overhead of calculating the density, we compute it only before expensive contractions. Moreover, once the average density exceeds 95%, we stop calculating In Figure 5] we compare the hybrid algorithm against the standard algorithm using either sparse or dense tensor formats.
Hybrid vs. sparse. The results of the comparison are shown in Figure 5a The hybrid algorithm is significantly faster on the vast majority of dense instances. It is more than T00x faster on 17 instances, and even more than 1000x faster on three instances. Overall, we achieve median speedup of 22x. We observe a small slowdown on small instances and on integer instances, which, as discussed earlier, are exceptional cases.
For sparse instances, the overhead of the hybrid over the sparse algorithm appears to be minimal. The median slowdown is approximately 1.08x. The hybrid algorithm takes twice as long for just four instances, primarily because it either does not switch to sparse or switches to sparse too late. However, switching too late only results in slowdown due to performance differences observed on integer instances. There are also three instances that benefit from the hybrid algorithm because they initially dense and thus benefit from the initial dense phase of the hybrid algorithm.
Hybrid vs. dense. The results of the comparison are shown Figure 5b) contrast to the comparison with the sparse algorithm, cannot expect speedup over the dense algorithm on dense instances, as additional operations are performed for the density test. However, in practice we observe small speedup on the integer instances from the third group, where sparse format because of their mixed density characteristics that can be seen in FFiuure 3 (right). The largest slowdowns incorrectly switch to the sparse format. On the other instances, the slowdown, and thus the overhead, on dense instances is indeed very small. The median slowdown is approximately 1.01x.
incorrectly switch to the
experienced on the remaining instances of the third group, where also
sparse instances, the hybrid algorithm clearly outperforms the dense algorithm, reaching a speedup of nearly 1000x very sparse instance. The median speedup is approximately 6.2x. The hybrid algorithm is notably slower than the dense algorithm on only one sparse instance, because the dense implementation performs better even though the instance is sparse.
Hybrid vs. virtual best. Overall, the hybrid algorithm closely matches the runtime of the faster static backend (the virtual best between dense and sparse) and occasionally outperforms it. The censity-estimation overhead negligible practice and provides clear gains on sparse instances. Measurable slowdowns occur mainly for mixed--ensiin cases and integer-typed inputs. As shown in Figure 6) the hybrid algorithm is more than twice as slow as the virtual best for only 8 of 158 instances, and the median slowdown is about 1.02x ( 2%). Conversely, the hybrid method can beat both static choices: our experiments it was substantially faster than the virtual best on of the 28 sparse instances.

---

2x
2x
2x
Slow Down Speed Up
10
10
Slowdown Speedup
10a
10
109
26
Dense Instances
Sparse Instances
Figure 6: Slowdown and speedup of hybrid algorithm vs. the virtual best algorithm.
Limitations
While our hybrid approach dynamically exploiting sparsity in einsum expressions demonstrates strong empirical performance across a range of benchmark instances, it also presents several limita- tions. First, the effectiveness of the hybrid strategy sensitive to the choice of sparsity threshold used to trigger format switching. Although our 5% threshold performs well across the dataset, this value may not generalize across all domains or applications. more adaptive or learned switching criterion could further improve performance but was not explored in this work.
Second, our current implementation only supports switching from dense sparse formats, and not vice versa. As demonstrated by mixed-density instances, some tensors become sparse mid-evaluation but later become dense again. For such cases, the optimal strategy may involve switching back and forth between sparse and dense formats based on the evolving average density. Our hybrid approach does not support this bidirectional switching, which limits its flexibility. Moreover, some contractions might benefit from mixing sparse and dense tensors. Supporting mixed-forma operations could further improve performance but beyond the scope of this work.
Third, our experiments concentrate on CPU performance because sparse einsum support on accel- erators (such as GPUs and TPUs) is currently device-agnostic in design and can be adapted to other hardware. In fact, different devices could be used for the dense and sparse phases, depending which device has better support for the respective format. However, we chose not to mix hardware in our experiments, as we consider the lack of sparse GPU einsum support a temporary limitation.
limited. Nevertheless, the hybrid algorithm is
Finally, our method assumes a static contraction order, which may not be optimal when sparsity evolves dynamically. Integrating contraction order optimization that is aware of evolving sparsity could offer further improvements but is non-trivial and left for future work.
6 Conclusions
Einsum is a key component of modern machine learning frameworks and also has applications in classical areas of artificial intelligence, such as probabilistic models, probabilistic inference, and constraint satisfaction. The range of tensor expressions that einsum implementations must support is extensive. Sparsity in einsum has traditionally been handled statically, meaning dense or sparse data structure is chosen before the evaluation begins. However, it is difficult to predict the most efficient data structure for given einsum expression, especially when it involves many input tensors. Our experiments show that sparsity often arises dynamically in weal-world einsum instances, with intermediate tensors becoming sparse during the evaluation, even if the input tensors are initially dense. To address this, we introduced hybrid algorithm that dynamically switches the data structure from dense to sparse based on the sparsity of the remaining tensors. The runtime of density checks negligible compared to the overall time required to evaluate the expression. Our hybrid algorithm is up three orders of magnitude faster than using an inefficient data structure and performs comparably to or faster than the optimal static choice.
10

---

Acknowledgments
This work was supported by the Carl Zeiss Foundation within the project Interactive Inference.
References
[1] Martín Abadi, Paul Barham, Jianmin Chen, Zhifeng Chen, Andy Davis, Jeffrey Dean, Matthieu Devin, Sanjay Ghemawat, Geoffrey Irving, Michael Isard, et al. TensorFlow: system for Large-Scale machine learning. In OSDI, 2016.
[2] Fahiem Bacchus, Shannon Dalmao, and Toniann Pitassi. Algorithms and complexity results for #SAT and Bayesian inference. In FOCS, 2003.
[3] Teodora Baluta, Zheng Leong Chua, Kuldeep S. Meel, and Prateek Saxena. Scalable quantitative verification for deep neural networks. In ICSE, 2021.
[4] Afrad Basheer, Yuan Feng. Christopher Ferrie, and Sanjiang Li. Ansatz-Agnostic Exponential Resource Saving in Variational Quantum Algorithms Using Shallow Shadows. In IJCAI, 2024.
[5] Harun Bayraktar, Ali Charara, and et al. Clark. cuQuantum SDK: A High-Performance Library for Accelerating Quantum Science. In QCE, 2023.
[6] Jacob D. Biamonte, Jason Morton, and Jacob W. Turner. Tensor Network Contractions for #SAT. Journal of Statistical Physics, 2015.
[7] Jacob D. Biamonte, Peter Wittek, Nicola Pancotti, Patrick Rebentrost, Nathan Wiebe, and Seth Lloyd. Quantum machine learning. Nat., 2017.
[8] Mark Blacher, Julien Klaus, Christoph Staudt, Sören Laue, Viktor Leis, and Joachim Giesen. Efficient and Portable Einstein Summation SQL. PACMMOD, 2023.
[9] Mark Blacher, Christoph Staudt, Julien Klaus, Maurice Wenig, Niklas Merk, Alexander Breuer, Max Engel, Sören Laue, and Joachim Giesen. Einsum Benehmark: Enabling the Development of Next-Generation Tensor Execution Engines. In NeurlPS, 2024.
[10] Diego Calanzone, Stefano Teso, and Antonio Vergari. Logically Consistent Language Models via Seuuo-Symboolic Integration. In ICLR, 2025.
[11] Lam Chi-Chung, Sadayappan, and Rephael Wenger. On optimizing a class of multi- dimensional loops with reduction for parallel execution. PPL, 1997.
[12] Daniel Cunnington, Mark Law, Jorge Lobo, and Alessandra Russo. Neuro-Symbolic Learning of Answer Set Programs from Raw Data. In IJCAI, 2023.
[13] G Daniel, Johnnie Gray, et al. Opt_ ein_einumm-ai python package for optimizing contraction order for einsum-like expressions. JOSS, 2018.
[14] Leonardo Dueñas-Osorin Kuldeep S. Meel, Roger Paredes, and Moshe Y. Vardi. Counting- based reliability estimation for power-transmission grids. AAAI, 2017.
[15] Gennaro Gala, Cassio P. de Campos, Antonio Vergari, and Erik Quaeghebeur. Scaling continuous latent variable models as probabilistic integral circuits. In NeurIPS, 2024.
[16] Artur García-Sáez and José Ignacio Latorre. An exact tensor network for the 3SAT problem. Quantum Inf. Comput., 2012.
[17] Zhangxiaowen Gong, Houxiang Ji, Christopher W. Fletcher, Christopher J. Hughes, and Josep Torrellas. SparseTraint Leveraging Dynamic Sparsity in Software for Training DNNs on General-Purpose SIMD Processors. In PACT, 2020.
[18] Johnnie Gray and Stefanos Kourtis. Hyper-optimized tensor network contraction. Quantum, 2021.

---

[19] Charles R. Harris, K. Jarrod Millman, Stéfan J. van der Walt, Ralf Gommers, Pauli Virtanen, David Couraapeau, Eric Wieser, Julian Taylor, Sebastian Berg, Nathaniel Smith, Robert Kern. et al. Array programming with NumPy. Nature, 2020.
[20] Anna L. D. Latour, Behrouz Babaki, Anton Dries, Angelika Kimmig, Guy Van den Broeck, and Siegfried Nijssen. Combining stochastic constraint optimization and probabilistic programming from knowledge compilation to constraint solving. In CP, 2017.
[21] Shiwei Liu, Tianlong Chen, Zahra Atashgahi, Xiaohan Chen, Ghada Sokar, Elena Mocanu, Mykola Pechenizkiy, Zhangyang Wang, and Decebal Constantin Mocanu. Deep Ensembling with No Overhead for either Training or Testing: The All-Round Blessings of Dynamic Sparsity. In ICLR, 2022.
[22] Igor L Markov and Yaoyun Shi. Simulating quantum computation by contracting tensor networks. SIAM Journal on Computing, 2008
[23] André Martins, Edoardo Ponti, Duarte Alves, Piotr Nawrot, and Saul Santos. Dynamic Sparsity in Machine Learning: Routing Information through Neural Pathways. NeurIPS Tutorial, 2024.
[24] Calarina Muslimani, Bram Grooten, Deepak Ranganatha Sastry Mamillapalli, Mykola Pech- enizkiy, Decebal Constantin Mocanu, and Matthew E. Taylor. Boosting Robustness in Preference-Based Reinforcement Learning with Dynamic Sparsity. CoRR, abs/2406.06495, 2024.
[25] Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, et al. Pytorch: An imperative style, performefoe deep learning library. In NeurlPS, 2019. 2019.
[26] PyData. Sparse einsum, 2018.
[27] Weijia Shi, Andy Shih, Adnan Darwiche, and Arthur Choi. On tractable representations of binary neural networks. In KR, 2020.
[28] Sahil Sidheckh and Sriraam Natarajan. Building Expressive and Tractable Probabilistic Genera- tive Models: A Review. In IJCAI, 2024.
[29] Daniel Silver, Tirthak Patel, and Devesh Tiwari. QUILT: Effective Multi-Class Classification on Quantum Computers Using an Ensemble of Diverse Quantum Classifiers. AAAI, 2022.
[30] Christoph Staudt, Mark Blacher, Julien Klaus, Farin Lippmann, and Joachim Giesen. Improved Cut Strategy for Tensor Network Contraction Orders. SEA, 2024.
[31] Nasib Ullah, Erik Schultheis, Mike Lasby, Yani Ioannou, and Rohit Babbar. Navigating Extremes: Dynamic Sparsity in Large Output Space. In NeurlPS, 2024.
[32] Antonio Vergari, YooJung Choi, Anji Liu, Stefano Teso, and Guy Van den Broeck. A composi- tional atlas of tractable circuit operations for probabilistic inference. In NeurlPS, 2021.
[33] Zhaofeng Wu, Ding Zhao, Qiao Liang, Jiahui Yu, Anmol Gulati, and Ruoming Pang. Dynamic Sparsity Neural Networks for Automatic Speech Recognition. In ICASSP, 2021.
[34] Lu Yin, Gen Li, Meng Fang, Li Shen, Tianjin Huang, Zhangyang Wang, Vlado Menkovski, Xiaolong Ma, Mykola Pechenizkiy, and Shiwei Liu. Dynamic Sparsity Is Channel-Level Sparsity Learner. In NeurlPS, 2023.
[35] Zhan Yu, Qiuhao Chen, Yuling Jiao, Yinan Li, Xiliang Lu, Xin Wang, and Jerry Zhijian Yang. Non-asymptotic approximation error bounds of parameterized quantum circuits. In NeurIPS, 2024.
[36] Ennan Zhai, Ang Chen, Ruzica Piskac, Mahesh Balakrishnan, Bingchuan Tian, Bo Song, and Haoliang Zhang. Check before you change: Preventing correlated failures in service updates. In NSDI, 2020.
[37] Bingyi Zhang and Viktor K. Prasanna. Dynasparse: Accelerating GNN Inference through Dynamic Sparsity Exploitation. In IPDPS, 2023.
[38] Pedro Zuidberg Dos Martires. Probabilistic neural circuits. In AAAI, 2024.
12

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
The abstract and introduction clearly state that the paper introduces a hybrid algorithm to exploit dynamic sparsity in einsum expressions, and supports this with both theoretical and empirical evidence.
Justification:
2. Limitations
Question:
Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: A dedicated Limitations section is included, discussing threshold sensitivity, lack of bidirectional format switching. limitations in backend optimization, and generaliza- tion to large-scale or mixed-density workloads.
3. Theory assumptions and proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [Yes]
Justification: All theoretical results (e.g., Lemma Lemma 2, Proposition 1) are stated with assumptions and justified by sketches or full proofs. A minor proof is deferred to the appendix.
4. Experimental result reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Answer: [Yes]
Justification: The paper describes the hardware, software, dataset used, and how experiments were run. All code and data are made available at the provided GitHub repository.
5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described in supplemental material?
Answer: [Yes]
Justification: The paper provides a GitHub link to a repository that includes both the code and instructions to run the experiments.
6. Experimental setting/details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
The experimental section specifies hardware, datasets, number of runs, median reporting, and other relevant settings. No learning parameters are involved, as the method is non-parametric.
Justification:
7. Experiment statistical significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [No]
Justification: Error bars are not included. Instead, results are reported using medians over five runs. Variability is low and the results are dominated by algorithmic complexity, not stochasticity.
13

---

8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed reproduce the experiments?
Question:
Answer: [Yes]
Justification: Section 4.1 describes the compute environment: 64-core Xeon server with 1.5TB RAM, and mentions that runtimes are in the order of seconds to hours, with outlier handling discussed.
9. Code of ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https: httttp:://nuuris..co http://ppp....pp..b../aaaacaaaaa
Answer: [Yes]
Justification: The research is entirely algorithmic and does not involve humans, personal data, or decision-making systems with ethical implications.
10. Broader impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [NA]
Justification: The work is foundational and methodological, focused on tensor expression evaluation. It does not directly interface with users society in way that warrants broader impact discussion.
11. Safeguards
Does the paper deseribe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, scraped datasets)?
Question:
Answer: [NA]
Justification: No such high-risk assets are released. The code is low-risk algorithmic implementation for tensor contraction.
12. Licenses for existing assets
Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Answer: [Yes]
Justification: All datasets and libraries (c.g., PyTorch, einsum benchmark) used in the paper are publicly available and properly cited.
13. New assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [Yes]
Justification: The provided GitHuu repository includes our sparse einsum implementation and instructions for reproducing all benchmarks and results.
14. Crowdsourcing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
Answer: [NA]
Justification: The paper does not involve any human subjects or crowdsourcing.
14

---

15. Institutional review board (IRB) approvals or equivalent for research with human subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/review based the requirements of your country institution) were obtained?
Answer: [NA]
Justification: The research does not involve human subjects and thus does not require IRB approval.
16. Declaration of LLM usage
Question: Does the paper describe the usage of LLMs if it is an important, original, or non-standard component of the core methods in this research?
Answer: [NA]
Justification: No large language models were used in the research methodology or con- tributed to any part of the core method development.
15


Reviews:
[{"review_id": "qxZrUA0yvs", "reviewer": "Reviewer_pjM9", "summary": "This paper talked about the limitations and constrains in tensor computing, where they modify tensor hypernetworks to speed up the computing process. The proposed method showed a great performance increase across tensors with different densities.", "questions": "Please see previous section.", "limitations": "Yes", "rating": 4, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "Strengths:\n1. This paper discussed the constrains and limitations in (sparse) tensor computations.\n2. They modified the tensor hypernetworks to speed up the computation. Such modification is suitable for tensors with varied densities.\n3. They provide a small C++ module plugs into PyTorch/opt_einsum for CPU capability.\n\nWeaknesses:\n1. One limitation I noticed and the author pointed out in the limitation section is how would this method work on GPU. Most tensor computation is utilizing GPU. If the proposed method can only work on CPU. It could greatly limit the applicability of the proposed method.\n2. Can the proposed method support sparse to dense fallback?\n3. How sufficient is the 5% threshold. How sensitive are results to this value?", "was_revised": true, "final_justification": "This paper brings efficient solutions in computing. I will keep the score as 4.", "rebuttal": "Thank you for reviewing our paper. We appreciate the opportunity to address the weaknesses that you have mentioned.\n\n**Weakness 1:** *One limitation I noticed and the author pointed out in the limitation section is how would this method work on GPU. Most tensor computation is utilizing GPU. If the proposed method can only work on CPU. It could greatly limit the applicability of the proposed method.*\n\nWe want to distinguish between the phenomenon of dynamic sparsity in *einsum* and the implementations that exploit it. Dynamic sparsity arises independently of the compute backend, that is, CPU or GPU. The proposed method would work on GPUs as well, but so far there are no established sparse *einsum* implementations for GPUs.\n\n**Weakness 2:** *Can the proposed method support sparse to dense fallback?*\n\nYes, this is possible. In this work, we kept the setup simple and focused on analyzing the phenomenon of dynamic sparsity in *einsum*. Nevertheless, our motivation is practical, and we will explore all options as we further engineer an efficient sparse *einsum* implementation that exploits dynamic sparsity.\n\n**Weakness 3:** *How sufficient is the 5% threshold. How sensitive are results to this value?*\n\nIn general, the threshold should not be chosen too large, because, as can be seen in Figure 3, sparsity can first decrease and then increase again, which can result in large performance overheads depending on the performance of the sparse *einsum* implementation. Switching later than necessary because of a threshold that is too small incurs only small overheads as long as the switch happens before the large, sparse contractions that typically dominate the end of the computations. On the benchmark, the optimal value is around 5%. Setting it as low as 1% would affect fewer than ten instances, all with minor sparse speed-ups."}, {"review_id": "4Vk7JRorZc", "reviewer": "Reviewer_FHEA", "summary": "This paper proposed a hybrid dense-sparse Einsum for computing large Eisnum expressions. The authors show that there exist non-trivial examples of Einsum programs, in the form of ${\\rm Grid}_n$, where the popular dense implementation is much slower than a coordinate-list-based sparse implementation. Through experiments, the authors show that their proposed hybrid dense-sparse Einsum is much performant that the baseline dense implementation.", "questions": "- It'd be better to spell out, in pseudocode, the algorithms for contraction for both sparse and dense tensors.\n- Specify what the ${\\rm Grid}$ model is used for. E.g., for modeling magnetic crystal structures? Given more context here would aid understanding of the setting of the method.\n- L147: ${\\rm Grid}_n$: Provide a large illustration / figure of this tensor network.\n- L148: Please add the definition $[n] = \\\\{i \\in \\mathbb{N} | 0 \\le i < n\\\\}$.", "limitations": "yes", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 3, "significance": 2, "originality": 2, "strengths_and_weaknesses": "Strengths:\n- Simple method that keeps track of the proportion of non-zero elements in tensors; and switch to sparse when below a threshold. This method is simple and performs well.\n- Detailed exposition on the ${\\rm Grid}_n$ tensor network.\n\nWeaknesses:\n- The background is too long in the paper, where the proposed solution appeared at Section 4.4, page 8. The hybrid algorithm is very simple (if non-zero entries >5%, dense). The substance of this paper does not seem to merit a long NeurIPS submission.\n- The proposed method works for ${\\rm Grid}$-type tensor networks, but actual tensor Einsum programs take many different forms (in deep learning; belief propagation; quantum computing; etc). The performance might not generalize to other forms of tensor networks. A good benchmark for this is the Einsum benchmark (https://benchmark.einsum.org/): https://proceedings.neurips.cc/paper_files/paper/2024/file/b1bbfdb9197bfc819a52c34dce493f85-Paper-Datasets_and_Benchmarks_Track.pdf, which I believe that the authors should do experiments on.\n- Tensor operations are now predominantly performed in a parallel way on specialized hardwares such as GPUs or TPUs. Analysis should be done over these heterogeneous hardware.", "was_revised": true, "final_justification": "I've read the comment and had a better understanding of the main contributions of the paper. I'm happy to increase the scores of the paper.", "rebuttal": "Thank you for reviewing our paper. Unfortunately, there are some inaccuracies in your review that we address below. We hope you take this into account and update your score accordingly.\n\n**Weakness 1:** *The proposed method works for Grid-type tensor networks, but actual tensor Einsum programs take many different forms (in deep learning; belief propagation; quantum computing; etc). The performance might not generalize to other forms of tensor networks. A good benchmark for this is the Einsum benchmark, which I believe that the authors should do experiments on.*\n\nAs we write in the second sentence of the experimental section (Section 4), our experiments were performed on this benchmark, which covers several application domains. \n\n**Weakness 2:** *The background is too long in the paper, where the proposed solution appeared at Section 4.4, page 8. The hybrid algorithm is very simple (if non-zero entries >5%, dense). The substance of this paper does not seem to merit a long NeurIPS submission.*\n\nWe do not consider the hybrid algorithm the main contribution of our paper, but rather the theoretical and practical analysis of dynamic sparsity in *einsum*. Therefore, we do not consider Section 3 as background material but as a contribution. The Grid example in Section 3 shows that there exist *einsum* expressions where the two data structures, dense and sparse, are exponentially separated, that is, the dense data structure dynamically becomes exponentially larger than the sparse one, whereas the overhead of sparse over dense is always polynomially bounded. Then, in Section 4, we show on the *einsum* benchmark that a separation can also be observed on practical *einsum* problems. Moreover, we show that dynamic sparsity can be exploited to speed up the evaluation of *einsum* expressions.     \n\n**Weakness 3:** *Tensor operations are now predominantly performed in a parallel way on specialized hardwares such as GPUs or TPUs. Analysis should be done over these heterogeneous hardware.*\n\nWe want to distinguish between the phenomenon of dynamic sparsity in *einsum* and implementations that exploit this phenomenon. Dynamic sparsity arises independently of the compute backend, that is, CPU, GPU, or TPU. The proposed method would work on TPUs and GPUs as well, but so far there are no established sparse *einsum* implementations for GPUs and TPUs. \n\n**Question/Comment 1:** *It'd be better to spell out, in pseudocode, the algorithms for contraction for both sparse and dense tensors.*\n\nContraction orders break down an *einsum* expression into a sequence of binary contractions, that is, contractions of two tensors. Implementationd of binary contraction algorithms in either case, dense and sparse, are highly advanced and still an active area of research. They are not in the focus of our paper.\n\n**Question/Comment 2:**  *Specify what the model is used for. E.g., for modeling magnetic crystal structures? Given more context here would aid understanding of the setting of the method.*\n\nThe Grid model is an instance of a #SAT (model counting) problem. Model counting has applications, among others, in probabilistic inference. For instance, the inference engine of the probabilistic programming language DICE uses model counting in its backend. However, here, the Grid model only serves as an example of a tensor contraction problem that exhibits dynamic sparsity. It serves as a formal proof that the difference between sparse and dense representations can be exponential. \n\n**Question/Comment 3:** *L147: Provide a large illustration / figure of this tensor network*\n\nThe Grid$_3$ model is illustrated in Figure 1.  \n\n**Question/Comment 4:** *L148: Please add the definition of $[n]$.*\n\nWe will do so."}, {"review_id": "tXdRLGokXe", "reviewer": "Reviewer_hrqu", "summary": "This paper evaluates the benefits of sparse tensor storage for einsum computations.  \nOn the theoretical side, the authors prove that for a well-chosen family of instances, dense storage can require exponentially more FLOPs for tensor contraction than sparse storage.  \nOn the practical side, they perform an analysis of sparsity evolution as contraction progresses, for a recently-introduced benchmark ([https://benchmark.einsum.org/](https://benchmark.einsum.org/)). This analysis shows that sparsity can sometimes appear halfway through the process, even when the input tensors are dense. Therefore, a new storage strategy is proposed, which switches from dense to sparse when the average density of the tensor sequence falls below a pre-specified threshold. Numerical experiments demonstrate that this strategy can lead to efficiency gains.", "questions": "- To justify a dynamic approach, the authors state that “in einsum expressions involving a large number of tensors, it is not known in advance whether intermediate tensors will become sparse.” Is that obvious? Are there heuristics that could guide a static storage choice?\n- Describe the contraction orders used in the experiments and their impact on the storage and FLOPs requirements.\n- Explain more precisely how the “not-yet contracted tensors” are defined and how their average density is calculated.\n- Add a comparison between the hybrid algorithm and the best of both dense and sparse algorithms, in order to back the claim from the abstract.\n- Defend your choice of sesum as a sparse einsum library, as opposed to the Python sparse library you mention in the introduction.\n\n---\n\nSatisfactory answers to most of these questions can lead me to increase my rating to accept.", "limitations": "Yes.", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 3, "significance": 2, "originality": 3, "strengths_and_weaknesses": "## Strengths\n\n### Quality\n\n- The separation theorem on model counting for logical formulas is a strong theoretical motivation for sparse storage.\n- The measurement of dynamic sparsity in the Einsum Benchmark provides empirical justification for an adaptive approach.\n\n### Clarity\n\n- The paper is well-written.\n- The experiments are reproducible.\n\n### Significance\n\n- Tensor contraction is a thriving topic in probabilistic inference, and this study opens up several promising research avenues, especially in Section 5.\n\n### Originality\n\n- The focus on CPU operations is unexpected in machine learning but justified for handling dynamic sparsity.\n- The adaptive storage policy is rather naive but it is possible that no one has benchmarked it before.\n\n## Weaknesses\n\n### Quality\n\n- The contraction order used in the experiments is never discussed in the main paper. In the appendix, we only read that “The benchmark includes one contraction order optimized for minimal intermediate size and one for minimal floating-point operations (flops).”\n- The abstract claims that “In our experiments on established benchmark einsum expressions, the hybrid solution is consistently competitive with or outperforms the better of the two static representations.”, but the hybrid algorithm is only compared to its purely dense and sparse counterparts *separately* (Figure 5). It should also be compared to the best of both for each instance.\n\n### Clarity\n\n- It is not clear to me how the hybrid algorithm defines “the average density of *the not-yet contracted tensors*”. After all, a tensor could be involved at several steps of the contraction order (which is defined on axes, if I understand correctly). To put it differently, contractions after step 1 can be defined on intermediate tensors which do not yet exist at the beginning of the process.\n\n### Significance\n\n- In the experiment code, a sparse tensor library called sesum (which I wasn’t able to find elsewhere) is apparently used as a binary wheel, which makes the experiments less reusable and generalizable.", "was_revised": true, "final_justification": "Raising my score to 5 following a convincing rebuttal, to make it clear that this paper deserves acceptance despite the concerns raised about GPU compatibility.", "rebuttal": "Thank you for reviewing our paper and your questions. We appreciate the opportunity to address them.\n\n**Question 1:** *To justify a dynamic approach, the authors state that “in einsum expressions involving a large number of tensors, it is not known in advance whether intermediate tensors will become sparse.” Is that obvious? Are there heuristics that could guide a static storage choice?*\n\nUnfortunately, this is not obvious. Sparsity in *einsum* is a truly dynamic property, and all the following are possible: 1. Sparse input tensors can lead to dense intermediate tensors. 2. Dense input tensors with only a few non-zero entries can lead to intermediate and result tensors that are hypersparse. 3. Sparse intermediate tensors can lead to dense tensors later in the contraction order. \n\nSparsity estimation techniques such as cardinality estimation, which are used in the context of databases, tend to be expensive and unreliable for large expressions. To the best of our knowledge, there are no sparsity estimation techniques specifically for *einsum*. Therefore, sparsity estimation for einsum is a promising avenue for future research. We will address this in the paper.\n\n**Question 2:** *Describe the contraction orders used in the experiments and their impact on the storage and FLOPs requirements.*\n\nWe use the contraction orders that are provided with the *einsum* benchmark, Reference [9] in the paper. We summarize key metrics of the contraction paths in Table 2 of the supplement. Regarding impact, dynamic sparsity is contraction order dependent. From our observations, flop-optimized contraction orders tend to result in sparser tensors than size-optimized orders. Therefore, we used flop optimized paths, except for problems that become too large to fit into RAM in the dense evaluation case. \n\n**Question 3:** *Explain more precisely how the “not-yet contracted tensors” are defined and how their average density is calculated.*\n\nAny tensor, input or intermediary, is contracted only once, because a contraction replaces the participating tensors by a new intermediary tensor. This can be explained using the example that is shown in Figure 1 in the paper. Before the evaluation of the *einsum* expression starts, twelve tensors (black circles) that are shown in Figure 1(left) are all not-yet contracted. Then, after the first contraction, shown in Figure 1(right), the four tensors at the center are contracted into an intermediate tensor resulting in nine tensors: eight original tensors and one new intermediate tensor. These nine tensors are not-yet contracted after the first contraction. Therefore, the density of the not-yet contracted tensors before the first contraction is the average density of the input tensors and the average density of the not-yet contracted tensors after the first contraction is the average density of the remaining nine tensors. Thus, the only new density that has to be computed after the first contraction is the density of the new intermediate tensor. The average density can then be updated accordingly. We will clarify this in the paper.  \n\n**Question 4:** *Add a comparison between the hybrid algorithm and the best of both dense and sparse algorithms, in order to back the claim from the abstract.* \n\nBesides the outlier instances, the comparison of hybrid vs the best of sparse and dense boils down to the plots of sparse instances using the sparse backend and dense instances using the dense backend. We will include a figure that shows this explicitly in supplemental material. Overall, the hybrid algorithm is only more than twice as slow as the best in eight out of 158 instances. In the median, we see a slowdown of roughly 1.02x, so about 2 percent slower, than the best of both. Notably, the hybrid approach can also be significantly faster than the best of sparse and dense. In our experiments, that happened on three out of 28 sparse benchmark problems. We will add a figure to the supplement to highlight this.\n\n**Question 5:** *Defend your choice of sesum as a sparse einsum library, as opposed to the Python sparse library you mention in the introduction.*\n\nIn the Python sparse library, each index is encoded in a separate integer. It turns out that this slows the evaluation of *einsum* expressions with tensors that have many indices. In the *sesum* implementation of *einsum*, the indices of a tensor are encoded in a single number, which leads to a smaller memory footprint and faster execution. \n\nEven on sparse instances with 1% non-zeros, such as *mc_2022_029* or *mc_rw_blasted_case1_b14_even3* from the benchmark, the Python sparse library can be significantly slower than PyTorch's dense *einsum* implementation.\n\nFor the reproducibility of our experiments, *sesum* is included in the project repository. Running `pip install sesum-0.3.8-py3-none-any.whl` or `uv sync` will make it available on any Linux or macOS machine."}, {"review_id": "KTcJzeUDCq", "reviewer": "Reviewer_DRa5", "summary": "The paper introduces the problem of choosing the correct data structure when evaluating einsum expressions on sparse data, particularly between tensors and coordinate arrays. They show on a family of grid formulae that the latter uses exponentially more flops than the former when contracting the representation of a grid. Then, experiments are performed to empirically support this.\nFurthermore, they propose an hybrid algorithm that switches the representation on the fly, showing positive results.", "questions": "- Is it possible to apply the hybrid approach to the actual training of a neural network? For example to a simple CNN on MNIST to see the benefits regarding training time and resource use.\n- Are there any related works that discuss the efficiency of structures on neural networks? It would be good to see if other approaches were considered before.", "limitations": "Yes, in Section 5.", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 2, "strengths_and_weaknesses": "Strengths\n\n- The paper addresses an important point about the efficiency of the structures we use during long routines of computation, which can lead to dramatic improvement in time and energy consumption.\n- The empirical results are founded in a proved proposition and align with them.The benchmark of einsum instances seems broad and has been recently published.\n- Well written and easy to read.\n\nWeaknesses\n\n- I find it strange that there is no related work. Although the approach seems completely novel, there must be research that speaks about transitions of structures with respect to some threshold. Adding some about this would improve to contextualize the originality of the approach and also to improve my confidence on the score.\n- Some of the limitations could be addressed in the main article such as switching back to a dense format, which should be relatively easy to implement, or experiments on GPUs, which are essential for high performance computation in the context of AI.\n\n\nMinor comments\n\n- Figure 2. Grid Site Length (n) -> Size. Also, why is dense not present on the right subfigure? I think the subfigures should be better explained.\n- Figure 3. Legend missing.", "was_revised": true, "final_justification": "The paper offers a principled analysis of dynamic sparsity in einsum evaluation, combining theoretical insight with a practical hybrid implementation. The method is well-motivated, and the experimental results align with the theoretical claims. I appreciate the clarifications provided in the rebuttal, especially regarding backend generality and use cases in neuro-symbolic reasoning. I maintain my positive score and increase my confidence.", "rebuttal": "Thank you for reviewing our paper and your questions. We appreciate the opportunity to address them.\n\n**Question 1:** *Is it possible to apply the hybrid approach to the actual training of a neural network? For example to a simple CNN on MNIST to see the benefits regarding training time and resource use.*\n\nThe dynamic switching between dense and sparse data structures based on average density that we introduce here is potentially also useful for multi-tensor operations used in CNNs. However, the benefits of dynamic switching become more pronounced for *einsum* expressions that involve many tensors. \n\n**Question 2:** *Are there any related works that discuss the efficiency of structures on neural networks? It would be good to see if other approaches were considered before.* \n\nDynamic sparsity is a known phenomenon in deep learning, see References  [17, 21, 24, 31, 33, 34, 37] in the paper. Traditionally, however, deep learning employs *einsum* expressions with only a few tensors. Only recently have large *einsum* expressions been used in machine learning areas such as probabilistic circuits, probabilistic neural circuits, and neuro-symbolic models that combine deep learning with logical constraints. Therefore, to the best of our knowledge, we are the first to study dynamic sparsity in large *einsum* expressions."}]

Author rebuttal:

Meta-review: The paper investigates dynamic sparsity in einsum computations, demonstrating a potential exponential gap between dense and sparse representations, and proposes a simple hybrid algorithm that switches between the two representations when density falls below a threshold. Experiments demonstrate that the hybrid solution is consistently competitive with, or outperforms, the better of the two static approaches.

The work offers a novel and well-founded treatment of dynamic sparsity as a hardware-agnostic property, combining theory with reproducible experiments. The authors convincingly addressed critiques on scope, writing, and related work during the rebuttal. Overall, the contribution is technically sound, relevant across domains, and is therefore recommended for acceptance.
Decision: Accept (poster)