So what is your solution, I'm curious.

Regarding the number of solvable boards. It is actually possible to calculate the exact number. The solvable boards correspond to the image of the matrix of moves. It will be a vector (sub)space over Z/2Z, so it will always be a power of 2, and its size will be determined by the rank of the matrix. For example for 5x5, there are 2^22 valid boards.

And regarding the number of solutions to a given board, (assuming you ignore the order of moves, because it doesn't matter) it will always be the exact same number, and in the case of 5x5 it will be 2^8. 8 is the nullity of the matrix of moves.

Note that 22+8=30, which is 44+33+22+11, which is the number of possible moves.

We might be in pretty similar spots with respect to a solver, then. For one with no conscious linear algebra or deep programming background, as long as XOR or modular addition is familiar (which I suppose narrows it), I think a possible sequence of realizations is:

1. The flippable regions need be flipped at most once, in any order. So, one could think of a solution as a collection of boolean variables, one for each region, indicating whether to flip it (true/1 if so, false/0 otherwise, feels natural).

2. A tile needs to be flipped an even number of times if it starts in the desired state, and an odd number of times otherwise. So: Give each tile an equation by forming the XOR (or addition modulo 2) of the variables for all the regions it inhabits, and then setting that expression equal to 0 if the tile starts out in the desired state, otherwise 1.

3. Observe that such equations are always easy to rearrange so as to solve for any desired variable.

4. Observe that we can therefore run through a familiar substitution process, and good things will happen.

What I find elegant here is that the resulting algorithm is known to many from early algebra and has a much more clearly bounded runtime than the backtracking searches that suit other logic puzzles. Of course, as you've evinced, the fuller representation with linear algebra opens many doors to both implementing this well and exploring its other properties.

Finding par, on the other hand, I'll have to consider further. I suspect that, after representing the solution space succinctly, you could get a fair way with a backtracking search therein (where to make a forward step is to concretely assign a value to a variable), with the search tree trimmed via the best par observed so far. If you have something drastically better than that, I'm sure that would be interesting!

Looking forward to seeing what you do next, also perhaps the potential blog post you mentioned elsewhere in this thread, etc.