To plot the absolute percent error

6 ビュー (過去 30 日間)
Naveen Krish
Naveen Krish 2022 年 3 月 18 日
コメント済み: Naveen Krish 2022 年 3 月 18 日
Hi, I'm doing a project on formation of disinfection byproducts during water treatment to solve this problem i have to plot the absolute percent error in the concentration A at each timestep of 3600 secs from ti = 0 to tf=12 hours. I need to compare the results between Euler's method and 4th order RK. I have wrriten the code but need to determine the absolute percent error. I have tried searching in MATLAB Answers pertaining to my topics of study, but to no avail. The MATLAB code and the error message are shown below:
%%Eulers method
nsteps = 12;
t = zeros (nsteps,1);
A = zeros (nsteps,1);
B = zeros(nsteps, 1);
P = zeros(nsteps,1);
A(1) = 1;
B(1) = 3;
C(1) = 0;
K = 5*10^-5
for k = 2:13
t(k) = t(k-1)+3600
A(k) = A(k-1)+(-K*A(k-1)*B(k-1))*3600;
B(k) = B(k-1)+(-Yb*(K*A(k-1)*B(k-1)))*3600;
P(k)= P(k-1)+ Yp*(K*A(k-1)*B(k-1))*3600;
end
%%RK Method
h=3600;
A0=1;
B0=3;
P0=0;
K=5*10^-5;
Yb=1;
Yp=0.15;
t = 0:h:43200;
fyt = @(t,y) [(-K*y(1)*y(2));
(-Yb*(K*y(1)*y(2)));
(Yp*(K*y(1)*y(2)))];
Y = zeros(3,numel(t))
Y(1,1) = 1.0;
Y(2,1) = 3.0;
Y(3,1) = 0;
for i=1 : numel(t)-1
k1 = fyt(t(i),Y(:,i));
k2 = fyt(t(i)+0.5*h,Y(:,i)+0.5*h*k1);
k3 = fyt(t(i)+0.5*h,Y(:,i)+0.5*h*k2);
k4 = fyt(t(i)+h,Y(:,i)+h*k3);
Y(:,i+1) = Y(:,i) + (h/6)*(k1+2*k2+2*k3+k4);
end
For the following parameter values: 𝑌𝐵 = 1, 𝑌𝑃 = 0.15, 𝐴0 = 1 mg/L, 𝐵0 = 3 mg/L, 𝑃0 = 0, 𝐾 = 5E − 5
Use the analytical solution to plot the absolute percent error in the concentration of A at each timestep for both Euler’s method and 4th order RK. (Notice how the error grows with each time step).
I have also attached the code/data for your convenience. I appreciate your help with troubleshooting the problem. Thanks for considering my request.
  2 件のコメント
Naveen Krish
Naveen Krish 2022 年 3 月 18 日
Naveen Krish
Naveen Krish 2022 年 3 月 18 日

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

回答 (0 件)

カテゴリ

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

製品


リリース

R2021b

Community Treasure Hunt

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

Start Hunting!

Translated by