В этом примере показано, как обучить Физику информированную нейронную сеть (PINN) численно вычислять решение уравнения Бургера при помощи ограниченной памяти BFGS (LBFGS) алгоритм.
Уравнение Бургера является дифференциальным уравнением с частными производными (PDE), которое возникает в различных областях прикладной математики. В частности, гидроаэромеханика, нелинейная акустика, газовая динамика и потоки трафика.
Учитывая вычислительную область, это примеры используют физику сообщила нейронной сети (PINN) [1] и обучают многоуровневую perceptron нейронную сеть, которая берет выборки как введено, где пространственная переменная, и переменная времени и возвращается , где u является решением уравнения Бургера:
с как начальное условие, и и как граничные условия.
Вместо обучения сеть с помощью trainNetwork
функция или использование пользовательского учебного цикла, который обновляет параграницы с помощью sgdmupdate
или подобные функции, этот пример оценивает настраиваемые параметры при помощи fmincon
функция (требует Optimization Toolbox™). fmincon
функция находит минимум ограниченных нелинейных многомерных функций.
Пример обучает модель путем осуществления этого, учитывая вход , выход сети выполняет уравнение Бургера, граничные условия и начальное условие. Чтобы обучить модель, экс-клен использует ограниченную память BFGS (LBFGS) алгоритм, который является приближенным методом ньютона, который аппроксимирует алгоритм Бройдена Флетчера Голдфарба Шэнно.
Обучение эта модель не требует данных о сборе заранее. Можно сгенерировать данные с помощью определения УЧП и ограничений.
Обучение модель требует набора данных узлов коллокации, которые осуществляют граничные условия, осуществите начальные условия и выполните уравнение Бургера.
Выберите 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 имеет ту выход .
Задайте параметры для каждой из операций и включайте их в структуру. Используйте формат parameters.OperationName_ParameterName
где parameters
структура, OperationName
имя операции (например, "fc1") и ParameterName
имя параметра (например, "Веса").
Алгоритм в этом примере требует, чтобы настраиваемые параметры были на первом уровне stucture, не используйте вложенные структуры на этом шаге. fmincon
функция требует, чтобы learnable, чтобы быть удвоилось.
Задайте количество слоев и количество нейронов для каждого слоя.
numLayers = 9; numNeurons = 20;
Инициализируйте параметры для первого полностью операция connect. Первое полностью операция connect имеет два входных канала.
parameters = struct; sz = [numNeurons 2]; parameters.fc1_Weights = initializeHe(sz,2,'double'); parameters.fc1_Bias = initializeZeros([numNeurons 1],'double');
Инициализируйте параметры для каждого остающегося промежуточного звена, полностью соединяют операции.
for layerNumber=2:numLayers-1 name = "fc"+layerNumber; sz = [numNeurons numNeurons]; numIn = numNeurons; parameters.(name + "_Weights") = initializeHe(sz,numIn,'double'); parameters.(name + "_Bias") = initializeZeros([numNeurons 1],'double'); end
Инициализируйте параметры для финала полностью операция connect. Финал полностью операция connect имеет ту выходной канал.
sz = [1 numNeurons]; numIn = numNeurons; parameters.("fc" + numLayers + "_Weights") = initializeHe(sz,numIn,'double'); parameters.("fc" + numLayers + "_Bias") = initializeZeros([1 1],'double');
Просмотрите сетевые параметры.
parameters
parameters = struct with fields:
fc1_Weights: [20×2 dlarray]
fc1_Bias: [20×1 dlarray]
fc2_Weights: [20×20 dlarray]
fc2_Bias: [20×1 dlarray]
fc3_Weights: [20×20 dlarray]
fc3_Bias: [20×1 dlarray]
fc4_Weights: [20×20 dlarray]
fc4_Bias: [20×1 dlarray]
fc5_Weights: [20×20 dlarray]
fc5_Bias: [20×1 dlarray]
fc6_Weights: [20×20 dlarray]
fc6_Bias: [20×1 dlarray]
fc7_Weights: [20×20 dlarray]
fc7_Bias: [20×1 dlarray]
fc8_Weights: [20×20 dlarray]
fc8_Bias: [20×1 dlarray]
fc9_Weights: [1×20 dlarray]
fc9_Bias: [1×1 dlarray]
Создайте функциональный model
, перечисленный в разделе Model Function в конце примера, который вычисляет выходные параметры модели глубокого обучения. Функциональный model
берет в качестве входа параметры модели и сетевые входные параметры, и возвращает выходной параметр модели.
Создайте функциональный modelGradients
, перечисленный в разделе Model Gradients Function в конце примера, который берет в качестве входа параметры модели, сетевые входные параметры и начальные и граничные условия, и возвращает градиенты потери относительно настраиваемых параметров и соответствующей потери.
fmincon
Целевая функцияСоздайте функциональный objectiveFunction
, перечисленный в fmincon
Раздел Objective Function примера, который возвращает потерю и градиенты модели. Функциональный objectiveFunction
берет в качестве входа, вектора из настраиваемых параметров, сетевых входных параметров, начальных условий, и имен и размеров настраиваемых параметров, и возвращает потерю, которая будет минимизирована fmincon
функционируйте и градиенты потери относительно настраиваемых параметров.
Задайте опции оптимизации:
Оптимизируйте использование fmincon
optmizer с алгоритмом LBFGS для не больше, чем 7 500 итераций и вычислений функции.
Оцените с допуском оптимальности 1e-5
.
Предоставьте градиенты для алгоритма.
options = optimoptions('fmincon', ... 'HessianApproximation','lbfgs', ... 'MaxIterations',7500, ... 'MaxFunctionEvaluations',7500, ... 'OptimalityTolerance',1e-5, ... 'SpecifyObjectiveGradient',true);
fmincon
Обучите сеть с помощью fmincon
функция.
fmincon
функция требует, чтобы настраиваемые параметры были заданы как вектор. Преобразуйте параметры в вектор с помощью paramsStructToVector
функция (присоединенный к этому примеру как вспомогательный файл). Чтобы преобразовать назад в структуру параметров, также возвратите названия параметра и размеры.
[parametersV,parameterNames,parameterSizes] = parameterStructToVector(parameters); parametersV = extractdata(parametersV);
Преобразуйте обучающие данные в dlarray
объекты с форматом 'CB'
(образуйте канал, пакет).
dlX = dlarray(dataX','CB'); dlT = dlarray(dataT','CB'); dlX0 = dlarray(X0,'CB'); dlT0 = dlarray(T0,'CB'); dlU0 = dlarray(U0,'CB');
Создайте указатель на функцию с одним входом, который задает целевую функцию.
objFun = @(parameters) objectiveFunction(parameters,dlX,dlT,dlX0,dlT0,dlU0,parameterNames,parameterSizes);
Обновите настраиваемые параметры с помощью fmincon
функция. В зависимости от количества итераций это может требовать времени к запущенному. Чтобы включить подробный многословный выход, установите 'Display'
опция оптимизации к 'iter-detailed'
.
parametersV = fmincon(objFun,parametersV,[],[],[],[],[],[],[],options);
Solver stopped prematurely. fmincon stopped because it exceeded the function evaluation limit, options.MaxFunctionEvaluations = 7.500000e+03.
Для предсказания преобразуйте вектор из параметров к структуре с помощью parameterVectorToStruct
функция (присоединенный к этому примеру как вспомогательный файл).
parameters = parameterVectorToStruct(parametersV,parameterNames,parameterSizes);
Для значений в 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
fmincon
Целевая функцияobjectiveFunction
функция задает целевую функцию, которая будет использоваться алгоритмом LBFGS.
Этот objectiveFunction
funtion берет в качестве входа, вектора из настраиваемых параметров parametersV
, сетевые входные параметры, dlX
и dlT
, начальные и граничные условия dlX0, dlT0
, и dlU0
, и имена и размеры настраиваемых параметров parameterNames
и parameterSizes
, соответственно, и возвращает потерю, которая будет минимизирована fmincon
функциональный loss
и вектор, содержащий градиенты потери относительно настраиваемых параметров gradientsV
.
function [loss,gradientsV] = objectiveFunction(parametersV,dlX,dlT,dlX0,dlT0,dlU0,parameterNames,parameterSizes) % Convert parameters to structure of dlarray objects. parametersV = dlarray(parametersV); parameters = parameterVectorToStruct(parametersV,parameterNames,parameterSizes); % Evaluate model gradients and loss. [gradients,loss] = dlfeval(@modelGradients,parameters,dlX,dlT,dlX0,dlT0,dlU0); % Return loss and gradients for fmincon. gradientsV = parameterStructToVector(gradients); gradientsV = extractdata(gradientsV); loss = extractdata(loss); 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))/2; % 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