-7
$\begingroup$
ε = .1; A = 0; ω0 = 1; ωf = 0; data = NDSolve[ {x''[t] == ε (1 - x[t]^2) x'[t] - ω0^2*x[t] + A*Cos[ωf*t], x[0] == 1, x'[0] == 0}, {t, 0, 100}]; Periodogram[data, ScalingFunctions -> "Absolute", PlotRange -> All] 

How would I create a suitable Power spectrum of the frequencies that this system has. For these set of values (at the top), there should ideally be only one unique frequency (ω0) since the other one is not in play.

Mathematica's example of power spectra utilizes the periodogram funcition, although I don't really understand what it does mathematically, so there is room for using another type of graph.

$\endgroup$
0

4 Answers 4

7
$\begingroup$
ε = 1/10; A = 0; ω0 = 1; ωf = 0; sol = NDSolve[{x''[ t] == ε*(1 - x[t]^2) x'[t] - ω0^2*x[t] + A*Cos[ωf*t], x[0] == 1, x'[0] == 0}, x[t], {t, -100, 100}]; Plot[x[t] /. sol, {t, -100, 100}] 

function plot

Periodogram[Flatten[Table[x[t] /. sol, {t, -100, 100, 0.5}]],PlotRange -> All] 

periodogram plot

$\endgroup$
6
$\begingroup$

Here's my shot at it:

ε = .1; A = 0; ω0 = 1; ωf = 0; sol = First[NDSolve[ {x''[t] == ε (1 - x[t]^2) x'[t] - ω0^2*x[t] + A*Cos[ωf*t], x[0] == 1, x'[0] == 0}, x, {t, 0, 100}]]; f = x /. sol; data = Table[f[t], {t, 0, 100}]; Periodogram[data, ScalingFunctions -> "Absolute"] 

enter image description here

$\endgroup$
13
  • $\begingroup$ I've worked with this one, but it does not give me the correct results. $\endgroup$ Commented Jul 3, 2013 at 13:54
  • 1
    $\begingroup$ A-ha... And we should guess: Are there any error messages ? What do they say ? Any attempts to fix the "problems" ? $\endgroup$ Commented Jul 3, 2013 at 14:06
  • $\begingroup$ There are no error messages, the problem is that I am not getting a frequency on the x-axis, and the situations where I expect chaos, the power spectrum should be more messy. $\endgroup$ Commented Jul 3, 2013 at 19:56
  • $\begingroup$ Are you using the same code ? If not, please, do provide the code you are using. $\endgroup$ Commented Jul 4, 2013 at 17:07
  • $\begingroup$ Well, I've been trying it with yours, but there is no luck, and I tried with the other answer, but everything should be above the x axis. $\endgroup$ Commented Jul 4, 2013 at 21:33
5
$\begingroup$

Consider the equation: $$\ddot{x} = (x,\dot{x},t)$$ The numerical solution generates a sequence of values of x at discrete values of t. This sequence is known as time series. Plotting x pk t might seem random or noisy, but often it contains harmonics at certain frequencies. Dominant frequencies in time series data can be investigated using the Discrete Fourier transform. For example we take N equally spaced values from the time series - $$x = x_{0},x_{1}, ..., x_{N}$$ The Discrete Fourier transform is defined as follows: $$X_{k}=\frac{1}{\sqrt{N}}\sum_{m=0}^{N}x_{m}e^{\frac{-2\pi ikm}{N}}, (k=0,1,2,...,N-1)$$ Generally $X_{k}$ is a complex number. To analyse the frequency structure we analyse the power spectrum $P(\omega _{k})$ defined by: $$P(\omega _{k})=X_{k}\overline{X_{k}}=\left | X_{k} \right |^{2}$$

ε = .1; A = 0; ω0 = 1; ωf = 0; ff = x[t] /. First[NDSolve[{x''[ t] == ε (1 - x[t]^2) x'[t] - ω0^2*x[t] + A*Cos[ωf*t], x[0] == 1, x'[0] == 0}, x, {t, 0, 300}]][[1]] Table[ff, {t, 0, 100}] // Fourier // Abs // ListLogPlot[#, PlotRange -> All] & 

Plot

$\endgroup$
4
  • $\begingroup$ I feel as if we are getting closer to the answer as I see you are using Fourier, although this isn't getting me the correct results. Here, there should only be one frequency at around .1. Can you tell me what you are plotting on the x and y axes? $\endgroup$ Commented Jul 19, 2013 at 15:32
  • $\begingroup$ @Slightly Correct me if I am wrong, but this is the answer to your question. Nothing more, nothing less. $\endgroup$ Commented Jul 22, 2013 at 13:19
  • $\begingroup$ This is getting closer to what I am looking for, although the values are not correct. Since I have ω0 = 1, I should only have one peak at 1 (on the x axis) $\endgroup$ Commented Jul 23, 2013 at 2:10
  • $\begingroup$ I don't agree. Thanks for the help, though. $\endgroup$ Commented Jul 23, 2013 at 15:48
-4
$\begingroup$

The code that gives me the best results is this.

ε = 3; A = 5; ω0 = 1; ωf = 1.788; sol = NDSolve[{x''[t] == ε*(1 - x[t]^2) x'[t] - ω0^2*x[t] + A*Sin[ωf*t], x[0] == 1, x'[0] == 0}, x[t], {t, 0, 500}, MaxSteps -> Infinity] Periodogram[Flatten[Table[x[t] /. sol, {t, 0, 500, 1}]], PlotRange -> All, ScalingFunctions -> "Absolute"] 
$\endgroup$
4
  • 1
    $\begingroup$ Well, you could have mentioned the parameters are different ... $\endgroup$ Commented Jul 23, 2013 at 19:35
  • $\begingroup$ It's a completely different graph. $\endgroup$ Commented Jul 23, 2013 at 22:44
  • 1
    $\begingroup$ Yeah, because the conditions are different :D :D + There is a Sin[$\omega ft$] term in your answer whereas there's a Cos[$\omega ft$] in our code, based on the equation you provided in the first place. Obviously you made a lot of mistakes and you did not acknowledge making them. $\endgroup$ Commented Jul 24, 2013 at 7:40
  • $\begingroup$ Oh yes, I do agree that the forcing function is different, but it does not change anything in the output. $\endgroup$ Commented Jul 24, 2013 at 13:23

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.