Why Mathematica gives me different results depending on the order I set the parameters when I solve this differential equation?

Question

I am asking Mathematica to solve a differential equation using DSolve, with two different parameters, q and m, and then evaluate at a point.

So I write the following lines of code:

x[y_] = w[y] /. 
    First[DSolve[{w'[y] == (1 - (2 m y^2)/(y^2 + q^2)^(3/2))^-1}, 
      w[y], y]] /. C[1] -> 0;
m = 1.; q = 1/2;

x[2]

Where the output is (I just chose a random point to evaluate this function, it does not really matter where):

-4.60489 - 1.14725*10^-19i .

However, If I first set the parameters (using the same values as before) and then let Mathematica solve it:

m = 1.; q = 1/2;

x[y_] = w[y] /. 
   First[DSolve[{w'[y] == (1 - (2 m y^2)/(y^2 + q^2)^(3/2))^-1}, w[y],
      y]] /. C[1] -> 0

x[2]

It yields:

Inactive[Integrate][(
 0.25 Sqrt[1.  + 4. K[1]^2] + 
  1. K[1]^2 Sqrt[1.  + 4. K[1]^2])/(-4. K[1]^2 + 
  0.25 Sqrt[1.  + 4. K[1]^2] + 
  1. K[1]^2 Sqrt[1.  + 4. K[1]^2]), {K[1], 1, 2}]

which is an integral that cannot be analytically solved.

How come this happens? Why would the order that I define my constants change the outcome of how Mathematica handles this differential equation? What I find even weirder, is that (at least to me) it seems like it would be simpler to set the parameters first and then solve it.

Answer 1

Both solutions are correct, and the behavior of Mathematica is (to me) not entirely surprising. To begin, plot the right side of the ODE.

ComplexPlot3D[(1 - (2 m y^2)/(y^2 + q^2)^(3/2))^-1 /. {m -> 1, q -> 1/2}, {y, 3}, 
    AxesLabel -> {Re[y], Im[y], None}, ImageSize -> Large, LabelStyle -> {15, Black, Bold}, 
    PlotRange -> {0, 6}]

Four poles, on the real axis, and two branch cuts, on the imaginary axis, are visible. On the other hand, for m -> 1/2, the poles move off the real axis.

With m and q unspecified, as in the first case in the Question, Mathematica cannot know the location of the poles or where in the complex plane the desired values of y lie. So, it appears to make unknown assumptions and proceed. With m and q specified, it knows were the poles are but not on which side of the branch cuts the desired values of y lie. Again, it makes assumptions and proceeds. I do find it surprising, however, that Mathematica makes different assumptions for the second case in the Question and the case considered in Ulrich Neumann's Answer, which differ only by inserting m -> 1. and m -> 1.

Now, address which solutions are correct.

s1f = DSolveValue[w'[y] == (1 - (2 m y^2)/(y^2 + q^2)^(3/2))^-1, w, y];
Plot[Evaluate[{w'[y] /. w -> s1f, (1 - (2 m y^2)/(y^2 + q^2)^(3/2))^-1} 
    /. {m -> 1, q -> 1/2} // Chop], {y, -5, 5}, ImageSize -> Large,
    PlotStyle->{Automatic, Dashed}, LabelStyle -> {15, Black, Bold}]

Visibly, the left and right sides of the ODE are equal. (Comparing the two with FullSimplify would be preferable, but I got tired of waiting for it to complete.) Doing the same for the second case in the Question and for Ulrich Neumann's Answer show that they too satisfy the ODE.

Answer 2

It seems to be dangerous to use y one both sides of the function definition x[y].

Perhaps this approach is what you're looking for

W[m_, q_] :=DSolveValue[{w'[y] == (1 - (2 m y^2)/(y^2 + q^2)^(3/2))^-1}, w , y] 
W[1,1/2][2]
(*(-3.21859 + 0. I) + C[1]*)