I solve particular boundary-value-problem for ODE with NDSolve “Shooting” method. The case is that solution is attained very slow, that seems to mean that boundary-value-problem which is supposed to be well-defined enough, in fact, is ill at some place. So i try figure out. First step is to see concrete values of a produced approximations to a solution. What language constructs should i use for that?

Simple example. Suppose we consider particle motion in vertical plane. We throw our particle from initial point with coordinates {x0,y0} := {0,0} and initial trajectory angle 45 degrees. And try to achieve point with coordinates {x1,y1} := {1,0} by varying particle initial velocity. We don’t know two things here: initial velocity and a duration of a motion. Here is how that toy problem can be presented and solved in *mathematica*:

`gravity = 10; bvpsol = NDSolve[ { { (* ODE system (5th order) *) {x''[u] / c[u]^2 == 0, y''[u] / c[u]^2 == -gravity, c'[u] == 0}, (* boundary conditions (5 items) *) {x[0] == y[0] == 0, x[1] == 1, y[1] == 0, x'[0] == y'[0]} } }, {x[u], y[u], c[u]}, u, Method -> {"Shooting", "StartingInitialConditions" -> {x[0] == y[0] == 0, x'[0] == y'[0] == 1, c[0] == 1}}] // Flatten; {dxdu, dydu} = D[{x[u], y[u]} /. bvpsol, u]; {vx0, vy0} = ({dxdu, dydu} / c[u] /. bvpsol) /. {u -> 0}; duration = c[u] /. bvpsol /. {u -> RandomReal[]}; ivpsol = NDSolve[{ (* ODE system *) {x''[t] == 0, y''[t] == -gravity}, (* initial values *) {x[0] == y[0] == 0, x'[0] == vx0, y'[0] == vy0} }, {x[t], y[t]}, {t, 0, duration}]; ParametricPlot[{x[t], y[t]} /. ivpsol, {t, 0, duration}, GridLines -> Automatic, AspectRatio -> 1/2] `

**Question:** Now what options or language construct should i use to see approximations which are produced NDSolve while solving boundary-value-problem?