# Solving second order ODEs with ANN

10 次查看（过去 30 天）
Dimitris Kastoris 2023-5-16

Can i use this method for second order ODEs and is this modification on modelGradients correct?
x = linspace(0,1,10000)';
inputSize = 1;
layers = [
featureInputLayer(inputSize,Normalization="none")
fullyConnectedLayer(10)
sigmoidLayer
fullyConnectedLayer(1)
sigmoidLayer];
dlnet = dlnetwork(layers);
numEpochs = 15;
miniBatchSize =100;
initialLearnRate = 0.1;
learnRateDropFactor = 0.3;
learnRateDropPeriod =5 ;
momentum = 0.9;
icCoeff = 7;
figure
set(gca,YScale="log")
lineLossTrain = animatedline(Color=[0.85 0.325 0.098]);
ylim([0 inf])
xlabel("Iteration")
ylabel("Loss (log scale)")
grid on
velocity = [];
iteration = 0;
learnRate = initialLearnRate;
start = tic;
% Loop over epochs.
for epoch = 1:numEpochs
% Shuffle data.
mbq.shuffle
% Loop over mini-batches.
while hasdata(mbq)
iteration = iteration + 1;
dlX = next(mbq);
% Evaluate the model gradients and loss using dlfeval and the modelGradients function.
% Update network parameters using the SGDM optimizer.
% To plot, convert the loss to double.
loss = double(gather(extractdata(loss)));
% Display the training progress.
D = duration(0,0,toc(start),Format="mm:ss.SS");
title("Epoch: " + epoch + " of " + numEpochs + ", Elapsed: " + string(D))
drawnow
end
% Reduce the learning rate.
if mod(epoch,learnRateDropPeriod)==0
learnRate = learnRate*learnRateDropFactor;
end
end
y = forward(dlnet,dlX);
% Evaluate the gradient of y with respect to x.
% Since another derivative will be taken, set EnableHigherDerivatives to true.
% Define ODE loss.
eq = dy + y/5 - exp(-(dlX / 5)) .* cos(dlX);
% Define initial condition loss.
ic = forward(dlnet,dlarray(0,"CB")) - 0 ;
% Specify the loss as a weighted sum of the ODE loss and the initial condition loss.
loss = mean(eq.^2,"all") + icCoeff * ic.^2;
end
y = forward(dlnet, dlX);
% Evaluate the gradient of y with respect to x.
dy = dlgradient(sum(y, "all"), dlX, 'EnableHigherDerivatives', true);
d2y = dlgradient(sum(dy, "all"), dlX, 'EnableHigherDerivatives', true);
% Define ODE loss.
eq = d2y + 1/5*dy + y + 1/5*exp(-dlX/5).*cos(dlX);
% Compute the initial conditions directly within this function
yAtZero = forward(dlnet, dlarray(0,"CB"));
% Use the gradient dy and evaluate it at x = 0
dyAtZero = forward(dlnet, dlarray(0,"CB"));
% Initial condition errors
icErrorY = yAtZero ; % because y(0) = 0
icErrorDY = dyAtZero - 1; % because y'(0) = 1
% Combine the ODE loss and the initial condition errors.
loss = mean(eq.^2, "all") + icCoeff * (icErrorY.^2 + icErrorDY.^2);
end
##### 1 个评论显示 -1更早的评论隐藏 -1更早的评论
Murat Balc? 2023-11-7，12:24
Hi Dimitris, aren't the following pieces of code identical to each other?
% Compute the initial conditions directly within this function
yAtZero = forward(dlnet, dlarray(0,"CB"));
% Use the gradient dy and evaluate it at x = 0
dyAtZero = forward(dlnet, dlarray(0,"CB"));

### 采纳的回答

Ranjeet 2023-6-8
Hi Dimitris,
I assume that the code you provided has been adopted from the following resource Solve Ordinary Differential Equation Using Neural Network.
The resource guides on training a neural network to estimate solution of a first order differential equation. The important point to be noted is the trained network estimates the analytical solution of the taken example ().
For the case of second order differential equation, we get the analytical solution in many forms including , etc. Here, are found from characteristic equation and are found from the initial conditions (Analytical solution exist in other forms as well).
Changes suggested to adopt the example for second order differential equation are:
1. Modify the modelGradients function to adopt the loss calculation based on second order differential equation taken to solve.
2. Test by increasing number of hidden layers as the analytical solution is now more complex.
3. Try changing other hyperparameters based on the obtained test results.
##### 1 个评论显示 -1更早的评论隐藏 -1更早的评论
thank you very much!

### 类别

Help CenterFile Exchange 中查找有关 Custom Training Loops 的更多信息

R2023a

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by