Foad S. Farimani
Foad S. Farimani

Reputation: 14016

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

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:

enter image description here

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.

Upvotes: 1

Views: 172

Answers (1)

Robert Dodier
Robert Dodier

Reputation: 17576

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

Upvotes: 1

Related Questions