# ODE simulation: issue with the size of arrays

2 Ansichten (letzte 30 Tage)
Alessandro Bellocchi am 13 Apr. 2023
Bearbeitet: Davide Masiello am 13 Apr. 2023
Hello, everyone. I am trying to simulate the evolution of the following model. However, I get an error with the size of the arrays (Unable to perform assignment because the size of the left side is 1-by-2 and the size of the right side is 1-by-330201). Can anyone help me solve it?
% Define the parameters
m = 0.25;
beta = 1.25;
rho = 0.04;
delta = 0.03;
alpha = 0.2;
pphi = 0.5;
mu = 0.2;
gamma = 1;
% Define the initial condition
phi0 = 10;
theta0 = 1;
% Define the values of sigma for the two situations
sigma_values = [0.1, 0.5];
% Define the time span for the simulation
tspan = [0 500];
% Define the Wiener process
dW = makedist('Normal', 'mu', 0, 'sigma', sqrt(diff(tspan)));
% Define the number of simulations
num_simulations = 5;
% Initialize arrays to store the phi values for each simulation
phi_values = zeros(length(sigma_values), length(tspan), num_simulations);
average_phi_values = zeros(length(sigma_values), length(tspan));
% Perform multiple simulations and store the phi values
for sim = 1:num_simulations
for i = 1:length(sigma_values)
% Define the current value of sigma
sigma = sigma_values(i);
% Define the function to be solved
finance = @(t, y) [
(1/gamma)*((m*y(2)^beta)/(rho+delta-(alpha*(beta*sigma)^2)/2)-(pphi-mu)-delta*y(1));
sigma*y(2)*dW.random()];
% Solve the differential equation numerically
[t, y] = ode45(finance, tspan, [phi0, theta0]);
% Store the phi values for the current simulation and sigma value
phi_values(i, :, sim) = y(:, 1)';
end
end
% Calculate the average phi value over all simulations for each time point and sigma value
for i = 1:length(sigma_values)
for j = 1:length(tspan)
average_phi_values(i, j) = mean(phi_values(i, j, :));
end
end
% Plot the average phi values against time for each sigma value
hold on
for i = 1:length(sigma_values)
plot(tspan, average_phi_values(i, :), 'LineWidth', 2);
end
% Add legend, labels, and title
legend(['\sigma = ', num2str(sigma_values(1))], ['\sigma = ', num2str(sigma_values(2))]);
xlabel('Time');
ylabel('Average \phi_t');
title('Evolution of Finance over Time for Different Values of \sigma and Wiener Process');
##### 1 Kommentar-1 ältere Kommentare anzeigen-1 ältere Kommentare ausblenden
Davide Masiello am 13 Apr. 2023
Unfortunately, you do not know the size of the t vector before the ODE integration itself.
Hence, you cannot initialize phi_values the way you did.
Even if you didn't initialize it, you'd probably end up with a size error, because each simulation will require a different number of steps to complete the integration.
This is due to the fact that MatLab ode solvers are adaptive solvers, with variable step size.

Melden Sie sich an, um zu kommentieren.

### Antworten (1)

Davide Masiello am 13 Apr. 2023
Bearbeitet: Davide Masiello am 13 Apr. 2023
In addition to my comment above, here's a solution to your problem.
It's based on the fact that you can specify tspan as a vector of points at which the solution must be reported.
This provides a homogeneous output regardless of the choice of the integration time steps.
clear,clc
% Define the parameters
m = 0.25;
beta = 1.25;
rho = 0.04;
delta = 0.03;
alpha = 0.2;
pphi = 0.5;
mu = 0.2;
gamma = 1;
% Define the initial condition
phi0 = 10;
theta0 = 1;
% Define the values of sigma for the two situations
sigma_values = [0.1, 0.5];
% Define the time span for the simulation
tspan = linspace(0,500,1000);
% Define the Wiener process
dW = makedist('Normal', 'mu', 0, 'sigma', sqrt(tspan(end)-tspan(1)));
% Define the number of simulations
num_simulations = 5;
% Initialize arrays to store the phi values for each simulation
phi_values = zeros(length(sigma_values), length(tspan), num_simulations);
average_phi_values = zeros(length(sigma_values), length(tspan));
% Perform multiple simulations and store the phi values
for sim = 1:num_simulations
for i = 1:length(sigma_values)
% Define the current value of sigma
sigma = sigma_values(i);
% Define the function to be solved
finance = @(t, y) [
(1/gamma)*((m*y(2)^beta)/(rho+delta-(alpha*(beta*sigma)^2)/2)-(pphi-mu)-delta*y(1));
sigma*y(2)*dW.random()];
% Solve the differential equation numerically
[t, y] = ode45(finance, tspan, [phi0, theta0]);
% Store the phi values for the current simulation and sigma value
phi_values(i, :, sim) = y(:, 1)';
end
end
% Calculate the average phi value over all simulations for each time point and sigma value
for i = 1:length(sigma_values)
for j = 1:length(tspan)
average_phi_values(i, j) = mean(phi_values(i, j, :));
end
end
% Plot the average phi values against time for each sigma value
plot(tspan, average_phi_values, 'LineWidth', 2);
Warning: Imaginary parts of complex X and/or Y arguments ignored.
% Add legend, labels, and title
legend(['\sigma = ', num2str(sigma_values(1))], ['\sigma = ', num2str(sigma_values(2))]);
xlabel('Time');
ylabel('Average \phi_t');
title('Evolution of Finance over Time for Different Values of \sigma and Wiener Process')
##### 0 Kommentare-2 ältere Kommentare anzeigen-2 ältere Kommentare ausblenden

Melden Sie sich an, um zu kommentieren.

### Kategorien

Mehr zu Financial Toolbox finden Sie in Help Center und File Exchange

R2022b

### Community Treasure Hunt

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

Start Hunting!

Translated by