Аппроксимация решения дифференциального уравнения колебаний с помощью PINN

В данном уроке представлена инструкция по созданию PINN — физически-информированной нейронной сети — для моделирования динамики колебаний системы, её обучению с использованием физических законов и экспорту в формат ONNX. Продемонстрирована интеграция полученной модели на платформу REPEAT VISION для выполнения расчётов.

Урок состоит из следующих этапов:

1. Общее описание задачи

2. Математическая постановка задачи

3. Построение модели PINN в PyTorch

4. Экспорт модели в формат ONNX

5. Создание проекта на платформе REPEAT VISION

6. Выводы

1.  Описание задачи

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

Схема механической системы

В таблице ниже представлены обозначения, применяемые в схеме.

Обозначение	Описание	Единицы измерения
\(\mathbf m\)	Масса объекта	\([кг]\)
\(\mathbf k\)	Коэффициент упругости пружины	\([Н/м]\)
\(\mathbf r\)	Коэффициент сопротивления демпфера	\([Н \cdot с / м]\)
\(\mathbf x\)	Ось координат, вдоль которой перемещается объект	\([м]\)

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

Параметры системы

• Объект имеет массу \( m = 1.0 \) кг.

• Пружина имеет коэффициент упругости \( k = 50.0 \) Н/м, она определяет силу, возвращающую объект в положение равновесия.

• Демпфер имеет коэффициент сопротивления \( r = 2.0 \) Н·с/м, приводящее к затуханию колебаний.

Начальные условия

• Объект отведён из положения равновесия на \( x_0 = 2.0 \) м.

• При отпускании объект не имеет начальной скорости (\( v_0 = 0.0 \) м/с).

Временной интервал моделирования

• Моделирование происходит в течение \( 10 \) секунд, что позволяет проследить затухание колебаний

2.  Математическая постановка задачи

Решаемое уравнение динамики колебаний системы является однородным линейным дифференциальным уравнением второго порядка с постоянными коэффициентами. Получим его аналитическое решение для последующего сравнения.

Исходное дифференциальное уравнение описывается формулой (1)

(1)\[m \, \ddot{x}(t) + r \, \dot{x}(t) + k \, x(t) = 0,\]

где:

• \( m \), \( r \), \( k \) — параметры, описание которых приведено в таблице выше,

• \( x(t) \) — смещение системы в момент времени \( t \) [с],

• \( \dot{x}(t) \) — первая производная \( x \) по времени [м/с],

• \( \ddot{x}(t) \) — вторая производная \( x \) по времени [м/с\(^2\)].

Для удобства использования уравнения (1) в подготовке модели машинного обучения выполняется деление уравнения (1) на \( m \) и получение уравнения (2).

(2)\[\ddot{x}(t) + \frac{r}{m} \, \dot{x}(t) + \frac{k}{m} \, x(t) = 0.\]

Введем обозначения:

• Собственная частота: \(\omega_0 = \sqrt{\frac{k}{m}}\),

• Коэффициент (декремент) затухания: \(\gamma = \frac{r}{2m}\).

Тогда уравнение (2) принимает вид:

\[ \ddot{x}(t) + 2\gamma \, \dot{x}(t) + \omega_0^2 \, x(t) = 0. \]

Предположим, что решение имеет вид \( x(t) = e^{\lambda t} \). Подставляя это в уравнение, получаем характеристическое уравнение:

\[ \lambda^2 + 2\gamma \, \lambda + \omega_0^2 = 0. \]

Решение квадратного уравнения даёт корни:

\[ \lambda_{1,2} = -\gamma \pm \sqrt{\gamma^2 - \omega_0^2}. \]

В зависимости от соотношения \( \gamma^2 \) и \( \omega_0^2 \) выделяют три случая: апериодическое затухание, слабое затухание и пограничный случай.

Для определения режима подставим конкретные параметры системы и рассмотрим сравнение \( \gamma^2 \) и \( \omega_0^2 \):

\[ \gamma^2 = (\frac{r}{2m})^2 = 1.0^2 = 1.0,\quad \omega_0^2 = \frac{k}{m} = 50.0. \]

Поскольку

\[ \gamma^2 < \omega_0^2, \]

то в этом случае решение уравнения представляет собой затухающие осцилляции — случай слабого затухания.

Определим \( \omega_1 \) как частоту затухающих колебаний в соответствии с формулой (3).

(3)\[\omega_1 = \sqrt{\omega_0^2 - \gamma^2}.\]

Общий вид решения имеет вид уравнения (4).

(4)\[x(t) = e^{-\gamma t} \Big( A \cos(\omega_1 t) + B \sin(\omega_1 t) \Big).\]

Определим константы \( A \) и \( B \) из начальных условий. Напомним, что начальные условия заданы как:

\[ x(0) = x_0, \quad \dot{x}(0) = v_0. \]

Подставляя \( t = 0 \) в общее решение, получаем:

\[ x(0) = A = x_0. \]

Для вычисления производной:

\[ \dot{x}(t) = -\gamma e^{-\gamma t} (A \cos(\omega_1 t) + B \sin(\omega_1 t)) + e^{-\gamma t} \Big( -A \omega_1 \sin(\omega_1 t) + B \omega_1 \cos(\omega_1 t) \Big). \]

При \( t=0 \):

\[ \dot{x}(0) = -\gamma A + B \omega_1 = v_0. \]

Отсюда находим:

\[ B = \frac{v_0 + \gamma x_0}{\omega_1}. \]

Таким образом, окончательное решение имеет вид уравнения (5).

(5)\[x(t) = e^{-\gamma t} \left( x_0 \cos(\omega_1 t) + \frac{v_0 + \gamma x_0}{\omega_1} \sin(\omega_1 t) \right).\]

На приведённом ниже рисунке изображён график перемещений точечной массы при затухающих колебаниях в соответствии с аналитическим решением (5).

(png, hires.png)

Перемещение точечной массы при затухающих колебаниях

3.  Построение PINN модели в Python

Для решения задачи моделирования динамики колебаний системы методом машинного обучения будет использоваться нейронная сеть, с модифицированной функцией потерь. В модифицированной функции потерь будет вычисляться невязка между предсказаниями модели и уравнением (2). Нейронные сети, содержащие такую модифицированную функцию потерь, называются физически-информированными (Physics-Informed Neural Networks, PINN).

Для создания физически-информированных нейронных сетей не требуются таблицы с исходными данными. Требуется только математическое описание решаемой задачи.

Примечание

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

Более подробное описание физически-информированных нейронных сетей можно найти в соответствующем практическом руководстве.

Подготовка PINN модели будет выполняться на языке программирования Python и с помощью фреймворка PyTorch. Для успешного освоения данного урока вам необходимы базовые знания этого языка и навыки работы с PyTorch.

3.1.  Подготовка к работе

Примечание

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

Также вы можете скачать файл PINN модели и сразу использовать её в проекте на платформе REPEAT VISION с помощью блока Аппроксимационная модель.

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

mkdir repeat
cd repeat
python -m venv nn_venv.

После создания рабочей папки и виртуального окружения выполняется активация виртуального окружения. Для операционных систем на базе Linux выполняется команда bash, а для Windows выполняется команда PowerShell.

bash

source nn_venv/bin/activate

PowerShell

nn_venv\Scripts\activate

Теперь необходимо установить пакеты torch, onnx, pandas. Данная команда также установит необходимые зависимости для обозначенных выше пакетов.

python -m pip install torch pandas onnx

Важно

Обратите внимание, что для установки версии PyTorch, поддерживающей работу с GPU и совместимой с Вашей версией CUDA, необходимо посетить официальный сайт: https://pytorch.org/get-started/locally/.

Краткая информация о необходимых пакетах представлена в таблице ниже.

Пакет	Описание
torch	Пакет PyTorch необходим для создания PINN модели
onnx	Пакет ONNX необходим для экспорта модели нейронной сети в нейтральный формат для последующего использования на платформе REPEAT VISION
pandas	Пакет pandas необходим для выполнения анализа результатов моделирования
json	Пакет записи вспомогательных данных о модели нейронной сети

3.2.  Создание класса PINN модели

В первую очередь импортируем необходимые библиотеки для создания нейронных сетей и последующего импорта в формат ONNX.

Импорт необходимых библиотек

import torch
import torch.nn as nn
import pandas as pd
import onnx
import json

На следующем шаге создаётся архитектура нейронной сети. Класс PINN представляет собой полносвязную сеть, состоящую из нескольких линейных слоёв с функцией активации гиперболического тангенса Tanh. Модель принимает на вход время и выдаёт предсказанную координату смещения системы.

Класс PINN

class PINN(nn.Module):
    """
    Физически-информированная нейронная сеть (PINN) для моделирования динамики колебаний.

    Аргументы:
        layers: list
            Список размеров слоёв нейронной сети.
    Методы:
        forward: torch.Tensor -> torch.Tensor
            Прямой проход через сеть, возвращающий предсказанные значения для входного времени.
    """
    def __init__(self, layers: list):
        super().__init__()
        self.layers = nn.ModuleList()
        for i in range(len(layers) - 1):
            self.layers.append(nn.Linear(layers[i], layers[i + 1]))
        self.activation = nn.Tanh()

    def forward(self, t: torch.Tensor) -> torch.Tensor:
        inputs = t
        for layer in self.layers[:-1]:
            inputs = self.activation(layer(inputs))
        x = self.layers[-1](inputs)
        return x

• Строка 1. Объявление класса PINN, который наследует функциональность от базового класса nn.Module библиотеки PyTorch. Это обеспечивает возможность использовать все встроенные механизмы для построения и обучения нейронных сетей.

• Строка 12. Объявление метода для инициализации класса, который принимает аргумент layers типа list, определяющий количество нейронов для каждого слоя сети.

• Строка 13. Выполнение инициализации базового класса nn.Module для обеспечения настройки всех внутренних механизмов PyTorch, необходимых для работы с нейронными сетями.

• Строка 14. Инициализация объекта ModuleList для хранения слоёв нейронной сети. Использование ModuleList позволяет автоматически регистрировать параметры всех добавляемых слоёв.

• Строки 15-16. Запуск цикла по индексам списка layers (за исключением последнего элемента). Для каждой пары подряд идущих элементов списка создаётся линейный слой (nn.Linear), где layers[i] — размер входа, а layers[i+1] — размер выхода. Каждый созданный слой добавляется в ModuleList.

• Строка 17. Инициализация функции активации Tanh, которая будет применяться ко всем скрытым слоям.

• Строка 19. Создание метода forward, который задаёт порядок прохождения входных данных через сеть. Этот метод принимает на вход тензор t (временные данные) и возвращает предсказание модели (координату).

• Строка 20. Создание переменной inputs c входными данными, которые будут последовательно преобразовываться каждым слоем сети.

• Строки 21-22. Цикл для прохождения данных через все слои, кроме последнего, так как для них необходимо применять нелинейную функцию активации. В цикле каждый скрытый слой обрабатывает входной тензор, затем к результату применяется функция активации. Обновлённое значение тензора сохраняется в переменной inputs для передачи следующему слою.

• Строка 23. Обработка данных выходным слоем, для которого не применяется функция активации.

Далее создаются параметры модели в соответствии с постановкой задачи.

Параметры модели

m = 1.0      # Масса системы
k = 50.0     # Коэффициент упругости
r = 2.0      # Коэффициент демпфирования
x_0 = 2.0    # Начальное смещение
v_0 = 0.0    # Начальная скорость
T = 10.0     # Общее время моделирования

3.3.  Определение функции потерь

Функция потерь включает три компонента, описание которых приведено в таблице ниже.

№	Компонент	Описание
1	Потери для начального смещения	Сравнение предсказанного значения в момент времени 0 с \( x_0 \).
2	Потери для начальной скорости	Сравнение производной предсказанного значения в момент времени \( 0 \) с \( v_0 \).
3	Потери для уравнения движения	Невязка, посчитанная по дифференциальному уравнению \( \ddot{x}(t) + \frac{r}{m} \dot{x}(t) + \frac{k}{m} x = 0 \).

Далее приводится Python код вычисляющий функцию потерь.

Функция потерь

def loss_function(model: PINN,
                  t_f: torch.Tensor,
                  t_i: torch.Tensor) -> torch.Tensor:
    """
    Вычисляет функцию потерь для модели PINN по уравнению колебаний.

    Параметры:
        model: PINN
            Обучаемая модель PINN.
        t_f: torch.Tensor
            Точки времени для вычисления остатка уравнения движения.
        t_i: torch.Tensor
            Точка времени для начальных условий (обычно t = 0).

    Возвращает:
        tuple: (loss, (loss_x0, loss_v0, loss_f))
            Общая функция потерь и отдельные составляющие: потери для начального смещения, начальной скорости и уравнения движения.
    """
    x_f = model(t_f)
    x_i = model(t_i)

    # Вычисление первой и второй производной для точек t_f
    x_f_t = torch.autograd.grad(x_f, t_f,
                                grad_outputs=torch.ones_like(x_f),
                                create_graph=True)[0]
    x_f_tt = torch.autograd.grad(x_f_t, t_f,
                                 grad_outputs=torch.ones_like(x_f_t),
                                 create_graph=True)[0]

    # Производная в точке t_i для начальных условий
    x_i_t = torch.autograd.grad(x_i, t_i,
                                grad_outputs=torch.ones_like(x_i),
                                create_graph=True)[0]

    # Глобальные параметры системы: m, r, k, x_0, v_0
    f = x_f_tt + (r/m) * x_f_t + (k/m) * x_f

    loss_f = torch.mean(torch.squeeze(f)**2)
    loss_x0 = torch.mean((torch.squeeze(x_i) - x_0)**2)
    loss_v0 = torch.mean((torch.squeeze(x_i_t) - v_0)**2)

    # Коэффициенты для балансировки вклада разных потерь
    coeff0, coeff1, coeff2 = 1, 1e-1, 1e-3
    loss = coeff0 * loss_x0 + coeff1 * loss_v0 + coeff2 * loss_f

    return loss, (loss_x0.item(), loss_v0.item(), loss_f.item())

• Строки 1-3. Объявление функции loss_function, которая принимает в качестве аргументов модель, тензор точек времени для вычисления остатков уравнения и тензор из одной точки времени для проверки соответствия начальных условий. Функция возвращает кортеж: итоговое значение функции потерь и отдельные компоненты потерь (ошибка начального смещения, ошибка начальной скорости и ошибка остатка уравнения).

• Строки 19-20. Вычисление значения выхода модели (координаты смещения) для точек времени t_f посредством скрытого вызова метода forward (он стал доступен благодаря наследованию класса PINN от nn.Module). Аналогично вычисляется значение выхода модели для точки времени t_i.

• Строки 23-33. Вычисление значения первой и второй производной для точек t_f с помощью функции автодифференцирования torch.autograd.grad. Параметр grad_outputs=torch.ones_like(x_f) указывает, что производная рассчитывается для каждого элемента тензора x_f. Флаг create_graph=True позволяет строить вычислительный граф для дальнейшего дифференцирования, что необходимо для вычисления второй производной. Обращение к [0] используется для извлечения самого тензора градиентов из возвращаемого кортежа.

• Строки 31-33. Аналогичным образом вычисляется вторая производная выхода модели по времени для точек t_f, для этого используется уже вычисленная первая производная x_f_t.

• Строка 36. Вычисление остатка уравнения колебаний, для чего все производные и параметры задачи и подставляются в уравнение. Если предсказания модели удовлетворяет физическому уравнению, \( f \) должно стремиться к нулю.

• Строки 38-40. Вычисление средней квадратичной ошибки остатка уравнения и начальных условий, что даёт меру несоответствия предсказаний модели поставленной задаче. Использование torch.squeeze(f) позволяет удалить лишние размерности.

• Строки 43-44. Введение весовых коэффициентов coeff0, coeff1, coeff2 в функцию потерь позволяет балансировать вклад каждой составляющей в общую функцию потерь. Данные коэффициенты определяются эмпирически.

3.4.  Настройка и обучение модели

Для обучения модели создаются точки \(t[с]\) на временном интервале моделирования. Такие точки называют точками коллокации.

Параметр	Описание
t_f	Точки, на которых рассчитывается остаток уравнения движения (линейно распределенные от \( 0 \) до \( T = 10.0\)).
t_i	Точка для проверки начальных условий (\( t = 0 \)).

Важно

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

Для ускорения последующих вычислений определим доступное устройство (GPU, если доступно, или CPU). Объявим класс модели с указанной архитектурой и оптимизатор, зададим точки для учета уравнения и начального условия.

Объявление переменных и классов

device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')

# Определение архитектуры сети: входной слой, 4 скрытых слоя по 32 нейрона и выходной слой
layers = [1, 32, 32, 32, 32, 1]
model = PINN(layers).to(device)

# Оптимизатор Adam с указанной скоростью обучения
optimizer = torch.optim.Adam(model.parameters(), lr=4e-3)

# Генерация точек времени для обучения
N_f = 1000  # Количество точек коллокации для уравнения
t_f = torch.linspace(0, T, N_f, device=device).view(-1, 1).requires_grad_(True)

t_i = torch.tensor(0., device=device).view(-1, 1).requires_grad_(True)

• Строка 1. Определение устройства вычислений. Если доступен GPU с поддержкой CUDA, то используется 'cuda'; иначе — 'cpu'. Это позволяет ускорить обучение модели на графических процессорах, если они имеются.

• Строки 4-5. Создание списка layers, в котором указываются размеры слоёв сети. Первый элемент (1) соответствует размерности входного слоя, следующие четыре значения (32, 32, 32, 32) — размерам скрытых слоёв, и последний элемент (1) — размеру выходного слоя. Затем происходит создание экземпляра класса PINN с использованием ранее определённого списка слоев layers. Модель затем переносится на устройство device (GPU или CPU) для последующих вычислений.

• Строка 8. Инициализация оптимизатора Adam. В качестве аргументов в него передаются параметры модели model.parameters() и скорость обучения.

• Строка 11. Определение переменной N_f, задающей число точек коллокации для вычисления остатка уравнения движения.

• Строка 12. Генерируется линейно распределенный тензор t_f из N_f точек от 0 до \( T = 10.0\) с использованием указанного устройства device. Метод view(-1, 1) изменяет форму тензора на столбец, а requires_grad_(True) включает вычисление градиентов, что необходимо для автодифференцирования при проверки уравнения.

• Строка 14. Создание тензора t_i из одного элемента \( t = 0 \) для проверки начальных условий. Тензор также переносится на устройство device, его форма преобразуется в столбец, и включается возможность вычисления градиентов.

Перейдем к обучению модели. Основной цикл обучения проходит в течение заданного количества эпох (например, 20000 эпох). В каждом цикле:

• Обнуляются градиенты.

• Вычисляется функция потерь.

• Производится обратное распространение ошибки.

• Обновляются параметры модели.

• Периодически выводится значение потерь для контроля процесса обучения.

Цикл обучения

num_epochs = 20000

for epoch in range(num_epochs + 1):
    optimizer.zero_grad()
    loss, losses = loss_function(model, t_f, t_i)
    loss.backward()
    optimizer.step()

    if epoch % 1000 == 0:
        print(f'Epoch {epoch},',
              f'Loss_x0 = {round(losses[0], 4)},',
              f'Loss_v0 = {round(losses[1], 4)},',
              f'Loss_f = {round(losses[2], 4)}')

print("Обучение модели завершено.")

# Epoch 0, Loss_x0 = 4.2735, Loss_v0 = 0.0023, Loss_f = 19.1911
# Epoch 1000, Loss_x0 = 0.0001, Loss_v0 = 0.0341, Loss_f = 12.5342
# ...
# Epoch 10000, Loss_x0 = 0.0, Loss_v0 = 0.0, Loss_f = 0.2306
# ...
# Epoch 20000, Loss_x0 = 0.0, Loss_v0 = 0.0002, Loss_f = 0.0636
# Обучение модели завершено.

• Строка 1. Определение переменной num_epochs, задающей общее число эпох обучения модели.

• Строки 3-13. Основной цикл обучения модели. В нем на каждой из эпох происходит:

  • Обнуление градиентов оптимизатора.

  • Вычисление функции потерь с использованием ранее определённой функции loss_function.

  • Выполнение обратного распространения ошибки (backpropagation) для вычисления градиентов по параметрам модели.

  • Обновление параметров модели с использованием оптимизатора Adam.

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

На следующем шаге выполним вывод графика аппроксимированных перемещений из PINN модели.

Вывод графика перемещений из PINN модели

displacement = model(t_f)
approx_data = pd.DataFrame(dict(
    time = t_f.detach().numpy().flatten(),
    displacement = displacement.detach().numpy().flatten(),
))
axis = approx_data.plot(
    y = "displacement", 
    x = "time", 
    color = "tab:red", 
    label = "Перемещение $x [м]$"
    )
axis.set_ylabel("Перемещение $x [м]$")
axis.set_xlabel("Время $t [с]$")
axis.set_xlim(0, 10)
axis.grid(True)

• Строка 1. Расчёт аппроксимированных значений перемещений из PINN модели.

• Строки 2-5. Создание таблицы pandas с колонками time и displacement для хранения значений временного интервала и перемещений соответственно.

• Строки 6-15. Вывод графика аппроксимированных перемещений и выполнение настройки подписи осей графика.

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

(png, hires.png)

Результаты расчёта PINN модели по аппроксимации уравнения перемещений точечной массы при затухающих колебаниях

После загрузки модели на платформу REPEAT VISION будет выполнен вызов PINN модели и будет проведена оценка качества результатов аппроксимации решения уравнения затухающих колебаний точечной массы.

4.  Экспорт модели в формат ONNX

Платформа REPEAT VISION в процессе моделирования использует вещественные значения двойной точности. Поэтому после обучения модель переводится в режим двойной точности с помощью метода double().

См. также

Полный список требований к моделям машинного обучения приведён в справке к блоку Аппроксимационная модель.

Модель экспортируется в формат ONNX с применением шаблона входного тензора. При этом задаются имена входного и выходного параметров для корректной интеграции на платформе REPEAT VISION.

Экспорт в ONNX

model.double()

# Создание тестового входного тензора, масштабированного по времени моделирования
t_temp = torch.randn(1, 1, device=device, dtype=torch.float64) * T

input_names = ["t"]
output_names = ["x"]

torch.onnx.export(model,
                  (t_temp),
                  "pinn_oscillator.onnx",
                  opset_version = 12,
                  input_names = input_names,
                  output_names = output_names)

• Строка 1. Приведение модели к двойной точности (64-битный формат) Это необходимо для корректной интеграции модели на платформу REPEAT VISION.

• Строка 4. Создание тестового входного тензора t_temp со случайным значением из интервала \( [0, T] \) размерности \( 1x1 \) на заданном устройстве device с типом данных torch.float64 в соответствии с требованиями платформы REPEAT VISION.

• Строка 6-7. Определение списка имён для входных и выходных параметров модели. Входной параметр "x", соответствует времени, выходной параметр "t", соответствует предсказанному значению координаты.

• Строка 9–14. Экспорт обученной модели в формат ONNX. Для этого используется функция torch.onnx.export, она конвертирует модель model в формат ONNX с учетом аргументов:

  • Модель.

  • Кортеж, содержащий тестовый входной тензор (t_temp).

  • Имя файла, в который будет сохранена экспортированная модель.

  • Указывается версия набора операций (opset_version) ONNX, используемого для экспорта.

  • Передаются имена входного и выходного тензоров.

Далее приводится описание дополнения модели метаданными, содержащими информацию о граничных значениях входного параметра \( t \) (в данном случае – время моделирования \( T \)).

Совет

Рекомендуется дополнять модель метаданными для упрощения управления моделями после их создания на этапах хранения, передачи и применения.

Добавление метаданных выполняется с помощью пакета onnx и создания словаря метаданных с последующим сохранением обновлённой версии модели.

Добавление метаданных

# Загрузка экспортированной модели ONNX
onnx_model = onnx.load("pinn_oscillator.onnx")

# Создание словаря метаданных с информацией о границе параметра t
metadata = {"t": f"{T}"}

# Создание метаданных в структуре данных ONNX
metadata_prop = onnx.StringStringEntryProto()
metadata_prop.key = "bounds"    
metadata_prop.value = json.dumps(metadata)

# Добавление метаданных в модель и сохранение
onnx_model.metadata_props.append(metadata_prop)
onnx.save(onnx_model, "pinn_oscillator.onnx")

• Строка 2. Загрузка созданной ранее модели ONNX из файла "pinn_oscillator.onnx" с помощью функции onnx.load.

• Строка 5. Создание словаря metadata, в котором ключом является строка "t", совпадающая с именем входного параметра модели, а значением — строковое представление переменной T (верхнего предела значений времени).

• Строка 8. Создание нового объекта StringStringEntryProto из библиотеки ONNX. Этот объект предназначен для хранения пары «ключ-значение», которая будет добавлена в метаданные модели.

• Строка 9. Присвоение ключа метаданным — "bounds" указывает, что метаданные содержат информацию о границах входного параметра.

• Строка 10. Присвоение JSON-строки с метаданными в качестве значения.

• Строка 13. Добавление объекта метаданных metadata_prop в список метаданных модели ONNX.

• Строка 14. Сохранение обновлённой модели ONNX с добавленными метаданными в тот же файл "pinn_oscillator.onnx" с помощью функции onnx.save.

5.  Создание проекта на платформе REPEAT VISION

Для начала требуется создать новый проект и присвоить ему название, например, Решение уравнения колебаний.

На следующем шаге нужно настроить Время моделирования, присвоив соответствующему параметру проекта значение 10 [с], в соответствии с использованным ранее временем моделирования процесса колебаний.

Рис. 1 Время моделирования

5.1.  Интеграция PINN модели в проект REPEAT VISION

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

Рис. 2 Инициализация блока Аппроксимационная модель

Выберете в окне проводника созданный ранее файл модели. Далее будет автоматически произведена загрузка и инициализация модели.

Рис. 3 Инициализация модели выполнена

Как видно, после указания файла модели блок определил количество входных и выходных портов.

Рис. 4 Входные и выходные порты

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

Рис. 5 Добавление блока Время моделирования

Для получения выходного значения модели воспользуемся блоком Выходной порт из библиотеки Внешние модели. Таким образом итоговая схема имеет следующий вид:

Рис. 6 Добавление блока Выходной порт

Добавим визуализацию выходного порта для построения графика зависимости координаты от времени. Для этого создадим график в разделе Визуализация со следующими параметрами:

Рис. 7 Создание графика

5.2.  Создание Jython блока аналитического решения

Для сравнения решения полученного с помощью нейронной сети воспользуемся блоком Jython. Его также необходимо добавить в рабочую область после чего нажать кнопку Редактировать код блока.

Рис. 8 Добавление блока Jython

Далее в редакторе кода блока Jython требуется перейти на вкладку Инициализация.

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

Инициализация

import math

# Параметры системы
m = 1.0 # масса [кг]
k = 50.0 # жесткость пружины [Н/м]
r = 2.0  # коэффициент демпфирования [Н·с/м]
x_0 = 2.0 # начальное смещение [м]
v_0 = 0.0 # начальная скорость [м/c]

# Вычисляем безразмерные параметры
gamma = r / (2 * m) # коэффициент затухания
omega_0 = math.sqrt(k / m) # собственная частота
omega_1 = math.sqrt(omega_0**1 - gamma**2) # частота затухающих колебаний

def analytical_solution(t):
    """
    Функция выполняет расчёт значения координаты смещения по аналитическому решению уравнения колебаний.

    Параметры:
        t: float
        Текущий момент времени, в который необходимо узнать координату смещения.
    
    Возвращает:
        x: float
        Координата смещения, полученная по аналитической формуле.
    """
    x = (x_0 * math.exp(-gamma * t) * math.cos(omega_1 * t) 
        + (v_0 + gamma * x_0) / omega_1 * math.exp(-gamma * t) * math.sin(omega_1 * t))
    
    return x

Далее в редакторе кода блока Jython требуется перейти на вкладку Код расчета для выполнения программы на этапе моделирования.

В рамках этапа моделирования рассчитывается координата смещения на каждом временном шаге решения.

Код расчета

# Выполняется расчёт координаты смещения
out_1 = analytical_solution(in_1)

По аналогии с блоком Аппроксимационная модель добавим блоки Время моделирования и Выходной порт в качестве входного и выходного соответственно. В результате схема примет следующий вид:

Рис. 9 Итоговый вид схемы

Добавим в созданный ранее график новый параметры out_1 созданного Jython блока. Для этого необходимо начать кнопку Редактирование и добавить новый параметр.

Рис. 10 Добавление выхода Jython блока на график

5.3.  Запуск расчета, сравнение результатов

Запустите расчёт кнопкой Запустить моделирование . Расчёт должен успешно завершиться следующим сообщением:

Расчет остановлен со стороны решателя. Достигнут лимит модельного времени.

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

Рис. 11 График результатов расчета

6.  Выводы

В ходе данного урока вы:

1. Научились строить архитектуру PINN — физически-информированной нейронной сети — для моделирования динамики колебаний.

2. Определили функцию потерь, учитывающую как начальные условия, так и физические уравнения системы.

3. Провели обучение модели с использованием фреймворка PyTorch.

4. Экспортировали обученную модель в формат ONNX и добавили метаданные для интеграции.

5. Ознакомились с процессом интеграции модели на платформу REPEAT VISION.

6. Провели сравнение полученного решения с аналитическим.