0

I have a system of three ODEs. This system oscillates. I need to determine the frequency of each of the three solutions.

My file. You can download and modify this.

I'm new to Mathematica and I know this question has been asked multiple times, but I can't get it working.

{xSol, ySol, zSol} = NDSolveValue[{ode1, ode2, ode3, ic2}, {x, y, z}, {t, 0, 200}]

Plotting xSol, ySol and zSol (three interpolating functions) gives:

Solutions

Can anyone show me how to determine the frequency of this thing using Reap and Sow or FFT. I've been instructed already how to do this using EcoEvo package, but I'm also wondering how to do it without such a package.

Both the Reap/Sow solution and the FFT come pretty close to what I have and need. But I fail in modifying the mentioned solutions so that it works with my code.

How do I individually determine the frequency of xSol, ySol and zSol?

My implementation of WhenEvent

For example, when I try to use WhenEvent, I used it like this:

{xSol, ySol, zSol} = NDSolveValue[{ode1, ode2, ode3, ic2, WhenEvent[x'[t] == 0, Sow[t]]}, {x, y, z}, {t, 0, 200}]

This gives me an error (which I presume to be about the WhenEvent) and my regular solutions. Because the solutions are correct, I think that the problem here is in the WhenEvent.

Errors

Error: Event location failed to converge to the requested accuracy or precision within 100 iterations between t = ... and t = ....

LATER UPDATE

{{xSol, ySol, zSol}, reapresult} = Reap[NDSolveValue[{ode1, ode2, ode3, ic2, WhenEvent[x[t] > .4 && x'[t] == 0,Sow[t]]}, {x, y, z}, {t, 0, 200} ]]

This outputs and empty list {}

I first forgot to add the Reap[]

All suggestions and solutions are highly appreciated! Many thanks

ralphjsmit
  • 187
  • 7
  • Perhaps it would be useful if you explained exactly how you are failing to adapt the solution you linked to. – Natas Oct 20 '20 at 12:04
  • 1
    Yes of course, thanks Natas! – ralphjsmit Oct 20 '20 at 12:09
  • OK, I think this has to do with the fact that x is (almost) zero for a long time. You can use (for instance) the condition x[t] == .8 && x'[t] > 0 instead, although this of course has to be modified for similar problems. I suggest you have a look at the documentation of WhenEvent. – Natas Oct 20 '20 at 12:16
  • You can look at my source code at (https://github.com/cklausme/EcoEvo/blob/master/EcoEvo/EcoEvo.nb) for inspiration. – Chris K Oct 20 '20 at 13:08
  • Hi Chris, thanks ;-) but as a Mathematica novice it's still one level too high – ralphjsmit Oct 20 '20 at 13:16

1 Answers1

1

I found the solution. My initial code

{{xSol, ySol, zSol}, reapresult} = Reap[NDSolveValue[{ode1, ode2, ode3, ic2, WhenEvent[x[t] > .4 && x'[t] == 0,Sow[t]]}, {x, y, z}, {t, 0, 200} ]]

gives an empty list for reapresult. But when I switch the order in the WhenEvent, it does give me the correct list.

{{xSol, ySol, zSol}, reapresult} = Reap[NDSolveValue[{ode1, ode2, ode3, ic2, WhenEvent[x'[t] == 0 && x[t] > 0.1,Sow[t]]}, {x, y, z}, {t, 0, 200} ]]

Thanks everyone for thinking along!

ralphjsmit
  • 187
  • 7