thanks Timm

2D rectangular channel, I set parabolic velocity for initial condition,

for cycle = 1:maxT

% MACROSCOPIC VARIABLES

rho = sum(fIn);

ux = reshape ( (cx * reshape(fIn,9,lx*ly)), 1,lx,ly) ./rho;*

uy = reshape ( (cy * reshape(fIn,9,lxly)), 1,lx,ly) ./rho;

```
% MACROSCOPIC (DIRICHLET) BOUNDARY CONDITIONS
% Inlet: pressure BC
rho(:,in,col)=1+0.001*sin(OMEGA*cycle);
ux(:,in,col)=1 - 1 ./ (rho(:,in,col)) .* ( ...
sum(fIn([1,3,5],in,col)) + 2*sum(fIn([4,8,7],in,col));
uy(:,in,col)=0;
```

% Outlet: pressure BC

rho(:,out,col) = 1;

ux(:,out,col) = -1 + 1 ./ (rho(:,out,col)) .* ( …

sum(fIn([1,3,5],out,col)) + 2*sum(fIn([2,6,9],out,col)) );

uy(:,out,col) = 0;

% MICROSCOPIC BOUNDARY CONDITIONS: INLET (Zou/He BC)

fIn(2,in,col) = fIn(4,in,col) + 2/3*rho(:,in,col).**ux(:,in,col);*

fIn(6,in,col) = fIn(8,in,col) + 1/2(fIn(5,in,col)-fIn(3,in,col)) …

+ 1/2rho(:,in,col).*uy(:,in,col) …*

+ 1/6rho(:,in,col).*ux(:,in,col);*

fIn(9,in,col) = fIn(7,in,col) + 1/2(fIn(3,in,col)-fIn(5,in,col)) …

- 1/2*rho(:,in,col).*uy(:,in,col) …*

+ 1/6rho(:,in,col).*ux(:,in,col);

```
% MICROSCOPIC BOUNDARY CONDITIONS: OUTLET (Zou/He BC)
fIn(4,out,col) = fIn(2,out,col) - 2/3*rho(:,out,col).*ux(:,out,col);
fIn(8,out,col) = fIn(6,out,col) + 1/2*(fIn(3,out,col)-fIn(5,out,col)) ...
- 1/2*rho(:,out,col).*uy(:,out,col) ...
- 1/6*rho(:,out,col).*ux(:,out,col);
fIn(7,out,col) = fIn(9,out,col) + 1/2*(fIn(5,out,col)-fIn(3,out,col)) ...
+ 1/2*rho(:,out,col).*uy(:,out,col) ...
- 1/6*rho(:,out,col).*ux(:,out,col);
```

% COLLISION STEP(incomresible model for feq)

for i=1:9

cu = 3*(cx(i)*ux+cy(i)**uy);*

fEq(i,:, =rho.*t(i)+ rho0.* t(i). …

( cu + 1/2(cu.*cu) - 3/2*(ux.^2+uy.^2) );

fOut(i,:, = fIn(i,:, - omega .* (fIn(i,:,:)-fEq(i,:,:));

end

OBSTACLE (BOUNCE-BACK)

for i=1:9

fOut(i,bbRegion) = fIn(opp(i),bbRegion);

end

```
% STREAMING STEP
for i=1:9
fIn(i,:,:) = circshift(fOut(i,:,:), [0,cx(i),cy(i)]);
end
```

thanks a lot