I am defining a following simplification for a functional E and a list
of atoms. E should pull out everything that does not depend on any of
the atoms in the list, but also it should pull out E(<...>). One can
think of E as a time averaging functional, and about the atoms in the
list as functions of time.
Here is an example ([x, y] is the list of the atoms)
(%i44) [E(y*E(x)*a), E(E(b*y*x)*a), E(a+b)];
(%o44) [a E(x) E(y), E(1) a b E(x y), E(1) (b + a)]
(%i45) [E(y*E(x)*a), E(E(b*y*x)*a), E(a+b)];
(%o45) [a E(x) E(y), E(1) a b E(x y), E(1) (b + a)]
My attempt is below. I have some questions. Am I on the right track or
it can be express in more simple way? Does it work in all situations? Is
there any better solution which allows for extensions? (I can think of
defining several such functionals and introducing a predicate which
determents if `mapatom' is a constant with respect to the functional.)
load("partition.mac")$
load("opsubst")$
/* determines if `mapatom' is a constant for an operator `oper'
(%i14) put('E, [x, y], 'listofnonconstants);
(%o14) [x, y]
(%i15) constantforop_expr(x, E);
(%o15) false
(%i16) constantforop_expr(a, E);
(%o16) true */
constantforop_mapatom(ma, oper):= block([loc: get(oper, 'listofnonconstants)],
if loc=false then true else
lfreeof(get(oper, 'listofnonconstants), ma))$
/* determines if `expr' is a constant for an operator `oper'
(%i14) put('E, [x, y], 'listofnonconstants);
(%o14) [x, y]
(%i15) constantforop_expr(x, E);
(%o15) false
(%i16) constantforop_expr(a, E);
(%o16) true
(%i17) constantforop_expr(x+y, E);
(%o17) false
(%i18) constantforop_expr(2*x*y, E);
(%o18) false
(%i19) constantforop_expr(a+b, E);
(%o19) true */
constantforop_expr(expr, oper):= block([],
catch(for ma in listofvars(hide_op(expr, oper)) do if not
constantforop_mapatom(ma, oper) then throw(false),
true))$
/* replaces E(<..>) by _HIDEFUN_ */
hide_op(expr, oper):=opsubst(oper=lambda([[arg]], _HIDEFUN_), expr)$
/* separate components of a product into two parts: constant and
non-constant with respect to operator E. Lists are placed into `ANSss'.
true is return if a constant part is not empty and not equal to 1
(%i20) put('E, [x, y], 'listofnonconstants);
(%o20) [x, y]
(%i21) partition_withconst(x*a);
(%o21) true
(%i22) [partition_withconst(x*a), 'ANSss];
(%o22) [true, ANSss]
(%i23) [partition_withconst(x*a), ANSss];
(%o23) [true, [[a], [x]]]
(%i24) partition_withconst(x*a);
(%o24) true
(%i25) [partition_withconst(x*a), ANSss];
(%o25) [true, [[a], [x]]]
(%i26) [partition_withconst(x+a), ANSss];
(%o26) [false, [[a], [x]]]
(%i27) [partition_withconst(a), ANSss];
(%o27) [true, [[a], [1]]]
(%i28) [partition_withconst(1), ANSss];
(%o28) [false, [[a], [1]]] */
partition_withconst(expr):= block([ret, prederror: true],
if expr#1 and constantforop_expr(expr, 'E) then ( ANSss:
[[expr], [1]], return(true)),
ret: partition_expression("*", lambda([v], constantforop_expr(v,
'E)), [], cons, "[", 'ANSss, expr),
ret and first(ANSss)#[1] and not emptyp(second(ANSss)))$
matchdeclare(pp, partition_withconst)$
tellsimpafter(E(pp), xreduce("*", first(ANSss))*E(xreduce("*", second(ANSss))))$
/* examples */
put('E, [x, y], 'listofnonconstants)$
[E(a), E(x*a), E(E(x)*a)];
[E(y*E(x)*a), E(E(b*y*x)*a), E(a+b)];