# Numerically solve a non-differential equation with InterpolatingFunction as output

I have some a function that I need to solve numerically over an entire domain (Ideally I would like an InterpolatingFunction object as output), which does not include any derivatives. I have found a way to do it using NDSolve, however I feel like there should be a better solution. Here is a simple example to show my NDSolve method:

Solving: $$y^2+y=x$$

The solutions: $$y(x) = \frac{1}{2}\left(-1 \pm \sqrt{1+4x}\right)$$

Numerically Solving using NDSolve:

NDSolve[D[y[x]^2 + y[x], x] == D[x, x] && y^2 + y == 0, y[x], {x, -5, 5}] //Quiet

This operation gives me two InterpolatingFunctions that, when plotted, match the analytical solutions, but feels very hacky. I have also noticed that for some functions it does not produce solutions for the entire function domain. Is there a better way to do this?

• I do what you did. Or use a DAE: Block[{F, x0 = 0}, F = Function[{x, y}, y^2 + y - x]; NDSolve[{F[x, y[x]] == 0, u'[x] == 1, u[x0] == x0, F[x0, y[x0]] == 0}, y, {x, -5, 5}] ] Mar 1 at 3:47
• @MichaelE2 Good to know! Do you know of any way to force NDSolve to cover the entire given domain? Sometimes when it hits large discontinuities, its output will not cover the entire requested domain. Mar 1 at 3:51
• I don't have a general way to deal with discontinuities. My first thought would be to use WhenEvent for finite jumps -- NDSolve does this automatically for known functions like Piecewise or UnitStep etc. Poles/asymptotes are more difficult. Mar 1 at 4:28
• Is your example oversimplified? Theequation is an ordinary equation not an differential equation. It can easily be solves using: Solve: Solve[y^2 + y == x, y] Mar 1 at 8:34