Contents

KETA01 Tenta 2010-03-11 CARMEN AREVALO

clear all
close all

Question 1

Ordered as B, T, X, C, BP, CH4, H2, n1, n2

A=zeros(9,9);
A(1,2)=1; A(1,3)=2; A(1,4)=3; A(1,6)=1; A(1,8)=-2.1;
A(2,1)=6; A(2,2)=8; A(2,3)=10; A(2,4)=12; A(2,5)=10;
  A(2,6)=4; A(2,7)=2; A(2,8)=-10.2; A(2,9)=-2;
A(3,1:4)=1; A(3,5)=2; A(3,8)=-1;
A(4,3)=-1; A(4,8)=0.35*(1-0.74);
A(5,2)=-1; A(5,8)=0.2*(1-.8);
A(6,4)=-1; A(6,8)=0.4*(1-.7);
A(7,8)=-5; A(7,9)=1;
A(8,1:7)=-1; A(8,5)=1000-1;
A(9,8)=1;
format short g
A
b=zeros(9,1);
b(9)=1000;
x=A\b
A =
            0            1            2            3            0            1            0         -2.1            0
            6            8           10           12           10            4            2        -10.2           -2
            1            1            1            1            2            0            0           -1            0
            0            0           -1            0            0            0            0        0.091            0
            0           -1            0            0            0            0            0         0.04            0
            0            0            0           -1            0            0            0         0.12            0
            0            0            0            0            0            0            0           -5            1
           -1           -1           -1           -1          999           -1           -1            0            0
            0            0            0            0            0            0            0            1            0
x =
          737
           40
           91
          120
            6
         1518
         3488
         1000
         5000

Sammansättning

faktor=100/sum(x(1:7));
B=x(1)*faktor
T=x(2)*faktor
X=x(3)*faktor
C=x(4)*faktor
BP=x(5)*faktor
CH4=x(6)*faktor
H2=x(7)*faktor
B =
       12.283
T =
      0.66667
X =
       1.5167
C =
     2
BP =
          0.1
CH4 =
         25.3
H2 =
       58.133

Question 2

R=1;
V=0.75;

Alternative solution 1 We need to find the zero of the function:

f=@(h)V-(pi*h.^2.*(3*R-h))/3;

Plot to find initial guess.

hm=linspace(0,2*R);
y=f(hm);
plot(hm,y)
grid
h=fzero(f,0.5)
h =
      0.53952

Alternative solution 2 Find roots of cubic polynomial in interval [0,2]

p=[pi/3 -pi*R 0 V];
h=roots(p)
h =
       2.9158
      0.53952
     -0.45528

The only root in [0,2] is

h(2)
ans =
      0.53952

Question 3

g=9.81;
C=0.55;
r=0.005;
h0=1.75;
R=1;
% Area of the hole
A=pi*r^2;

Right-hand side of diff eq:

hprime=@(t,h)-C*A*sqrt(2*g*h)./(pi*h.*(2*R-h));

We take a span of one day:

tspan=[0,60*60*24];
type myevents10.m
options=odeset('events',@myevents10);
[t,h,te,he,ie] = ode15s(hprime,tspan,h0,options);
function [value,isterminal,direction] = myevents10(t,h)
%MYEVENTS10 Events function for prob 3 in Tenta KETA01 2010-3-11
   % Find when the value of h is zero
value=h;
isterminal=1;
direction=0;

end


The time is given in seconds

sol=[te,he]
% time in minutes
Tmin=te/60
% time in hours
Thr=Tmin/60
sol =
        24023                4.5857e-007 + 8.277e-013i
Tmin =
       400.38
Thr =
        6.673