# How to find all multi-indices satisfying given contraints?

Posted 9 years ago
5701 Views
|
10 Replies
|
5 Total Likes
|
 Hey there,First, in case someone isn't familiar with multi-indices, they are essentially lists of non-negative integers, e.g. $\alpha = (\alpha_1,\ldots,\alpha_n).$ Further, denote $\left| \alpha \right|_1 = \sum_{i=1}^n \left|\alpha_i\right|$. What I would like is some efficient way of finding all multi-indices for a given $n$ satisfying $$a \leq \left| \alpha \right|_1 \leq b,$$with the additional constraint that $\alpha_i > 0$ for each $1 \leq i \leq n$. The integers $a$ and $b$ are of course given.More specifically, I need to sum over all such multi-indices (passing each as an argument to a function), but that's simple enough once I have a list of the multi-indices.-Jonatan
10 Replies
Sort By:
Posted 9 years ago
 You could just have a1 go from 1 to b-n+1, a2 from 1 to (b-n+1)-a1, a3 from 1 to (b-n+1)-(a1+a2), etc.One can use Reduce to find all actual multiindex sets, but that 's overkill. Likewise, could look at exponents in the expansion of (1+a1+a2+...+an)^(b-a)*(a1+...+an)^a. Again, overkill.
Posted 9 years ago
 Just table it out With[{b= ...., n = ...}, allowedInd = Table[{a1,a2,a3,...}, {a1,1,b-n+1},{a2, 1, b-n+2-a1},{a3, 1, b-n+3-a1-a2},...] ] and use it afterwards. Be careful to formulate correct Table-iterators.Another simple way is to include If-statements into the generic summand term, giving 0 if the multi-index is not allowed to contribute.
Posted 9 years ago
 Thank you, Daniel and Udo, that seems to work. I do have a few issues with this solution though, I hope there's some easy fix for these.Firstly, the constraint $a \leq \left|\alpha\right|_1$ is missing from these solutions. As you say, this is fixable with a simple If-statement, and that's fine, however I can't see a way of doing this without leaving nulls/zeros/empty sets in the output. Secondly, I think I failed to mention this earlier but ideally I'd like the output to be as a list of lists, such as the output from David's suggested method.I'm sure both of these problems are easily fixable, I'm just new to Mathematica so I don't yet know how to do it easily. My first thought was Flatten, but then I would need to know for sure how many levels the output list has; also, that still leaves the task of weeding out the nulls/whatnots.Finally, there's the issue that the suggested method only works if $n$ is fixed -- I would prefer to implement a function which simply takes $a, b$ and $n$ as input. Again, this is mostly a case of my question being unclear, sorry for that.Thank you for the help so far. =)-JonatanP.S. As a side note, in case someone else ends up needing this for something: the appropriate upper limits for the iterators are b-n+1, b-n+2-a1, b-n+3-a1-a2, and so forth. Otherwise, if we have $b = 5$ and $n = 3$, say, then the last limit becomes at most $5 - 3 + 1 - 1 - 1 = 1$, when it should be 3 in the case where $\alpha_1 = \alpha_2 = 1$.
Posted 9 years ago
 Lower conditions as well as the whole job can be done by brute force up to a certain n In[41]:= With[{b = 7, n = 5, a = 1}, If[a < 0 || n < 0 || b < 1 || b - n < a, Print["You will get the empty set!"]]; If[a < n, Print["Your conditions run empty on lower boundary!"]]; Union[Flatten[ DeleteMissing[ Table[If[a <= Plus @@ {a1, a2, a3, a4, a5} <= b, X[a1, a2, a3, a4, a5], Missing[]], {a1, 1, b - 1}, {a2, 1, b - 1}, {a3, 1, b - 1}, {a4, 1, b - 1}, {a5, 1, b - 1}], n]]] /. X -> List ] During evaluation of In[41]:= Your conditions run empty on lower boundary! Out[41]= {{1, 1, 1, 1, 1}, {1, 1, 1, 1, 2}, {1, 1, 1, 1, 3}, {1, 1, 1, 2, 1}, {1, 1, 1, 2, 2}, {1, 1, 1, 3, 1}, {1, 1, 2, 1, 1}, {1, 1, 2, 1, 2}, {1, 1, 2, 2, 1}, {1, 1, 3, 1, 1}, {1, 2, 1, 1, 1}, {1, 2, 1, 1, 2}, {1, 2, 1, 2, 1}, {1, 2, 2, 1, 1}, {1, 3, 1, 1, 1}, {2, 1, 1, 1, 1}, {2, 1, 1, 1, 2}, {2, 1, 1, 2, 1}, {2, 1, 2, 1, 1}, {2, 2, 1, 1, 1}, {3, 1, 1, 1, 1}} if you want multi-indices disregarding the order, it is simply In[42]:= With[{b = 7, n = 5, a = 1}, If[a < 0 || n < 0 || b < 1 || b - n < a, Print["You will get the empty set!"]]; If[a < n, Print["Your conditions run empty on lower boundary!"]]; Union[Sort /@ Flatten[DeleteMissing[ Table[If[a <= Plus @@ {a1, a2, a3, a4, a5} <= b, X[a1, a2, a3, a4, a5], Missing[]], {a1, 1, b - 1}, {a2, 1, b - 1}, {a3, 1, b - 1}, {a4, 1, b - 1}, {a5, 1, b - 1}], n]]] /. X -> List ] During evaluation of In[42]:= Your conditions run empty on lower boundary! Out[42]= {{1, 1, 1, 1, 1}, {1, 1, 1, 1, 2}, {1, 1, 1, 1, 3}, {1, 1, 1, 2, 2}} this last thing brings you to an efficient way to do it: Start with the trivial allowed lower bound multiindex {1,1,1,..., Max[1,a-n+1]}, raise the last position until the limit b is reached ({1,1,1,...,b-n+1}) now distribute the IntegerPartitions/@ Range[Max[1,a-n+1],b-n+1-Max[1,a-n+1]] over the trivial allowed lower bound multiindex, with other words In[51]:= With[{b = 7, n = 5, a = 1}, If[a < 0 || n < 0 || b < 1 || b - n < a, Print["You will get the empty set!"]]; If[a < n, Print["Your conditions run empty on lower boundary!"]]; Join[{{1, 1, 1, 1, Max[1, a - n + 1]}}, {1, 1, 1, 1, Max[1, a - n + 1]} + PadLeft[#, n] & /@ Flatten[IntegerPartitions /@ Range[Max[1, a - n + 1], b - n + 1 - Max[1, a - n + 1]], 1]] ] During evaluation of In[51]:= Your conditions run empty on lower boundary! Out[51]= {{1, 1, 1, 1, 1}, {1, 1, 1, 1, 2}, {1, 1, 1, 1, 3}, {1, 1, 1, 2, 2}} if they are needed with respect to order, the result must be shuffled back to un-order.Check for errors please, and find possibly other, even more efficient ways to do it!
Posted 9 years ago
 Thank you, that got me just what I needed! There was some bug fixing needed, but that was quite trivial once I had the right idea from your code. This is what I finally ended up with: With[{b=10,n=5,a=8}, If[a < 0 || n < 0 || b < 1 || b < a, Print["You will get the empty set!"]]; If[a < n, Print["Your conditions run empty on lower boundary!"]]; Flatten[Permutations /@ (ConstantArray[1,n] + PadLeft[#,n]& /@ Flatten[IntegerPartitions[#,n]& /@ Range[Max[0, a-n], b-n], 1]), 1] ] Basically, your code worked fine when $a < n$, so all I had to do was fix the Range-call and replace the "initial vector" with all ones, and everything worked perfectly. As a nice side-effect, replacing the initial vector with a call to ConstantArray removed the only place in the function with a hard-coded value for $n$, allowing it to work for any $n$ given as input.Thanks again, your help was invaluable. =)-JonatanP.S. Since the solution ended up taking quite a different form, Daniel's solution wasn't needed in the end, but thank you for that as well nonetheless, it was nice to learn how I might go about doing it if I find myself needing it later on.P.P.S. I tested the efficiency, and with $(a,b,n) = (10,20,5)$, your "brute force" solution took 15.7 seconds to run, the code above took 0.02 seconds. Very nice improvement, if I had gone about solving this myself I doubt I would've found anything better than the brute force solution with my limited knowledge.
Posted 9 years ago
 Sorry if I come back to this, but the Parma Polyhedron Library might be of interest in cases where IntegerPartitions[] might take too big memory requirements (usual around b - n = 80 in the language of the current problem on a 8 GB RAM computer). I tried to fit a Henrich package Multidimensional Convex Polyhedra and Searching for Magic Squares to this multi-index problem, but all I fiddled out was the null space, possibly because I'm unable to define an equation in this case, only inequalities are available.
Posted 9 years ago
 An interesting addition, but in my work computing the multi-indices will never actually be the bottleneck, since both the computation time and memory usage for the remainder of my program is exponential in $b$.
Posted 9 years ago
 This might be what you have in mind. indicesFromLowToHigh[l_, h_, n_, var_] := Module[ {vars = Array[var, n]}, Table[{var[j], 1, h - n + j - Sum[var[k], {k, j - 1}]}, {j, n}] ] Example:indicesFromLowToHigh[2, 7, 3, x] Out[361]= {{x[1], 1, 5}, {x[2], 1, 6 - x[1]}, {x[3], 1, 7 - x[1] - x[2]}} For actual use as a sequence of iterators one would have something like Sequence@@indicesFromLowToHigh[...].
Posted 9 years ago
 This disrespects the lower bound if a > n In[9]:= (* Lichtblau *) Remove[indicesFromLowToHigh] indicesFromLowToHigh[l_, h_, n_, var_] := Module[{vars = Array[var, n]}, Table[{var[j], 1, h - n + j - Sum[var[k], {k, j - 1}]}, {j, n}]] In[12]:= indicesFromLowToHigh[8, 10, 5, x] Out[12]= {{x[1], 1, 6}, {x[2], 1, 7 - x[1]}, {x[3], 1, 8 - x[1] - x[2]}, {x[4], 1, 9 - x[1] - x[2] - x[3]}, {x[5], 1, 10 - x[1] - x[2] - x[3] - x[4]}} In[18]:= Select[ Flatten[Table[X[x1, x2, x3, x4, x5], {x1, 1, 6}, {x2, 1, 7 - x1}, {x3, 1, 8 - x1 - x2}, {x4, 1, 9 - x1 - x2 - x3}, {x5, 1, 10 - x1 - x2 - x3 - x4}]] /. X -> List, (Plus @@ # < 8 || Plus @@ # > 10) &] Out[18]= {{1, 1, 1, 1, 1}, {1, 1, 1, 1, 2}, {1, 1, 1, 1, 3}, {1, 1, 1, 2, 1}, {1, 1, 1, 2, 2}, {1, 1, 1, 3, 1}, {1, 1, 2, 1, 1}, {1, 1, 2, 1, 2}, {1, 1, 2, 2, 1}, {1, 1, 3, 1, 1}, {1, 2, 1, 1, 1}, {1, 2, 1, 1, 2}, {1, 2, 1, 2, 1}, {1, 2, 2, 1, 1}, {1, 3, 1, 1, 1}, {2, 1, 1, 1, 1}, {2, 1, 1, 1, 2}, {2, 1, 1, 2, 1}, {2, 1, 2, 1, 1}, {2, 2, 1, 1, 1}, {3, 1, 1, 1, 1}} but if they become deleted afterwards, it does fit In[19]:= (* Lehtonen *) lehto = With[{b = 10, n = 5, a = 8}, If[a < 0 || n < 0 || b < 1 || b < a, Print["You will get the empty set!"]]; If[a < n, Print["Your conditions run empty on lower boundary!"]]; Flatten[Permutations /@ (ConstantArray[1, n] + PadLeft[#, n] & /@ Flatten[IntegerPartitions[#, n] & /@ Range[Max[0, a - n], b - n], 1]), 1] ]; In[26]:= lichtb = DeleteCases[ Flatten[Table[ If[8 <= Plus @@ {x1, x2, x3, x4, x5} <= 10, X[x1, x2, x3, x4, x5], Y[x1, x2, x3, x4, x5]], {x1, 1, 6}, {x2, 1, 7 - x1}, {x3, 1, 8 - x1 - x2}, {x4, 1, 9 - x1 - x2 - x3}, {x5, 1, 10 - x1 - x2 - x3 - x4}]], _Y] /. X -> List; In[28]:= Sort[lehto] == Sort[lichtb] Out[28]= True so one could call this subtle brute force.
Posted 9 years ago
 Have you tried Subsets (if you want the indices in order) or Tuples (for the more general case)? You may then have to add a shift.The most common case is to pick k indices from n in ascending order. Say 3 indices from 5: Subsets[Range[5], {3}] Giving: {{1, 2, 3}, {1, 2, 4}, {1, 2, 5}, {1, 3, 4}, {1, 3, 5}, {1, 4, 5}, {2, 3, 4}, {2, 3, 5}, {2, 4, 5}, {3, 4, 5}}