76 Computational Modelling in Hydraulic and Coastal Engineering
le=Dx*(nx-1);
% Initial values of u and h;
uo= sqrt(g*r*ho/2/le/fc);
for i=1:nx
u(i)=uo;
h(i)=ho-ho*(i-1)/(nx-1);
hmax(i)=0;
u1(i)=uo;
h1(i)=ho-ho*(i-1)/(nx-1);
h1max(i)=0;
u2(i)=uo;
h2(i)=ho-ho*(i-1)/(nx-1);
h2max(i)=0;
end
for n=1:tn
t=Dt*n;
% Downstream condition for three vane closure times;
if t u(nx)=uo*(tv-t)/tv;
end
if t u1(nx)=uo*(tv1-t)/tv1;
end
if t u2(nx)=uo*(tv2-t)/tv2;
end
% Calculation of new water velocities;
for k=2:nx-1
ff=2*fc*u(k)*abs(u(k))/r;
u(k)=u(k)-Dt/Dx*g*(h(k)-h(k-1))-Dt*ff;
ff1=2*fc*u1(k)*abs(u1(k))/r;
u1(k)=u1(k)-Dt/Dx*g*(h1(k)-h1(k-1))-Dt*ff1;
ff2=2*fc*u2(k)*abs(u2(k))/r;
u2(k)=u2(k)-Dt/Dx*g*(h2(k)-h2(k-1))-Dt*ff2;
end
% Calculation of piezometric head;
for k=2:nx-1
h(k)=h(k)-Dt/Dx*c^2/g*(u(k+1)-u(k));
h1(k)=h1(k)-Dt/Dx*c^2/g*(u1(k+1)-u1(k));
h2(k)=h2(k)-Dt/Dx*c^2/g*(u2(k+1)-u2(k));
% Max piezometric head value for each pipe section;
if h(k)>hmax(k)
hmax(k)=h(k);
end
if h1(k)>h1max(k)
h1max(k)=h1(k);
end
if h2(k)>h2max(k)
h2max(k)=h2(k);
end