a test of defrules expansion code



On Nov. 21, Robert Dodier wrote:

> OK, here's my attempt at a solution.
>
> I used defrule instead of tellsimpafter to define rules.
> I chose defrule in order to work around a tellsimpafter
> idiosyncrasy; namely that multiple tellsimpafter rules are not all
> executed -- there is a flag AFTERFLAG* or something like that
> which, I guess, is supposed to prevent infinite loops, but also prevents
> multiple rules for the same operator from firing.
>
> I'd be happy to explain the rules in more detail (assuming they work
> as intended). I think the built-in rule machinery is the best way to
> go here -- doing stuff like factoring out scalars by hand is essentially
> just reinventing stuff that already exists in the rule code.
--------------------------------

I have found a disagreement with the correct answer for the
case of expanding the commutator

  [b^^3, c^^5]

by using b = Sig[2] and c = Sig[3], two of the
well known Pauli matrices (which are in fact
idempotent: the square of each reduces to the unit
2 x 2 matrix), so the correct answer is
just [b, c] once the Pauli matrices are inserted.

However, we first do the general expansion for
arbitrary b and c, first using the defrule code,
and then using the known formula for the
expansion of [B^^n,  C^^m]  (well, I don't really
know where it is referenced, but it is fairly obvious
if one first fully expands the first arg, and then
fully expands the second arg (see formula
in batch file listing below).

The defrule code produces a zero matrix result, while
the formal expansion formula produces the
correct answer ( 2*%i * Sig[1] ).

Perhaps the defrule code needs some slight
tweaking for the expansion part.

Ted Woollett

------------------------------
batch file run:

Maxima 5.22.1 http://maxima.sourceforge.net
using Lisp GNU Common Lisp (GCL) GCL 2.6.8 (a.k.a. GCL)
Distributed under the GNU Public License. See the file COPYING.
Dedicated to the memory of William Schelter.
The function bug_report() provides bug reporting information.
"2010-11-25"

(%i1) batch ("comm.mac")$
read and interpret file: #pc:/work5/comm.mac
(%i2) " test expansion with a pair of Pauli matrices "
(%i3) Sig[2]:matrix([0,-%i],[%i,0])
(%i4) Sig[3]:matrix([1,0],[0,-1])
(%i5) comm_replace:Sig[2] . Sig[3]-Sig[3] . Sig[2]
(%o5) matrix([0,2*%i],[2*%i,0])
(%i6) " required value for [b^^3, c^^5] is "
(%i7) comm_val:Sig[2]^^3 . Sig[3]^^5-Sig[3]^^5 . Sig[2]^^3
(%o7) matrix([0,2*%i],[2*%i,0])
(%i8) " ----------------------------------------"
(%i9) " defrule code method of expansion "
(%i10) load(commutators)
(%o10) "c:/work5/commutators.mac"
(%i11) declare([s1,s2,s3,s4],scalar)
(%i12) declare([a,b,c,d],nonscalar)
(%i13) ee:comm(b^^3,c^^5)
(%o13) comm(b^^3,c^^5)
(%i14) ee_reduce:comm_simp(ee)
(%o14) c . c . c . c . comm(b,c) . b . b+4*c . c . c . comm(b,c) . c . b . b
                                        +6*c . c . comm(b,c) . c . c . b . b
                                        +4*c . comm(b,c) . c . c . c . b . b
                                        +comm(b,c) . c . c . c . c . b . b
                                        +2*b . c . c . c . c . comm(b,c) . b
                                        +8*b . c . c . c . comm(b,c) . c . b
                                        +12*b . c . c . comm(b,c) . c . c . 
b
                                        +8*b . c . comm(b,c) . c . c . c . b
                                        +2*b . comm(b,c) . c . c . c . c . b
                                        +b . b . c . c . c . c . comm(b,c)
                                        +4*b . b . c . c . c . comm(b,c) . c
                                        +6*b . b . c . c . comm(b,c) . c . c
                                        +4*b . b . c . comm(b,c) . c . c . c
                                        +b . b . comm(b,c) . c . c . c . c
(%i15) ee_r1:ratsubst(comm_replace,comm(b,c),ee_reduce)
(%i16) ee_r2:subst([b = Sig[2],c = Sig[3]],ee_r1)
(%o16) matrix([0,0],[0,0])
(%i17) " which does not agree with comm_val  "
(%i18) comm_val
(%o18) matrix([0,2*%i],[2*%i,0])
(%i19) " ----------------------------------------------"
(%i20) "The general expansion of [B^^n, C^^m] should be "
(%i21) comm_exp(B,n,C,m):=sum(
                sum(B^^(j-1) . C^^(k-1) . comm(B,C) . C^^(m-k) . 
B^^(n-j),k,1,
                    m),j,1,n)
(%i22) " test for our case "
(%i23) dd_reduce:comm_exp(b,3,c,5)
(%o23) c^^4 . comm(b,c) . b^^2+c^^3 . comm(b,c) . c . b^^2
                              +c^^2 . comm(b,c) . c^^2 . b^^2
                              +c . comm(b,c) . c^^3 . b^^2
                              +comm(b,c) . c^^4 . b^^2+b^^2 . c^^4 . 
comm(b,c)
                              +b^^2 . c^^3 . comm(b,c) . c
                              +b^^2 . c^^2 . comm(b,c) . c^^2
                              +b^^2 . c . comm(b,c) . c^^3
                              +b^^2 . comm(b,c) . c^^4
                              +b . c^^4 . comm(b,c) . b
                              +b . c^^3 . comm(b,c) . c . b
                              +b . c^^2 . comm(b,c) . c^^2 . b
                              +b . c . comm(b,c) . c^^3 . b
                              +b . comm(b,c) . c^^4 . b
(%i24) dd_r1:ratsubst(comm_replace,comm(b,c),dd_reduce)
(%i25) dd_r2:subst([b = Sig[2],c = Sig[3]],dd_r1)
(%o25) matrix([0,2*%i],[2*%i,0])
(%i26) " which agrees with comm_val "
(%i27) comm_val
(%o27) matrix([0,2*%i],[2*%i,0])
(%i28) " ---------------------------"

/*
=============================
batch file comm.mac


" test expansion with a pair of Pauli matrices "$

Sig[2] : matrix ([0,-%i],[%i,0])$

Sig[3] : matrix ([1,0],[0,-1])$

comm_replace : Sig[2] . Sig[3] - Sig[3] . Sig[2];

 " required value for [b^^3, c^^5] is "$

comm_val : Sig[2]^^3 . Sig[3]^^5  - Sig[3]^^5 . Sig[2]^^3 ;

" ----------------------------------------"$

" defrule code method of expansion "$


load(commutators);

 declare ([s1, s2, s3, s4], scalar)$

 declare ([a, b, c, d], nonscalar)$

 ee : comm (b^^3, c^^5);

 ee_reduce : comm_simp (ee);

 ee_r1 : ratsubst (comm_replace, comm (b,c), ee_reduce )$


 /* ee_r1 : subst (comm = comm_replace, ee_reduce); */

 ee_r2 : subst ([b = Sig[2], c = Sig[3]], ee_r1 );

 " which does not agree with comm_val  "$

 comm_val;

 " ----------------------------------------------"$

 "The general expansion of [B^^n, C^^m] should be "$

 comm_exp (B,n,C,m) :=

  sum (sum (B^^(j-1) . C^^(k-1) . comm (B,C) . C^^(m-k) . 
B^^(n-j),k,1,m),j,1,n)$

  " test for our case "$

 dd_reduce : comm_exp (b,3,c,5);

 dd_r1 : ratsubst (comm_replace, comm(b,c), dd_reduce)$

 dd_r2 : subst ( [b = Sig[2], c = Sig[3]], dd_r1 );

 " which agrees with comm_val "$

 comm_val;

 " ---------------------------"$

 */

=================================
Robert Dodier's defrule code is the file
commutator.mac, in which the replacement
FOO --->  comm_simp had been made,
and we have included the final expansion
inside the definition of comm_simp.

 /*  commutators.mac



/* elw: 11-23-10, FOO --> comm_simp
 * and include expand in def.
 */


declare (bilinear, feature);
bilinearp (e) := featurep (e, bilinear);
declare (comm, bilinear);

matchdeclare (xx, bilinearp);
matchdeclare ([aa, bb, cc], all);
defrule (rbilineara1, xx (aa + bb, cc), map (lambda ([u], xx (u, cc)), bb));
defrule (rbilineara2, xx (cc, aa + bb), map (lambda ([u], xx (cc, u)), bb));

matchdeclare (xx, bilinearp);
matchdeclare (aa, scalarp);
matchdeclare (bb, nonscalarp);
matchdeclare (cc, all);
defrule (rbilinearm1, xx (aa * bb, cc), aa * xx (bb, cc));
defrule (rbilinearm2, xx (cc, aa * bb), aa * xx (cc, bb));

nctimesp (e) := not atom(e) and op(e) = ".";
matchdeclare (aa, nctimesp);
matchdeclare (cc, all);
defrule (rcommp1, comm (aa, cc), first (aa) . comm (rest (aa), cc) +
comm (rest (aa, -1), cc) . last (aa));
defrule (rcommp2, comm (cc, aa), first (aa) . comm (cc, rest (aa)) +
comm (cc, rest (aa, -1)) . last (aa));

matchdeclare (aa, nonscalarp);
matchdeclare (nn, integerp);
defrule (rncexpt, aa^^nn, apply (".", makelist (aa, i, 1, nn)));

/* put antisymmetric declaration here so built-in rules don't
interfere w/ rule definitions */
declare (comm, antisymmetric);

comm_simp (e) :=
   (e : apply1 (e, rbilineara1, rbilineara2, rbilinearm1, rbilinearm2),
    block ([dotexptsimp : false],
        apply1 (applyb1 (e, rncexpt), rcommp1, rcommp2),
        expand (%%)));

/* example:
 *
    load (commutators);
    declare ([s1, s2, s3, s4], scalar);
    declare ([a, b, c], nonscalar);
    e : comm (3*s1*a^^2 + 4*s2*b^^3 ,  (s3*a^^4)/5 - 6*s4*c^^5);
    comm_simp (e);
 *
 */

*/