The y_RKF(0) = 1 indicates the initial conditions of the ode or something entirely else? cause if i change it to y_RKF(1) = 1 it works but it doesn't output the right solution cause the initial conditions are wrong. Thanks in advance.
RKF array indices error
1 次查看(过去 30 天)
显示 更早的评论
Hello and Merry Christmas to all,
i am trying to solve an ode using RKF with 6th order accuracy.
The form of the ode is y'(t) = f(t,y(t)) and specifically is :
y' = sqrt(y), with initial values y(0) = 1 and the exact solution is y(t) = 1/4(t + 2)^2
i am using this code that i get from mathworks documentation about RKF method
clc; close all; clear all;
%% Inputs and Declaration
h = 0.1; % solution stepsize
x = 0:h:1; % input vector
% Function declaration
f = @(x, y) sqrt(y);
% Initial conditions
y_RKF(0) = 1;
gamma = [16/135, 0, 6656/12825, 28561/56430, -9/50, 2/55]; % vector related to RKF procedure
%% Core: Runge Kutta 6th order procedure
for i = 1:length(x)-1
k1 = h*f(x(i), y_RKF(i)); k1 = double(k1);
k2 = h*f(x(i)+h/4, y_RKF(i)+k1/4); k2 = double(k2);
k3 = h*f(x(i)+3/8*h, y_RKF(i)+3/32*k1+9/32*k2); k3 = double(k3);
k4 = h*f(x(i)+12/13*h, y_RKF(i)+1932/2197*k1-7200/2197*k2+7296/2197*k3); k4 = double(k4);
k5 = h*f(x(i)+h, y_RKF(i)+439/216*k1-8*k2+3680/513*k3-845/4104*k4); k5 = double(k5);
k6 = h*f(x(i)+h/2, y_RKF(i)-8/27*k1+2*k2-3544/2565*k3+1859/4104*k4-11/40*k5); k6 = double(k6);
K = [k1, k2, k3, k4, k5, k6];
y_RKF(i+1) = y_RKF(i) + sum(K.*gamma); % new solution
end
%% Analytical Exact Solution
y_exact = 1/4 .* (x + 2) .^ 2;
error_RKF = y_exact - y_RKF; % error between
But i get this error :
Array indices must be positive integers or logical values.
Error in rkf4_5 (line 13)
y_RKF(0) = 1;
Any help why this happens? Thanks in advance
回答(1 个)
Cris LaPierre
2020-12-30
编辑:Cris LaPierre
2020-12-30
y_RKF(0) = 1 creates an error because MATLAB does not use 0-based indexing. All indexing must use positive integers (1, 2, 3...) or logical indexing.
What your initial condition y(0)=1 means is y=1 when t=0. This corresponds to and index of 1 because x(1)=0.
Last thought - your error is currently the elementwise difference between corresponding values in y_exact and y_RKF. You probably want the total error, so sum the result.
error_RKF = sum(y_exact - y_RKF)
3 个评论
Cris LaPierre
2020-12-30
编辑:Cris LaPierre
2020-12-30
You are confusing function input and indexing. Here, y is a function of t. In your code, you define t as the variable x.
h = 0.1; % solution stepsize
x = 0:h:1; % input vector
So when you say you want y(0)=1, you mean the value of y at t=0 to be 1, or

x and y are variables, each containing values. You want to set the initial condition, or the value of y when t=0, to 1. The first value of x is 0, so the corresponding value in y (the first value) should be set to 1.
x(1)
% this works
y(1)=1
% This does not
y(0)=1
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Partial Differential Equation Toolbox 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!