Free surface flows 97
jr=nth-1;
end
Q(1)=Qin(jr)+(Qin(jr+1)-Qin(jr))*(k*Dt-(jr-1)*dt)/dt;
% Calculation of the velocity and discharge;
for i=2:nx-1
if h(i-1)<0.001
u(i)=0;
else u(i)=u(i)+Dt*g*(S(i)-(h(i)-h(i-1))/
Dx-u(i)*abs(u(i))/C(i)^2*2/(h(i)+h(i-1)));
end
Q(i)=u(i)*b(i)*(h(i)+h(i-1))/2;
end
% Downstream boundary condition;
Q(nx)=b(nx)*h(nx-1)*sqrt(g*h(nx-1));
% Estimation of the water depths:
for i=1:nx-1
h(i)=h(i)-Dt/Dx/(b(i)+b(i+1))*2*(Q(i+1)-Q(i));
end
% Data stored at different time intervals for plotting;
if k == round(nt/9)
for m=1:nx
h1t(m)=h(m);
Q1t(m)=Q(m);
end
end
if k == round(nt/6)
for m=1:nx
h2t(m)=h(m);
Q2t(m)=Q(m);
end
end
if k == round(nt/4.5)
for m=1:nx
h3t(m)=h(m);
Q3t(m)=Q(m);
end
end
if k == round(nt/3)
for m=1:nx
h4t(m)=h(m);
Q4t(m)=Q(m);
end
end
if k == round(nt/2)
for m=1:nx
h5t(m)=h(m);
Q5t(m)=Q(m);
end
end
if k == nt