Two Step Adam Bashford Method

104 ビュー (過去 30 日間)
B
B 2021 年 9 月 21 日
回答済み: Veena 2023 年 8 月 6 日
I am trying to make a function that implements the two step Adam Bashford Method to solve an ODE
function [t, w, h] = abs2(f, a, b, alpha, n)
%AB2 Two-step Adams Bashforth method
% [t, w, h] = ab2(f, a, b, alpha, n) performs the two-step Adams Bashforth
% method for solving the IVP y' = f(t,y) with initial condition y(a) = alpha
% taking n steps from t = a to t = b. The first step from t = a to t = a + h
% is performed using the modified Euler method.
h=(b-a)/n;
t=a:h:b;
w=zeros(1,length(t));
w(1)=alpha;
%modified euler method
for i=2
k1=h*f(t(i),w(i));
k2=h*f(t(i)+h,w(i)+k1);
w(i+1)=w(i)+1/2*(k1+k2);
end
for i=3:length(t)
w(i+1)=w(i)+(3/2)*h*f(t(i),w(i))-.5*h*f(t(i-1),w(i-1));
end
end
Code to call the function
f=@(t,y) 3*t+y/t;
alpha=5;
a=1;
b=2;
n=3;
[t, w, h] = ab2(f, a, b, alpha, n)
%% Output %%
t = 1×4
1.0000 1.3333 1.6667 2.0000
w = 1×5
5.0000 0 1.6333 3.9567 6.9492
h = 0.3333
As you can the output doesn't look right at all. Any help would be appreciated

採用された回答

Alan Stevens
Alan Stevens 2021 年 9 月 21 日
As follows
f=@(t,y) 3*t+y/t;
alpha=5;
a=1;
b=2;
n=3;
[t, w, h] = abs2(f, a, b, alpha, n);
plot(t,w,'-o'),grid
xlabel('t'),ylabel('w')
disp(['h = ' num2str(h)])
h = 0.33333
function [t, w, h] = abs2(f, a, b, alpha, n)
%AB2 Two-step Adams Bashforth method
% [t, w, h] = ab2(f, a, b, alpha, n) performs the two-step Adams Bashforth
% method for solving the IVP y' = f(t,y) with initial condition y(a) = alpha
% taking n steps from t = a to t = b. The first step from t = a to t = a + h
% is performed using the modified Euler method.
h=(b-a)/n;
t=a:h:b;
w=zeros(1,length(t));
w(1)=alpha;
%modified euler method
k1=h*f(t(1),w(1));
k2=h*f(t(1)+h,w(1)+k1);
w(2)=w(1)+1/2*(k1+k2);
for i=2:length(t)-1 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
w(i+1)=w(i)+(3/2)*h*f(t(i),w(i))-.5*h*f(t(i-1),w(i-1));
end
end
  1 件のコメント
B
B 2021 年 9 月 22 日
thanks

サインインしてコメントする。

その他の回答 (1 件)

Veena
Veena 2023 年 8 月 6 日

Source code for lagrange inverse interpolation

カテゴリ

Help Center および File ExchangeMathematics についてさらに検索

製品


リリース

R2021a

Community Treasure Hunt

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

Start Hunting!

Translated by