Mathematica Mathematica - Find max value NDSolve Plot

AI Thread Summary
The discussion focuses on finding the maximum value of a function derived from a numerical solution to a differential equation using Mathematica's NDSolve. A user shares their code for solving and plotting the equation but struggles to extract the maximum value from the plotted results. Another user suggests deconstructing the interpolating function to access the y-values directly and provides a method for looping through time values to find the maximum. Additionally, there is clarification on the use of the ReplaceAll operator (/.), which is necessary for plotting the results accurately. The conversation highlights the importance of understanding Mathematica's syntax and functions for effective data analysis.
carey1000
Messages
2
Reaction score
0
Dear forum members:

After numerically solving a differential equation and plotting the results I would like to determine the single maximum value in the plotted range but do not know how.

The code below works for numerically solving the differential equation and plotting the results.

NDSolve[{x''[t] + x[t] - 0.167 x[t]^3 == 0.005 Cos[t + -0.0000977162*t^2/2], x[0] == 0, x'[0] == 0}, x, {t, 0, 1000}]

Plot[Evaluate[x[t] /. s], {t, 0, 1000},
Frame -> {True, True, False, False}, FrameLabel -> {"t", "x"}, FrameStyle -> Directive[FontSize -> 15], Axes -> False]

Thank you!
Carey
 
Physics news on Phys.org
I can't find an easy way to do it, with mapping or anything, but deconstructing the interpolating function can get it:

Max[s[[1]][[1]][[2]][[4]][[3]]]

This is basically the y-values of the array. Just found it by looking at what creates an interpolating function. I'm not sure if its always at this location, but I don't see why it wouldn't be.
 
Tabulate the data and use Max@
 
Maxdata = Evaluate[x[1]/.s];

Then, find a way to loop t = 2 thru 1000:
(Newmax = If[Evaluate[x[t]/.s]>Maxdata, Evaluate[x[t]/.s, Maxdata];
Maxdata = Newmax;)
 
carey1000 said:
Dear forum members:

..

Plot[Evaluate[x[t] /. s], {t, 0, 1000},
Frame -> {True, True, False, False}, FrameLabel -> {"t", "x"}, FrameStyle -> Directive[FontSize -> 15], Axes -> False]

Carey

To plot the results, why you need to divide by " .s "??
 
panzer7910 said:

To plot the results, why you need to divide by " .s "??

In the help manual of Mathematica, it gives the example:
Code:
s = NDSolve[{y'[x] == y[x] Cos[x + y[x]], y[0] == 1}, y, {x, 0, 30}]
Plot[Evaluate[y[x] /. s], {x, 0, 30}, PlotRange -> All]

Also you could find the symbol /. means ReplaceAll (/.) in its help manual.
Code:
In[1]:= {x, x^2, y, z} /. x -> a

Out[1]= {a, a^2, y, z}

In[2]:= {x, x^2, y, z} /. x -> {a, b}

Out[2]= {{a, b}, {a^2, b^2}, y, z}

In[3]:= Sin[x] /. Sin -> Cos

Out[3]= Cos[x]
Combine these two, you could find that s represents the result of this numerical solution, and
Code:
Evaluate[y[x] /. s]
means using the numerical results s to represent the value of y[x] as y[x] itself does not have a definition in the ordinary form so the plot function cannot plot it unless you use something(here is the numerical results s) to represent it
 
Last edited:
Will this do?

In[1]:= f = x/.NDSolve[{x''[t]+x[t]-0.167x[t]^3==0.005 Cos[t+ -0.0000977162*t^2/2], x[0]==0, x'[0]==0}, x, {t,0,1000}][[1,1]];
Plot[f[x], {x, 0, 1000}]

Out[2]= ...PlotSnipped...

In[3]:= {*Looks like max is between 775 and 825, turn off warnings and hunt for it*)
Off[FindMaximum::eit];
Sort[Table[FindMaximum[{f[x], 0<x<1000}, {x,y,y-1,y+1}], {y,775,825}], OrderedQ[{First[#1], First[#2]}] &]

Out[4]= {
{-0.699441, {x -> 796.}}, {-0.695645, {x -> 803.}}, {-0.688966, {x -> 789.}},
{-0.681066, {x -> 810.}}, {-0.661475, {x -> 782.}}, {-0.66001, {x -> 817.}},
{-0.637412, {x -> 824.}}, {-0.614744, {x -> 775.}}, {0.094972, {x -> 820.}},
{0.108052, {x -> 778.}}, {0.121419, {x -> 813.}}, {0.144966, {x -> 785.}},
{0.145051, {x -> 806.}}, {0.16006, {x -> 799.}}, {0.161367, {x -> 792.}},
{0.412197, {x -> 795.}}, {0.414114, {x -> 788.}}, {0.423985, {x -> 802.}},
{0.432996, {x -> 781.}}, {0.444894, {x -> 809.}}, {0.469508, {x -> 816.}},
{0.470473, {x -> 774.}}, {0.492212, {x -> 823.}}, {1.02978, {x -> 779.}},
{1.05629, {x -> 821.}}, {1.05931, {x -> 786.}}, {1.0706, {x -> 814.}},
{1.07659, {x -> 793.}}, {1.08044, {x -> 807.}}, {1.08311, {x -> 800.}},
{1.16416, {x -> 780.}}, {1.17012, {x -> 787.}}, {1.17933, {x -> 794.}},
{1.19063, {x -> 801.}}, {1.20206, {x -> 808.}}, {1.21145, {x -> 815.}},
{1.21692, {x -> 822.}}, {1.22094, {x -> 779.646}}, {1.22094, {x -> 779.647}},
{1.23615, {x -> 786.619}}, {1.23615, {x -> 786.621}}, {1.24892, {x -> 793.61}},
{1.24892, {x -> 793.612}}, {1.25869, {x -> 800.615}}, {1.25869, {x -> 800.617}},
{1.26499, {x -> 807.63}}, {1.26499, {x -> 807.632}}, {1.26624, {x -> 821.672}},
{1.26624, {x -> 821.673}}, {1.26753, {x -> 814.65}}, {1.26753, {x -> 814.652}}}
 

Similar threads

Replies
1
Views
2K
Replies
1
Views
2K
Replies
1
Views
2K
Replies
1
Views
5K
Replies
1
Views
2K
Replies
4
Views
3K
Replies
6
Views
2K
Back
Top