Free surface flows 133
for j=3:nx-2
uu(j)=uun(j);
ub(j)=ubn(j);
Hbfinal(j)=hbn(j);
Htotal(j)=hun(j)+hbn(j);
end
Eu=0;
El=0;
for j=3:nx-2
Eu=Eu+(hu(j)+hu(j-1))/2*uu(j)^2;
El=El+(hb(j)+hb(j-1))/2*ub(j)^2;
end
Hupper(i)=hu(nx-2);
Hlower(i)=hb(nx-2);
end
plot(1:j,Hbfinal,'r','Linewidth',1.5)
hold on
plot(1:j,Htotal,'b','Linewidth',1.5)
axis([3,48,0,21])
xlabel('Distance x 200 [m]'), ylabel('Water surface and
interface elevation [m]')
text(10,18.9,'Water free-surface')
text(20,11.5,'Interface at the end of simulation (t =
nt*Dt)')
text(25,18,'(Wind direction is from left to right)')
figure, plot(1:i,Hupper,'b','Linewidth',1.5)
hold on
plot(1:i,Hlower,'r','Linewidth',1.5)
axis([0,100000,6.5,13.5])
xlabel('Time x 10 [s]'), ylabel('Surface oscillation at
right boundary [m]')
text(40000,12,'Upper layer')
text(40000,8,'Lower layer')
PROBLEM 5.6
Solve the same exercise by making the following suggested changes while
keeping the rest of the data constant:
1. Change the relative density difference to 0.01, 0.001 and 0.0001. Run
the simulations until steady-state conditions are achieved. Compare
and comment on your results. Note that full upwelling may destabilize the solution.
2. Change the wind speed to 1 m/s, 5 m/s, 25 m/s and 50 m/s. Run the
simulations until steady-state conditions are achieved. Compare and
comment on your results.
3. Change the code as needed, so that the program simulates a left-toright wind of 10 m/s for 6 hours and then a right-to-left wind of