Stream function in Channel flow

I am trying to simulate channel flow with D2Q9 LBM.

bbbbbbbb
i o
i o
i o
bbbbbbbb

where i is uniform velocity inlet
b -> wall boundary
o ->open boundary condition

I took this code to calculate stream function from one program, but not getting proper result.


strf(0,0)=0.0
do i=0,n
rhoav=0.5*(rho(i-1,0)+rho(i,0))
if(i.ne.0) strf(i,0)=strf(i-1,0)-rhoav0.5(v(i-1,0)+v(i,0))
do j=1,m
rhom=0.5*(rho(i,j)+rho(i,j-1))
strf(i,j)=strf(i,j-1)+rhom0.5(u(i,j-1)+u(i,j))
end do
end do


Please help me to figure out mistake.