2D Lid driven cavity

Lattice Boltzmann Method LBM: Theory
#1

Hi
I was just a beginner. I was trying a fortran code for 2D lid driven cavity flow. I have used bounce bance back rule for walls. Initially I have assigned zero to the boundary nodes 0(left),nx+1(right) and the bottom wall(0) and assigned 1 to top node.

I have applied collision to the boundary nodes(0,nx+1:0,ny+1). (is it recommended for half-way bounce back rule?)
I have not applied streaming to boundary nodes. Should I need to apply streaming to the boundary nodes also?

But i still didnt get the desired velocity plot.

my algorithm is as follows:
initialisation of u,v,f
boundary treatment
calculation of macrovariables
calculation of feq
collision
streaming and again boundary treatment.

I have attached my code for reference.

subroutine inithydro

    implicit double precision(a-h,o-z)
    include'lbe.par'
       do j = 1, ny
          do i = 1, nx
             rho1(i,j)=rho
             u(i,j) = u0
             v(i,j) = 0.0
           enddo
        endif
    end

subroutine hydrovar

    implicit double precision(a-h,o-z)
    include'lbe.par'
    do j = 1,ny
       do i = 1,nx
       rho1(i,j) = ( f(0,i,j) + f(1,i,j) + f(2,i,j) +
 .              f(3,i,j) + f(4,i,j) + f(5,i,j) +
 .              f(6,i,j) + f(7,i,j) + f(8,i,j) )

       u(i,j) = ( f(1,i,j) - f(3,i,j) + f(5,i,j) -
 .              f(6,i,j) - f(7,i,j) + f(8,i,j)) / rho1(i,j)

       v(i,j) = ( f(5,i,j) + f(2,i,j) + f(6,i,j)
 .              - f(7,i,j) - f(4,i,j) - f(8,i,j)) / rho1(i,j)
       rho_boun = ( f(0,i,j) + f(1,i,j) + f(3,i,j)
 .              +2*( f(2,i,j) + f(5,i,j) + f(6,i,j)))
       enddo
       end

subroutine collis

    implicit double precision(a-h,o-z)
    include'lbe.par'
       do j = 0,ny+1
          do i = 0,nx+1
             do k = 0,npop-1
             f(k,i,j) = f(k,i,j) * ( 1.0d0 - omega)
 .                      + omega * feq(k,i,j)
          enddo
       enddo
    enddo
    return
    end

subroutine streaming
implicit double precision(a-h,o-z)
include’lbe.par’

    do j = ny,1,-1
       do i =1, nx
          f(2,i,j) = f(2,i,j-1)
          f(6,i,j) = f(6,i+1,j-1)
       enddo
    enddo

    do j = ny,1,-1
       do i = nx,1,-1
          f(1,i,j) = f(1,i-1,j)
          f(5,i,j) = f(5,i-1,j-1)
       enddo
    enddo

    do j = 1,ny
       do i = nx,1,-1
       f(4,i,j) = f(4,i,j+1)
       f(8,i,j) = f(8,i-1,j+1)
       enddo
    enddo

    do j = 1,ny
       do i = 1,nx
         f(3,i,j) = f(3,i+1,j)
         f(7,i,j) = f(7,i+1,j+1)
      enddo
    enddo

    return
    end

    subroutine bc
    implicit double precision(a-h,o-z)
    include'lbe.par'

c left wall - bounce back
do j = 1,ny
f(1,0,j) = f(3,1,j)
f(5,0,j) = f(7,1,j+1)-(1/2)(f(2,1,j+1)-f(4,1,j+1))
f(8,0,j) = f(6,1,j-1)+(1/2)(f(2,1,j-1)-f(4,1,j-1))
enddo
c right wall - bounce back
do j = 1,ny
f(3,nx+1,j) = f(1,nx,j)
f(7,nx+1,j) = f(5,nx,j-1)+(1/2)(f(2,nx,j-1)-f(4,nx,j-1))
f(6,nx+1,j) = f(8,nx,j+1)-(1/2)(f(2,nx,j+1)-f(4,nx,j+1))
enddo
c North wall - moving

    do i = 1,nx
       f(4,i,ny+1) =f(2,i,ny)
       f(7,i,ny+1) =f(5,i-1,ny)+(1/2)*(f(1,i-1,ny)-f(3,i-1,ny))
 .                          -((1/2)*rho_boun*ux)
       f(8,i,ny+1) =f(6,i+1,ny)+(1/2)*(f(1,i+1,ny)-f(3,i+1,ny))
 .                          +((1/2)*rho_boun*ux)
    enddo

c South wall - bounce back
do i = 1,nx
f(2,i,0) = f(4,i,1)
f(5,i,0) = f(7,i+1,1)-(1/2)(f(1,i+1,1)-f(3,i+1,1))
f(6,i,0) = f(8,i-1,1)+(1/2)(f(1,i-1,1)-f(3,i-1,1))
enddo
return
end

Can anyone help me.

Many thanks in advance.
Ramki.

#2

If you are using half way bounce back, you will have fluid nodes from i=1…nx, j=1…ny, On all of these nodes you do collision and streaming, which requires calculating the moments and equilibrium function. For the boundary conditions you need the f_i that are entering these nodes. For half way bounce back you would say

f_b(fluid,t+1)=f*_a(fluid, t)

where a and b are opposite directions, ‘fluid’ means the fluid nodes next to a boundary and * denotes the post collision distribution. For this method it may be easier to do streaming first and then collisions. Let’s take the south wall as an example. You apply the LB evolution equation everywhere, including j=1. Now fill the “buffer” nodes outside of the fluid domain at j=0 withe the distribution functions you need to find. That is, f_2(j=0)=f_4*(j=1), etc etc. Note that f_4 here is post collision because you have just done the collisions to all nodes in the domain (they are simply the value of f_4 at the end of a timestep). Now you move to the next time step (t+1) and the first thing yo do is streaming. So the f_2 you have defined via bounce back is now streamed into the fluid domain,f_2(fluid)=f_2(buffer). This is precisely what the equation above says.

Basically, you apply streaming and collisions to all fluid nodes, and then fill some ‘buffer’ nodes with the particles that will be streamed into the fluid domain at the next time step. (For very accurate results (especially with low Re numbers), you should used an MRT or TRT model to get rid of the slip error, but BGK should still give reasonable results)

Note also that since the wall is no located half way between nodes (eg between j=0 and j=1) then, if you have divide your unit lengthscale into N fluid nodes then your grid spacing is dx=1/N. For example, if you have Poiesuille flow with a channel height L=1 in the non-dimensional system and lattice fluid nodes form j=1,…N with walls at j=0.5 and j=N+0.5 then dx=1/N.

#3

Thank you pleb.