Factoring a separable integral with a product of independent integrands

Question

I would like to input

someFunction[Integrate[p[x] p[y], {x, -1, 1}, {y, -1, 1}]]

and to get the following output

(Integrate[p[x], {x, -1, 1}])^2

Is this possible? (Function p is undefined, I am interested in a symbolic expression only.)

UPD: The question is not that how to define someFunction for this particular case. I am interested in a method which will work for any multiple integral: if there is a product function inside such integral then the answer has to be factorized as well.

Best wishes,

Dmitri.

Answer 1

This will work with any number of independent integrants. Define:

repl[l_] := # /. Thread[#[[All, 1]] -> Table[x, {Length[#]}]] &[l]

inTfaC[int_] := Times @@ MapThread[Integrate[#1, #2] &, 
   repl /@ {First[#], Rest[#]} &[
    int /. {Integrate -> List, Times -> List}]]

Now verify:

test = Integrate[p[x] p[y] q[z] r[s] r[u], 
{x, -1, 1}, {y, -1, 1}, {z, a, 2}, {s, 3, b}, {u, 3, b}]

inTfaC[test]

Answer 2

Adapting linearExpand from my answer to How to do algebra on unevaluated integrals?, we can come up with some transformations to factor separable multiple integrals. The function someFunction internally deals with and returns Inactive integrals, which can be evaluated with Activate, if appropriate or desired.

Examples

someFunction[Integrate[p[x] p[y], {x, -1, 1}, {y, -1, 1}]]

someFunction[Integrate[p[x] p[y], {x, -1, 1}, {y, -1, 100}]]

someFunction[Integrate[p[x] p[y], {x, -1, 1}, {y, -1, x}]]

someFunction[Integrate[p[x] (p[y] + q[y]), {x, -1, 1}, {y, -1, 1}]^2]

(* Vitaliy Kaurov's example *)
test = Integrate[p[x]*p[y]*q[z]*r[s]*r[u],
 {x, -1, 1}, {y, -1, 1}, {z, 0, 2}, {s, 3, 10}, {u, 3, 10}];
someFunction[test]
% /. changeVar[x]

someFunction[Integrate[q[x], {x, 0, 2}] + test]

Outline

We need some auxilliary functions.

• iterated converts a multiple integral into an iterated integral.
• linearExpand applies linearity properties to an integral, distributing the integral over a sum and factoring out constants.
• factorConstants recursively factors constants out of nested integrals using linearExpand.
• changeVar changes the variable of integration.

These can be applied by Simplify via the TransformationFunctions option. The trick, and it can be tricky, is devising a ComplexityFunction that will prefer a result in a factored out form. Note that we use changeVar only to combine two integrals that are equivalent. At the end we can change the variable in the independent integrals to the same one.

Code dump

ClearAll[linearExpand, iterated, factorConstants, changeVar, someFunction];

linearExpand[e_, x_, head_] := 
  e //. {op : head[arg_Plus, __] :> Distribute[op], 
    head[arg1_Times, var_, rest___] /; ! FreeQ[var, x] :> 
     With[{dependencies = Internal`DependsOnQ[#, x] & /@ List @@ arg1}, 
      Pick[arg1, dependencies, False] *
       head[Pick[arg1, dependencies, True], var, rest]]};

iterated[Integrate[f_, dom : {_Symbol, _, _} ..]] := 
  Fold[Inactive[Integrate], f, Reverse@{dom}];
iterated[Inactive[Integrate][f_, dom : {_Symbol, _, _} ..]] := 
  Fold[Inactive[Integrate], f, Reverse@{dom}];

factorConstants[Inactive[Integrate][f_, {v_Symbol, a_, b_}]] := 
  Inactive[Integrate][
    f /. j : Inactive[Integrate][_, _] :> factorConstants[j], {v, a, b}] /.
      i : Inactive[Integrate][_, {v, _, _}] :> 
    linearExpand[i, v, Inactive[Integrate]];
factorConstants[x_] := x;

changeVar[Inactive[Integrate][f_, {v_, a_, b_}]] := 
  i : Inactive[Integrate][g_, {w_, a, b}] /; 
    Simplify[f == g /. v -> w] :> (i /. w -> v);
changeVar[v_] := 
  i : Inactive[Integrate][g_, {w_Symbol, _, _}] :> (i /. w /; FreeQ[g, Integrate] -> v);

someFunction[integral_] := 
  With[{v = 
     FirstCase[Hold[integral], 
      HoldPattern[(Integrate | Inactive[Integrate])[_, {x_, _, _}, ___]] :> x,
      Missing[], 
      Infinity]},
   (Simplify[
      integral /.
       {i : (Integrate | Inactive[Integrate])[_, {_Symbol, _, _}, {_Symbol, _, _} ..] :> 
          iterated@Inactivate[i, Integrate],
        i : (Integrate | Inactive[Integrate])[_, {_Symbol, _, _}] :>
          Inactivate[i, Integrate]},
      TransformationFunctions -> {Automatic,
        # /. i : Inactive[Integrate][f_, {_, _, _}] :> factorConstants[i] &,
        # /. changeVar /@ 
           DeleteDuplicates@
            Cases[#, Inactive[Integrate][f_, {v_, a_, b_}], Infinity] &},
      ComplexityFunction -> (LeafCount[#] +
          10 Count[{#}, 
            Inactive[Integrate][Inactive[Integrate][__], __], Infinity] +
          5 Count[{#}, 
            Inactive[Integrate][f_, __] /; ! FreeQ[f, Integrate], Infinity] -
          Count[{#}, Power[Inactive[Integrate][f_, __], _], Infinity] +
          Length@DeleteDuplicates@
            Cases[{#}, Inactive[Integrate][f_, {x_, _, _}] :> x, Infinity] &)
    ] /. changeVar[v]) /; FreeQ[v, Missing]
   ];
someFunction[expr_] := expr;

The first two terms after LeafCount of the ComplexityFunction penalize nested integrals, the first one penalizing nested integrals in which the constants have not been factored out. The next term, being subtracted, rewards gathering powers of integrals. The last term penalizes extra variables of integration, which drives Simplify to prefer changing variables to the same thing.

Answer 3

Maybe for this case use:

someFunction[int_] := (int /. _[y] :> 1/2)^2