How to do abstract matrix operation?

Question

In this simple example, I have 2 lists, each of which contains 3 matrices。 $J$ here is also a matrix with the same dimension. What I want is $aJc+bJ(b+c)$, but if I use dot product, the result is as follows,

x = {a, b, b}; y = {c, b, c}; 
Sum[x[[i]] .J. y[[i]], {i, 1, 3}] // FullSimplify
(* return is a.J.c + b.J.b + b.J.c *)

When I use multiplication, then Mathematica goes too far:

x = {a, b, b}; y = {c, b, c}; 
Sum[x[[i]] *J* y[[i]], {i, 1, 3}] // FullSimplify
(* return is (b^2 + (a + b) c) J*)

So how to get what I want ?

Answer 1

We can surely generalize this (to include scalar multiples, for instance), but here's an implementation that factors the expression written in terms of Dots:

dotFactor[expr_] := expr //. {
  Dot[x__, z1_] + Dot[x__, z2_] :> Dot[x, z1 + z2], 
  Dot[z1_, x__] + Dot[z2_, x__] :> Dot[z1 + z2, x], 
  Dot[1, x__] | Dot[x__, 1] :> Dot[x]
 }

Then, using

x = {a, b, b};
y = {c, b, c};

we get

expr = Inner[#1.J.#2 &, x, y]
expr // dotFactor
(* a.J.c + b.J.b + b.J.c *)
(* a.J.c + b.J.(b + c) *)

Also,

x = {a, b, b, d, a};
y = {c, b, c, d, b};
expr = Inner[#1.J.#2 &, x, y]
expr // dotFactor
(* a.J.b + a.J.c + b.J.b + b.J.c + d.J.d *)
(* (a + b).J.(b + c) + d.J.d *)