5
$\begingroup$

First-order ODEs of the form $$P(x,y)dx+Q(x,y)dy=0$$ are said to be exact if $$\frac{\partial P}{\partial y}=\frac{\partial Q}{\partial x}.$$ However, if the partial derivates are different, it can be found an integrating factor $\mu(x,y)$ such that, when multiplying the original ODE with it, it becomes exact: $$\frac{\partial (\mu P)}{\partial y}=\frac{\partial (\mu Q)}{\partial x}.$$

Is there any way to automatically finding an integrating factor for first-order ODEs with Mathematica?

$\endgroup$
1
  • $\begingroup$ Can you get one from the solution? $\endgroup$
    – Michael E2
    Jul 9, 2020 at 13:38

2 Answers 2

6
$\begingroup$

I did not know there was a build-in function to find integrating factor, but here is a basic implementation from scratch.

In this, the form $m(x,y) dx + n(x,y) dy = 0$ is used, as I am more used to it. So in place of our $P$ it is now $m(x,y)$ and in place of your $Q$ it is now $n(x,y)$

You call it as

m = -2*Exp[2*x]*x^3 - 2*Exp[y];
n = Exp[y] x;
getIntegratingFactor[m, n, x, y]

Mathematica graphics

m = Exp[y] - x; n = Exp[y]*(Exp[y] + x);
getIntegratingFactor[m, n , x, y]

Mathematica graphics

m = 2*x*y; n = -2*x^2 + y^2;
getIntegratingFactor[m, n , x, y]

Mathematica graphics

m = -(-x y - 1);
n = (4 x^3 y - 2 x^2);
getIntegratingFactor[m, n , x, y]

Mathematica graphics

m = x^2 + y^2 + 2 x; n = 2 y;
getIntegratingFactor[m, n , x, y]

Mathematica graphics

code

getPatterns[expr_, pat_] := 
  Last@Reap[expr /. a : pat :> Sow[a], _, Sequence @@ #2 &];

getIntegratingFactor[m_, n_, x_, y_] := Module[{a, b, r, s, mu, t},
   (*find integrating factor for m*dx+n*dy=0*)
   (*version 1.0 alpha, July 9, 2020 10 AM*)
   
   If[Simplify[D[m, y] - D[n, x] == 0],
    Return["It is allready exact, no integrating factor needed",Module]];
   
   a = Simplify[(D[m, y] - D[n, x])/n];
   If[Length[getPatterns[a, y]] == 0,
    Return[Row[{"Integrating factor is mu=",Exp[Integrate[a, x]]}], Module]];
   
   b = Simplify[(D[n, x] - D[m, y])/m];
   If[Length[getPatterns[b, x]] == 0,
    Return[Row[{"Integrating factor is mu=", Exp[Integrate[b, y]]}],Module]];

   r = (D[n, x] - D[m, y])/(x*m - y*n);
   r = Simplify[r];
   r = r /. (x^s_.*y^s_.) -> t^s;
   If[Length[getPatterns[r, x]] == 0 && Length[getPatterns[r, y]] == 0,
     mu = Simplify[Exp[Integrate[r, t]]];
     mu = mu /. t -> (x*y);
     Return[Row[{"Integrating factor is mu=", mu}], Module]
    ,
     Print["Unable to find integrating factor"];
   ]
   ];

The helper function getPatterns used above is thanks to Carl Woll.

To help clarify the code above, here is Fortran-like flow chart of the algorithm. Drawing done using ipe Latex drawing diagram.

enter image description here

Bug reports are always welcome.

$\endgroup$
2
  • $\begingroup$ It is a good code. Its very compact but difficult to read for novice. For example, I could not see step where t=xy is defined and then testing A and B for t. I got it later. I think it can be expanded for other values of t. Good work. $\endgroup$
    – Aschoolar
    May 26 at 20:40
  • $\begingroup$ @Aschoolar I could not see step where t=xy is defined and then testing A and B for t it is done in this code r = r /. (x^s_.*y^s_.) -> t^s; near 70% down the code with the code that follows that. That handle the last case in the chart (the R case) $\endgroup$
    – Nasser
    May 26 at 22:44
5
$\begingroup$

There's an internal function used to construct the solution. If you don't want to reverse-engineer an integrating factor from the solution, you can use DSolve`DSolveFirstOrderODEDump`IntegratingFactor:

Block[{P, Q},
 P = (Cos[x] - Sin[x]) Sin[y];
 Q = Cos[x] Cos[y];
 mu = DSolve`DSolveFirstOrderODEDump`IntegratingFactor[
   Q, P, -D[P, y] + D[Q, x], x, y];
 mu -> D[mu*P, y] - D[mu*Q, x] // Simplify
 ]
(*  E^x -> 0  *)

Returns $Failed when it fails:

Block[{P, Q},
 P = E^(x y);
 Q = Cos[x] Cos[y];
 mu = DSolve`DSolveFirstOrderODEDump`IntegratingFactor[
   Q, P, -D[P, y] + D[Q, x], x, y];
 mu -> D[mu*P, y] - D[mu*Q, x] // Simplify
 ]
(*  $Failed -> $Failed (E^(x y) x + Cos[y] Sin[x])  *)
$\endgroup$
2
  • $\begingroup$ thanks for pointing out this function. I am still trying to figure how to use it to reproduce the results I have. But will try to figure it out. $\endgroup$
    – Nasser
    Jul 9, 2020 at 15:07
  • 1
    $\begingroup$ @Nasser If I substitute P=m, Q=n, I get all your results. $\endgroup$
    – Michael E2
    Jul 9, 2020 at 17:22

Your Answer

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

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