3
$\begingroup$

I'm an absolute beginner using Mathematica just to solve a pair of non-linear ODEs and I'm struggling a bit with the lingo. The solutions come out as x -> InterpolatingFunction[-], y -> InterpolatingFunction[-]. Both are functions purely of time.

All I want to do now is to extract a few values from these (or the corresponding plots if easier), namely gradients at certain points, and points with certain gradients. Things like Grad and FindMaximum don't seem to work.

Can anyone help? Even if its just a bit of code I can copy to get them?

Here is what I have so far:

s = NDSolve[{x'[t] == -(y[t]^-1) x[t]^(7/2) (0.99*10^-9) - y[t] x[t]^(-1/2) (0.24*10^-3), 
             y'[t] == x[t]^(5/2) (0.99*10^-9) - (y[t]^2) (x[t]^(-3/2)) (1.21*10^-3), 
             x[0] == 1.5*10^7, y[0] == 10^10},
            {x, y}, {t, 0, 4000}]

Two ODEs with initial conditions and time boundaries.

And I plot it

TP = Plot[x[t] /. s, {t, 0, 400}]
NP = Plot[y[t] /. s, {t, 0, 400}]
$\endgroup$
  • 1
    $\begingroup$ I'm missing the actual question, perhaps, but NDSolve[] returns a function you can differentiate. For instance: yf = y /. First[NDSolve[{y'[x] == -y[x], y[0] == 1}, y, {x, 0, 1}]]; Table[yf'[x], {x, 0, 1, 1/10}] $\endgroup$ – J. M.'s discontentment Jan 25 '17 at 14:25
  • 1
    $\begingroup$ Can you provide some code of your specific problem? $\endgroup$ – mgamer Jan 25 '17 at 14:46
  • $\begingroup$ code now included $\endgroup$ – Duncan Mark Horne Jan 25 '17 at 17:45
2
$\begingroup$

I may be misunderstanding your question here, but it appears to me that FindMaximum etc work quite well with your interpolating functions.

For instance, using the solutions s as defined by your code, we can plot e.g. y[t] that has a clear maximum:

Plot[y[t] /. s, {t, 0, 400}]

Mathematica graphics

Clearly the function has a maximum for $t$ around 100. We can use FindMaximum, FindArgMax, or NArgMax to find the value of t for which the maximum arises:

FindMaximum[y[t] /. s, t]
(* Out: {2.35387*10^10, {t -> 119.186}} *)

FindArgMax[y[t] /. s, t]
(* Out: {119.186} *)

NArgMax needs the most hand-holding: first, it doesn't like its argument to be a list, so we use y[t] /. First@s to extract only the interpolation function, instead of a list of one item. Second, its initial guess leads it to evaluate the interpolation function outside of its domain; we prevent that by imposing constraints on the value of $t$ that are well within the interpolation domain and close to the max.

NArgMax[{y[t] /. First@s, 100 < t < 400}, t]
(* Out: 119.186 *)

As you can see all results agree.


You can also obtain derivatives directly from the interpolation functions as y'[t] /. s. For instance:

Plot[y'[t] /. s, {t, 0, 400}, PlotRange -> All]

Mathematica graphics

You could also look for that maximum directly by finding the zeroes of the first derivative:

FindRoot[(y'[t] /. s) == 0, {t, 100}]
(* Out: {t -> 119.186} *)
| improve this answer | |
$\endgroup$

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.