## Numerical Analysis

Problem. Frequency power spectrum

Consider the following nonlinear 3rd-order system:

dx1

dt

= −(x2 + x3),

dx2

dt

= x1 +

1

5

x2,

dx3

dt

=

1

5

+ (x1 − ϵ)x3.

This equation is called the R¨ossler’s equation, and it has been applied to study chemical reactions

and electrical circuits. Integrate this equation numerically using the initial conditions x1(0) = 2.0,

x2(0) = 2.0, x3(0) = 0 for 4 different choices of parameter ϵ: ϵ = 2.5, 3.0, 4.1, 5.0. Provide the following

outputs for each values of ϵ.

1. Present a phase portrait of x1, x2 and x3 in a three dimensional phase space (i.e., take x1, x2

x3 as x−, y−, z−axis, respectively) by using plot3 function (see Matlab help for detail). Note

that the time series of a solution may contain some initial transient. Therefore when plotting the

trajectories in the phase space, omit (do not show) the initial transient trajectory. Assume the

initial transient time as roughly 15 basic periods of the solution x1(t) (basic period is the period

of the most dominant oscillation, i.e., the Fourier mode that has the largest energy).

Comments: To find the period of the dominant Fourier mode in the time series of x1, you compute

a frequency spectrum of the x1(t) first, and then from the frequency spectrum find the frequency at

which largest spectral energy is attained. From that frequency you can compute the basic period,

and by using the period you can remove the transient part in the phase portrait. The following

script is one of the examples for identifying the time array index it at which the transient part

ends.

% Estimate transient limit

[maxpower, ip]=max(power);

T1=2*pi/omega(ip);

[tmp, it]=min(abs(t-15*T1));

Here the power is a spectral energy density array and omega is a discrete frequency which must

be provided from the time series of x1 beforehand. Then, first, we find the maximum spectral

energy and find the index ip at which the energy is maximum. Then find the basic period T1

from the discrete frequency at the maximum energy, i.e., omega(ip). Then we find the time index

that is nearest to 15 times the dominant period. It can be found by searching a minimum value

of abs(t-15*T1), and in this search the corresponding index (it) of array t is also obtained. By

using the index it, you can plot a phase portrait after the transient part.

2. Present a frequency spectrum of the signal x1(t), and find the most dominant frequency and its

corresponding period.

Comments: To obtain a reasonable frequency spectrum, you must obtain a sufficiently long

time signal. You must explore and choose sampling interval until the resulting spectrum is not

significantly affected by the aliasing error.