Re: [Maxima-discuss] Z-Transforms

[Mike Miller]

On Fri, Jul 08, 2016 at 12:12:15 -0400, Przemek Klosowski wrote:
> On 07/07/2016 03:11 PM, Stavros Macrakis (Σταῦρος Μακράκης) wrote:
> > a dozen-line Octave program to do what can be done in two lines in
> > Maxima. [...]
> > 
> > T[n,k]:=if k < 1 or k > n then 0 elseif k = 1 or k = n then 1
> >          else k*T[n-1,k]+(n-k+1)*T[n-1,k-1]$
> > 
> > Then
> > 
> > genmatrix(T,11,10) =>
> > 
> >  [ 1   0      0        0        0        0       0        0       0
> > 0   ]
> >  [                                   ]
> >  [ 1   1      0        0        0        0       0        0       0
> > 0   ]
> >  [                                   ]
> >  [ 1   4      1        0        0        0       0        0       0
> > 0   ]
> >  [                                   ]
> >  [ 1   11     11       1        0        0       0        0       0
> > 0   ]
> >  [                                   ]
> >  [ 1   26     66      26        1        0       0        0       0
> > 0   ]
> >  [                                   ]
> >  [ 1   57    302      302      57        1       0        0       0
> > 0   ]
> >  [                                   ]
> >  [ 1  120    1191    2416     1191      120        1        0       0
> > 0   ]
> >  [                                   ]
> >  [ 1  247    4293    15619    15619     4293      247       1       0
> > 0   ]
> >  [                                   ]
> >  [ 1  502   14608    88234   156190    88234     14608     502      1
> > 0   ]
> >  [                                   ]
> >  [ 1  1013  47840   455192   1310354  1310354   455192    47840    1013
> > 1   ]
> >  [                                   ]
> >  [ 1  2036  152637  2203488  9738114  15724248  9738114  2203488  152637
> > 2036 ]
> 
> 
> I thought "challenge accepted" and tried to write this in a vectorized
> shortcut, but I couldn't hack it. I start by
> 
> T=eye(11,10);T(:,1)=1; 
> repmat(1:10,11,1).*shift(T,1)+(repmat((1:11)',1,10)-repmat(1:10,11,1)+1).*shift(shift(T,1),1,2)
> 
> but I couldn't see how to nicely write the rest without explicit loops. Can
> anyone suggest an approach?

I don't have an expression, but I think it would involve rewriting that
nicely formulated recursive definition into a closed-form non-recursive
definition for each element of the array, because that's the way the
Octave language works.

So the first step would be to redefine the above so that one can compute
the value of T(n,k) in the lower triangle independently of any other
element of T.

As a non-mathematician, that looks harder than it's worth thinking
about, so I would just write the function in C++ and it becomes a fast
function call in Octave.

-- 
mike