Why does InverseLaplaceTranform give an incorrect solution?

Question

Bug introduced in 9.0.1, fixed in 11.1.

---

I was trying to answer this question, but I ended up asking this.

Method 1

A workaround is to differentiate the ODE once,

y'[x] == x + Integrate[y[r], {r, 1, x}]
Newode=D[%, x]

y''[x] == 1 + y[x]

Now use DSolve to solve the new ODE.

DSolve[{Newode, y'[0] == 0, y[0] == 1}, y[x], x]//FullSimplify

{{y[x] -> -1 + 2 Cosh[x]}}

Method 2

We can also use LaplaceTransform to solve the integro-differential equation:

ode = y'[x] == x + Integrate[y[r], {r, 1, x}]
LaplaceTransform[ode, x, s]
Solve[%, LaplaceTransform[y[x], x, s]]
InverseLaplaceTransform[%, s, x]

{{y[x] -> -1}}

Here are my questions:

1. Why do the two solutions differ?
2. Why is the LaplaceTranform solution not correct?
3. The above two solutions are incorrect. Why?

Answer 1

When you take derivative of $y^{\prime}\left( x\right) =x+\int_{1} ^{x}y\left( r\right) dr$, then using Leibniz rule of integration

\begin{align*} y^{\prime\prime}\left( x\right) & =1+\overset{0}{\overbrace{\int_{1} ^{x}\frac{\partial y\left( r\right) }{\partial x}dr}}+\left. y\left( r\right) \right\vert _{r=x}\frac{dx}{dx}-\overset{0}{\overbrace{\left. y\left( r\right) \right\vert _{r=1}\frac{d1}{dx}}}\\ & =1+y\left( x\right) \end{align*}

This ODE has the solution $-1+2\cosh\left( x\right) $ for the initial conditions given. But the same solution is obtained if the lower limit of integration is zero. i.e. using $y^{\prime}\left( x\right) =x+\int_{0} ^{x}y\left( r\right) dr$. This is due to Leibniz rule. Any lower limit constant will give same solution. So using this form, with zero lower limit, allows Mathematica to give the correct answer, since now it know about the rule $\mathcal{L} \left( \int_{0}^{t}y\left( \tau\right) d\tau\right) =\frac{Y\left( s\right) }{s}$. But when you had $1$ as lower limit, it could not use this rule. Hence using lower limit of zero, gives the correct answer

ode = y'[x] == x + Integrate[y[r], {r, 0, x}];
LaplaceTransform[ode, x, s];
Solve[%, LaplaceTransform[y[x], x, s]];
InverseLaplaceTransform[%, s, x] /. {y[0] -> 1, y'[x] -> 0} // ExpToTrig

Answer 2

First of all, as mentioned in this comment, the answer for the 3rd question is, the new added boundary for Newode should be

newbc = ode /. x -> 1
(* y'[1] == 1 *)

Then we can get the desired result:

ode = y'[x] == x + Integrate[y[r], {r, 1, x}];
bc = y[0] == 1;

sol = First@DSolve[{Newode, newbc, bc}, y, x]
(* {y -> Function[{x}, (E^-x (-E + 2 E^2 - E^x + 2 E^(2 x) - E^(2 + x) + E^(1 + 2 x)))/(
   1 + E^2)]} *)

There seems to be no direct way to deduce the value of y'[0], so LaplaceTransform isn't suitable for the original problem linked in the question .

Then let's deal with the first 2 questions. I think the initial value problem formed by the fake initial condition y'[0] == 0 is interesting, because it seems to expose a bug of InverseLaplaceTransform. It fails to correctly handle unevaluated Integrate with an assumption related to s in this case:

teqn = LaplaceTransform[ode, x, s];
tsol = Solve[teqn, LaplaceTransform[y[x], x, s]][[1, 1, -1]] /. y@0 -> 1 // Apart
(* (1 + s^2)/(s (-1 + s^2))- Integrate[y[x], {x, 0, 1}, Assumptions-> s > 0]/(-1 + s^2) *)

If the assumption inside Integrate is removed, then InverseLaplaceTransform will evaluate correctly:

solgeneral = InverseLaplaceTransform[(1 + s^2)/(s (-1 + s^2)) - 
  Integrate[y[x], {x, 0, 1}]/(-1 + s^2), s, x]
(* -1 + (1/2)*E^x*(2 - Integrate[y[x], {x, 0, 1}]) + 
   ((1/2)*(2 + Integrate[y[x], {x, 0, 1}]))/E^x *)

Notice

ode /. x -> 0 /. y'[0] -> 0
(* 0 == Integrate[y[r], {r, 1, 0}] *)

We find the final solution:

solgeneral /. HoldPattern@Integrate[__] :> 0 // Simplify
(* -1 + E^-x + E^x *)

It's the same as the result of Method 1.

---

Update

In v8.0.4 InverseLaplaceTransform will return unevaluated, also, here's the response from WRI:

…… It does appear that InverseLaplaceTransform is not behaving as expected with the above function. I have forwarded an incident report to our developers with the information you provided.……

So I think it's safe to call the current behavior of InverseLaplaceTransform a bug.

---

Update 2

In v11.1 InverseLaplaceTransform returns unevaluated again. Since the incorrect result is no longer returned, I think we can say the bug is fixed, though in a not that perfect way…