poiseuille-boundary condition Zou/He

Hi,

again I have a question to the Zou/He boundary conditions for the d2q9 model. I applied poiseuille flow to a small gap between 2 parallel walls.It’s working quite good, but at the outlet I have an increase in speed of about 10 percent. So I tried to give the top and bottom nodes of the outlet (actually wall nodes) a special treat, but still the results don’t improve much. Following is the code.


    do z=1,nz_lb
        if (spalt_2D(x,z) == 2) then                        ! gap
            rho_outlet = f_2D(0,x,z) + f_2D(1,x,z) + f_2D(2,x,z) + f_2D(3,x,z) &
                       + f_2D(4,x,z) + f_2D(5,x,z) + f_2D(6,x,z) + f_2D(7,x,z) + f_2D(8,x,z)

            u_outlet = (f_2D(4,x,z) + f_2D(0,x,z) + f_2D(2,x,z) + 2.d0 *(f_2D(1,x,z) 
                        + f_2D(5,x,z) + f_2D(8,x,z)))/rho_outlet - 1.d0
            
            f_2D(3,x,z) = f_2D(1,x,z) - 2.d0/3.d0 * rho_outlet * u_outlet
            f_2D(6,x,z) = f_2D(8,x,z) + 1.d0/2.d0 * (f_2D(4,x,z) - f_2D(2,x,z)) 
                             - 1.d0/6.d0 * u_outlet * rho_outlet
            f_2D(7,x,z) = f_2D(5,x,z) + 1.d0/2.d0 * (f_2D(2,x,z) - f_2D(4,x,z)) 
                             - 1.d0/6.d0 * u_outlet * rho_outlet
        end if
        if (z==1) then                                      ! lower wall
            rho_outlet = f_2D(0,x,z) + f_2D(1,x,z) + f_2D(2,x,z) + f_2D(3,x,z) &
                       + f_2D(4,x,z) + f_2D(5,x,z) + f_2D(6,x,z) + f_2D(7,x,z) + f_2D(8,x,z)
            
            f_2D(3,x,z) = f_2D(1,x,z)
            f_2D(2,x,z) = f_2D(4,x,z)
            f_2D(6,x,z) = f_2D(8,x,z)
            f_2D(5,x,z) = 1.d0/2.d0 * (rho_outlet - (f_2D(0,x,z) + f_2D(1,x,z) + f_2D(2,x,z) 
                             + f_2D(3,x,z) + f_2D(4,x,z) + f_2D(6,x,z) + f_2D(8,x,z)))
            f_2D(7,x,z) = 1.d0/2.d0 * (rho_outlet - (f_2D(0,x,z) + f_2D(1,x,z) + f_2D(2,x,z) 
                             + f_2D(3,x,z) + f_2D(4,x,z) + f_2D(6,x,z) + f_2D(8,x,z)))
        else if (z==nz_lb) then                            ! upper wall
            rho_outlet = f_2D(0,x,z) + f_2D(1,x,z) + f_2D(2,x,z) + f_2D(3,x,z) &
                       + f_2D(4,x,z) + f_2D(5,x,z) + f_2D(6,x,z) + f_2D(7,x,z) + f_2D(8,x,z)
                                                             
            f_2D(3,x,z) = f_2D(1,x,z)
            f_2D(7,x,z) = f_2D(5,x,z)
            f_2D(4,x,z) = f_2D(2,x,z)
            f_2D(6,x,z) = 1.d0/2.d0 * (rho_outlet - (f_2D(0,x,z) + f_2D(1,x,z) + f_2D(2,x,z) 
                             + f_2D(3,x,z) + f_2D(4,x,z) + f_2D(5,x,z) + f_2D(7,x,z)))
            f_2D(8,x,z) = 1.d0/2.d0 * (rho_outlet - (f_2D(0,x,z) + f_2D(1,x,z) + f_2D(2,x,z) 
                             + f_2D(3,x,z) + f_2D(4,x,z) + f_2D(5,x,z) + f_2D(7,x,z)))
        end if
    end do

Does somebody have an idea about the reason for the increase of velocity? The easiest way would be to apply symmetric boundary conditions and the trouble is gone, but for latter calculations it could be a difficult choice.

it is difficult to read the code and find the error. but most matlab example codes on lbmethod.org implement zou/he. could you use them to debug your boundary condition?