kill(all)$

maxima : programming: multinomial version 1.0
peter.vlasschaert@gmail.com,27/07/2022

math person :
------------------
The Trinomial theorem | Arithmetic and Geometry Math Foundations 56 | N J Wildberger
→ https://www.youtube.com/watch?v=k4z2Y8Y_r-M
--------------------------------------------------------------------------------------------------------------------

--------------------------------------------------------------------------------------------------------------------
see: website
maxima : programming :Binomial_theorem Newton,version 1.0
--------------------------------------------------------------------------------------------------------------------

load packages ,option:
-----------------------------

load(functs)$/*combination*/
load(lrats)$
ratprint:false$

multinomial
----------------

ordergreat(a,b,c)$

properties(sum);

declare(m,integer)$

askinteger(m);

p1:(a+b)^m = ' sum(combination(m,j)*a^(m-j)*b^j,j,0,m);

p1a:subst([b=b+c],p1);

p1b:(b+c)^j = ' sum(combination(j,k)*b^(j-k)*c^k,k,0,j);

p1c:lhs(p1a)=subst([p1b],rhs(p1a));

p1d:subst([m=2],p1c);

result : see Binomium van Newton

p1e:lhs(p1d)=expand(ev(rhs(p1d),nouns));

check :

p1f:lhs(p1d) =expand(lhs(p1d));

p1g:is(rhs(p1e)=rhs(p1f));

load("sym/compile")$ /*needed to use :multinomial*/

p2a:multinomial(3,[1,1,1]);/* 3!/(1!*1!*1!) , 3! =6 */

3!;

my definition : mymultinomial
-----------------------------------------
pp : number.
qq : list of numbers.

mymultinomial(pp,qq):=block([],
s1:pp,
s2:makelist(qq[i],i,1,length(qq)),
s3:s1!/product(s2[i]!,i,1,length(qq)),
   return(s3)     
)$

syntax : function(args) : = block([local var.],code,return(result))

→ optimize code :

grind(mymultinomial)$

p2b:mymultinomial(3,[1,1,1]);

check :

is(p2a=p2b);

theory :
----------

multinomial(n,[n1,n2]) =(n!/(n1!*n2!)) =C(n,n1)*C(n-n1,n2)
------------------------------------------------------------------------------
constraint : n1+n2 = n
-------------------------------------------------------------------------------
rem : combination(n,k) = C(n,k)

declare([n,n1,n2],integer)$

p3:n!/(n1!*n2!) =combination(n,n1)*combination(n-n1,n2);

assume(equal(n,n1+n2))$/* n=n1+n2*/

map(askinteger,[n,n1,n2]);

p3a:lhs(p3) =  ev(rhs(p3),nouns);

see: Menu: Simply → Factorials and Gamma →Convert factorials

p3b:lhs(p3)=makefact(rhs(p3a));

p3b1:subst([n=n1+n2],rhs(p3b));

p3b2:ev(p3b1,n1+n2=n);

check :

is(lhs(p3b)=p3b2);

? coef of  a^2*b^1*c^2 for (a+b+c)^5
-----------------------------------------------
------------------------------------------------

p4:a^2*b^1*c^2;

p4a:args(p4);

a → select:2 out 5 > rest = 3

p4b1:'combination(5,2)=combination(5,2);

b → select:1 out 5-2 > rest = 2

p4b2:'combination(5-2,1)=combination(5-2,1);

c → select:2 out 5-2-1

p4b3:'combination(5-2-1,2)=combination(5-2-1,2);

result :

p4b:p4b1*p4b2*p4b3;

→ simple way : computation.

p4c1:multinomial(5,[2,1,2]);

p4c2:multinomial(5,[2,1,2]);

rem : p4b=p4c1=p4c2

theory ; number of terms
(x1+x2+x3+...+xr)^n
------------------------------------------
answer: combination(n+r-1,r-1)

p5:(a+b+c)^2;/*example :multinomial*/

→ answer :

p5a:'combination(2+3-1,3-1) =combination(2+3-1,3-1);

p5b:(a+b+c)^2=expand((a+b+c)^2);

p5c:length(rhs(p5b));

check : number of terms "p5"

is(rhs(p5a)=p5c);

------------------------------------------------------------------------------

theory :
----------

(x1+x2+...+xr)^n = sum(multinomial(n,[i1,i2,...,ir])*product(xj^ij,j,1,r))
--------------------------------------------------------------------------------------------
constraint : i1+i2+...+ir = n
---------------------------------------------------------------------------------------------

example : (x+y)^3
-----------------------

l:[]$
ll:[]$
for i:0 thru 3 do(
    for j:0 thru 3 do(
    if i+j =3 then   
    (
     l:endcons([i,j],l) ,
     ll:endcons(mymultinomial(3,[i,j])*x^i*y^j,ll)
     )       
    )
);

l;

ll;

lll:sum(ll[i],i,1,length(ll));

maxima result :

p6:(x+y)^3 = expand((x+y)^3);

check :lll  ?  rhs(p6)

is(lll=rhs(p6));

kill(l,ll,lll);

lll;

example : (x+y+z)^v, v=2
---------------------------------
term :x^i*y^j*z^k,see : "p4"

l:[]$
ll:[]$

v:2;/* power:multinomial*/

for i:0 thru v do(
    for j:0 thru v do(
        for k:0 thru v do(
    if i+j +k=v then   
    (
     l:endcons([i,j],l) ,
     ll:endcons(mymultinomial(v,[i,j,k])*x^i*y^j*z^k,ll)
     )       
    )
  )
)$

ll;

lll:sum(ll[i],i,1,length(ll));

maxima result :

p7:(x+y+z)^2 = expand((x+y+z)^2);

check :lll  ?  rhs(p7)

is(lll=rhs(p7));

------------------------------------------------------------------------------------------