2
$\begingroup$

I am searching for a while now, but I don't seem to be able to find an Answer for my Problem - if I am just not able to search properly, I am really sorry.

I Simplified my Problem to the following example: 

eqnA := {a''[x] + a[x] == b*Sin[a[x]]};
condA := {a[0] == a'[0] == 1};
stopCond := {WhenEvent[Evaluate[Re[a[x]] <= 0], xMax = x; "StopIntegration"]}
system := Join[eqnA, condA, stopCond];
sol = ParametricNDSolve[system, a, {x, 0, \[Infinity]}, {b}]
Plot[Evaluate@Table[a[b][x] /. sol, {b, .8, 1.2, .05}], {x, 0, xMax},PlotRange -> All]

Now, the Plot I want to create is a Parametric Plot of a over xMax. (Yes, I understand, that in this example it would be just a line - if it is any help: I am trying to solve the TOV equations and create the Mass over Radius diagram with central density as parameter)

I think, my understanding of WhenEvent and the uses of it might be flawed.

Any Tips are appreciated!

Edit: Since I think I did not formulate my problem clearly enough I try to explain a bit further:

With the post from george2079 I was able to get a bit further, but It still does not work. 

f[b_] := (a[b] /. sol)["Domain"][[1, 2]]
g[b_] := a[b][f[b]] /. sol
ParametricPlot[{f[b], g[b]}, {b, 0, 2}, PlotRange -> All]

I would expect that to produce something close to a line along the x-axis, but instead I get an error message. Plotting f[b] and g[b] individually works just fine.

Improve this question
$\endgroup$
2
  • $\begingroup$ By "Parametric Plot of a over xMax", do you mean a plot of the point {xMax, a[b][xMax]} as b varies? Note that a[b][xMax] would approximately be zero. $\endgroup$
    – Michael E2
    Commented May 17, 2016 at 14:09
  • $\begingroup$ If I understand you correctly, yes. xMax[b] on the x-axis, on the y-axis a[b][xMax], and I am interested in the resulting line when varying b (in my simple example a straight line along the x-axis) $\endgroup$
    – Zoldor
    Commented May 17, 2016 at 14:16

1 Answer 1

Reset to default
3
$\begingroup$

Your WhenEvent is just fine. If you look at an individual solution you can see where it stopped:

            (a[.8] /. sol)["Domain"]

{{0., 3.6351}}

The trouble is Plot extrapolates* past the solution domain, so you need to specify each plot domain individually:

Show[Table[
  Plot[ #[x] /. sol ,
    Evaluate[{x, Sequence @@ (First@#["Domain"])}],
    PlotStyle -> Hue[RandomReal[]]] &@(a[b] /. sol) ,
     {b, .8, 1.2, .05} ], PlotRange -> All]

enter image description here

Of course knowing whats going on you might just specify PlotRange->{0,Automatic} as well.

If you want to plot the stop point vs. b :

Plot[ (a[b] /. sol)["Domain"][[1, 2]] , {b, .8, 1.2}]

enter image description here

*and it really ought to throw a warning about it...

Edit: here is a parametric plot of the end value of each solution:

ParametricPlot[
  (Flatten@Last@Transpose[(a[b] /. sol)[{"Grid", "ValuesOnGrid"}]]),
       {b, .8, 1.2}, PlotRange -> {Automatic, {-1, 1}}, 
       AxesOrigin -> {Automatic, -1}, AspectRatio -> 1/GoldenRatio]

enter image description here

letting ParametricPlot autorange we can see the accuracy of the solver meeting your zero condition.

enter image description here

Improve this answer
$\endgroup$
4
  • $\begingroup$ Thank you, this already looks really useful! Maybe I can solve my Problem playing with "Domain". But I think my question formulated a bit confusing (Sorry for that, english is not my native tongue). I tried to explain, that the Plot I really want to create would be something like: ParametricPlot[{xMax[b],a[b]},{b, .8, 1.2}] $\endgroup$
    – Zoldor
    Commented May 18, 2016 at 5:41
  • $\begingroup$ responding to comment, a[b] at what x? a[b][xMax[b]] is zero due to your WhenEvent $\endgroup$
    – george2079
    Commented May 18, 2016 at 15:25
  • $\begingroup$ yes, a[b][xMax[b]] is zero in my example code, but in the Parametric Plot it should be a Line on the x-axis (since a is on the y-axis and because of my condition). Like I said, this is because i used a simplified ODE as example. Since this causes a lot of confusion: My full system are the coupled Tolman-Oppenheimer-Volkoff equations, with the WhenEvent triggered by the pressure reaching approx. 0 and the Plot should be the Mass at r=rMax over the max. of the Radius (which I get via the WhenEvent). But the problem is the same as in the example ODE that I used. $\endgroup$
    – Zoldor
    Commented May 19, 2016 at 4:40
  • $\begingroup$ Thank you, this does exactly what I need! $\endgroup$
    – Zoldor
    Commented May 20, 2016 at 16:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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