7.109 Bug in multilinear in Maple V.5 (25.5.99)

7.109.1 Tullia Dymarz

Why doesn’t this work...

> define(h,multilinear); 
> h(h(2*a,b),c); 
Error, (in h) too many levels of recursion
 

while this works...

> define(n,linear); 
> n(n(2*a,b),c); 
 
                           2 n(n(a, b), c)
 

It is corrected with Maple 6. (U. Klein)

7.109.2 Robert Israel

It’s a bug. The programmers evidently didn’t expect that the arguments to h might contain h. The result is that whenever x or y contains h, h(x,y) attempts to evaluate itself, resulting in an infinite recursion.

The "linear" code doesn’t contain this bug, but it has another one:

> n(a*b); 
 Error, (in n) too many levels of recursion
 

Here is a way to patch h so that it will work:

h := readlib(procmake)(subs('`&function`(op,`&expseq`(`&local`[2]))' = 
 '`&function`(op,`&expseq`(`&function`(subs,`&expseq`(_FUNCNAME = 
 `&args`[-2], `&local`[2]))))' , readlib(procbody)(h)));
 

If h is a binary function, you could also define it to be multilinear "by hand" as follows:

> define(h, 
  h(a::nonunit(algebraic)+b::nonunit(algebraic),c::anything) = h(a,c)+h(b,c), 
  h(a::anything,b::nonunit(algebraic)+c::nonunit(algebraic)) = h(a,b)+h(a,c), 
  h(a::nonunit(constant)*b::nonunit(algebraic),c::anything) = a*h(b,c), 
  h(a::anything,b::nonunit(constant)*c::nonunit(algebraic)) = b*h(a,c));
 

7.109.3 Andre Heck (9.6.99)

| From: Robert Israel <[e-mail]>

As the author of the original code for defining multilinear mappings (1991), I can only confirm that it is a bug and that it was introduced later (1996) by WMI programmers. The original code was still present in Maple V Version Release 4:

 
    |\^/|     Maple V Release 4 (University of Amsterdam) 
._|\|   |/|_. Copyright (c) 1981-1996 by Waterloo Maple Inc. All rights 
 \  MAPLE  /  reserved. Maple and Maple V are registered trademarks of 
 <____ ____>  Waterloo Maple Inc. 
      |       Type ? for help. 
> define(multilinear(h)); 
> h(h(2*a,b),c); 
                                2 h(h(a, b), c)
 

Interested in the old code?

> interface(verboseproc=3): 
> print(`define/multilinear`); 
proc(OpName) 
local t; 
global _X; 
option 
`Copyright (c) 1991 by the University of Waterloo. All rights reserved.`; 
    if nargs <> 1 or not type(OpName, 'string') then 
        ERROR(`invalid arguments`) 
    fi; 
    if assigned(_X) then _X := '_X' fi; 
    t := proc() 
        local func, k, term, terms, cf, coefficient; 
        option remember, `Copyright (c) 1991 by the University of Waterl\ 
        oo. All rights reserved.`; 
            if nargs < 1 then ERROR(`no element`) 
            elif nargs = 1 then 
                if type(args[1], 'constant') then RETURN(args[1]*'_X(1)') 
                elif type(args[1], `+`) then RETURN(map(_X, args[1])) 
                elif type(args[1], `*`) and hastype(args[1], 'constant') 
                then 
                    cf := select(type, args[1], 'constant'); 
                    term := args[1]/cf; 
                    RETURN(cf*_X(term)) 
                fi 
            elif 2 <= nargs then 
                if type(args[1], `+`) then 
                    RETURN(map(_X, args[1 .. nargs])) 
                else for k from 2 to nargs do 
                        if type(args[k], `+`) then 
                            func := subs({_KK = k, 'F' = _X}, proc() 
                                    F(args[2 .. _KK], args[1], 
                                    args[_KK + 1 .. nargs]) 
                                end); 
                            RETURN(map(func, args[k], args[1 .. k - 1], 
                                args[k + 1 .. nargs])) 
                        fi 
                    od 
                fi; 
                terms := NULL; 
                coefficient := 1; 
                for k to nargs do 
                    if type(args[k], 'constant') and args[k] <> 1 then 
                        coefficient := coefficient*args[k]; term := 1 
                    elif 
                    type(args[k], `*`) and hastype(args[k], 'constant') 
                    then 
                        cf := select(type, args[k], 'constant'); 
                        coefficient := coefficient*cf; 
                        term := args[k]/cf 
                    else term := args[k] 
                    fi; 
                    terms := terms, term 
                od; 
                if coefficient <> 1 then RETURN(coefficient*_X(terms)) fi 
            fi; 
            RETURN('_X(args)') 
        end; 
    if member(OpName, {'func', 'k', 'term', 'terms', 'cf', 'coefficient'}) 
    then t := subs(OpName = 'b', eval(t)) 
    fi; 
    if type(OpName, protected) then 
        unprotect('OpName'); lprint(`Warning: new definition for`, OpName) 
    fi; 
    OpName := subs(_X = OpName, eval(t)); 
    NULL 
end
 

This piece of code will work if you change the name of the procedure, say to ‘mlinear‘, and change the first type checking into if nargs <> 1 or not type(OpName, 'symbol')

 
    |\^/|     Maple V Release 5 (WMI Campus Wide License) 
._|\|   |/|_. Copyright (c) 1981-1997 by Waterloo Maple Inc. All rights 
 \  MAPLE  /  reserved. Maple and Maple V are registered trademarks of 
 <____ ____>  Waterloo Maple Inc. 
      |       Type ? for help. 
> read mlinear; 
> mlinear(h); 
> h(h(2*a,b),c); 
                                2 h(h(a, b), c)
 

Note that I have not tested whether other problems arise with this code in Maple V Release 5.