230 Computational Modelling in Hydraulic and Coastal Engineering
end
% Main program;
for j=2:ny-2
for i=2:nx-2
h=fm(i,j,fDep);
if h>0.1
u=fm(i,j,fCoU);
tbx=1000*fb*u^2;
tbxs=tbx/(1000*up2);
qsx(i,j)=0.05*tbxs^2.5/fb*qss;
v=fm(i,j,fCoV);
tby=1000*fb*v^2;
tbys=tby/(1000*up2);
qsy(i,j)=0.05*tbys^2.5/fb*qss;
end
end
for j=1:ny
qsx(nx-1,j)=qsx(nx-2,j);
qsy(nx-1,j)=qsy(nx-2,j);
end
end
for j=2:ny-2
for i=2:nx-2
dh(i,j)=(qsx(i+1,j)-qsx(i,j))/Dd+(qsy(i,j+1)qsy(i,j))*Dt/Dd;
end
end
for j=1:ny
for i=1:nx
DD(i,j)=dh(i,j)*nt;
% Avoid plotting trivial sediment erosion/deposition
effects;
if abs(DD(i,j))<0.01;
DD(i,j)=0;
end
end
end
[i,j]=meshgrid(1:1:nx,1:1:ny);
df1=fDep(:,3);
df2=reshape(df1,70,40);
df2=df2';
An=ones(size(df2));
idx=find(df2== 0);
An(idx)=0;
figure;
pcolor(An);
colormap(jet);
shading flat;
hold on;
i=1:nx;