0
$\begingroup$

How can i extract some value from NDSolve. We know that gives us a InterpolatingFunction. But i need to use some specific value. For example this is the NDSolve: enter image description here

This is the equation used

un2 := -\[Alpha] y[t]/((x[t] - \[Beta]*Sign[Sin[2 \[Pi] t]])^2 + y[t]^2);
vn2 := \[Alpha] (x[t] - \[Beta]*Sign[Sin[2 \[Pi] t]])/((x[t] - \[Beta]*Sign[Sin[2 \[Pi] t]])^2 + y[t]^2);

The code of NDSolve:

r2[t_] := {x[t], y[t]}
solv2[x0_, y0_, tmax_] := NDSolve[{r2'[t] == {un2, vn2}, r2[0] == {x0, y0}}/. {\[Alpha] -> 1, \[Beta] -> 0.4}, r2[t], {t, 0, tmax}, StartingStepSize -> 0.01, Method -> {"FixedStep", Method -> {"ExplicitRungeKutta", DifferenceOrder -> 4}}, MaxStepFraction -> Infinity]

And i hope to get, for example, the value of x and y when t=0 (which is set in the initial condition). I was trying to do something like this: enter image description here

This is the code i used to extract value:

x01 = 0.5; y01 = 0.5;
sopx1 = x[0] /. solv2[x01, y01, 10][[1]]
sopy1 = y[0] /. solv2[x01, y01, 10][[1]]

but as you notice, it doesn't really work.

$\endgroup$
  • 1
    $\begingroup$ Welcome on MMA.SE! Could you write copy-and-pastable code so that it's easier to help you? $\endgroup$ – anderstood Oct 19 at 8:34
  • 1
    $\begingroup$ Thanks for letting me know. I have added all the parts that is necessary. If you need more information, plz let me know. And thanks for the help $\endgroup$ – Noah Ren Oct 19 at 9:24
  • $\begingroup$ Strongly related, if not duplicate: mathematica.stackexchange.com/a/67036/1871 $\endgroup$ – xzczd Oct 19 at 10:02
  • $\begingroup$ Thanks, I am going check on this. $\endgroup$ – Noah Ren Oct 19 at 10:21
  • $\begingroup$ I have seen the link. It is helpful. However, i can plot the graph or manipulate it. It is just, if i want to use the value from the list of the solution, individually. I actually can do this with the similar code i posted, but i can't do the same thing here $\endgroup$ – Noah Ren Oct 19 at 10:25
0
$\begingroup$

The quick fix is this:

x01 = 0.5; y01 = 0.5;
sopx1 = x[t] /. solv2[x01, y01, 10][[1]] /. t -> 0
sopy1 = y[t] /. solv2[x01, y01, 10][[1]] /. t -> 0

and the reason why this is necessary becomes clear when you look at the result of

solv2[x01, y01, 10]

which returns two replacement rules for x[t] and y[t]. These only work when the argument is exactly t, they do not match for arbitrary arguments as your code implicitly assumes. This is why you need to first replace x[t] and only after that replacement can insert the value for t.

There are some recommendations (mostly in accordance with Bill's answer) which make handling such situation somewhat easier, here is how I would write your code:

un2[α_,β_][x_, y_, t_] := -α*y[t]/((x[t] - β*Sign[Sin[2*Pi*t]])^2 + y[t]^2);
vn2[α_,β_][x_, y_, t_] := (
  α*(x[t] - β*Sign[Sin[2*Pi*t]])/((x[t] - β*Sign[Sin[2*Pi*t]])^2 + y[t]^2
);

solv2[α_,β_, x0_, y0_, tmax_] := Module[{x, y, t},
  NDSolveValue[{
      x'[t] == un2[α,β][x, y, t],
      y'[t] == vn2[α,β][x, y, t],
      x[0] == x0,
      y[0] == y0
    },
    {x, y}, {t, 0, tmax},
    Method -> {"FixedStep", Method -> {"ExplicitRungeKutta", DifferenceOrder -> 4}},
    StartingStepSize -> 0.01,
    MaxStepFraction -> Infinity
  ]
]

{xsol, ysol} = solv2[1., 0.4, 0.5, 0.5, 10]

xsol[0]
ysol[0]

and these are the recommended changes (depending on your situation, not all of them might be applicable or advantageous, but you should know about them and probably only deviate when having a good reason to do so):

  • make usage of parameters explicit instead of implicit. Wether or not to use subvalues to differentiate between "variables" and "parameters" is a matter of taste but somewhat common practive, I think. Avoiding implicit paramters helps to prevent your results becoming dependent on the order you evaluate your cells and avoids errors in general.

  • use NDSolveValue, it is newer than NDSolve and basically does the same thing. The difference is that it returns the solution as interpolating functions directly instead of rules. That is a bit less general and less conform to other similar functions but almost always what you want when solving differential equations.

  • do not use arguments in the definitions of the dependent variables, that will make NDSolve and NDSolveValue return the pure functions without arguments and helps to avoid the problems you have seen.

Think of the result as a mathematical function stored in a variable. You can use this variable just as if you would have defined a function with that name in almost every function that accepts functions as arguments. That is much more convenient for further analysis in almost all cases. Plotting functions, evaluating for given values or even getting derivatives is very convenient then:

Table[{xsol[t], ysol[t]}, {t, {0, 1, 2}}]
Plot[{xsol[t], ysol[t]}, {t, 0, 10}]

xsolderiv = Derivative[1][xsol]
Plot[xsolderiv[t], {t, 0, 10}]

or something a bit more advanced:

Show[
  Plot[{xsol[t], ysol[t]}, {t, 0, 10}],
  Graphics[{
    PointSize[Large], Red, 
    Point[
      Extract[NMinimize[{xsol[t], 0 < t < 9}, t], {{2, 1, 2}, {1}}]
    ]
  }]
]
$\endgroup$
  • $\begingroup$ Thanks for the explanation! $\endgroup$ – Noah Ren Oct 20 at 9:06
0
$\begingroup$

Try this more conventional form of writing a differential equation.

x0=0.5; y0=0.5; tmax=10; α=1.; β=0.4;
un2 = -α y[t]/((x[t]-β*Sign[Sin[2 Pi t]])^2+y[t]^2);
vn2 = α (x[t] - β*Sign[Sin[2 Pi t]])/((x[t] - β*Sign[Sin[2 Pi t]])^2 + y[t]^2);
sol={x[t],y[t]}/.
  NDSolve[{x'[t]==un2,y'[t]==vn2,x[0]==x0,y[0]==y0},{x[t],y[t]},{t,0,tmax}][[1]];
Plot[sol,{t,.9,1.1}]
sol/.t->1

which instantly gives you

{-0.294470, -0.287727}

and the plot seems to agree with the displayed numeric values.

Then you can try changing it back into your style, one tiny step at a time, and see at what point it breaks.

$\endgroup$
  • $\begingroup$ Thanks for point it out. My previous work was doing with such form, and it works really well, as mention in the comments in the question, but i couldn't do it this way. And this time, i had more than two hundred lines of code already written, it will be hash to change all of them. But you solution is absolutely correct. Thank you very much $\endgroup$ – Noah Ren Oct 20 at 8:37
  • $\begingroup$ @NoahRen Someone once wrote that some languages have one way of solving each problem. Mathematica is not like that, there are 6 or 12 ways of solving every problem. Are each really different? Likely no, they differ only in simple details that make no real difference. But many answers here will be followed by 2 or 4 other answers saying "no no, use my slight differences because that is what I use." For problems you will solve again and again pick one simple method that you can memorize and use without making too many mistakes and use that method every time. If you don't you will have problems. $\endgroup$ – Bill Oct 20 at 18:27
0
$\begingroup$

NDSolveValue does the trick:

solv2[x0_, y0_, tt_, tmax_] := 
  NDSolveValue[{r2'[t] == {un2, vn2}, 
    r2[0] == {x0, y0}} /. {\[Alpha] -> 1, \[Beta] -> 0.4}, 
    r2[tt],
    {t, 0, tmax},
    StartingStepSize -> 0.01, 
    Method -> {"FixedStep", Method -> {"ExplicitRungeKutta", DifferenceOrder -> 4}}, 
    MaxStepFraction -> Infinity]

solv2[0.5, 0.5, 0, 10]
(*    {0.5, 0.5}    *)

solv2[0.5, 0.5, 1, 10]
(*    {-0.29447, -0.287727}    *)
$\endgroup$
  • $\begingroup$ Thank you very much. I thought i could tick for accepting all of the all answers, because all of yours solution are indeed correct. But, if i modify with NDSolveValue, then all the plot i had before has to be modified. And sadly, i don't have enough time. If i could accept the solution for all of you. I will definitely tick the box. Thank you the help $\endgroup$ – Noah Ren Oct 20 at 8:33
  • $\begingroup$ @NoahRen Still, you can upvote if you think the answer is useful (to you or someone else having a similar question). $\endgroup$ – anderstood Oct 20 at 10:04
  • $\begingroup$ Yes, i have upvote the all the solution. $\endgroup$ – Noah Ren Oct 21 at 3:13

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.