Surface gravity water waves 143
Dx=2;
Dt=0.1;
nx=200;
nt=3040;
r=g*(Dt/Dx)^2/2;
nx2=nx/2;
% Option for variable bed friction;
for i=1:nx
fb(i)=fbo;
end
% Bathymetry data for sloping bed;
for i=1:nx
% For horizontal bed active next line;
% h(i)=ho;
% For sloping bed activate next line;
h(i)=ho-(ho-.5)*i/nx;
bed(i)=5-h(i);
end
for i=1:nx
z(i)=0;
zo(i)=0;
zm(i)=0;
end
% Estimation of the eddy viscosity;
for i=1:nx
ev(i)=0;
if a>h(i)/2.5
ev(i)=ec*h(i)*sqrt(g*h(i));
end
end
c=sqrt(g*h(1));
L=c*T;
zr1=0;
% Main program;
for k=1:nt
% Open sea boundary conditions: Incident and back
radiated wave;
zr2=z(2)-a*sin(2*pi*((k-1)*Dt/T-Dx/L));
zr1=zr1+Dt/Dx*c*(zr2-zr1);
zn(1)=a*sin(2*pi*k*Dt/T)+zr1;
% Solution within the main domain;
for i=2:nx-1
br=ev(i)*Dt/
Dx^2*(z(i+1)-2*z(i)+z(i-1)-zo(i+1)+2*zo(i)-zo(i-1));
tt=r*((h(i+1)+h(i))*(z(i+1)-z(i))-(h(i)+h(i-1))*
(z(i)-z(i-1)));
zn(i)=2*z(i)-zo(i)+tt-Dt*fb(i)*(z(i)-zo(i))+br;
end