$\begingroup$I am writing a 3 Step AB Method code to output the computed solutions, time points and step sizes and a log-log plot of the errors. Whilst I input the Time Domain, RHS of the function, Initial Conditions and step sizes.
I am getting the following error, when trying to get my error term for the Log-Log plot and I don't know why:
Undefined operator '-' for input arguments of type 'function_handle'. Error in AB32 (line 39) errAB3 = abs(y(i+2) - G);
I am putting the following in to the command window:
[y, t, h] = AB32([-1,2], @(t, y) -3*y+6*t+5, 2*exp(3)-1 , 39.1592, 39.1474, 10.^-4, @(t) 2*exp(-3*t)+2*t+1)
My code is the following:
% To solve y' = -3y+6t+5 s.t. y(-1) = 2exp(3)-1 for -1 <= t <= 2 using % 3-step Adams-Bashforth method. %Inputs Below: % Dom = Time Domain i.e. 2-(-1)=3 % f = RHS Function % F = Initial Solution % F2 = 2nd Step found using RK2 % F3 = 3rd Step found using RK2 % h = step size. % G = Exact Solution function [y, t, h] = AB32(dom, f, F, F2, F3, h, G) f = @(t, y) -3*y + 6*t + 5; % RHS T = diff(dom); % This is length of the time interval for which you're solving for, i.e. 2-0 = 2. N = ceil(T/h); % total number of times steps h = T/N; % Preallocations: t = zeros(N+1, 1); y = zeros(N+1, 1); % Initializations: y(1) = F; % Initial value of the ODE IVP % Compute the values at point 2 & 3 using RK2 Method: t(2) = -1+h; y(2) = F2; t(3) = -1+(2*h); y(3) = F3; % Main loop for marching N steps: for i = 3:N t(i+3) = -1+(i-1)*h; % time points y(i+3) = y(i+2) + (h/12)*(23*f(t(i+2), y(i+2)) - (16*f(t(i+1), y(i+1))) + 5*f(t(i), y(i))); % 3-step Adams-Bashforth method!!! end %Find Error in Estimation errAB3 = abs(y(i+2) - G); %Plotting of Log-Log Plot loglog(h, errAB3, '.-r', 'MarkerSize', 15) legend('Error in Adams-Bashforth 3 Step Method') xlabel('Step Size h') ylabel('error at point x') title('Error in Adams-Bashforth 3 Step Method') 2 Answers
$\begingroup$Disclaimer: I don't know Matlab.
It appears that G is a function, and that you're supposed to pass a t value into it when using its results as a number. y(i+2) - G(i+2) or similar should work.
You program a function integrate(..) getting the relevant arguments, doing the start steps and the multi-step integration. This you have more or less done. Next you need a main loop over the different 'h', something like
M = 300 h = 0.1^linspace(2,6,M) err = zeros(M,1) for k =1:M [y,t] = integrate(...,h(k)) err(k)=max(abs(y-G(t))) end loglog(h,err,...) If you put everything in one file, do not forget to close subroutines with end.
