# Julia and phase portraits of dynamic systems

3r33385. 3r33394.

3r33333. using Pkg

pkg "add Compose # master"

pkg "add Gadfly # master"

3r33385. 3r33394.

Not the most convenient option, but while waiting for the update 3r33333. PlotlyJS

can be for the sake of comparison and try. 3r33333. 3r33385. 3r33394.

You also need a powerful tool for solving differential equations 3r-338. DifferentialEquations

3r33333. 3r33385. 3r33394.

3r33333. ]add DifferentialEquations # An alternative way to add a package.

3r33385. 3r33394.

This is the most extensive and supported package, providing a lot of methods for solving difures. Now we proceed to the phase portraits. 3r33333. 3r33385. 3r33394.

The solutions of ordinary differential equations are often more conveniently portrayed in a different way than 3r-3319. 3r3355.

, and in

phase space

, along whose axes the values of each of the functions found are deposited. The argument t is included in the graphs only parametrically. In the case of two ODEs, such a graph, the phase portrait of the system, is a curve on the phase plane and is therefore particularly visual. 3r33333. 3r33385. 3r33394. 3r362. Oscillator

3r33385. 3r33394.

Dynamics of a harmonic oscillator with attenuation

described by the following system of equations:

3r33385. 3r33394.

3r33333.

3r33333.

3r33333. 3r33385. 3r33394.

and regardless of the initial conditions (x? y0), it comes to equilibrium, a point (?0) with a zero angle of deviation and zero speed. 3r33333. 3r33385. 3r33394.

At

3r388.

The decision takes the form peculiar to the classical oscillator. 3r33333. 3r33385. 3r33394.

3r33333. using DifferentialEquations, Gadfly # we include packages in the project - a terribly long process

3r33394. function garmosc (ω = pi, γ = ? x0 = 0.? y0 = 0.1)

3r33394. A =[0 1

-ω^2 γ]3r33394. 3r33394. syst (u, p, t) = A * u # ODE system

3r33394. u0 =[x0, y0]# start cond-ns

tspan = (0.? 4pi) # time period

3r33394. prob = ODEProblem (syst, u? tspan) # problem to solve

sol = solve (prob, RK4 (), reltol = 1e-? timeseries_steps = 4)

3r33394. N = length (sol.u)

J = length (sol. U[1])

3r33394. U = zeros (N, J)

3r33394. for i in 1: N, j in 1: J

U[i,j]= sol.u[i] [j]# rewrite the answers in a matrix

more comfortably. end

U

end

3r33385. 3r33394.

Let's sort the code. The function takes the frequency and attenuation parameter values as well as the initial conditions. Nested function

syst ()

sets the system. In order to come out single-line, they used matrix multiplication. Function

solve ()

taking a huge number of parameters is very flexible to adjust the problem, but we only indicated the method of solution - Runge-Kutta 4 (3r3128. there are many other 3r-3129.), relative error, and the fact that not all points need to be saved, but only every fourth . In the variable

sol

the answer matrix is preserved, and

sol.t

contains a vector of time values, and

sol.u

system solutions at these times. All this is quietly drawn in

Plots

, and for

Gadfly

will have to complete the answer in a matrix more comfortably. We do not need times, so we only return solutions. 3r33333. 3r33385. 3r33394.

Construct the phase portrait:

3r33385. 3r33394.

3r33333. Ans0 = garmosc ()

plot (x = Ans0[:,1], y = Ans0[:,2])

3r33385. 3r33394.

Among the functions of a high order for us is particularly noticeable 3r3333330. broadcast (f,[1, 2, 3]) which is in function

f

alternately substitutes the values of the array[1, 2, 3]. Short record

` f. ([1, 2, 3]) 3r33376. . This is useful to vary the frequency, attenuation parameter, initial coordinate, and initial velocity. 3r33333. 3r33385. 3r33394. `

` `** Hid under the spoiler [/b] **

** **

** 3r33333. L = Array {Any} (undef, 3) # created an array of anything (for graph layers) 3r3r94. clr =["red" "green" "yellow"]3r33394. 3r33394. function ploter (Answers)**

for i = 1: 3

L[i]= layer (x = Answ[i] [:,1], y = Answers[i] [:,2], Geom.path, Theme (default_color = clr[i]))

end

p = plot (L[1], L[2], L[3])

end

3r33394. Ans1 = garmosc. (W2w2w235.) # Frequency

p1 = ploter (Ans1)

3r33394. Ans2 = garmosc. (Pi,[0.1 0.3 0.5]) # Attenuation

p2 = ploter (Ans2)

3r33394. Ans3 = garmosc. (Pi, 0.?[0.2 0.5 0.8])

p3 = ploter (Ans3)

3r33394. Ans4 = garmosc. (Pi, 0.? 0.?[0.2 0.5 0.8])

p4 = ploter (Ans4)

3r33394. set_default_plot_size (16cm, 16cm)

3r33394. vstack (hstack (p? p2), hstack (p? p4)) # blind horizontally and vertically

** 3r33385. 3r33394. **

** 3r3204. 3r33333. 3r33395. 3r33395. 3r33385. 3r33394. **

** Now we will consider not small oscillations, but oscillations of an arbitrary amplitude: 3r-3388. 3r33385. 3r33394. **

** 3r33333. **

** 3r33333. 3r33333. 3r33385. 3r33394. **

** The graph of the solution on the phase plane is not an ellipse (which indicates the nonlinearity of the oscillator). The less we choose the amplitude of oscillations (i.e., the initial conditions), the less non-linearity will appear (therefore, the small oscillations of the pendulum can be considered harmonic). 3r33333. 3r33385. 3r33394. **

** **** Hid under the spoiler [/b] **

** **

** 3r33333. function ungarmosc (ω = pi, γ = ? x0 = 0.? y0 = 0.1)**

3r33394. function syst (du, u, p, t)

du[1]= u[2]3r33394. du[2]= -sin (ω * u[1]) + γ * u[2]3r33394. end

3r33394. u0 =[x0, y0]# start cond-ns

tspan = (0.? 2pi) # time period

3r33394. prob = ODEProblem (syst, u? tspan) # problem to solve

sol = solve (prob, RK4 (), reltol = 1e-? timeseries_steps = 4)

3r33394. N = length (sol.u)

J = length (sol. U[1])

3r33394. U = zeros (N, J)

3r33394. for i in 1: N, j in 1: J

U[i,j]= sol.u[i] [j]# rewrite the answers in a matrix

more comfortably. end

U

end

3r33394. Ans1 = ungarmosc. (W2w2w235.)

p1 = ploter (Ans1)

3r33394. Ans2 = ungarmosc. (Pi,[0.1 0.4 0.7])

p2 = ploter (Ans2)

3r33394. Ans3 = ungarmosc. (Pi, 0.?[0.1 0.5 0.9])

p3 = ploter (Ans3)

3r33394. Ans4 = ungarmosc. (Pi, 0.? 0.?[0.1 0.5 0.9])

p4 = ploter (Ans4)

3r33394. vstack (hstack (p? p2), hstack (p? p4))

** 3r33385. 3r33394. **

** Brusselator 3r33385. 3r33394. **

** This model describes some autocatalytic chemical reaction in which diffusion plays a certain role. The model was proposed in 1968 by Lefebvre and Prigogine 3r-3388. 3r33385. 3r33394. **

** 3r33333. **

** 3r33333. 3r33333. 3r33385. 3r33394. **

** Unknown features 3r33232. reflect the dynamics of the concentration of intermediate products of a chemical reaction. Parameter model 3r33320. makes sense of the initial concentration of the catalyst (the third substance). 3r33333. 3r33385. 3r33394. **

** In more detail, the evolution of the phase portrait of the Brusselator can be observed by performing calculations with a different parameter 3r-3319. 3r33320. . When it is increased, the node will first gradually shift to a point with coordinates (? 3r33333. ) Until it reaches 3r3333314. bifurcation r3r3315. values. 3r33320. = 2. At this point, a qualitative reorganization of the portrait takes place, which is expressed in the birth of the limiting cycle. With a further increase of 3r33333. 3r33320. there is only a quantitative change in the parameters of this cycle. 3r33333. 3r33385. 3r33394. **

** **** Hid under the spoiler [/b] **

** **

** 3r33333. function brusselator (μ = 0.? u0 =[0.1; 0.1])**

3r33394. function syst (du, u, p, t)

du[1]= - (μ + 1) * u[1]+ u[1]* u[1]* u[2]+ 1

du[2]= μ * u[1]- u[1]* u[1]* u[2]3r33394. end

3r33394. tspan = (0.? 10) # time period

3r33394. prob = ODEProblem (syst, u? tspan) # problem to solve

sol = solve (prob, RK4 (), reltol = 1e-? timeseries_steps = 1)

3r33394. N = length (sol.u)

J = length (sol. U[1])

3r33394. U = zeros (N, J)

3r33394. for i in 1: N, j in 1: J

U[i,j]= sol.u[i] [j]# rewrite the answers in a matrix

comfortably. end

U

end

3r33394. L = Array {Any} (undef, 10)

3r33394. function ploter (Answers)

for i = 1: length (Answers)

L[i]= layer (x = Answ[i] [:,1], y = Answers[i] [:,2], Geom.path)

end

plot (L[1], L[2], L[3], L[4], L[5], L[6], L[7], L[8], L[9], L w222. L end

3r33394. SC =[[0 0.5],[0 1.5],[2.5 0],[1.5 0],[0.5 1],[1 0],[1 1],[1.5 2],[0.1 0.1],[0.5 0.2]]

3r33394. Ans1 = brusselator. (0.? SC)

Ans2 = brusselator. (0.? SC)

Ans3 = brusselator. (1.? SC)

Ans4 = brusselator. (2.? SC)

p1 = ploter (Ans1)

p2 = ploter (Ans2)

p3 = ploter (Ans3)

p4 = ploter (Ans4)

3r33394. set_default_plot_size (16cm, 16cm)

3r33394. vstack (hstack (p? p2), hstack (p? p4))

3r33385. 3r33394.

3r33381. 3r33333. 3r33395. 3r33395. 3r33385. 3r33394.

It's enough for today. Next time, we’ll try to learn how to use another graphic package on new tasks, at the same time getting a hand on Julia’s syntax. In the process of solving problems, they slowly begin to be traced not so directly to the pros and cons rather, to pleasantness and inconvenience — a separate conversation should be devoted to this. And of course, I would like to see something more serious in my games with functions - so I urge julists to tell more about serious projects, this would help many people and me in particular in learning this wonderful language. 3r33333. 3r33395. 3r33394. 3r33394.

! function (e) {function t (t, n) {if (! (n in e)) {for (var r, a = e.document, i = a.scripts, o = i.length; o-- ;) if (-1! == i[o].src.indexOf (t)) {r = i[o]; break} if (! r) {r = a.createElement ("script"), r.type = "text /jаvascript", r.async =! ? r.defer =! ? r.src = t, r.charset = "UTF-8"; var d = function () {var e = a.getElementsByTagName ("script")[0]; e.parentNode.insertBefore (r, e)}; "[object Opera]" == e.opera? a.addEventListener? a.addEventListener ("DOMContentLoaded", d,! 1): e.attachEvent ("onload", d ): d ()}}} t ("//mediator.mail.ru/script/2820404/"""_mediator") () (); 3r33393. 3r33394. 3r33395.