Численное решение ОДУ второго порядка с odeint

@Zauberbiest · Регистрация: 27.08.2016

Студворк — интернет-сервис помощи студентам

Здравствуйте! Решаю диффур $https://www.cyberforum.ru/cgi-bin/latex.cgi?\ddot{y}+4y=0$ при начальных условиях $https://www.cyberforum.ru/cgi-bin/latex.cgi?y(0) = 1, \dot{y}(0) = 1$ , аналитическим решением которого будет $https://www.cyberforum.ru/cgi-bin/latex.cgi?y = \cos (2x)+\frac{1}{2}\sin(2x)$ .
Для численного решения с odeint понижаю порядок путем замены y' = z(y), y'' = z'(y)*y' = z'(y)*z(y) = z'z.
Получившееся уравнение с разделяющимися переменными интегрирую и провожу обратную замену. Получаем слеудующее
$https://www.cyberforum.ru/cgi-bin/latex.cgi?{\dot{y}}^{2}+4{y}^{2}+const = 0$
вспоминаем начальные условия:
$https://www.cyberforum.ru/cgi-bin/latex.cgi?{\dot{y}}^{2}+4{y}^{2}-5 = 0$
разделяем диффур на "левую" и "правую" части:
$https://www.cyberforum.ru/cgi-bin/latex.cgi?\frac{dy}{dx}=\sqrt{5-4{y}^{2}}$

Полученное уравнение вроде должно решаться численно с помощью odeint.
Однако, я что-то не то делаю в коде:

Python
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
from scipy.integrate import odeint   # Для численного решения ОДУ.
import numpy as np                   # Для массивов и математики.
import matplotlib.pyplot as plt      # Для графиков.
 
# Правая часть диффура.
def f(y, x):
    return np.sqrt(5.0-4.0*y**2)
# Аналитическое решение
def analytical(x):
    return np.cos(2.0*x)+0.5*np.sin(2.0*x)
 
x = np.arange(0.0, 1.0, 0.001)    # Массив точек интегрирования.
y0 = [1.0, 1.0]    # Список начальных значений ( y'(0) = 1, y(0) = 1 ).
 
solution = odeint(f, y0, x)
print solution[0]    # Решение в точке x[0] = 0
 
print 'numerical:'
plt.plot( x, solution )    # График численного решения.
plt.show()                 # Отобразить.
print 'analytical:'
plt.plot( x, analytical(x) )    # График аналитического решения
plt.show()                      # Отобразить.

Потому что выдает он мне вместо численного решения кракозябры (меняющиеся при разных отрезках интегрирования).
Подскажите, где я накосячил? В возвращаемых функцией значениях? Или где-то еще?

@Zauberbiest · 27.11.2017, 09:59 **[ТС]**

Судя по всему, дело в том, что в функции f ("правой" части диффура), стоит выражение под корнем, которое может принимать (и принимает) отрицательные значения. Если сузить отрезок интегрирования до [0; 0.1], то численное решение получается (и даже хорошо "совпадает" с аналитическим).
Но я же хочу получить решение на "нормальном" отрезке, так чтобы был виден характер функции. Что делать?

@Zauberbiest · 27.11.2017, 10:35 **[ТС]**

Оказалось, что можно не понижать порядок исходного уравнения, а записать в функции f() его "правую часть" (выразить y''=-4y):

Python
1
2
3
4
def f(u, x):
    y = u[0]
    z = u[1]
    return [z, -4.0*y]

@Zauberbiest · 27.11.2017, 11:50 **[ТС]**

т.е. получается, порядок-то мы конечно понизили... Заменой y' = p, y'' = p'.
И в такой записи мы численно решаем систему ОДУ первого порядка:
$https://www.cyberforum.ru/cgi-bin/latex.cgi?\left\{\begin{matrix}\\ \frac{dp}{dx} = -4y\\ \frac{dy}{dx} = p\end{matrix}\right.$

Новые блоги и статьи Все статьи Все блоги /
Midnight Chicago Blues kumehtar 24.03.2026 Такой Midnight Chicago Blues, знаешь?. . Когда вечерние улицы становятся ночными, а ты не можешь уснуть. Ты идёшь в любимый старый бар, и бармен наливает тебе виски. Ты смотришь на пролетающие. . .	Контроль уникальности заводского номера - вариант №2 Maks 24.03.2026 В отличие от предыдущего варианта добавлено прерывание циклов, также добавлены новые переменные для сохранения контекста ошибки перед прерыванием цикла: Процедура ПередЗаписью(Отказ, РежимЗаписи,. . .	SDL3 для Desktop (MinGW): Вывод текста со шрифтом TTF с помощью библиотеки SDL3_ttf на Си и C++ 8Observer8 24.03.2026 Содержание блога Финальные проекты на Си и на C++: finish-text-sdl3-c. zip finish-text-sdl3-cpp. zip	Жизнь в неопределённости kumehtar 23.03.2026 Жизнь — это постоянное существование в неопределённости. Например, даже если у тебя есть список дел, невозможно дойти до точки, где всё окончательно завершено и больше ничего не осталось. В принципе,. . .
Модель здравоСохранения: работники работают быстрее после её введения. anaschu 23.03.2026 geJalZw1fLo Корпорация до введения программа здравоохранения имела много невыполненных работниками заданий, после введения программы количество заданий выросло. Но на выплатах по больничным это. . .	Контроль уникальности заводского номера - вариант №1 Maks 23.03.2026 Алгоритм контроля уникальности заводского (или серийного) номера на примере документа выдачи шин для спецтехники с табличной частью в КА2. Данные берутся из регистра сведений, по которому настроено. . .	Хочу заставить корпорации вкладываться в здоровье сотрудников: делаю мат модель здравосохранения anaschu 22.03.2026 e7EYtONaj8Y Z4Tv2zpXVVo https:/ / github. com/ shumilovas/ med2. git	Программный отбор элементов справочника по группе Maks 22.03.2026 Установка программного отбора элементов справочника "Номенклатура" из модуля формы документа в КА2. В качестве фильтра для отбора справочника служит группа номенклатуры. Отбор по наименованию. . .

@Zauberbiest 99 / 2 / 0 Регистрация: 27.08.2016 Сообщений: 39
	27.11.2017, 09:59 [ТС]
	Судя по всему, дело в том, что в функции f ("правой" части диффура), стоит выражение под корнем, которое может принимать (и принимает) отрицательные значения. Если сузить отрезок интегрирования до [0; 0.1], то численное решение получается (и даже хорошо "совпадает" с аналитическим). Но я же хочу получить решение на "нормальном" отрезке, так чтобы был виден характер функции. Что делать? Миниатюры 0

@Zauberbiest 99 / 2 / 0 Регистрация: 27.08.2016 Сообщений: 39
	27.11.2017, 11:50 [ТС]
	т.е. получается, порядок-то мы конечно понизили... Заменой y' = p, y'' = p'. И в такой записи мы численно решаем систему ОДУ первого порядка: $https://www.cyberforum.ru/cgi-bin/latex.cgi?\left\{\begin{matrix}\\ \frac{dp}{dx} = -4y\\ \frac{dy}{dx} = p\end{matrix}\right.$ 0