Skip to main content

Simple Euler Method

##Usage:Call Octave from terminal
##and then call EulerMethodUmitAlkus.m
##from octave and finally
##press enter. That's all.
##Simple Euler Method
##Constants and initializations
x=[]; ## initial empty vector for x
y=[]; ## initial empty vector for y
x(1)=1; ## initial value of x
y(1)=1; ## initial value of y
h=1E-3; ## increment in x
dery=[]; ## 1st derivative of y wrt x
n=1; ## inital loop index for while
## enter the while loop for the interval x=[1,2]
while (x(n)<=2)
x(n+1)=x(n)+h;
dery(n+1)=x(n)*x(n)-2*y(n)/x(n); ##given
y(n+1)=y(n)+h*dery(n); ##Euler method
n++;
endwhile
##exit from the 1st while loop
##Modified Euler Method
##Constant and initializations
ymid=[]; ## empty vector function evaluated at x midpoint
xmid=[]; ## empty vector func. of midpoints in x
ymid(1)=1; ## inital value for ymid.
derymid=[]; ## derivative of y at midpoints
##Enter the 2nd while loop
n=1;
while (x(n)<=2)
xmid(n)=x(n)+h/2;
derymid(n)=xmid(n)*xmid(n)-2*ymid(n)/xmid(n);
ymid(n+1)=ymid(n)+derymid(n)*h/2;
n++;
endwhile
## Plot for Simple Euler Method
subplot(2,1,1)
hold on
plot(x,y,'r-');
title('Simple Euler Method');
xlabel('x');
plot(x,dery,'b-');
legend('y','dy/dx')
hold off
##Plot for Modified Euler Method
subplot(2,1,2)
plot(x,ymid,'g-');
title('Modified Euler Method');
xlabel('x');
ylabel('y-modified');
print('-dpsc','EulerMethodUmitAlkus.ps')
save -text EulerMethodUmitAlkus.dat

Comments

Popular posts from this blog

NEWTON-RAPSON METHOD FOR HEAT FLOW

##Constants and initializations
a=5.67E-8; ## Stefan-Boltzman constant[Watt/meter^2Kelvin^4]
e=0.8; ## Rod surface emissivity [Dimensionless]
h=20; ## Heat transfer coefficient of air flow [W/m^2-K]
Tinf=Ts=25; ## Temperature of air and the walls of the close[Celcius]
D=0.1; ## Diameter of the rod[meter]
I2R=100; ## Electric power dissipated in rod (Ohmic Heat)[W]
T=[]; ## Temperature of the rod[*C]
T(1)=25; ## Initial guess of the temperature of the rod[*C]
Q=[]; ## Heat function [W]
Qp=[]; ## First derivative of Q wrt T [W/C*].
for i=1:100
Q(i)=pi*D*(h*(T(i)-Tinf)+e*a*(T(i)^4-Ts^4))-I2R;
Qp(i)=pi*D*(h+4*e*a*T(i)^3);
T(i+1)=T(i)-Q(i)/Qp(i); ## Newton-Rapson Method
endfor
printf('The steady state temperature is %f\n',T(i+1))
save -text HeatFlowTemp.dat
## The plot
t=1:100; ##temperature
for n=1:100
H(n)=pi*D*(h*(t(n)-Tinf)+e*a*(t(n)^4-Ts^4))-I2R;
endfor
plot(t,H)
xlabel('T(Celcius)');
ylabel('Q(Watt)');
legend('Q(T)');
title('Heat flow vs Temperatu…