Simplify expression: overcome a barrier in the complexity function

Question

I would like to simplify

x=n[a] f[b, a] + f[a, b]

using 2 rules

r1[e_] := e /. {n[x_] :> 1 - h[x]}
r2[e_] := e /. {f[x_, y_] :> -f[y, x]}

where the expected answer is h[a]f[a,b]. Using

Simplify[n[a] f[b, a] + f[a, b], TransformationFunctions -> {Automatic, r1, r2}]

returns unaltered expression. In contrast, the transformation functions are operative at the intermediate steps, and Simplify works as expected:

Simplify[1 - n[a], TransformationFunctions -> {Automatic, r1, r2}]
Simplify[-f[b, a], TransformationFunctions -> {Automatic, r1, r2}]
Simplify[-n[a] f[a, b] + f[a, b], TransformationFunctions -> {Automatic, r1, r2}]
(* h[a] *)
(* f[a, b] *)
(* f[a, b] h[a] *)

I was trying to guide Simplify by ComplexityFunction, however, no success. I suspect, all intermediate forms have higher LeafCount due to extra minuses. But how to overcome that? I appreciate your help.

Comment

I prefer that Simplify works selectively and applies the given transformation rules only if the final complexity is lower in the standard definition---therefore, I am talking about the barrier here. This can be best illustrated by an example. Besides xdefined above, consider one more expression that does not need simplification:

y=n[c] f[c, d]

The results should be as follows:

OverBarrierSimplify[x]
(* f[a, b] h[a] *)
OverBarrierSimplify[y]
(* n[c] f[c, d] *)
OverBarrierSimplify[x+y]
(* f[a, b] h[a] + n[c] f[c, d]*)

In this example, the desired simplification function OverBarrierSimplify is expected to figure out that in n[a] f[b, a] + f[a, b]+ n[c] f[c, d] 1st and 2nd terms can be combined together, whereas the 3rd term does not need any simplification. This leads to the reduction of the LeafCount from 16 to 13.

Afterthoughts

I am eager to learn about a good approach for this case. However, I am aware that this is a difficult conceptual problem. This becomes obvious while doing searches for similar questions. Here is incomplete list:

• Rules for the Gamma function
• Transformation rules are not applied
• A need to separately consider symbolic and numeric cases for the transformation
• Using assumptions along with TransformationFunctions
• Custom complexity function

These are useful discussions, but how do you approach this kind of problems in general? Do you try to "help" Mathematica and do some intermediate simplifications by hands (when they are obvious)? Or do you completely avoid Simplify and use purpose-tuned rules for ReplaceAll? Any feedback is greatly appreciated!

Answer 1

Snooping into what Simplify is doing (using definitions of x, y, r1, and r2 from OP):

snoop = Function[f,
  Function[, Echo[f@##, HoldForm@HoldForm@f@##], HoldAllComplete],
  HoldAllComplete
];

Simplify[x + y,
  TransformationFunctions -> {Automatic, snoop@r1, snoop@r2},
  ComplexityFunction -> snoop@Simplify`SimplifyCount
]

seems to suggest that indeed it's backing up immediately after transformation function increases complexity. I don't know how to change that behavior, so in similar situations I try to use "manual" replacements.

First thing we need to do is to gather all replacement rules, for given relation, that can replace smallest possible relevant expressions, found in relation, that are explicitly expressible by other expressions from that relation. With relation from OP $\forall_x \; n(x) + h(x) = 1$ we can explicitly express: $n(x) = 1 - h(x)$ and $h(x) = 1 - n(x)$. With $\forall_{x,y} \; f(x,y) = -f(y,x)$ relation there's nothing we need to change, although for this relation one might prefer to choose canonical ordering of arguments instead of trying to Simplify it, for the sake of example let's include it. So our full set of rules will be:

rules = {n[x_] :> 1 - h[x], h[x_] :> 1 - n[x], f[x_, y_] :> -f[y, x]};

Now we could take an expression, find all occurrences of left hand side patterns from our rules, and create list of possible expressions obtained by replacing only single occurrence. Then we could simplify them, choose simplest results and repeat this procedure on them until simplest results stop changing. This procedure can be automated by following function:

singleReplaceList = Function[{expr, rules, levelspec, n},
  Join @@ (
    Function[r,
      ReplacePart[expr, # -> Replace[Extract[expr, #], r]]& /@
        Position[expr, First@r, levelspec, n]
    ] /@
      rules
  )
];
selectGenerator // ClearAll
selectGenerator[Infinity, _] = Identity;
selectGenerator[n_, f_] := f[#, UpTo@n]&;

replaceSimplify // ClearAll
replaceSimplify // Options = {
  TimeConstraint -> 300,
  MaxIterations -> Infinity,
  ComplexityFunction -> Automatic,
  Simplify -> Automatic,
  "SimplifiedSelect" -> Automatic,
  "MaxSimplified" -> 2,
  "ReplacementsSelect" -> Automatic,
  "MaxReplacements" -> 10^3,
  "MaxReplacementsPerRule" -> Automatic
};
replaceSimplify[expr_, rules : {(_Rule | _RuleDelayed)...}, opts : OptionsPattern[]] :=
  Module[{pool = {{None, expr}}, maxRepl = OptionValue@"MaxReplacements", replSimpl},
    TimeConstrained[
      With[
        {
          complSimpl = With[{cf = Replace[OptionValue@ComplexityFunction, Automatic :> Simplify`SimplifyCount]},
            {cf@#, #}& @* Replace[OptionValue@Simplify, Automatic :>
              With[{simplOpts = FilterRules[{opts, Options@replaceSimplify}, Options@Simplify]},
                Simplify[#, simplOpts]&
              ]
            ]
          ],
          replSelect = Replace[OptionValue@"ReplacementsSelect", Automatic :> selectGenerator[maxRepl, RandomSample]],
          maxReplPerRule = Replace[OptionValue@"MaxReplacementsPerRule", Automatic :> maxRepl],
          simplSelect = Replace[OptionValue@"SimplifiedSelect", Automatic :> selectGenerator[OptionValue@"MaxSimplified", Take]]
        }
        ,
        call : replSimpl@{_, x_} := (
          call = {};
          complSimpl /@ replSelect@singleReplaceList[x, rules, {0, Infinity}, maxReplPerRule]
        );
        pool = {complSimpl@expr};
        FixedPoint[
          (pool = simplSelect@Union[Sequence @@ (replSimpl /@ #), pool])&,
          pool,
          OptionValue@MaxIterations
        ]
      ]
      ,
      Replace[OptionValue@TimeConstraint, {tot_, _} :> tot]
    ];
    pool[[1, 2]]
  ]

Let's see how replaceSimplify works, on some test expressions, compared to Simplify with appropriate TransformationFunctions:

{
  #,
  Simplify[#, TransformationFunctions -> Prepend[ReplaceAll /@ rules, Automatic]] & /@ #,
  replaceSimplify[#, rules] & /@ #
} &@
  {1 - n[a], 1 - h[a], -f[b, a], f[a, b] - f[a, b] n[a], f[a, b] + f[b, a] n[a] + f[c, d] n[c], 2 f[a, b] + 2 f[b, a] n[a] + 3 f[c, d] n[c] - f[a, b] f[c, d] n[c] - f[b, a] f[c, d] n[a] n[c]} // TableForm
(* 1 - n[a]    1 - h[a]    -f[b, a]    f[a, b] - f[a, b] n[a]    f[a, b] + f[b, a] n[a] + f[c, d] n[c]    2 f[a, b] + 2 f[b, a] n[a] + 3 f[c, d] n[c] - f[a, b] f[c, d] n[c] - f[b, a] f[c, d] n[a] n[c]
   h[a]        n[a]        f[a, b]     f[a, b] h[a]              f[a, b] + f[b, a] n[a] + f[c, d] n[c]    3 f[c, d] n[c] + f[a, b] (2 + f[d, c] n[c]) + f[b, a] n[a] (2 + f[d, c] n[c])
   h[a]        n[a]        f[a, b]     f[a, b] h[a]              f[a, b] h[a] + f[c, d] n[c]              3 f[c, d] n[c] + f[a, b] h[a] (2 + f[d, c] n[c]) *)

By default replaceSimplify uses only two simplest simplified replacement results in next step. If this doesn't give desired result we can take obtained expression and run replaceSimplify on it with increased "MaxSimplified" option:

(res = replaceSimplify[2 x + 3 y - x y + x^2 - y^3 + x y^2, rules]) // AbsoluteTiming
res // Simplify`SimplifyCount
(res = replaceSimplify[res, rules, "MaxSimplified" -> 8]) // AbsoluteTiming
res // Simplify`SimplifyCount
(* {0.553314, f[a, b]^2 h[a]^2 + 3 f[c, d] n[c] + f[d, c]^3 n[c]^3 + f[a, b] h[a] (2 + f[d, c] n[c] + f[c, d]^2 n[c]^2)} *)
(* 52 *)
(* {0.526003, 3 f[c, d] n[c] + f[d, c]^3 n[c]^3 + f[a, b] h[a] (2 + f[a, b] h[a] + f[d, c] n[c] + f[c, d]^2 n[c]^2)} *)
(* 48 *)

replaceSimplify is expected to be much slower than vanilla Simplify on larger expressions. If it's using too much memory reducing value of "MaxReplacements" may help.

Answer 2

Try this:

((n[a] f[b, a] + f[a, b] /. {n[u_] -> 1 - h[u]}) // Expand) /. 
 f[b, a] -> -f[a, b]

(*   f[a, b] h[a]   *)

Have fun!