Can operators that act over indices be defined in Maxima?

24 Views Asked by At

I need to define a differentiation operator that acts over the indices of the elements x[t] of a time series. It should act this way: B x[t] --> x[t-1].

I have seen here that one can establish the dependency of a symbol over a variable by depends (for example, depends(b, [R1, R2]);). This would allow to pass a function as an argument to another function. However, this doesn't work for my purpose, since if I want an expression like (1+B)^n involving the operator acting over a term x[t] it won't give the expected result.

Can this kind of operator be defined in maxima?

1

There are 1 best solutions below

0
Robert Dodier On

There isn't anything built-in, but it is possible to get something working via pattern matching rules. (Pattern matching in the sense of expression pattern matching -- stuff like 4*a - b matches x + y with x equal to 4*a and y equal to - b.)

Here is what I came up with. This defines B as a "backup" operator. announce_rules_firing is an undocumented Easter egg which can help make it clearer what's going on.

announce_rules_firing: true;

matchdeclare (cc, lambda ([e], freeof (B, e)));
matchdeclare (aa, all);
matchdeclare (nn, all);
matchdeclare (xx, all);
matchdeclare (bb, lambda ([e], not freeof (B, e)));

tellsimpafter (B (xx), subst ('t = 't - 1, xx));
tellsimpafter ((B^nn)(xx), subst ('t = 't - nn, xx));
tellsimpafter ((cc*B)(xx), cc*subst ('t = 't - 1, xx));
tellsimpafter ((cc*B^nn)(xx), cc*subst ('t = 't - nn, xx));

simp: false;
tellsimpafter ((aa + bb)(xx), aa + if not atom (bb) and op (bb) = "+" then map (lambda ([e], e(xx)), bb) else bb (xx));
simp: true;

Examples. It is assumed shift_operator.mac contains the preceding stuff.

(%i2) load ("shift_operator.mac");
(%o2)                  shift_operator.mac
(%i3) B(y[t]);
By Brule1 , B(y[t]) --> y[t-1] 
(%o3)                        y
                              t - 1
(%i4) B(y[t] - x[t]);
By Brule1 , B(y[t]-x[t]) --> y[t-1]-x[t-1] 
(%o4)                    y      - x
                          t - 1    t - 1
(%i5) (B^2)(u[t]);
By subvarrule1 , (B^2)(u[t]) --> u[t-2] 
(%o5)                        u
                              t - 2
(%i6) (k*B)(w[t]);
By subvarrule2 , (B*k)(w[t]) --> k*w[t-1] 
(%o6)                       k w
                               t - 1
(%i7) (k*B^n)(y[t]);
By subvarrule3 , (B^n*k)(y[t]) --> k*y[t-n] 
(%o7)                       k y
                               t - n
(%i8) ((1 + B)^4)(x[t]);
By subvarrule4 , ((B+1)^4)(x[t]) --> ((B+1)^4)(x[t]) 
                                 4
(%o8)                     (B + 1) (x )
                                    t
(%i9) expand (((1 + B)^4)(x[t]));
By subvarrule4 , ((B+1)^4)(x[t]) --> ((B+1)^4)(x[t]) 
By subvarrule1 , (B^4)(x[t]) --> x[t-4] 
By subvarrule3 , (4*B^3)(x[t]) --> 4*x[t-3] 
By subvarrule3 , (6*B^2)(x[t]) --> 6*x[t-2] 
By subvarrule2 , (4*B)(x[t]) --> 4*x[t-1] 
By subvarrule4 , (B^4+4*B^3+6*B^2+4*B+1)(x[t]) --> 
  4*x[t-1]+6*x[t-2]+4*x[t-3]+x[t-4]+1 
(%o9)      4 x      + 6 x      + 4 x      + x      + 1
              t - 1      t - 2      t - 3    t - 4

This much is enough to get going -- maybe you can say more about what-all you need to do, and we can think about ways to extend this stuff.