Matlab/Mupad symbolic simplification for pretty output

Question

I need to automatically simplify some symbolic expressions but the simplify function of matlab can't do everything i need.

Example:

simplify(expand((ax + bx)^2 + (ay + by)^2))

Which results in the output

ax^2 + 2*ax*bx + ay^2 + 2*ay*by + bx^2 + by^2

So i've tried to create my own rules to make the output prettier

function [ result ] = simplify_pretty( term )
    read(symengine, 'simplify_pretty_rules.mu');
    result = feval(symengine, 'Simplify', term, 'SelectRules = simplify_pretty_rules')
end

with "simplify_pretty_rules.mu"

simplify_pretty_rules := proc()
begin
    [
        Rule(#X^2 - 2 * #X * #Y + #Y^2, (#X - #Y)^2),
        Rule(#X^2 + 2 * #X * #Y + #Y^2, (#X + #Y)^2)
    ]
end_proc:

This works for input like...

simplify_pretty(expand((ax + bx)^2 + (ay + by)^2))

...but when changed to...

simplify_pretty(expand(-(ax + bx)^2 - (ay + by)^2))

...it doesn't

So the questions is: Is it possible to create rules which work in (nearly) all situations? What am I doing wrong?

---

Edit: This seems to be a general problem. If an expression has some kind of factor/is disarranged it seems to be quite hard for matlab to figure this out. I need a robust way to do this kind of stuff even if it takes some time to process.

---

Edit2:

What type of expressions I'm talking about?
The expressions contain only mathematical simple operations (add/sub/mul/div).

What do I want to get simplified?
One of the most common things I want to simplify are binomials which got expanded (see the original question).

Another common thing contained in the expressions expanded dot/cross products which could get simplified to cos(a)/sin(a) expressions which might get merged together in the next steps. But this is more advanced so I'm not sure if there are useful automatic solutions for it.

What is "pretty" in my case?
I guess most important is the length of the expressions (as shorter as better). But it would also be neat to optimize expressions like

a/2 + b/2 + c/2

to

(a + b + c) / 2

Answer 1

For me it works if you add a #n element to your rules, so that your simplify_pretty_rules file looks like this:

simplify_pretty_rules := proc()
begin
    [
    Rule(#n * #X^2 - #n * 2 * #X * #Y + #n * #Y^2, #n * (#X - #Y)^2),
    Rule(#n * #X^2 + #n * 2 * #X * #Y + #n * #Y^2, #n * (#X + #Y)^2)
    ]
end_proc:

I believe in this case the #n refers to a constant term. This gives me:

term = 'expand((ax - bx)^2 - 3 *(ay + by)^2)';
feval(symengine, 'Simplify', term, 'SelectRules = simplify_pretty_rules')

ans =  (ax - bx)^2 - 3*(ay + by)^2

So it seems to work.