Borrar filtros
Borrar filtros

End a Calculation if a condition is satisfied.

6 visualizaciones (últimos 30 días)
DIP
DIP el 9 de Feb. de 2017
Respondida: DIP el 20 de Feb. de 2017
I have a calculation to perform:
Part A (Number of iterations are known)
ODE is solved using Runge Kutta 4 and I get an answer
Part B (I dont know the iterations to solve this)
ODE is solved using euler implicit
Condition (Part A-Part B <= 10e-5)
The program should end
Can you help me out with the code ?
[EDITED, Jan, moved from comment section:]
clear all;close all; clc;
format long;
u=1.5;
S=-33;
L=12/100;
N=11;
delx=L/(N-1);
K(1)=0.414;
y(1)=0;
C(1)=0.414;
x(1)=0;
%4th Order RK Method
for i=1:N-1
x(i+1)=x(i)+delx;
k1 = delx*(S/u*C(i));
k2 = delx*(S/u*C(i)+S/u*k1/2);
k3 = delx*(S/u*C(i)+S/u*k2/2);
k4 = delx*(S/u*C(i)+S/u*k3);
C(i+1) = C(i) + (1/6)*(k1+2*k2+2*k3+k4);
end
rk4=C(i+1);
N_euler=7839;
for j=1:7839 % i am not supposed to know that the iteration ends here,
it should be dynamic. (When N_euler = 7839, rk4-K(j+1)=10e-5)
%Euler Implicit
delxeuler=L/(N_euler-1);
y(j+1)=y(j)+delxeuler;
K(j+1)=(K(j)/delxeuler)*(1/((1/delxeuler)-(S/u)));
if abs(rk4-K(j+1))<0.00001
break
end
end
plot(y,K,'-.+','color','g')
title('Concentration of CO vs. Distance');
xlabel('Axial(x) Direction [m]');
ylabel('Concentration of CO[mol/m3]');
hold on
plot(x,C,'-o','color','r')
legend('Euler Implicit: N=7389','Runge Kutta 4th Order: N=11');

Respuesta aceptada

DIP
DIP el 20 de Feb. de 2017
while abs(rk4-imp)>=0.00001
N_euler = N_euler+1;
delxeuler=L/(N_euler-1);
for j=1:N_euler-1
%Euler Implicit
y(j+1)=y(j)+delxeuler;
K(j+1)=(K(j)/delxeuler)*(1/((1/delxeuler)-(S/u)));
imp=K(j+1);
err=rk4-imp;
end
end

Más respuestas (1)

Rik
Rik el 9 de Feb. de 2017
The command you are looking for is 'while' instead of 'for'.
So replace the "for j=1:7839" in your code with "while abs(rk4-K(j+1))>0.00001" and you should be golden.
Note: with while loops you have to increment your counter yourself, otherwise your code will loop forever, so add a "j=j+1;" below your while
  2 comentarios
DIP
DIP el 9 de Feb. de 2017
but K(j+1) is not defined. and will throw an error. Also, how do I define the stepsize delxeuler ?
Rik
Rik el 9 de Feb. de 2017
There are two options:
initialize K before the loop to a value start will result in starting the loop (e.g. K=10^6*rk4),
or you copy the content of the loop, so it runs once before it gets to the while

Iniciar sesión para comentar.

Categorías

Más información sobre Loops and Conditional Statements en Help Center y File Exchange.

Community Treasure Hunt

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

Start Hunting!

Translated by