The following two expressions are identical from a human point of view:

a = Integrate[f[x, y], x, y]

b = Integrate[Integrate[f[x, y], x], y]

The same goes here:

a = Integrate[f[x] + g[x], x]

b = Integrate[f[x], x] + Integrate[g[x], x]

However, Mathematica considers them different and FullSimplify[a - b] produces nasty no-zero results in both cases. This due to the fact that Mathematica "wraps" the inner integrate in b into () in the first case and then cannot thread + over Integrate in the second case.

Apparently there might be some edge cases when the order of integrations is important but I am not interested in those. The function under integral f[x, y] is a simplification of what's going on and so I cannot rely on its name. The reason why I need that is that I combine several such equations and together they are supposed to produce some conservation law - basically the sum of several of such integrals must come out as exact zero. However, some of them are defined as triple integrals, some as double integrals, some as integrals by x first and then by y and some as the other way around (or something like that). I cannot define all of them the same way because I use intermediate pieces for other calculations.

How can I explain to Mathematica that a == b or, better, convert b into a.

In a sense, I need standard associativity and distributivity algebraic rules applied to integrals and the ability to convert an integral over an integral into a multiple integral (up to triple integral), e.g.:

integrate[integrate[a * b[x, y] * c[x], x] - d * e[y], y]

should unwrap into:

a * integrate[b[x, y] * c[x], x, y] - d * integrate[e[y], y]

and so on and so forth. Of course, any of these a, b, c are not just a, b, c but could be some complicated expressions.

Following the advice below, I am fine introducing an integrate instead of using built in Integrate until the last moment, except that defining the transformation rules for this integrate seems elusive.



3 Answers 3


For you first pair of integrals, let us introduce the following rule:

rule = Integrate[Integrate[g_[x, y], {x, -Infinity, Infinity}], {y, -Infinity, 
    Infinity}] :>Integrate[g[x, y], {x, -Infinity, Infinity}, {y, -Infinity, Infinity}];


a = Integrate[
  f[x, y], {x, -Infinity, Infinity}, {y, -Infinity, Infinity}];

b = Integrate[
  Integrate[f[x, y], {x, -Infinity, Infinity}], {y, -Infinity,Infinity}];


enter image description here


b /. rule

returns the integral without the internal parentheses:

enter image description here


a - b /. rule

(*  0  *)

Concerning your second example, one can apply the function Distribute to the first integral:

a = Distribute@Integrate[f[x] + g[x], {x, -Infinity, Infinity}];
b = Integrate[f[x], {x, -Infinity, Infinity}] + 
  Integrate[g[x], {x, -Infinity, Infinity}];

returning the following:

enter image description here


a - b

(*  0  *)

Edit To address your question:

To transform the integrals of your last example the rule should be slightly modified to adopt it to the integrals in question:

rule = Integrate[Integrate[Times[g1_[x, y], g2_[x]], x], y] :> 
  Integrate[g1[x, y]*g2[x], x, y];

After that, everything can be done the same way as previously:

expr1=Integrate[Integrate[a*b[x, y]*c[x], x] - d*e[y], y]

    Distribute[expr1] /. rule

enter image description here

A minor note. If you did not previously assign the functions b[x,y], c[x] and so on, it is not necessary to use the lowercase integrate. The expression anyway will not be integrated.

Another story, if you have previously assigned these functions, but want to prevent the integration. It is possible, of course, to use the lowercase integrate as you did in the question.

There is, however, another possibility that I prefer. It is to inactivate the integrals:

expr2 = Inactivate[expr1, Integrate]

The advantage here is that the integrals are shown on the screen in the traditional way. To indicate that they are inactivated they are shown blended. Like this, it is easier to look at the expression. Try it! Then the application of the same procedure:

expr3=Distribute[expr2] /. rule

returns the same result as above. Later, you can activate it, if needed:


In general, my experience shows that the rule should be adapted to fit the form of integrals you deal with.

Have fun!

  • $\begingroup$ Both works. Thanks! $\endgroup$ Jan 16, 2022 at 13:34
  • $\begingroup$ I've updated the question with the result that I got after applying your transformations. $\endgroup$ Jan 16, 2022 at 13:52
  • $\begingroup$ @Konstantin Konstantinov There you need to work a bit further. Please post your equation in Mma form. $\endgroup$ Jan 16, 2022 at 18:50
  • $\begingroup$ It is far worse than I thought. Give me a few days, please. $\endgroup$ Jan 18, 2022 at 20:44
  • $\begingroup$ I've updated the question with the clarification. $\endgroup$ Jan 18, 2022 at 22:20

There's Flatten:

Flatten[Integrate[Integrate[f[x, y], x], y]]
(*  Integrate[f[x, y], x, y]  *)


  1. Expanding integrals is answered here (linearExpand): How to do algebra on unevaluated integrals?

  2. Combining integrals is (almost) answered here (linearCombine): Pushing Mathematica's FullSimplify to a global complexity minimum

linearCombine needs to be tweaked for Integrate, something like this, maybe:

  SetAttributes[f, NumericFunction];
  constantQ = 
   NumericQ[# /.
      s_Symbol /; MemberQ[Attributes[s], Constant] :> f[0]] &];
linearCombine[e_Plus, func_, constants_List : {}] := Block[constants,
   SetAttributes[#, Constant] & /@ constants;
      Replace[e, {c_?constantQ*func[a_] :> c*a, func[a_] :> a}, 1]
      ] /; MatchQ[e, _[((_?constantQ)*func[_] | func[_]) ..]]
linearCombine[e_Plus, func_, constants_List : {}] := e;

linearCombine[Integrate[f[x], x] + 2 a Integrate[g[x], x], 
 HoldPattern@Integrate[#, x] &, {a}]
(*  \[Integral](f[x] + 2 a g[x]) \[DifferentialD]x  *)

You could try to add definitions to "Integrate". However, the danger to upset something is big. Therefore, I would advice to define a new integration routine that then calls the real integration.

You can then define any definition you like for the new function without upsetting something.

Let's call the new routine "integrate", lowercase.:

integrate[a_[x_], {x_, x1_, x2_}] = Integrate[a[x], {x, x1, x2}];

Now e.g. if the intergand is a sume of 2 terms we can change it into 2 integrals like:

integrate[a_[x_] + b[x_], {x_, x1_, x2_}] = 
 Integrate[a[x], {x, x1, x2}] + Integrate[b[x], {x, x1, x2}];

integrate[a[x] + b[x], {x, 0, 1}]

enter image description here

Or if we have a nested integral, we may change it to a double integral like:

integrate[integrate[a_[x_, y_], {x_, x1_, x2_}], {y_, y1_, y2_}] = 
  Integrate[Integrate[a[x, y], {x, x1, x2}], {y, y1, y2}] ;

 Integrate[f[x, y], {x, -Infinity, Infinity}], {y, -Infinity, 

enter image description here

  • $\begingroup$ 1. What about double integral (the first example)? 2. a[x] and b[x] could be anything, not just a and b and so I cannot "hard code" them. The new integrate should work for any input and thread + and - over integrate. Also whatever does not depend from integration variable should be "extractable" out of integral. Otherwise Integrate[a * b[x], x] won't be equal to a * Integrate[b[x], x]. $\endgroup$ Jan 16, 2022 at 13:11
  • $\begingroup$ Note you may add any definition you like. I added a definition for a nested integral. $\endgroup$ Jan 16, 2022 at 13:32

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.