How to generate a multidimensional array of indexed variables in Maxima?

Question

I want to take a list of non-negative integers D=[d1,...,dm] and and generate a multidimensional array of indexed symbols A in the form of:

where 0<=i_j<=d_j. For example if D=[2,3] then A should be

[[a_[0,0],a_[0,1],a_[0,2]],
 [a_[1,0],a_[1,1],a_[1,2]]]

For this case I could nest two for loops to generate the said array, however D does not necessarily have a length of 2 and I don't know how to nest an arbitrary number of for loops!

I would appreciate if you could help me know how I can generate A from D.

P.S. What I want to finally achieve is to create a multivariate polynomial as explained here.

@Eelvex No I have not. Would you be so kind to share an example below in a post? — Foad S. Farimani, Aug 03 '18 at 20:22
@Eelvex I tested this. this is amazing. Is it also possible to use for generating the X tensor as I have explained [here](https://stackoverflow.com/questions/51628270/how-to-generate-a-symbolic-multivariate-polynomial-of-a-given-dimension-in-sympy)? — Foad S. Farimani, Aug 03 '18 at 20:33
I suppose you could use `genmatrix(lambda ([i,j], x[1]^i*x[2]^j), 2, 3, 0, 0) * genmatrix(a, 2, 3, 0, 0))` and then `list_matrix_entries(%)` and `apply("+", %)`. I don't know if that would be the best way to do it. — Eelvex, Aug 04 '18 at 14:34
@Eelvex Great. I was wondering if you would write this is a post below? — Foad S. Farimani, Aug 06 '18 at 07:01

Robert Dodier · Accepted Answer · 2018-08-02T16:49:23.407

Here's one way to do it. The essential part is that I called cartesian_product to construct the list of all combinations of indices, and then arrayapply to create the subscripted expressions.

(%i11) ii:setify(makelist(i, i, 0, n)), n=2;
(%o11)                      {0, 1, 2}
(%i12) apply (cartesian_product, makelist (ii, m)), m=3;
(%o12) {[0, 0, 0], [0, 0, 1], [0, 0, 2], [0, 1, 0], [0, 1, 1], 
[0, 1, 2], [0, 2, 0], [0, 2, 1], [0, 2, 2], [1, 0, 0], 
[1, 0, 1], [1, 0, 2], [1, 1, 0], [1, 1, 1], [1, 1, 2], 
[1, 2, 0], [1, 2, 1], [1, 2, 2], [2, 0, 0], [2, 0, 1], 
[2, 0, 2], [2, 1, 0], [2, 1, 1], [2, 1, 2], [2, 2, 0], 
[2, 2, 1], [2, 2, 2]}
(%i13) map (lambda ([l], arrayapply (_a, l)), %);
(%o13) {_a       , _a       , _a       , _a       , _a       , 
          0, 0, 0    0, 0, 1    0, 0, 2    0, 1, 0    0, 1, 1
_a       , _a       , _a       , _a       , _a       , 
  0, 1, 2    0, 2, 0    0, 2, 1    0, 2, 2    1, 0, 0
_a       , _a       , _a       , _a       , _a       , 
  1, 0, 1    1, 0, 2    1, 1, 0    1, 1, 1    1, 1, 2
_a       , _a       , _a       , _a       , _a       , 
  1, 2, 0    1, 2, 1    1, 2, 2    2, 0, 0    2, 0, 1
_a       , _a       , _a       , _a       , _a       , 
  2, 0, 2    2, 1, 0    2, 1, 1    2, 1, 2    2, 2, 0
_a       , _a       }
  2, 2, 1    2, 2, 2
(%i14) grind (%);

{_a[0,0,0],_a[0,0,1],_a[0,0,2],_a[0,1,0],_a[0,1,1],_a[0,1,2],
 _a[0,2,0],_a[0,2,1],_a[0,2,2],_a[1,0,0],_a[1,0,1],_a[1,0,2],
 _a[1,1,0],_a[1,1,1],_a[1,1,2],_a[1,2,0],_a[1,2,1],_a[1,2,2],
 _a[2,0,0],_a[2,0,1],_a[2,0,2],_a[2,1,0],_a[2,1,1],_a[2,1,2],
 _a[2,2,0],_a[2,2,1],_a[2,2,2]}$
(%o14)                        done

This is just working at the top-level interactive prompt; if you need to construct a function, I think you'll see how to do it.

EDIT: Here's a way to create the polynomial.

(%i16) S : {[0, 0, 0], [0, 0, 1], [0, 0, 2], [0, 1, 0], [0, 1, 1], 
[0, 1, 2], [0, 2, 0], [0, 2, 1], [0, 2, 2], [1, 0, 0], 
[1, 0, 1], [1, 0, 2], [1, 1, 0], [1, 1, 1], [1, 1, 2], 
[1, 2, 0], [1, 2, 1], [1, 2, 2], [2, 0, 0], [2, 0, 1], 
[2, 0, 2], [2, 1, 0], [2, 1, 1], [2, 1, 2], [2, 2, 0], 
[2, 2, 1], [2, 2, 2]} $

(%i17) L : listify (S) $ 

(%i18) A : map (lambda ([l], arrayapply (_a, l)), L);
(%o18) [_a       , _a       , _a       , _a       , _a       , 
          0, 0, 0    0, 0, 1    0, 0, 2    0, 1, 0    0, 1, 1
_a       , _a       , _a       , _a       , _a       , 
  0, 1, 2    0, 2, 0    0, 2, 1    0, 2, 2    1, 0, 0
_a       , _a       , _a       , _a       , _a       , 
  1, 0, 1    1, 0, 2    1, 1, 0    1, 1, 1    1, 1, 2
_a       , _a       , _a       , _a       , _a       , 
  1, 2, 0    1, 2, 1    1, 2, 2    2, 0, 0    2, 0, 1
_a       , _a       , _a       , _a       , _a       , 
  2, 0, 2    2, 1, 0    2, 1, 1    2, 1, 2    2, 2, 0
_a       , _a       ]
  2, 2, 1    2, 2, 2
(%i19) U : map (lambda ([l], product (u[i]^l[i], i, 1, length(l))), L);
                2                 2   2   2      2  2
(%o19) [1, u , u , u , u  u , u  u , u , u  u , u  u , u , 
            3   3   2   2  3   2  3   2   2  3   2  3   1
           2                          2      2      2
u  u , u  u , u  u , u  u  u , u  u  u , u  u , u  u  u , 
 1  3   1  3   1  2   1  2  3   1  2  3   1  2   1  2  3
    2  2   2   2      2  2   2      2         2     2   2  2
u  u  u , u , u  u , u  u , u  u , u  u  u , u  u  u , u  u , 
 1  2  3   1   1  3   1  3   1  2   1  2  3   1  2  3   1  2
 2  2      2  2  2
u  u  u , u  u  u ]
 1  2  3   1  2  3
(%i20) A.U;                        
        2  2            2    2               2    2            2
(%o20) u  u  _a        u  + u  u  _a        u  + u  _a        u
        1  2   2, 2, 2  3    1  2   2, 1, 2  3    1   2, 0, 2  3
                 2  2              2  2                    2
 + u  _a        u  u  + _a        u  u  + u  _a        u  u
    1   1, 2, 2  2  3     0, 2, 2  2  3    1   1, 1, 2  2  3
                 2                 2              2
 + _a        u  u  + u  _a        u  + _a        u
     0, 1, 2  2  3    1   1, 0, 2  3     0, 0, 2  3
    2  2                 2                    2
 + u  u  _a        u  + u  u  _a        u  + u  _a        u
    1  2   2, 2, 1  3    1  2   2, 1, 1  3    1   2, 0, 1  3
                 2                 2
 + u  _a        u  u  + _a        u  u  + u  _a        u  u
    1   1, 2, 1  2  3     0, 2, 1  2  3    1   1, 1, 1  2  3
 + _a        u  u  + u  _a        u  + _a        u
     0, 1, 1  2  3    1   1, 0, 1  3     0, 0, 1  3
    2  2              2                 2
 + u  u  _a        + u  u  _a        + u  _a
    1  2   2, 2, 0    1  2   2, 1, 0    1   2, 0, 0
                 2              2
 + u  _a        u  + _a        u  + u  _a        u
    1   1, 2, 0  2     0, 2, 0  2    1   1, 1, 0  2
 + _a        u  + u  _a        + _a
     0, 1, 0  2    1   1, 0, 0     0, 0, 0

Note that the ordering of terms within each product doesn't conform to what humans would consider the usual convention, e.g. [1]^2*u[2]^2*_a[2,2,2]*u[3]^2 is the first term. Maxima is ordering the terms according to the subscripts, therefore _a[2,2,2] comes after u[1] and before u[3]. In some contexts this coincides with what humans expect, but here it doesn't; in any event, Maxima is consistent in hope of making programmatic manipulation work better.

(%i21) grind (%);

u[1]^2*u[2]^2*_a[2,2,2]*u[3]^2+u[1]^2*u[2]*_a[2,1,2]*u[3]^2
                              +u[1]^2*_a[2,0,2]*u[3]^2
                              +u[1]*_a[1,2,2]*u[2]^2*u[3]^2
                              +_a[0,2,2]*u[2]^2*u[3]^2
                              +u[1]*_a[1,1,2]*u[2]*u[3]^2
                              +_a[0,1,2]*u[2]*u[3]^2
                              +u[1]*_a[1,0,2]*u[3]^2
                              +_a[0,0,2]*u[3]^2
                              +u[1]^2*u[2]^2*_a[2,2,1]*u[3]
                              +u[1]^2*u[2]*_a[2,1,1]*u[3]
                              +u[1]^2*_a[2,0,1]*u[3]
                              +u[1]*_a[1,2,1]*u[2]^2*u[3]
                              +_a[0,2,1]*u[2]^2*u[3]
                              +u[1]*_a[1,1,1]*u[2]*u[3]
                              +_a[0,1,1]*u[2]*u[3]
                              +u[1]*_a[1,0,1]*u[3]
                              +_a[0,0,1]*u[3]
                              +u[1]^2*u[2]^2*_a[2,2,0]
                              +u[1]^2*u[2]*_a[2,1,0]
                              +u[1]^2*_a[2,0,0]
                              +u[1]*_a[1,2,0]*u[2]^2
                              +_a[0,2,0]*u[2]^2
                              +u[1]*_a[1,1,0]*u[2]
                              +_a[0,1,0]*u[2]+u[1]*_a[1,0,0]
                              +_a[0,0,0]$
(%o21)                        done

Thanks a lot. this is very close to what I want, except it is not a multidimensional array/list. What I want to finally achieve is to create a multivariate polynomial as explained [here](https://stackoverflow.com/questions/51628270/how-to-generate-a-symbolic-multivariate-polynomial-of-a-given-dimension-in-sympy). and for that I think I need a multidimensional list to perform the Frobenius inner product later. — Foad S. Farimani, Aug 02 '18 at 06:54
@Foad I've edited my answer to include a way to create the polynomial. Let me know if I've misunderstood the goal. — Robert Dodier, Aug 02 '18 at 16:44
you are awesome. I'm gonna first try the code, and then edit my post to fit with the answer you gave me. I was also wondering if you could be so kind to take a look at my follow up question [here](https://cs.stackexchange.com/questions/95886/algorithm-for-using-power-series-to-numerically-solve-a-partial-differential-equ)? — Foad S. Farimani, Aug 02 '18 at 18:57

How to generate a multidimensional array of indexed variables in Maxima?

1 Answers1