Contaminant and sediment transport by advection and diffusion 215
jmm=(ny-1)*Dd;
% Initial position of particles;
for i=1:ipp
x(i)=xo;
y(i)=yo;
end
% Number of decaying particles on every time step;
depp=round(ipp*Dt*dec);
diff=sqrt(6*dif/Dt);
% Main program;
for k=1:nt
% Movement of particles by advection and diffusion;
for i=1:ipp
X=fix(x(i)/Dd)+1;
Y=fix(y(i)/Dd)+1;
if X>2 && X2 && Y
u=fm(X,Y,fCoU);
ur=fm(X+1,Y,fCoU);
v=fm(X,Y,fCoV);
vo=fm(X,Y+1,fCoV);
uu=(u+ur)/2;
vv=(v+vo)/2;
tempx=x(i);
tempy=y(i);
ax=rand;
x(i)=x(i)+(2*(ax-0.5)*diff+uu)*Dt;
ay=rand;
y(i)=y(i)+(2*(ay-0.5)*diff+vv)*Dt;
X=fix(x(i)/Dd)+1;
Y=fix(y(i)/Dd)+1;
% Deflection on solid boundarues;
h=fm(X,Y,fDep);
if h<0.1
x(i)=tempx;
y(i)=tempy;
end
end
end
% Decaying of particles;
if dec>0
for i=1:depp
ad=rand;
X=fix(ad*ipp)+1;
x(X)=0;
y(Y)=0;
end
end
% Computation of particle concentration in grid cells;
for i=1:nx
jmm=(ny-1)*Dd;
% Initial position of particles;
for i=1:ipp
x(i)=xo;
y(i)=yo;
end
% Number of decaying particles on every time step;
depp=round(ipp*Dt*dec);
diff=sqrt(6*dif/Dt);
% Main program;
for k=1:nt
% Movement of particles by advection and diffusion;
for i=1:ipp
X=fix(x(i)/Dd)+1;
Y=fix(y(i)/Dd)+1;
if X>2 && X
ur=fm(X+1,Y,fCoU);
v=fm(X,Y,fCoV);
vo=fm(X,Y+1,fCoV);
uu=(u+ur)/2;
vv=(v+vo)/2;
tempx=x(i);
tempy=y(i);
ax=rand;
x(i)=x(i)+(2*(ax-0.5)*diff+uu)*Dt;
ay=rand;
y(i)=y(i)+(2*(ay-0.5)*diff+vv)*Dt;
X=fix(x(i)/Dd)+1;
Y=fix(y(i)/Dd)+1;
% Deflection on solid boundarues;
h=fm(X,Y,fDep);
if h<0.1
x(i)=tempx;
y(i)=tempy;
end
end
end
% Decaying of particles;
if dec>0
for i=1:depp
ad=rand;
X=fix(ad*ipp)+1;
x(X)=0;
y(Y)=0;
end
end
% Computation of particle concentration in grid cells;
for i=1:nx
