(こちらもしばらくは工事中。
DifferentialEquations.jl: Scientific Machine Learning (SciML) Enabled Simulation and Estimation
Julia だけでなく、Python, R でも使えるらしい。
Julia で使えるようにするには
using Pkg Pkg.add("DifferentialEquations") using DifferentialEquations |
“Solving ODEs in Julia” by Chris Rackauckas
using DifferentialEquations f(u,p,t) = 1.01*u u0 = 1/2 tspan = (0.0,1.0) prob = ODEProblem(f,u0,tspan) sol = solve(prob, Tsit5(), reltol=1e-8, abstol=1e-8) using Plots plot(sol,linewidth=5,title="Solution to the linear ODE with a thick line", xaxis="Time (t)",yaxis="u(t) (in μm)",label="My Thick Line!") # legend=false plot!(sol.t, t->u0*exp(1.01t),lw=3,ls=:dash,label="True Solution!") savefig("diff.png") |
こちらは , reltol=1e-8, abstol=1e-8 としてある。 それも込めて、プログラムの書き方を見てほっとする。 最初の例からこういうのを見せるべきだ。
ソルバーの指定が出来るそうで、, Tsit5() というのがそれに相当する (Tsitouras 5/4 Runge-Kutta method)。 どういうソルバーが良いかは、 「ODE Solvers」 に詳しく書いてある。 MATLAB だと ode45 というのがあるけれど、 翻訳は DP5() だとか (どういう意味だろう)。 それよりは大抵の場合 Tsit5() の方が効率的だとか。 Tsitouras というのはギリシャの人か。 バリバリ論文を出しているんだ (名前が珍しいので、検索して論文探すの簡単)。 知らなかった。
(2022/1/9追記)
少しツッコミを入れる。
1.01 のようなマジック・ナンバーを入れる (2箇所別々) のは拙い。
パラメーターを渡す仕掛けがあり、それを使うべきだと思う。
using DifferentialEquations f(u,p,t) = p*u u0 = 1/2 tspan = (0.0,1.0) a=1.01 prob = ODEProblem(f,u0,tspan,a) sol = solve(prob, Tsit5(), reltol=1e-8, abstol=1e-8) using Plots plot(sol,linewidth=5,title="Solution to the linear ODE with a thick line", xaxis="Time (t)",yaxis="u(t) (in μm)",label="My Thick Line!") # legend=false plot!(sol.t, t->u0*exp(a*t),lw=3,ls=:dash,label="Exact Solution!") savefig("diff.png") |
多次元の場合の例として、Loenz 方程式があげられている。
using DifferentialEquations function lorenz!(du,u,p,t) σ,ρ,β = p du[1] = σ*(u[2]-u[1]) du[2] = u[1]*(ρ-u[3]) - u[2] du[3] = u[1]*u[2] - β*u[3] end u0 = [1.0,0.0,0.0] p = (10,28,8/3) tspan = (0.0,100.0) prob = ODEProblem(lorenz!,u0,tspan,p) sol = solve(prob) using Plots plot(sol,vars=(1,2,3)) |
何と MATLABDiffEq.jl というのもあるらしい (MATLAB をインストールしてあることが必要だとか)。
(某月某日) あ、面白いの見つけた。