9.1. Метод Эйлера#
Простейшим методом решения (9.1)
является метод Эйлера. В этом методе производная \(u'\) заменяется конечной разностью
В геометрической интерпретации это означает, что метод приближает \(y_{i+1}\) значением касательной через \((t_i, y_i)\) в точке \(t_i + \tau\).
Пусть дана задача Коши для функции \(u\). Значения \(u(t)\) при \(t>0\) находим из рекуррентной формулы
На примере метода Эйлера введём несколько понятий, применяющихся к разностным схемам [89].
Фиксируем точку \(t\) и построим последовательность сеток \(\omega_\tau\) (9.2) таких, что \(\tau \to 0\) и \(t_i = i \tau = t\). Говорят, что метод сходится в точке \(t\), если
Метод называют сходящимся на отрезке \((0, T]\), если он сходится в каждой точке отрезка.
Если при этом
то говорят, что метод имеет \(p\)-ый порядок точности.
Альтернативное определение сходимости [16] заключается в рассмотрении \(y_\tau\) и проекции \(u_\tau\) как сеточных функций (т.е. функций, определённых на \(\omega_\tau\)), тогда требуют \(\|y_\tau - u_\tau\|_{F_\tau} \to 0\), где \(F_\tau\) – пространство сеточных функций.
При подстановке решения \(u\) в численный метод получается выражение для погрешности аппроксимации разностного уравнения на решении исходного уравнения, или невязки. Например, для метода Эйлера (9.3) получим вид невязки \(\psi_i\)
Если невязка
то говорят, что разностный метод аппроксимирует исходное дифференциальное уравнение.
При этом, если невязка имеет вид
то говорят о \(p\)-ом порядке аппроксимации.
Несложно найти порядок аппроксимации метода Эйлера. Сама невязка имеет вид (9.4), для её получения разложим \(u_{i+1}-u_i\) по формуле Тейлора
Подставим полученное соотношение в невязку (9.4)
Таким образом, метод Эйлера имеет первый порядок аппроксимации.
Помимо сходимости и аппроксимации, вводят различные понятия устойчивости. Одной из обобщённых является устойчивость по правой части. Данная устойчивость [16] означает, что при малом возмущении правой части уравнения \(f + \epsilon\) разностный метод будет иметь (а) единственное решение \(z_\tau\), (б) решение \(z_\tau\) будет не сильно отличаться от невозмущённого, т.е. \(\|z_\tau - y_\tau\| \le c \|\epsilon_\tau\|\).
Показано [16], что из аппроксимации и устойчивости по правой части следует сходимость метода с порядком точности, совпадающем с порядком аппроксимации.
9.1.1. Реализация#
В этой главе для постановки задачи Коши (скалярной) мы будем пользоваться структурой данных CauchyODEProblem
.
struct CauchyODEProblem{T<:Real,F<:Function}
bound::Tuple{T,T} # отрезок интегрирования
u₀::T # начальное значение интегрируемой функции
f::F # правая часть ОДУ
function CauchyODEProblem(; f::Function, tstart::Real, tend::Real, u₀::Real)
new{Float64, typeof(f)}(
float.(tuple(tstart, tend)),
float(u₀),
f,
)
end
end
Явный метод Эйлера
"""
euler(problem; nsteps)
Решает задачу Коши `problem` явным методом Эйлера за `nsteps` шагов.
"""
function euler(problem::CauchyODEProblem; nsteps::Integer)
u = Vector{Float64}(undef, nsteps + 1)
u[1] = problem.u₀
tstart, tend = problem.bound
trange = range(tstart, tend; length=nsteps+1)
τ = step(trange)
@inbounds for i in 1:nsteps
tᵢ, uᵢ = trange[i], u[i]
u[i+1] = uᵢ + τ * problem.f(tᵢ, uᵢ)
end
return trange, u
end
Решим задачу Коши
явным методом Эйлера и сравним с точным решением \(u = 2 \cosh t\).
problem = CauchyODEProblem(
f=(t, u) -> -u + 2 * exp(t),
tstart=0,
tend=1,
u₀=2,
)
plt = plot(; xlabel="time", leg=:topleft)
plot!((t) -> 2 * cosh(t); label=L"2 \times \cosh{t}", lw=2, xlabel=L"t", ylabel=L"u(t)")
for n in (1, 2, 3, 10)
t, u = euler(problem; nsteps=n)
plot!(t, u; label="euler, nsteps = $n", marker=:o)
end
plt