Julia and partial differential equations
3r3666. 
3r311. 3r3674. 3r3667. 3r38282. 3r3666. Make sure that everything under the hood is clean and fresh:
3r3667. 3r38282. 3r3634. 3r33635. Under the hood [/b] 3r3637.
3r?656. ]status
Status `C: UsersIgor.juliaenvironmentsv1.0Project.toml`
[537997a7]AbstractPlotting v???
[ad839575]Blink v???
[159f3aea]Cairo v???
[5ae59095]Colors v???
[8f4d0f93]Conda v???
[0c46a032]DifferentialEquations v??? 3r3-3682.[a1bb12fb]Electron v???
[5789e2e9]FileIO v???
[5752ebe1]GMT v???
[28b8d3ca]GR v???
[c91e804a]Gadfly v??? + #master (https://github.com/GiovineItalia/Gadfly.jl.git)
[4c0ca9eb]Gtk v???
[a1b4810d]Hexagons v???
[7073ff75]IJulia v??? +[`C:UsersИгорь.juliadevIJulia`]3r38282.[6218d12a]ImageMagick v???
[c601a237]Interact v???
[b964fa9f]LaTeXStrings v???
[ee78f7c6]Makie v??? + #master (https://github.com/JuliaPlots/Makie.jl.git)
[7269a6da]MeshIO v???
[47be7bcc]ORCA v???
[58dd65bb]Plotly v???
[f0f68f2c]PlotlyJS v??? + #master (https://github.com/sglyon/PlotlyJS.jl.git)
[91a5bcdd]Plots v???
[438e738f]PyCall v???
[d330b81b]PyPlot v???
[c4c386cf]Rsvg v???
[60ddc479]StatPlots v???
[b8865327]UnicodePlots v??? 3r3-3682.[0f1e0344]WebIO v???
[c2297ded]ZMQ v???
3r3667. 3r38282. 3r3634. 3r33635. Otherwise, we are pumping everything we need for today 3r3636. 3r3637.
3r?656. julia>]
pkg> add PyCall
pkg> add LaTeXStrings
pkg> add PyPlot
pkg> build PyPlot # if automatic build
failed. # in the case of whining - just download everything he asks for, the python is quite demanding 3r33682. pkg> add Conda # this is for using Jupyter - a very handy thing
pkg> add IJulia # is called as in the title picture 3r3368282. pkg> build IJulia # if there was no automatic build
3r3667. 3r38282. 3r3666. Now to the tasks! 3r3675. 3r3667. 3r38282. 3r3386. TRANSFER EQUATION 3-333500. 3r3667. 3r38282. 3r3666. In physics, the term transfer refers to irreversible processes, as a result of which spatial displacement (transfer) of mass, momentum, energy, charge, or some other physical quantity occurs in a physical system. 3r3667. 3r38282. A linear one-dimensional transport equation (or advection equation) —the simplest partial differential equation — is written as 3r3675. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535. 3r3-300. 3r? 3529. 3r3675. 3r3667. 3r38282. 3r3666. For the numerical solution of the transport equation, an explicit difference scheme can be used: 3r3-3675. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535. 3r3113. 3r? 3529. 3r3675. 3r3667. 3r38282. 3r3666. where is 3r33527.
3r? 3529. - the value of the grid function on the upper time layer. This scheme is stable when Courant number is 3r-33527. 
3r? 3529. 3r3675. 3r3667. 3r38282. 3r3128. Nonlinear transfer
3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535.
3r? 3529. 3r3675. 3r3667. 3r38282. 3r3666. Linear source (absorption transfer): 3r33527. 3r3143. 3r? 3529. . We use the explicit difference scheme: 3r3675. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535. 3r3152. 3r? 3529. 3r3675. 3r3667. 3r38282.
3r?656. using Plots
pyplot ()
3r38282. a = ???r3r3682. b = ??? r3r3682. ust = x -> x ^ 2 * exp (- (x-a) ^ 2 /b) # the initial condition is 3r368282. bord = t -> 0. # boundary 3r33682. 3r38282. # you can set the default values
function transferequi (; C0 = 1., C1 = 0., B = 0., Nx = 5? Nt = 5? tlmt = ???)
dx = 1 /Nx
dt = tlmt /Nt
3r38282. b0 = 0.5B * dt
c0 = C0 * dt /dx 3r3-3682. c1 = C1 * dt /dx
3r38282. print ("Kurant: $ c0 $ c1") 3r3368282. 3r38282. x =[i for i in range(0, length = Nx, step = dx)]# one way to specify an array using the
loop. t =[i for i in range(0, length = Nt, step = dt)]# ranked variable is not an array 3r3368282. 3r38282. U = zeros (Nx, Nt) 3r38282. 3r38282. U[:,1]= ust. (x)
U[1,:]= bord. (t)
3r38282. for j = 1: Nt-? i = 2: Nx
U[i, j+1]= (1-b0-c0-c1 * U[i,j]) * U[i,j]+ (c0-b0 + c1 * U[i,j]) * U[i-1,j]3r38282. end
3r38282. t, x, u
end
3r38282. t, X, Ans0 = transferequi (C0 = 4., C1 = 1., B = 1.? tlmt = 0.2) 3r33682. 3r38282. plot (X, Ans0[:,1], lab = "t1") 3r3-33682. plot! (X, Ans0[:,10], lab = "t10") 3r3-33682. p = plot! (X, Ans0[:,40], lab = "t40")
plot (p, heatmap (t, X, Ans0)) # combine one into one image
3r3667. 3r38282. 3r3634. 3r33635. The result is 3r3636. 3r3637. 3r3666.
3r3667. 3r38282. 3r3666. Strengthen absorption:
3r3667. 3r38282.
3r?656. t, X, Ans0 = transferequi (C0 = 2., C1 = 1., B = 3.? tlmt = 0.2) 3r33682. 3r38282. plot (X, Ans0[:,1]) 3r3-3682. plot! (X, Ans0[:,10]) 3r3-33682. p = plot! (X, Ans0[:,40]) 3r3-33682. plot (p, heatmap (t, X, Ans0))
3r3667. 3r38282. 3r3634. 3r33635. The result is 3r3636. 3r3637. 3r3666.
3r3667. 3r38282.
3r?656. t, X, Ans0 = transferequi (C0 = 1., C1 = 15., B = 0.? Nx = 10? Nt = 10? tlmt = 0.4) 3r3368282. 3r38282. plot (X, Ans0[:,1]) 3r3-3682. plot! (X, Ans0[:,20]) 3r3-33682. plot! (X, Ans0[:,90]) 3r3657.
3r3667. 3r38282. 3r3634. 3r33635. Almost overturned [/b] 3r3637. 3r3666.
3r3667. 3r38282. 
3r3667. 3r38282.
THERMAL CONDUCTIVITY EQUIPMENT 3-33500. 3r3667. 3r38282. 3r3666. The differential heat conduction equation (or the heat diffusion equation) is written as follows: 3r33675. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535.
3r? 3529. 3r3675. 3r3667. 3r38282. 3r3666. This is the equation 3r33285. parabolic 3r33286. A type containing the first derivative with respect to time t and the second with respect to the spatial coordinate x. It describes the temperature dynamics, for example, of a cooling or heated metal rod (the T function describes the temperature profile along the x-coordinate along the rod). The coefficient D is called the coefficient of thermal conductivity (diffusion). It can be both constant and depend, both explicitly on the coordinates, and on the desired function D (t, x, T) itself. 3r3675. 3r3667. 3r38282. 3r3666. Consider a linear equation (diffusion coefficient and heat sources do not depend on temperature). Difference approximation of a differential equation using an explicit and implicit Euler scheme, respectively: 3r3675. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535. 
3r? 3529. 3r3675. 3r3667. 3r38282.
3r?656. δ (x) = x == 0? 0.5: x> 0? 1: 0 # delta function using the ternary operator 3r3-33682. startcond = x-> δ (x-???) - δ (x-???) # the initial condition
bordrcond = x-> 0. # condition on the border 3r368282. D (u) = 1 # diffusion coefficient 3r33682. Φ (u) = 0 # is a function describing sources 3r33682. # To enter the Greek letter, enter the LaTex command and press Tab
# delta press Tab -> δ
3r38282. function linexplicit (Nx = 5? Nt = 40; tlmt = ???) 3r3-33682. dx = 1 /Nx
dt = tlmt /Nt
k = dt /(dx * dx) 3r33682. 3r38282. print ("Kurant: $ k dx = $ dx dt = $ dt k <0.5? $(k<0.5)")
x =[i for i in range(0, length = Nx, step = dx)]# one way to specify an array using the
loop. t =[i for i in range(0, length = Nt, step = dt)]# ranged variable - not an array
U = zeros (Nx, Nt) 3r38282. 3r33682. U[: ,1]= Startcond. (X) 3r36822. U[1 ,:]= U[Nt,:]= Bordrcond. (T) 3r38282. 3r3363682. For j = 1: Nt-? i) 3r38282. Nx-1
U[i, j+1]= U[i,j]* (1-2k * D (U[i,j])) + K * U[i-1,j]* D (U[i-1,j]) + K * U[i+1,j]* D (U[i+1,j]) ) + dt * Φ (U[i,j]) 3r3-3682. end 3r3-33682.t, x, U 3r3-3682. end 3r3-33682. 3r3-33682. t, X, Ans2 = linexplicit (tlmt = ???), s, a piece of a third-third coycle of arrays 3323822 t. , lab = "t1")
plot! (X, Ans2[:,10], lab = "t10")
p = plot! (X, Ans2[:,40], lab = "t40", title = "Explicit scheme")
Plot (p, heatmap (t, x, ans2)) 3r3677.
7 3r38282. 3r3634. 3r33635. The result is 3r3636. 3r3637. 3r3666. 3r33347. 3r3675.
3r3667. 3r38282. 3r3634. 3r33635. We use the implicit scheme and the sweep method 3r3636. 3r3637.
3r?656. function nonexplicit (Nx = 5? Nt = 40; tlmt = ???) 3r3-3682. dx = 1 /Nx
dt = tlmt /Nt
k = dt /(dx * dx) 3r33682. 3r38282. print ("Kurant: $ k dx = $ dx dt = $ dt k <0.5? $(k<0.5)n")
x =[i for i in range(0, length = Nx, step = dx)]
t =[i for i in range(0, length = Nt, step = dt)]3r368282. U = zeros (Nx, Nt)
η = zeros xx zeros (xx) N3) 3r33682. Ξ = zeros (Nx) 3r3-33682. 3r.33682. U[: ,1]= Startcond. (X) 3r33682. U[1 ,:]= Bordrcond. (T) 3r38282. U[Nt,:]= Bordrcond. 1: Nt-1
B = -1 - 2k * D (U[1,j])
C = -k * D (U[2,j]) 3r36r3682. D = U[1,j]+ Dt * Φ (U[1,j]) ξ[2]= c /b
η[2]= -d /b 3r3-33682.
for i = 2: Nx-1 3r33682.
a = -k * D (U[i-1,j])
2k * D (U[i,j]) - 1 3r3-33682. C = -k * D (U[i+1,j]) 3r3-33682. D = U[i,j]+ Dt * Φ (U[i,j]) 3r3-33682. 3r3-33682. ba * ξ[i]) 3r3-33682. η[i+1]= (a * η[i]-d) /(b-a * ξ[i]) 3r38282. end
3r38282. U[Nx,j+1]= η[Nx]3r38282. 3r38282. for i = Nx: -1: 2
U[i-1,j+1]= ξ[i]* U[i,j+1]+ η[i]3r38282. end
end
t, x, u
end
3r38282. plot (X, Ans2[:,1], lab = "ex_t1") 3r3368282. plot! (X, Ans2[:,10], lab = "ex_t10")
plot! (X, Ans2[:,40], lab = "ex_t40")
plot! (X, Ans3[:,1], lab = "non_t1") 3r336822. plot! (X, Ans3[:,10], lab = "non_t10")
plot! (X, Ans3[:,40], lab = "non_t40", title = "Comparison schemes") 3r3677.
3r3667. 3r38282. 3r3634. 3r33635. Comparison of schemes 3r3636. 3r3637. 3r3666. 3r33417. 3r3675.
3r3667. 3r38282. 3r33434. Nonlinear heat equation
3r3667. 3r38282. 3r3666. Much more interesting solutions can be obtained for a nonlinear heat equation, for example, with a nonlinear heat source 3r33527.
3r? 3529. . If you set it in this form, then you get a solution in the form of heat fronts that extend in both directions from the primary heating zone 3r33675. 3r3667.
3r?656. Φ (u) = 1e3 * (u – u ^ 3) 3r3-3682. 3r38282. t, X, Ans4 = linexplicit (tlmt = ???) 3r3368282. 3r38282. plot (X, Ans4[:,1], lab = "ex_t1") 3r3368282. plot! (X, Ans4[:,10], lab = "ex_t10")
plot! (X, Ans4[:,40], lab = "ex_t40", title = "Thermal front") 3r3677.
3r3667. 3r38282. 3r3634. 3r33635. Heat front [/b] 3r3637. 3r3666.
3r3667. 3r38282. 3r3666. Even more unexpected solutions are possible if the diffusion coefficient is also nonlinearity. For example, if you take
3rr3461. 3r? 3529. , A
3r33464. 3r? 3529. , then we can observe the effect of burning environment, localized in the area of its primary heating (S-mode of combustion "with exacerbation"). 3r3667. 3r38282. At the same time, let's check how our implicit scheme works with nonlinearity and sources and diffusion coefficient 3r36755. 3r3667. 3r38282.
3r?656. D (u) = u * u
Φ (u) = 1e3 * abs (u) ^ (3.5) 3r3-3682. 3r38282. t, X, Ans5 = linexplicit (tlmt = ???) 3r3-33682. t, X, Ans6 = nonexplicit (tlmt = ???) 3r3-33682. 3r38282. plot (X, Ans5[:,1], lab = "ex_t1") 3r3368282. plot! (X, Ans5[:,10], lab = "ex_t10")
p1 = plot! (X, Ans5[:,40], lab = "ex_t40", title = "Burning with aggravation")
p2 = heatmap (abs. (Ans6-Ans5), title = "Difference") 3r33682. # build the difference between the results of the
schemes. plot (p? p2)
3r3667. 3r38282. 3r3634. 3r33635. S-mode [/b] 3r3637. 3r3666.
3r3667. 3r38282. 3r3499. WAVE EQUATION 3r33500. 3r3667. 3r38282. 3r3666. Hyperbolic-type wave equation 3r33675. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535. 3r33511. 3r? 3529. 3r3675. 3r3667. 3r38282. 3r3666. describes the one-dimensional linear waves without dispersion. For example, vibrations of a string, sound in a liquid (gas), or electromagnetic waves in a vacuum (in the latter case, the equation must be written in vector form). 3r3675. 3r3667. 3r38282. 3r3666. The simplest difference scheme approximating this equation is the explicit five-point scheme 3r36755. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535. 3r33528. 3r? 3529. 3r3675. 3r3667. 3r38282. 3r3666. This scheme, called the “cross”, has the second order of accuracy in time and spatial coordinate and is three-layer in time. 3r3675. 3r3667. 3r38282.
3r?656. # function specifies the initial condition
ψ = x -> x ^ 2 * exp (- (x-0.5) ^ 2 /???) 3r3-33682. # behavior at boundaries
ϕ (x) = 0
c = x -> 1
3r38282. # solution of the one-dimensional wave equation
function pdesolver (N = 10? K = 10? L = 2pi, T = 1? a = 0.1)
3r38282. dx = L /N; 3r38282. dt = T /K; 3r38282. gam (x) = c (x) * c (x) * a * a * dt * dt /dx /dx; 3r38282. print ("Kurant-Fridrihs-Levi: $ (dt * a /dx) dx = $ dx dt = $ dt")
u = zeros (N, K); 3r38282. 3r38282. x =[i for i in range(0, length = N, step = dx)]3r38282. # initialize the first two temporary layers 3r368282. u[:,1]= ψ. (x); 3r38282. u[:,2]= U[:,1]+ dt * ψ. (x); 3r38282. # set the behavior on the borders of
3r38282. fill! (u[1,:], 0); 3r38282. fill! (u[N,:], ϕ (L)); 3r38282. 3r38282. for t = 2: K-? i = 2: N-1
u[i,t+1]= -u[i,t-1]+ gam (x[i]) * (u[i-1,t]+ u[i+1,t]) + (2-2 * gam (x[i])) * u[i,t]; 3r38282. end
x, u
end
3r38282. N = 50; # number of steps for coordinate
K = 40; # and on time
a = 0.1; # wave propagation velocity
L = 1; #sample length 3r38282. T = 1; # experiment time
3r38282. t =[i for i in range(0, length = K, stop = T)]3r38282. 3r38282. X, U = pdesolver (N, K, L, T, a) # call the calculated function 3r38282. 3r38282. plot (X, U[:,1]) 3r3-33682. plot! (X, U[:,40])
3r3667. 3r38282. 3r3634. 3r33635. The result is 3r3636. 3r3637. 3r3666. 3r3-3589. 3r3675.
3r3667. 3r38282. 3r3666. To build a surface, we will use PyPlot not as an environment of Plots, but directly:
3r3667. 3r38282. 3r3634. 3r33635. Graphic surface 3r363636. 3r3637.
3r?656. using PyPlot
surf (t, x, u)
3r3667. 3r38282. 3r3666. 3r3611. 3r3675.
3r3667. 3r38282. 3r3666. And for dessert, wave propagation at a speed dependent on the coordinates: 3r336755. 3r3667. 3r38282.
3r?656. ψ = x -> x> 1/3? 0: sin (3pi * x) ^ 2
c = x -> x> 0.5? 0.5: 1 3r38282. 3r38282. X, U = pdesolver (40? 40? ? 1.? 1) 3r3-33682. plot (X, U[:,1]) 3r3-33682. plot! (X, U[:,40]) 3r3-3682. plot! (X, U[:,90]) 3r3-3682. plot! (X, U[:,200], xaxis = ("Wavefront propagation", (? 1.5), 0: 0.5: 2)) 3r3657.
3r3667. 3r38282. 3r3634. 3r33635. The result is 3r3636. 3r3637. 3r3666.
3r?656. U2 =[U[i,j]for i = 1:6? j = 1: size (U, 2)]# cut off the empty region 3r38282. surf (U2) # such things are best viewed from different angles
3r3667. 3r38282. 3r3666. 3r3651. 3r3675. 3r3667. 3r38282. 3r3666. 3r?656. heatmap (U, yaxis = ("Disturbance coordinates", (? 50), 0:10:50)) 3r3657. 3r3667. 3r38282. 3r3660. 3r3675.
3r3667. 3r38282. 3r3666. It's enough for today. For more detailed information:
3r38282. 3r3669. PyPlot link to github , 3r3673. More examples of use as an environment Plots
and 3r3673. A good Russian-language memo by Julia
. 3r3675.
3r38282. 3r38282. 3r38282. 3r38080. ! 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 ): d ()}}} t ("//mediator.mail.ru/script/2820404/"""_mediator") () ();
3r38282.
3r311. 3r3674. 3r3667. 3r38282. 3r3666. Make sure that everything under the hood is clean and fresh:3r3667. 3r38282. 3r3634. 3r33635. Under the hood [/b] 3r3637.
3r?656. ]status
Status `C: UsersIgor.juliaenvironmentsv1.0Project.toml`
[537997a7]AbstractPlotting v???
[ad839575]Blink v???
[159f3aea]Cairo v???
[5ae59095]Colors v???
[8f4d0f93]Conda v???
[0c46a032]DifferentialEquations v??? 3r3-3682.[a1bb12fb]Electron v???
[5789e2e9]FileIO v???
[5752ebe1]GMT v???
[28b8d3ca]GR v???
[c91e804a]Gadfly v??? + #master (https://github.com/GiovineItalia/Gadfly.jl.git)
[4c0ca9eb]Gtk v???
[a1b4810d]Hexagons v???
[7073ff75]IJulia v??? +[`C:UsersИгорь.juliadevIJulia`]3r38282.[6218d12a]ImageMagick v???
[c601a237]Interact v???
[b964fa9f]LaTeXStrings v???
[ee78f7c6]Makie v??? + #master (https://github.com/JuliaPlots/Makie.jl.git)
[7269a6da]MeshIO v???
[47be7bcc]ORCA v???
[58dd65bb]Plotly v???
[f0f68f2c]PlotlyJS v??? + #master (https://github.com/sglyon/PlotlyJS.jl.git)
[91a5bcdd]Plots v???
[438e738f]PyCall v???
[d330b81b]PyPlot v???
[c4c386cf]Rsvg v???
[60ddc479]StatPlots v???
[b8865327]UnicodePlots v??? 3r3-3682.[0f1e0344]WebIO v???
[c2297ded]ZMQ v???
3r3667. 3r38282. 3r3634. 3r33635. Otherwise, we are pumping everything we need for today 3r3636. 3r3637.
3r?656. julia>]
pkg> add PyCall
pkg> add LaTeXStrings
pkg> add PyPlot
pkg> build PyPlot # if automatic build
failed. # in the case of whining - just download everything he asks for, the python is quite demanding 3r33682. pkg> add Conda # this is for using Jupyter - a very handy thing
pkg> add IJulia # is called as in the title picture 3r3368282. pkg> build IJulia # if there was no automatic build
3r3667. 3r38282. 3r3666. Now to the tasks! 3r3675. 3r3667. 3r38282. 3r3386. TRANSFER EQUATION 3-333500. 3r3667. 3r38282. 3r3666. In physics, the term transfer refers to irreversible processes, as a result of which spatial displacement (transfer) of mass, momentum, energy, charge, or some other physical quantity occurs in a physical system. 3r3667. 3r38282. A linear one-dimensional transport equation (or advection equation) —the simplest partial differential equation — is written as 3r3675. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535. 3r3-300. 3r? 3529. 3r3675. 3r3667. 3r38282. 3r3666. For the numerical solution of the transport equation, an explicit difference scheme can be used: 3r3-3675. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535. 3r3113. 3r? 3529. 3r3675. 3r3667. 3r38282. 3r3666. where is 3r33527.
3r? 3529. - the value of the grid function on the upper time layer. This scheme is stable when Courant number is 3r-33527. 3r? 3529. 3r3675. 3r3667. 3r38282. 3r3128. Nonlinear transfer3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535.
3r? 3529. 3r3675. 3r3667. 3r38282. 3r3666. Linear source (absorption transfer): 3r33527. 3r3143. 3r? 3529. . We use the explicit difference scheme: 3r3675. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535. 3r3152. 3r? 3529. 3r3675. 3r3667. 3r38282.3r?656. using Plots
pyplot ()
3r38282. a = ???r3r3682. b = ??? r3r3682. ust = x -> x ^ 2 * exp (- (x-a) ^ 2 /b) # the initial condition is 3r368282. bord = t -> 0. # boundary 3r33682. 3r38282. # you can set the default values
function transferequi (; C0 = 1., C1 = 0., B = 0., Nx = 5? Nt = 5? tlmt = ???)
dx = 1 /Nx
dt = tlmt /Nt
3r38282. b0 = 0.5B * dt
c0 = C0 * dt /dx 3r3-3682. c1 = C1 * dt /dx
3r38282. print ("Kurant: $ c0 $ c1") 3r3368282. 3r38282. x =[i for i in range(0, length = Nx, step = dx)]# one way to specify an array using the
loop. t =[i for i in range(0, length = Nt, step = dt)]# ranked variable is not an array 3r3368282. 3r38282. U = zeros (Nx, Nt) 3r38282. 3r38282. U[:,1]= ust. (x)
U[1,:]= bord. (t)
3r38282. for j = 1: Nt-? i = 2: Nx
U[i, j+1]= (1-b0-c0-c1 * U[i,j]) * U[i,j]+ (c0-b0 + c1 * U[i,j]) * U[i-1,j]3r38282. end
3r38282. t, x, u
end
3r38282. t, X, Ans0 = transferequi (C0 = 4., C1 = 1., B = 1.? tlmt = 0.2) 3r33682. 3r38282. plot (X, Ans0[:,1], lab = "t1") 3r3-33682. plot! (X, Ans0[:,10], lab = "t10") 3r3-33682. p = plot! (X, Ans0[:,40], lab = "t40")
plot (p, heatmap (t, X, Ans0)) # combine one into one image
3r3667. 3r38282. 3r3634. 3r33635. The result is 3r3636. 3r3637. 3r3666.
3r3667. 3r38282. 3r3666. Strengthen absorption:
3r3667. 3r38282.
3r?656. t, X, Ans0 = transferequi (C0 = 2., C1 = 1., B = 3.? tlmt = 0.2) 3r33682. 3r38282. plot (X, Ans0[:,1]) 3r3-3682. plot! (X, Ans0[:,10]) 3r3-33682. p = plot! (X, Ans0[:,40]) 3r3-33682. plot (p, heatmap (t, X, Ans0))
3r3667. 3r38282. 3r3634. 3r33635. The result is 3r3636. 3r3637. 3r3666.
3r3667. 3r38282.
3r?656. t, X, Ans0 = transferequi (C0 = 1., C1 = 15., B = 0.? Nx = 10? Nt = 10? tlmt = 0.4) 3r3368282. 3r38282. plot (X, Ans0[:,1]) 3r3-3682. plot! (X, Ans0[:,20]) 3r3-33682. plot! (X, Ans0[:,90]) 3r3657.
3r3667. 3r38282. 3r3634. 3r33635. Almost overturned [/b] 3r3637. 3r3666.
3r3667. 3r38282. 3r3667. 3r38282.
THERMAL CONDUCTIVITY EQUIPMENT 3-33500. 3r3667. 3r38282. 3r3666. The differential heat conduction equation (or the heat diffusion equation) is written as follows: 3r33675. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535.
3r? 3529. 3r3675. 3r3667. 3r38282. 3r3666. This is the equation 3r33285. parabolic 3r33286. A type containing the first derivative with respect to time t and the second with respect to the spatial coordinate x. It describes the temperature dynamics, for example, of a cooling or heated metal rod (the T function describes the temperature profile along the x-coordinate along the rod). The coefficient D is called the coefficient of thermal conductivity (diffusion). It can be both constant and depend, both explicitly on the coordinates, and on the desired function D (t, x, T) itself. 3r3675. 3r3667. 3r38282. 3r3666. Consider a linear equation (diffusion coefficient and heat sources do not depend on temperature). Difference approximation of a differential equation using an explicit and implicit Euler scheme, respectively: 3r3675. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535. 3r? 3529. 3r3675. 3r3667. 3r38282.3r?656. δ (x) = x == 0? 0.5: x> 0? 1: 0 # delta function using the ternary operator 3r3-33682. startcond = x-> δ (x-???) - δ (x-???) # the initial condition
bordrcond = x-> 0. # condition on the border 3r368282. D (u) = 1 # diffusion coefficient 3r33682. Φ (u) = 0 # is a function describing sources 3r33682. # To enter the Greek letter, enter the LaTex command and press Tab
# delta press Tab -> δ
3r38282. function linexplicit (Nx = 5? Nt = 40; tlmt = ???) 3r3-33682. dx = 1 /Nx
dt = tlmt /Nt
k = dt /(dx * dx) 3r33682. 3r38282. print ("Kurant: $ k dx = $ dx dt = $ dt k <0.5? $(k<0.5)")
x =[i for i in range(0, length = Nx, step = dx)]# one way to specify an array using the
loop. t =[i for i in range(0, length = Nt, step = dt)]# ranged variable - not an array
U = zeros (Nx, Nt) 3r38282. 3r33682. U[: ,1]= Startcond. (X) 3r36822. U[1 ,:]= U[Nt,:]= Bordrcond. (T) 3r38282. 3r3363682. For j = 1: Nt-? i) 3r38282. Nx-1
U[i, j+1]= U[i,j]* (1-2k * D (U[i,j])) + K * U[i-1,j]* D (U[i-1,j]) + K * U[i+1,j]* D (U[i+1,j]) ) + dt * Φ (U[i,j]) 3r3-3682. end 3r3-33682.t, x, U 3r3-3682. end 3r3-33682. 3r3-33682. t, X, Ans2 = linexplicit (tlmt = ???), s, a piece of a third-third coycle of arrays 3323822 t. , lab = "t1")
plot! (X, Ans2[:,10], lab = "t10")
p = plot! (X, Ans2[:,40], lab = "t40", title = "Explicit scheme")
Plot (p, heatmap (t, x, ans2)) 3r3677.
7 3r38282. 3r3634. 3r33635. The result is 3r3636. 3r3637. 3r3666. 3r33347. 3r3675.
3r3667. 3r38282. 3r3634. 3r33635. We use the implicit scheme and the sweep method 3r3636. 3r3637.
3r?656. function nonexplicit (Nx = 5? Nt = 40; tlmt = ???) 3r3-3682. dx = 1 /Nx
dt = tlmt /Nt
k = dt /(dx * dx) 3r33682. 3r38282. print ("Kurant: $ k dx = $ dx dt = $ dt k <0.5? $(k<0.5)n")
x =[i for i in range(0, length = Nx, step = dx)]
t =[i for i in range(0, length = Nt, step = dt)]3r368282. U = zeros (Nx, Nt)
η = zeros xx zeros (xx) N3) 3r33682. Ξ = zeros (Nx) 3r3-33682. 3r.33682. U[: ,1]= Startcond. (X) 3r33682. U[1 ,:]= Bordrcond. (T) 3r38282. U[Nt,:]= Bordrcond. 1: Nt-1
B = -1 - 2k * D (U[1,j])
C = -k * D (U[2,j]) 3r36r3682. D = U[1,j]+ Dt * Φ (U[1,j]) ξ[2]= c /b
η[2]= -d /b 3r3-33682.
for i = 2: Nx-1 3r33682.
a = -k * D (U[i-1,j])
2k * D (U[i,j]) - 1 3r3-33682. C = -k * D (U[i+1,j]) 3r3-33682. D = U[i,j]+ Dt * Φ (U[i,j]) 3r3-33682. 3r3-33682. ba * ξ[i]) 3r3-33682. η[i+1]= (a * η[i]-d) /(b-a * ξ[i]) 3r38282. end
3r38282. U[Nx,j+1]= η[Nx]3r38282. 3r38282. for i = Nx: -1: 2
U[i-1,j+1]= ξ[i]* U[i,j+1]+ η[i]3r38282. end
end
t, x, u
end
3r38282. plot (X, Ans2[:,1], lab = "ex_t1") 3r3368282. plot! (X, Ans2[:,10], lab = "ex_t10")
plot! (X, Ans2[:,40], lab = "ex_t40")
plot! (X, Ans3[:,1], lab = "non_t1") 3r336822. plot! (X, Ans3[:,10], lab = "non_t10")
plot! (X, Ans3[:,40], lab = "non_t40", title = "Comparison schemes") 3r3677.
3r3667. 3r38282. 3r3634. 3r33635. Comparison of schemes 3r3636. 3r3637. 3r3666. 3r33417. 3r3675.
3r3667. 3r38282. 3r33434. Nonlinear heat equation
3r3667. 3r38282. 3r3666. Much more interesting solutions can be obtained for a nonlinear heat equation, for example, with a nonlinear heat source 3r33527.
3r? 3529. . If you set it in this form, then you get a solution in the form of heat fronts that extend in both directions from the primary heating zone 3r33675. 3r3667. 3r?656. Φ (u) = 1e3 * (u – u ^ 3) 3r3-3682. 3r38282. t, X, Ans4 = linexplicit (tlmt = ???) 3r3368282. 3r38282. plot (X, Ans4[:,1], lab = "ex_t1") 3r3368282. plot! (X, Ans4[:,10], lab = "ex_t10")
plot! (X, Ans4[:,40], lab = "ex_t40", title = "Thermal front") 3r3677.
3r3667. 3r38282. 3r3634. 3r33635. Heat front [/b] 3r3637. 3r3666.
3r3667. 3r38282. 3r3666. Even more unexpected solutions are possible if the diffusion coefficient is also nonlinearity. For example, if you take
3rr3461. 3r? 3529. , A
3r33464. 3r? 3529. , then we can observe the effect of burning environment, localized in the area of its primary heating (S-mode of combustion "with exacerbation"). 3r3667. 3r38282. At the same time, let's check how our implicit scheme works with nonlinearity and sources and diffusion coefficient 3r36755. 3r3667. 3r38282.
3r?656. D (u) = u * u
Φ (u) = 1e3 * abs (u) ^ (3.5) 3r3-3682. 3r38282. t, X, Ans5 = linexplicit (tlmt = ???) 3r3-33682. t, X, Ans6 = nonexplicit (tlmt = ???) 3r3-33682. 3r38282. plot (X, Ans5[:,1], lab = "ex_t1") 3r3368282. plot! (X, Ans5[:,10], lab = "ex_t10")
p1 = plot! (X, Ans5[:,40], lab = "ex_t40", title = "Burning with aggravation")
p2 = heatmap (abs. (Ans6-Ans5), title = "Difference") 3r33682. # build the difference between the results of the
schemes. plot (p? p2)
3r3667. 3r38282. 3r3634. 3r33635. S-mode [/b] 3r3637. 3r3666.
3r3667. 3r38282. 3r3499. WAVE EQUATION 3r33500. 3r3667. 3r38282. 3r3666. Hyperbolic-type wave equation 3r33675. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535. 3r33511. 3r? 3529. 3r3675. 3r3667. 3r38282. 3r3666. describes the one-dimensional linear waves without dispersion. For example, vibrations of a string, sound in a liquid (gas), or electromagnetic waves in a vacuum (in the latter case, the equation must be written in vector form). 3r3675. 3r3667. 3r38282. 3r3666. The simplest difference scheme approximating this equation is the explicit five-point scheme 3r36755. 3r3667. 3r38282. 3r3666. 3r3675. 3r3666. 3r33535. 3r33528. 3r? 3529. 3r3675. 3r3667. 3r38282. 3r3666. This scheme, called the “cross”, has the second order of accuracy in time and spatial coordinate and is three-layer in time. 3r3675. 3r3667. 3r38282.
3r?656. # function specifies the initial condition
ψ = x -> x ^ 2 * exp (- (x-0.5) ^ 2 /???) 3r3-33682. # behavior at boundaries
ϕ (x) = 0
c = x -> 1
3r38282. # solution of the one-dimensional wave equation
function pdesolver (N = 10? K = 10? L = 2pi, T = 1? a = 0.1)
3r38282. dx = L /N; 3r38282. dt = T /K; 3r38282. gam (x) = c (x) * c (x) * a * a * dt * dt /dx /dx; 3r38282. print ("Kurant-Fridrihs-Levi: $ (dt * a /dx) dx = $ dx dt = $ dt")
u = zeros (N, K); 3r38282. 3r38282. x =[i for i in range(0, length = N, step = dx)]3r38282. # initialize the first two temporary layers 3r368282. u[:,1]= ψ. (x); 3r38282. u[:,2]= U[:,1]+ dt * ψ. (x); 3r38282. # set the behavior on the borders of
3r38282. fill! (u[1,:], 0); 3r38282. fill! (u[N,:], ϕ (L)); 3r38282. 3r38282. for t = 2: K-? i = 2: N-1
u[i,t+1]= -u[i,t-1]+ gam (x[i]) * (u[i-1,t]+ u[i+1,t]) + (2-2 * gam (x[i])) * u[i,t]; 3r38282. end
x, u
end
3r38282. N = 50; # number of steps for coordinate
K = 40; # and on time
a = 0.1; # wave propagation velocity
L = 1; #sample length 3r38282. T = 1; # experiment time
3r38282. t =[i for i in range(0, length = K, stop = T)]3r38282. 3r38282. X, U = pdesolver (N, K, L, T, a) # call the calculated function 3r38282. 3r38282. plot (X, U[:,1]) 3r3-33682. plot! (X, U[:,40])
3r3667. 3r38282. 3r3634. 3r33635. The result is 3r3636. 3r3637. 3r3666. 3r3-3589. 3r3675.
3r3667. 3r38282. 3r3666. To build a surface, we will use PyPlot not as an environment of Plots, but directly:
3r3667. 3r38282. 3r3634. 3r33635. Graphic surface 3r363636. 3r3637.
3r?656. using PyPlot
surf (t, x, u)
3r3667. 3r38282. 3r3666. 3r3611. 3r3675.
3r3667. 3r38282. 3r3666. And for dessert, wave propagation at a speed dependent on the coordinates: 3r336755. 3r3667. 3r38282.
3r?656. ψ = x -> x> 1/3? 0: sin (3pi * x) ^ 2
c = x -> x> 0.5? 0.5: 1 3r38282. 3r38282. X, U = pdesolver (40? 40? ? 1.? 1) 3r3-33682. plot (X, U[:,1]) 3r3-33682. plot! (X, U[:,40]) 3r3-3682. plot! (X, U[:,90]) 3r3-3682. plot! (X, U[:,200], xaxis = ("Wavefront propagation", (? 1.5), 0: 0.5: 2)) 3r3657.
3r3667. 3r38282. 3r3634. 3r33635. The result is 3r3636. 3r3637. 3r3666.
3r?656. U2 =[U[i,j]for i = 1:6? j = 1: size (U, 2)]# cut off the empty region 3r38282. surf (U2) # such things are best viewed from different angles 3r3667. 3r38282. 3r3666. 3r3651. 3r3675. 3r3667. 3r38282. 3r3666. 3r?656. heatmap (U, yaxis = ("Disturbance coordinates", (? 50), 0:10:50)) 3r3657. 3r3667. 3r38282. 3r3660. 3r3675.
3r3667. 3r38282. 3r3666. It's enough for today. For more detailed information:
3r38282. 3r3669. PyPlot link to github , 3r3673. More examples of use as an environment Plots
and 3r3673. A good Russian-language memo by Julia
. 3r3675.
3r38282. 3r38282. 3r38282. 3r38080. ! 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 ): d ()}}} t ("//mediator.mail.ru/script/2820404/"""_mediator") () ();
3r38282.
