# Internal Shooting Method of NDSolve in combination with NDSolveReinitialize?

To explain my problem, I am trying to extend the BVP problem example from the help that illustrates how to use the shooting method of NDSolve:

sols = First[NDSolve[{x''[t] + Sin[x[t]] == 0, x[0] == x[10] == 0}, x, t,
Method -> {"Shooting", "StartingInitialConditions" -> {x[0] == 0, x'[0] == 1.75}}]];
Plot[Evaluate[x[t] /. sols], {t, 0, 10}]


I managed to split up the problem as following:

ndssdata = First[NDSolveProcessEquations[{x''[t] + Sin[x[t]] == 0, x[0] == 0, x[10] == 0},
x, t,
Method -> {"Shooting", "StartingInitialConditions" -> {x[0] == 0, x'[0] == 1.75}}]];

NDSolveIterate[ndssdata, 10];
ndsol1 = NDSolveProcessSolutions[ndssdata];
Plot[Evaluate[x[t] /. ndsol1], {t, 0, 10}]


The solution is the same and everything works as expected:

However, when I try to change the boundary conditions with NDSolveReinitialize it complains about a restriction that initial conditions should be specified at a single point, which does not make sense for a BVP:

newstate = First[NDSolveReinitialize[ndssdata, {x[0] == 0, x[10] == 1}]]


NDSolveReinitialize::ndcinit: Initial conditions should be specified at a single point.

I was hoping to speed up my problem, in which I plan solving the same set of equations with the internal shooting method many times with varying boundary conditions. I know how to do it with my own FindRoot shooting approach, but this would be more elegant. Does anyone know how to do this?

Unfortunately the "shooting" is done by NDSolveProcessEquations, which converts the BVP into an IVP. Thus NDSolveReinitialize[ndssdata] is basically operating on an IVP, and you won't be able to approach your problem in this way. If we examine the state data, we see that the initial conditions are already there.

{sdb, sdf} = ndssdata@"SolutionData"
(*
{{0., None, {0., 1.87817}, {1.87817, 0.}, {}, {}, {}, {}},
{0., None, {0., 1.87817}, {1.87817, 0.}, {}, {}, {}, {}}}
*)

NDSolveSolutionDataComponent[sdf, "TemporalDerivatives"]
(*  {1.87817, 0.}  *)


If we compare with the OP's computed solution ndsol1, we see they are the same, except for a tiny error that I presume is due to round-off in the InterpolatingFunction.

{x'[0], x''[0]} /. ndsol1
(*  {1.87817, -1.0842*10^-19}  *)


If you want further evidence, then we can trace the call to NDSolveShootingImplementationShootForInitialConditions that implements the shooting method. We see it is called when processing the equations:

Trace[
ndssdata =
First[NDSolveProcessEquations[
{x''[t] + Sin[x[t]] == 0, x[0] == 0, x[10] == 0},
x, t,
Method -> {"Shooting", "StartingInitialConditions" -> {x[0] == 0, x'[0] == 1.75}}]],
_NDSolveShootingImplementationShootForInitialConditions,
TraceInternal -> True
]
(* output of the call with arguments *)


But it is not called when the integration is iterated:

Trace[
NDSolveIterate[ndssdata, 10],
_NDSolveShootingImplementationShootForInitialConditions,
TraceInternal -> True
]
(*  {}  *)


One could check NDSolveProcessSolutions, too -- still nothing.