我正在尝试使用julia提供的DifferentialEquations.jl,它工作正常,直到我尝试在二阶常微分方程上使用它。例如,考虑二阶常微分方程
x''(t) = x'(t) + 2* x(t),带初始条件
x'(0) = 0,x(0) =1
它有一个解析解: x(t) =2/3exp(-t)+1/3exp(2t)。
为了从数值上解决这个问题,我运行以下代码:
using DifferentialEquations;
function f_simple(ddu, du, u, p, t)
ddu[1] = du[1] + 2*u[1]
end;
du0 = [0.]
u0 = [1.]
tspan = (0.0,5.0)
prob2 = SecondOrderODEProblem(f_simple, du0, u0, tspan)
sol = solve(prob2,reltol=1e-8, abstol=1e-8);有了这些,
sol(3)[2] = 122.57014434362732而解析解的结果是134.50945587649028,所以我在这里有点迷惑。
发布于 2021-05-07 07:54:29
根据the documentation for DifferentialEquations.jl的说法,Vern7()适合于非刚性方程的高精度解:
sol = solve(prob2, Vern7(), reltol=1e-8, abstol=1e-8)
julia> println(sol(3)[2])
134.5094558872943在我的机器上,这与解析解非常接近。我不太确定使用的默认方法是什么:文档表明solve有一些方法可以在没有指定求解器的情况下选择适当的求解器。
有关Vern7()的更多信息,请访问Jim Verner's page on Runge-Kutta algorithms。
https://stackoverflow.com/questions/67373015
复制相似问题