В этом примере показано, как решить уравнение Бургера с помощью глубокого обучения.
Уравнение Бургера является дифференциальным уравнением с частными производными (PDE), которое возникает в различных областях прикладной математики. В частности, гидроаэромеханика, нелинейная акустика, газовая динамика и потоки трафика.
Учитывая вычислительную область, это примеры используют физику сообщила нейронной сети (PINN) [1] и обучают многоуровневую perceptron нейронную сеть, которая берет выборки как введено, где пространственная переменная, и переменная времени и возвращается , где u является решением уравнения Бургера:
с как начальное условие, и и как граничные условия.
Пример обучает модель путем осуществления этого, учитывая вход , выход сети выполняет уравнение Бургера, граничные условия и начальное условие.
Обучение эта модель не требует данных о сборе заранее. Можно сгенерировать данные с помощью определения УЧП и ограничений.
Обучение модель требует набора данных узлов коллокации, которые осуществляют граничные условия, осуществите начальные условия и выполните уравнение Бургера.
Выберите 25 равномерно распределенных моментов времени, чтобы осуществить каждое из граничных условий и .
numBoundaryConditionPoints = [25 25]; x0BC1 = -1*ones(1,numBoundaryConditionPoints(1)); x0BC2 = ones(1,numBoundaryConditionPoints(2)); t0BC1 = linspace(0,1,numBoundaryConditionPoints(1)); t0BC2 = linspace(0,1,numBoundaryConditionPoints(2)); u0BC1 = zeros(1,numBoundaryConditionPoints(1)); u0BC2 = zeros(1,numBoundaryConditionPoints(2));
Выберите 50 равномерно распределенных пространственных точек, чтобы осуществить начальное условие .
numInitialConditionPoints = 50; x0IC = linspace(-1,1,numInitialConditionPoints); t0IC = zeros(1,numInitialConditionPoints); u0IC = -sin(pi*x0IC);
Группа вместе данные для начальных и граничных условий.
X0 = [x0IC x0BC1 x0BC2]; T0 = [t0IC t0BC1 t0BC2]; U0 = [u0IC u0BC1 u0BC2];
Выберите 10 000 точек, чтобы осуществить выход сети, чтобы выполнить уравнение Бургера.
numInternalCollocationPoints = 10000; pointSet = sobolset(2); points = net(pointSet,numInternalCollocationPoints); dataX = 2*points(:,1)-1; dataT = points(:,2);
Создайте datastore массивов, содержащий обучающие данные.
ds = arrayDatastore([dataX dataT]);
Задайте многоуровневую perceptron архитектуру с 9 полностью операции подключения с 20 скрытыми нейронами. Первое полностью операция connect имеет два входных канала, соответствующие входным параметрам и . Последнее полностью операция connect имеет ту выход .
Задайте параметры для каждой из операций и включайте их в struct. Используйте формат parameters.OperationName.ParameterName
где parameters
struct, OperationName
имя операции (например, "fc1") и ParameterName
имя параметра (например, "Веса").
Задайте количество слоев и количество нейронов для каждого слоя.
numLayers = 9; numNeurons = 20;
Инициализируйте параметры для первого полностью операция connect. Первое полностью операция connect имеет два входных канала.
parameters = struct; sz = [numNeurons 2]; parameters.fc1.Weights = initializeHe(sz,2); parameters.fc1.Bias = initializeZeros([numNeurons 1]);
Инициализируйте параметры для каждого остающегося промежуточного звена, полностью соединяют операции.
for layerNumber=2:numLayers-1 name = "fc"+layerNumber; sz = [numNeurons numNeurons]; numIn = numNeurons; parameters.(name).Weights = initializeHe(sz,numIn); parameters.(name).Bias = initializeZeros([numNeurons 1]); end
Инициализируйте параметры для финала полностью операция connect. Финал полностью операция connect имеет ту выходной канал.
sz = [1 numNeurons]; numIn = numNeurons; parameters.("fc" + numLayers).Weights = initializeHe(sz,numIn); parameters.("fc" + numLayers).Bias = initializeZeros([1 1]);
Просмотрите сетевые параметры.
parameters
parameters = struct with fields:
fc1: [1×1 struct]
fc2: [1×1 struct]
fc3: [1×1 struct]
fc4: [1×1 struct]
fc5: [1×1 struct]
fc6: [1×1 struct]
fc7: [1×1 struct]
fc8: [1×1 struct]
fc9: [1×1 struct]
Просмотрите параметры первого полносвязного слоя.
parameters.fc1
ans = struct with fields:
Weights: [20×2 dlarray]
Bias: [20×1 dlarray]
Создайте функциональный model
, перечисленный в разделе Model Function в конце примера, который вычисляет выходные параметры модели глубокого обучения. Функциональный model
берет в качестве входа параметры модели и сетевые входные параметры, и возвращает выходной параметр модели.
Создайте функциональный modelGradients
, перечисленный в разделе Model Gradients Function в конце примера, который берет в качестве входа параметры модели, сетевые входные параметры и начальные и граничные условия, и возвращает градиенты потери относительно настраиваемых параметров и соответствующей потери.
Обучите модель в течение 3 000 эпох с мини-пакетным размером 1 000.
numEpochs = 3000; miniBatchSize = 1000;
Чтобы обучаться на графическом процессоре, если вы доступны, задайте среду выполнения "auto"
. Используя графический процессор требует Parallel Computing Toolbox™ и поддерживаемого устройства графического процессора. Для получения информации о поддерживаемых устройствах смотрите Поддержку графического процессора Релизом (Parallel Computing Toolbox) (Parallel Computing Toolbox).
executionEnvironment = "auto";
Задайте опции оптимизации ADAM.
initialLearnRate = 0.01; decayRate = 0.005;
Обучите сеть с помощью пользовательского учебного цикла.
Создайте minibatchqueue
возразите, что процессы и управляют мини-пакетами данных во время обучения. Для каждого мини-пакета:
Формат данные с размерностью маркирует 'BC'
(пакет, канал). По умолчанию, minibatchqueue
объект преобразует данные в dlarray
объекты с базовым типом single
.
Обучайтесь на графическом процессоре согласно значению executionEnvironment
переменная. По умолчанию, minibatchqueue
объект преобразует каждый выход в gpuArray
если графический процессор доступен.
mbq = minibatchqueue(ds, ... 'MiniBatchSize',miniBatchSize, ... 'MiniBatchFormat','BC', ... 'OutputEnvironment',executionEnvironment);
Преобразуйте начальные и граничные условия в dlarray
. Для точек входных данных задайте формат с размерностями 'CB'
(образуйте канал, пакет).
dlX0 = dlarray(X0,'CB'); dlT0 = dlarray(T0,'CB'); dlU0 = dlarray(U0);
Если обучение с помощью графического процессора, преобразуйте начальную букву и условия к gpuArray
.
if (executionEnvironment == "auto" && canUseGPU) || (executionEnvironment == "gpu") dlX0 = gpuArray(dlX0); dlT0 = gpuArray(dlT0); dlU0 = gpuArray(dlU0); end
Инициализируйте параметры для решателя Адама.
averageGrad = []; averageSqGrad = [];
Ускорьте функцию градиентов модели использование dlaccelerate
функция. Чтобы узнать больше, смотрите, Ускоряют Пользовательские Учебные Функции Цикла.
accfun = dlaccelerate(@modelGradients);
Инициализируйте график процесса обучения.
figure C = colororder; lineLoss = animatedline('Color',C(2,:)); ylim([0 inf]) xlabel("Iteration") ylabel("Loss") grid on
Обучите сеть.
Для каждой итерации:
Считайте мини-пакет данных из мини-пакетной очереди
Оцените градиенты модели и потерю с помощью ускоренных градиентов модели и dlfeval
функции.
Обновите скорость обучения.
Обновите настраиваемые параметры с помощью adamupdate
функция.
В конце каждой эпохи обновите учебный график со значениями потерь.
start = tic; iteration = 0; for epoch = 1:numEpochs reset(mbq); while hasdata(mbq) iteration = iteration + 1; dlXT = next(mbq); dlX = dlXT(1,:); dlT = dlXT(2,:); % Evaluate the model gradients and loss using dlfeval and the % modelGradients function. [gradients,loss] = dlfeval(accfun,parameters,dlX,dlT,dlX0,dlT0,dlU0); % Update learning rate. learningRate = initialLearnRate / (1+decayRate*iteration); % Update the network parameters using the adamupdate function. [parameters,averageGrad,averageSqGrad] = adamupdate(parameters,gradients,averageGrad, ... averageSqGrad,iteration,learningRate); end % Plot training progress. loss = double(gather(extractdata(loss))); addpoints(lineLoss,iteration, loss); D = duration(0,0,toc(start),'Format','hh:mm:ss'); title("Epoch: " + epoch + ", Elapsed: " + string(D) + ", Loss: " + loss) drawnow end
Проверяйте эффективность ускоренной функции путем проверки частоты успешных обращений и уровня заполняемости.
accfun
accfun = AcceleratedFunction with properties: Function: @modelGradients Enabled: 1 CacheSize: 50 HitRate: 99.9984 Occupancy: 2 CheckMode: 'none' CheckTolerance: 1.0000e-04
Для значений в 0,25, 0.5, 0.75, и 1, сравнивают ожидаемые значения модели глубокого обучения с истинными решениями уравнения Бургера с помощью ошибка.
Установите целевые времена тестировать модель в. В течение каждого раза вычислите решение в 1 001 равномерно распределенной точке в области значений [-1,1].
tTest = [0.25 0.5 0.75 1]; numPredictions = 1001; XTest = linspace(-1,1,numPredictions); figure for i=1:numel(tTest) t = tTest(i); TTest = t*ones(1,numPredictions); % Make predictions. dlXTest = dlarray(XTest,'CB'); dlTTest = dlarray(TTest,'CB'); dlUPred = model(parameters,dlXTest,dlTTest); % Calcualte true values. UTest = solveBurgers(XTest,t,0.01/pi); % Calculate error. err = norm(extractdata(dlUPred) - UTest) / norm(UTest); % Plot predictions. subplot(2,2,i) plot(XTest,extractdata(dlUPred),'-','LineWidth',2); ylim([-1.1, 1.1]) % Plot true values. hold on plot(XTest, UTest, '--','LineWidth',2) hold off title("t = " + t + ", Error = " + gather(err)); end subplot(2,2,2) legend('Predicted','True')
Графики показывают, как близко предсказания к истинным значениям.
solveBurgers
функция возвращает истинное решение уравнения Бургера во времена t
как обрисовано в общих чертах в [2].
function U = solveBurgers(X,t,nu) % Define functions. f = @(y) exp(-cos(pi*y)/(2*pi*nu)); g = @(y) exp(-(y.^2)/(4*nu*t)); % Initialize solutions. U = zeros(size(X)); % Loop over x values. for i = 1:numel(X) x = X(i); % Calculate the solutions using the integral function. The boundary % conditions in x = -1 and x = 1 are known, so leave 0 as they are % given by initialization of U. if abs(x) ~= 1 fun = @(eta) sin(pi*(x-eta)) .* f(x-eta) .* g(eta); uxt = -integral(fun,-inf,inf); fun = @(eta) f(x-eta) .* g(eta); U(i) = uxt / integral(fun,-inf,inf); end end end
Модель обучена путем осуществления этого, учитывая вход выход сети выполняет уравнение Бургера, граничные условия и начальное условие. В частности, два количества способствуют потере, которая будет минимизирована:
,
где и .
Здесь, соответствуйте узлам коллокации на контуре вычислительной области и счета и граничное и начальное условие. точки во внутренней части области.
Вычисление требует производных из выхода из модели.
Функциональный modelGradients
берет в качестве входа, параметры модели parameters
, сеть вводит dlX
и dlT
, начальные и граничные условия dlX0
, dlT0
, и dlU0
, и возвращает градиенты потери относительно настраиваемых параметров и соответствующей потери.
function [gradients,loss] = modelGradients(parameters,dlX,dlT,dlX0,dlT0,dlU0) % Make predictions with the initial conditions. U = model(parameters,dlX,dlT); % Calculate derivatives with respect to X and T. gradientsU = dlgradient(sum(U,'all'),{dlX,dlT},'EnableHigherDerivatives',true); Ux = gradientsU{1}; Ut = gradientsU{2}; % Calculate second-order derivatives with respect to X. Uxx = dlgradient(sum(Ux,'all'),dlX,'EnableHigherDerivatives',true); % Calculate lossF. Enforce Burger's equation. f = Ut + U.*Ux - (0.01./pi).*Uxx; zeroTarget = zeros(size(f), 'like', f); lossF = mse(f, zeroTarget); % Calculate lossU. Enforce initial and boundary conditions. dlU0Pred = model(parameters,dlX0,dlT0); lossU = mse(dlU0Pred, dlU0); % Combine losses. loss = lossF + lossU; % Calculate gradients with respect to the learnable parameters. gradients = dlgradient(loss,parameters); end
Модель, обученная в этом примере, состоит из ряда полностью операций подключения с tanh операцией между каждым.
Функция модели берет в качестве входа параметры модели parameters
и сеть вводит dlX
и dlT
, и возвращает выход dlU
модели.
function dlU = model(parameters,dlX,dlT) dlXT = [dlX;dlT]; numLayers = numel(fieldnames(parameters)); % First fully connect operation. weights = parameters.fc1.Weights; bias = parameters.fc1.Bias; dlU = fullyconnect(dlXT,weights,bias); % tanh and fully connect operations for remaining layers. for i=2:numLayers name = "fc" + i; dlU = tanh(dlU); weights = parameters.(name).Weights; bias = parameters.(name).Bias; dlU = fullyconnect(dlU, weights, bias); end end
Maziar Raissi, Париж Пердикэрис и Джордж Эм Карниадакис, Физика Информированное Глубокое обучение (Первая часть): управляемые данными Решения Нелинейных Дифференциальных уравнений с частными производными https://arxiv.org/abs/1711.10561
К. Бэсдевэнт, М. Девилл, П. Холденванг, Ж. Лакруа, Дж. Оуэззэни, Р. Пеирет, П. Орланди, А. Патера, решения для Спектральной и конечной разности уравнения Burgers, Компьютеров & жидкостей 14 (1986) 23–41.
dlarray
| dlfeval
| dlgradient
| minibatchqueue