Автор статьи: Артем Михайлов

Дифференциальные уравнения — это мощный инструмент, который используется во многих областях науки и техники. От физики и химии до экономики и биологии, они предоставляют основу для моделирования различных процессов. Дифференциальные уравнения описывают динамику системы, позволяя нам понять, как система меняется со временем или в ответ на различные воздействия.

Python, являясь одним из самых популярных языков программирования в научном сообществе, предлагает отличные возможности для работы с дифференциальными уравнениями. Он обладает мощными библиотеками для символьных вычислений, численного решения уравнений и визуализации, что делает его замечательным инструментом для изучения и применения дифференциальных уравнений.

В этой статье мы рассмотрим, как использовать Python для решения различных типов дифференциальных уравнений, начиная от простейших обыкновенных дифференциальных уравнений и заканчивая более сложными системами и частными дифференциальными уравнениями. Мы будем использовать несколько различных методов, включая как аналитические, так и численные подходы, и рассмотрим различные библиотеки Python, которые можно использовать для этих целей.

Введение в Python для решения дифференциальных уравнений


Python, как уже было упомянуто, является отличным инструментом для работы с дифференциальными уравнениями. Этот язык программирования предлагает несколько мощных библиотек, которые специализируются на решении таких уравнений, позволяя нам проводить как аналитические, так и численные вычисления.

Одной из ключевых библиотек, которую мы будем использовать в этой статье, является SymPy. SymPy — это библиотека Python для символьных математических вычислений, которая позволяет нам проводить аналитическое решение дифференциальных уравнений. Это означает, что мы можем получить точное математическое решение наших уравнений, когда это возможно.

Для численного решения дифференциальных уравнений мы будем использовать библиотеку SciPy. SciPy — это основная библиотека для научных вычислений в Python, которая предоставляет множество функций для численного решения дифференциальных уравнений, включая различные методы, такие как метод Эйлера и метод Рунге-Кутты.

Наконец, для решения частных дифференциальных уравнений мы будем использовать библиотеку FiPy. FiPy является специализированной библиотекой Python для решения частных дифференциальных уравнений методом конечных объемов.

Чтобы начать работу с этими библиотеками, вам нужно их установить. Это можно сделать, используя менеджер пакетов Python pip:

pip install sympy scipy fipy matplotlib

После установки этих библиотек вы готовы начать решать дифференциальные уравнения с Python! В следующих разделах мы подробно рассмотрим, как использовать каждую из этих библиотек для решения различных типов дифференциальных уравнений.

Аналитическое решение дифференциальных уравнений


Аналитическое решение дифференциальных уравнений означает поиск точного математического формула, которая описывает решение. Это может быть особенно полезно, когда у нас есть необходимость в точности, или когда мы хотим получить глубокое понимание поведения нашей системы. Однако, стоит отметить, что не все дифференциальные уравнения могут быть решены аналитически.

Использование SymPy для аналитического решения


Для аналитического решения дифференциальных уравнений мы будем использовать библиотеку SymPy. SymPy предоставляет функции для манипулирования и решения символьных математических выражений, что делает её идеальным инструментом для этой задачи.

Пример 1: Обыкновенное дифференциальное уравнение первого порядка
Прежде всего, давайте рассмотрим пример аналитического решения простого обыкновенного дифференциального уравнения первого порядка с использованием SymPy.

import sympy as sp

# Определяем переменные
t = sp.symbols('t')
y = sp.Function('y')

# Определяем уравнение
equation = sp.Eq(y(t).diff(t), -2 * y(t))

# Решаем уравнение
solution = sp.dsolve(equation)

print(solution)

В этом примере мы решаем дифференциальное уравнение первого порядка dy/dt = -2y, где y — это функция, которую мы хотим найти, а t — это независимая переменная, по отношению к которой мы дифференцируем. Функция dsolve в SymPy используется для решения дифференциальных уравнений.

Пример 2: Обыкновенное дифференциальное уравнение второго порядка
Теперь давайте рассмотрим пример обыкновенного дифференциального уравнения второго порядка.

# Определяем уравнение
equation = sp.Eq(y(t).diff(t, t), -y(t))

# Решаем уравнение
solution = sp.dsolve(equation)

print(solution)

В этом случае, мы решаем уравнение вида d^2y/dt^2 = -y, где y — это функция, которую мы хотим найти, а t — это независимая переменная. diff(t, t) означает взятие второй производной по t.

Аналитическое решение дифференциальных уравнений может быть мощным инструментом для получения точных и информативных решений. Однако, как мы увидим в следующих разделах, не всегда возможно (или практично) получить аналитическое решение, и в этих случаях мы должны обратиться к численным методам.

Численное решение дифференциальных уравнений


В предыдущем разделе мы рассмотрели основы численного решения дифференциальных уравнений с использованием метода Эйлера. Однако, метод Эйлера — это всего лишь один из многих численных методов, доступных для решения дифференциальных уравнений. Другие методы, такие как методы Рунге-Кутты, могут быть более точными или стабильными для определенных типов дифференциальных уравнений.

Метод Рунге-Кутты второго порядка (RK2)


Метод Рунге-Кутты второго порядка (RK2), также известный как метод средней точки, предлагает более точное приближение, чем метод Эйлера, за счет добавления дополнительного шага, который учитывает скорость изменения в середине каждого временного интервала. Давайте рассмотрим пример решения того же дифференциального уравнения dy/dt = -2y с использованием метода RK2.

# Определяем начальное условие
y[0] = 1

# Решаем уравнение с помощью метода Рунге-Кутты второго порядка
for i in range(1, time.shape[0]):
    k1 = alpha * y[i-1]
    k2 = alpha * (y[i-1] + dt * k1 / 2)
    y[i] = y[i-1] + dt * k2

# Визуализируем результат
plt.plot(time, y)
plt.xlabel('Time')
plt.ylabel('y')
plt.title('Solution of dy/dt = -2y using RK2 method')
plt.show()

Метод Рунге-Кутты четвертого порядка (RK4)


Один из наиболее известных и широко используемых численных методов — это метод Рунге-Кутты четвертого порядка (RK4). RK4 улучшает точность по сравнению с RK2, используя четыре «предсказания» скорости изменения на разных моментах времени в течение каждого временного интервала.

# Определяем начальное условие
y[0] = 1

# Решаем уравнение с помощью метода Рунге-Кутты четвертого порядка
for i in range(1, time.shape[0]):
    k1 = alpha * y[i-1]
    k2 = alpha * (y[i-1] + dt/2 * k1)
    k3 = alpha * (y[i-1] + dt/2 * k2)
    k4 = alpha * (y[i-1] + dt * k3)
    y[i] = y[i-1] + dt/6 * (k1 + 2*k2 + 2*k3 + k4)

# Визуализируем результат
plt.plot(time, y)
plt.xlabel('Time')
plt.ylabel('y')
plt.title('Solution of dy/dt = -2y using RK4 method')
plt.show()

Численные методы, такие как методы Эйлера и Рунге-Кутты, являются мощными инструментами для решения дифференциальных уравнений, когда аналитические решения недоступны или слишком сложны. Однако, они также вводят численные ошибки, особенно при больших временных шагах или в долгосрочной перспективе. Поэтому всегда важно проверять качество численного решения и выбирать подходящий метод и размер шага для конкр# Численное решение дифференциальных уравнений

Решение систем дифференциальных уравнений


Системы дифференциальных уравнений встречаются во многих областях науки и инженерии, когда несколько величин взаимосвязаны и взаимодействуют друг с другом.

Давайте рассмотрим простую систему из двух дифференциальных уравнений первого порядка, которая описывает динамику хищник-жертва, известную как модель Лотки-Вольтерра:

dx/dt =  alpha*x -   beta*x*y
dy/dt = delta*x*y - gamma*y

где x — популяция жертв (например, кроликов), y — популяция хищников (например, лис), а alpha, beta, delta и gamma — параметры, которые описывают взаимодействие между двумя популяциями.

Мы можем решить эту систему уравнений с помощью метода Рунге-Кутты 4-го порядка следующим образом:

import numpy as np
import matplotlib.pyplot as plt

# Задаем параметры
alpha = 1.0
beta = 0.1
delta = 0.075
gamma = 1.5

# Задаем начальные условия
x = 10
y = 5

# Задаем временной интервал и шаг
t = 0
t_end = 100
dt = 0.01

# Создаем списки для сохранения значений популяции и времени
time = []
x_values = []
y_values = []

while t < t_end:
    # Подсчет инкрементов для x и y
    dx1 = dt * (alpha * x - beta * x * y)
    dy1 = dt * (delta * x * y - gamma * y)
    
    dx2 = dt * (alpha * (x + dx1/2) - beta * (x + dx1/2) * (y + dy1/2))
    dy2 = dt * (delta * (x + dx1/2) * (y + dy1/2) - gamma * (y + dy1/2))
    
    dx3 = dt * (alpha * (x + dx2/2) - beta * (x + dx2/2) * (y + dy2/2))
    dy3 = dt * (delta * (x + dx2/2) * (y + dy2/2) - gamma * (y + dy2/2))
    
    dx4 = dt * (alpha * (x + dx3) - beta * (x + dx3) * (y + dy3))
    dy4 = dt * (delta * (x + dx3) * (y + dy3) - gamma * (y + dy3))
    
    # Обновляем значения x и y
    x += (dx1 + 2*dx2 + 2*dx3 + dx4) / 6
    y += (dy1 + 2*dy2 + 2*dy3 + dy4) / 6
    
    # Сохраняем текущие значения
    time.append(t)
    x_values.append(x)
    y_values.append(y)
    
    # Обновляем время
    t += dt

# Рисуем график численности популяций
plt.plot(time, x_values, label='Preys')
plt.plot(time, y_values, label='Predators')
plt.xlabel('Time')
plt.ylabel('Population')
plt.legend()
plt.show()

В этом примере мы использовали метод Рунге-Кутты 4-го порядка для решения системы дифференциальных уравнений, описывающих динамику популяций хищников и жертв. Этот метод дает хорошую точность и стабильность для широкого диапазона дифференциальных уравнений и их систем.

Но не забывайте, что для более сложных систем, особенно тех, которые включают в себя стохастические или частичные дифференциальные уравнения, может потребоваться использование более сложных методов или специализированных программных библиотек, таких как scipy.integrate в Python или deSolve в R.

Решение частных дифференциальных уравнений


Частные дифференциальные уравнения (ЧДУ) являются основным инструментом описания многих естественных и технических явлений. Они находят широкое применение в физике, инженерии, экономике и других науках.

Пример — Уравнение теплопроводности.

Уравнение теплопроводности — это классический пример ЧДУ. Оно описывает, как температура меняется с течением времени и пространства в данной области:

где u — температура, t — время, x — расстояние, и k — коэффициент теплопроводности.

В Python мы можем использовать метод конечных разностей для решения этого уравнения. Сначала импортируем необходимые библиотеки:

import numpy as np
from scipy import sparse
from scipy.sparse.linalg import spsolve
import matplotlib.pyplot as plt

Создадим сетку:

N = 50 # Количество узлов в сетке
x = np.linspace(0, 1, N)

Определим коэффициент k и начальные условия:

k = 0.1
u0 = np.sin(np.pi * x)

Создадим матрицу системы уравнений:

diagonals = np.zeros((3, N)) # Создаем матрицу с тремя диагоналями
diagonals[0,:] = -1
diagonals[1,:] = 2
diagonals[2,:] = -1
A = sparse.spdiags(diagonals, [-1,0,1], N, N, format="csr") # Создаем разреженную матрицу

Решим систему уравнений и визуализируем результат:

for i in range(100): # 100 временных шагов
    b = k * u0
    u = spsolve(A, b) # Решаем систему уравнений
    plt.plot(x, u) # Визуализируем результат
    u0 = u
plt.show()

Результат — это график, показывающий, как температура меняется с течением времени и пространства.

Это всего лишь один из многих примеров, как можно решить ЧДУ в Python. Важно отметить, что сложность решения может значительно возрастать в зависимости от характеристик конкретного ЧДУ, включая его порядок, тип, и начальные и граничные условия.

Заключение


Надеемся, что эта статья помогла вам понять, как использовать Python для решения дифференциальных уравнений, и что она станет отправной точкой для дальнейшего изучения и использования Python в этой области. Python — мощный инструмент для проведения вычислений и анализа данных. Его применение в решении дифференциальных уравнений — всего лишь один из многих примеров его эффективности.

Напоследок хочу порекомендовать вам бесплатный урок, в рамках которого эксперты OTUS расскажут, что общего между программистом и математиком и в чем между ними разница.
Будет показано, что языки программирования строятся на тех же принципах, что доказательства математических теорем, а выразительные способности языков программирования ограничиваются понятиями полноты и непротиворечивости.

  1. Зарегистрироваться на бесплатный урок

Комментарии (17)


  1. Anton_Menshov
    17.07.2023 14:00
    +14

    PDE - partial differential equations. Но в русском языке устоявшаяся терминология "уравнения в частных производных". Аббревиатуры ЧДУ употребляется очень редко и является калькой с английского. Используются аббревиатуры ОДУ (обыкновенные дифференциальные уравнения) и УРЧП (уравнения в частных производных), иногда встречаются и другие вариации. Можно продолжать и дальше по статье.

    Если пишете\переводите статьи по мат физике, пожалуйста, уважайте терминологию, для математики это важно. Здесь нужны писатели (или копирайтеры\переводчики) с соответствующим опытом.


  1. belch84
    17.07.2023 14:00
    +4

    В предыдущем разделе мы рассмотрели основы численного решения дифференциальных уравнений с использованием метода Эйлера.
    На самом деле в предыдущем разделе вы рассматривали аналитическое решение дифференциальных уравнений, а вовсе не метод Эйлера, метод Эйлера вы вообще не рассматривали
    Численные методы, такие как методы Эйлера и Рунге-Кутты, являются мощными инструментами для решения дифференциальных уравнений, когда аналитические решения недоступны или слишком сложны. Однако, они также вводят численные ошибки, особенно при больших временных шагах или в долгосрочной перспективе. Поэтому всегда важно проверять качество численного решения и выбирать подходящий метод и размер шага для конкр# Численное решение дифференциальных уравнений
    На самом деле для реального решения уравнений методами Рунге-Кутты используется адаптивное изменение шага, которое позволяет на участках, где решение изменяется мало, двигаться с бОльшим шагом, а там, где изменения решения значительны, измельчать шаг до необходимого минимума.
    Вообще, все упомянутое в комментарии (вместе с термином «частные дифференциальные уравнения»), наводит на мысль, что автору в подготовке статьи помогает какой-то не очень умный ИИ


  1. Travisw
    17.07.2023 14:00
    +3

    Даже формулы уравнения нормальные не сделали и еще sympy умеет pdsolve


  1. dlinyj
    17.07.2023 14:00
    -1

    Если что, на хабре можно нормально формулы писать:

    В этом случае, мы решаем уравнение вида d^2y/dt^2 = -y, где y — это функция, которую мы хотим найти, а t — это независимая переменная. diff(t, t) означает взятие второй производной по t.

    Если правильно распарсил, то будет вот так:

    $\frac{\text{d}^{2y}}{\text{d}t^{2}}=-y$


    Сама формула:

    \frac{\text{d}^{2y}}{\text{d}t^{2}}=-y

    Набирал тут: primat.org/mathred/mathred.html


    1. dlinyj
      17.07.2023 14:00

      administrator поломались формулы на хабре, видимо.
      Починились, а то не отображались.


    1. dlinyj
      17.07.2023 14:00

      Вот пример, крутой статьи с формулами и питоном, чтобы понять уровень https://habr.com/ru/articles/748430/


    1. Travisw
      17.07.2023 14:00
      -1

      почему у тебя в числителе y в степени?


      1. dlinyj
        17.07.2023 14:00
        -1

        Опечатка. Если правильно парсить, то возможно автор имел в виду это:


        d^2y/dt^2 = -y


        тогда будет вот так


        $d^{\frac{2y}{dt^{2}}}=-y$


        Скобки, надо не забывать скобки! Поэтому хорошо бы авторам осваивать latex и делать нормальные формулы, чтобы не вводить в заблуждение.


        1. Travisw
          17.07.2023 14:00
          -1

          Левый d в числитель а 2 в степень, а "y" справа


          1. dlinyj
            17.07.2023 14:00
            -1

            Я записываю исключительно то, что они дали в своей формуле.


            1. Travisw
              17.07.2023 14:00

              ты не правильно распарсил короче


              1. dlinyj
                17.07.2023 14:00
                -1

                Вы можете исправить положение и написать как правильно. Написал ровно то, что там было записано.


                1. Travisw
                  17.07.2023 14:00

                  Нафиг надо


  1. Uranix
    17.07.2023 14:00
    +1

    Использовать spsolve для трехдиагональной матрицы - это из пушки по воробьям. Трехдиагональные системы решают прогонкой, в питоне для этого есть solve_banded.

    Ну и метод для уравнения теплопроводности никак не может решать уравнение теплопроводности просто потому что в нем даже размер шага по времени не участвует.


  1. kbtsiberkin
    17.07.2023 14:00
    +1

    Странно расписывать простейшие методы по методичкам. Идём в scipy и там лежит вполне сносный integrate.RK45, например. С адаптивным шагом и любой размерностью задачи.

    Летняя практика у студентов где-то идёт что ли?


    1. mikko_kukkanen
      17.07.2023 14:00

      Рунге-Кутт он же только для задачи Коши. Для краевых задач в scipy ничего нет кроме странного solve_bvp.


  1. mikko_kukkanen
    17.07.2023 14:00
    +1

    В библиотеках python по численным методам много подводных камней. Начиная с того, что в python нет нормальной процедуры даже нахождения первой производной функции, заданной именно численно. т.е. набором данных. Библиотека на эту тему удалена из scipy, видимо, из-за безысходности. Два предлагаемых взамен репозитория https://github.com/pbrod/numdifftools и https://github.com/maroba/findiff/ - не о том. Особой проблемы посчитать по 3- или 5- точечной схеме, конечно, нет, но никто не делает как полагается из соображений наилучшей аппроксимации.

    В scipy для диффуров первого порядка нет никакого упоминания про неединственность решения, ничего про условие Липшица.

    С задачей Коши все выглядит более-менее, а вот с краевыми задачами для диффуров второго порядка в scipy.integrate, как-то черезчур наворочено.

    Что касается уравнений в частных производных, особенно в 2 или 3 размерных координатах, то либо профессиональные CAD-средства типа solidworks, либо ничего. Про МКЭ вообще никто не вспоминает.

    Мое мнение - в python нормальных готовых решений нет, все самим придется писать. Учите численные методы и диффуры с функаном.