I want a program that calculates the drag and lift force

Lattice Boltzmann Method LBM: Algorithms
#1

I want a program that calculates the drag and lift force.
flow around a cylinder problem please help me.

1 Like
#2

Dear sofiane,

If you run a self made code, it is quite simple to determine the vertical and horizontal forces exerted on a body.

You can check the momentum exchange method (MEM) as a good starting point.
you can find many information resources about this algorithm on internet. Briefly, it is an explicit method that sums up all the distribution functions that hit the boundary with a straightforward formula.

I hope that this guides you to the right direction.

Cheers

#3

@Puigar : I read about the momentum exchange method and tried to apply it but it didn’t work. I think i missed a point . So, just as a verification:
-Should i place the force calculus after collision and before streamin?
-The distribution functions “f” appearing in the force expression are in the same node “x” (x node at the obstacle boundary) or no ?
-The expression of momentum exchange is computed at the same time step no ? Or it uses a previous value ?
Here is the article i used “https://www.researchgate.net/profile/Li-Shi_Luo/publication/230691844_Viscous_flow_computations_with_the_method_of_lattice_Boltzmann_equation/links/0fcfd50ae507934a70000000.pdf” page 336.

#4

-Should i place the force calculus after collision and before streamin?

probably this is the best way to do it properly.

-The distribution functions “f” appearing in the force expression are in the same node “x” (x node at the obstacle boundary) or no ?

it depends. You should be aware that the name of momentum exchange means that you should have an exchange (of F distributions AT the physical boundary). therefore, it will extremely depend on which algorithm you are using to define the wall, which finally defines at which exact position you have the physical wall.

if you use:

  • fullway BB: you sum up all the distributions that are going to hit the boundary wall PLUS all of the distributions that are outgoing right from the obstacle to the fluid.

  • half way BB: you then only need sum up the distributions that are bound to hit the obstacle boundary the next time step and multiply each distribution by 2, because as we define the wall at the middle of the 2 nodes (property of the halfway BB), the outgoing F distribution is exactly the same.

  • second order BB scheme: then it becomes trickier. Before streaming, you store and sum all the distributions that go to the boundary. After the streaming and BB are applied, you add the modified distributions that face the opposite direction to the stored distributions. (this one is better to do when you completely dominate the other aforementioned two BB algorithms)

-The expression of momentum exchange is computed at the same time step no ? Or it uses a previous value ?

Computationally depends on the type of bounce back.

  • fullway: you can do it in a single operation because you already have all the variables. you can do it both before and after streaming step. the normal thing however is to do it before streaming and after collision
  • halfway: you also can do it on a single operation before or after the streaming. Theoretically the MEM is calculated at t+(Dt)/2 on this case. You better do it before streaming because it is also the normal way.
  • second order BB: here there are two sum operations. one before and the other right after the streaming step.

cheers!

1 Like
#5

Thanks a lot for this detailed answer :slight_smile: I’ve been searching for long time for this but didn’t find a clear answer.
Just a little question, what do you mean by “second order bounce back”? I thought that the fullway is first order and the halfway is second order. So, the halfway is the second order type, or you mean another ( name? ) type of bounce back?

#6

for second order BB I meant the second order BB for curved geometries (Mei-Luo-Shy method, for example), which regards the physical walls in other places than the middle between two nodes.

#7

% *************************************** %
% ********** Cd/Cl Computation ********** %
% *************************************** %
% Determination of the drag and lift on the profil

function [Fx, Fy, Cd, Cl] = CdCl(f,ftemp,Nxlimg,Nylims,Nxlimd,Nylimn,Num,R,uref,rhoref,C)

global bin binC

% Computation of the forces applied to the profil
k = 1;
Fx=0;
Fy=0;
FX=zeros(Nxlimd+2,Nylimn+2);
FY=zeros(Nxlimd+2,Nylimn+2);
Balayage=0;
for i = Nxlimg:Nxlimd
for j = Nylims:Nylimn
if bin(i,j,Num)==1
if bin(i+1,j,Num)==0
FX(i,j)=FX(i,j)+(f(i+1,j,1)+ftemp(i+1,j,3));
Balayage=Balayage+1;
end
if bin(i+1,j+1,Num)==0
FX(i,j)=FX(i,j)+(f(i+1,j+1,5)+ftemp(i+1,j+1,7));
FY(i,j)=FY(i,j)+(f(i+1,j+1,5)+ftemp(i+1,j+1,7));
Balayage=Balayage+1;
end
if bin(i,j+1,Num)==0
FY(i,j)=FY(i,j)+(f(i,j+1,2)+ftemp(i,j+1,4));
Balayage=Balayage+1;
end
if bin(i-1,j+1,Num)==0
FX(i,j)=FX(i,j)-(f(i-1,j+1,6)+ftemp(i-1,j+1,8));
FY(i,j)=FY(i,j)+(f(i-1,j+1,6)+ftemp(i-1,j+1,8));
Balayage=Balayage+1;
end
if bin(i-1,j,Num)==0
FX(i,j)=FX(i,j)-(f(i-1,j,3)+ftemp(i-1,j,1));
Balayage=Balayage+1;
end
if bin(i-1,j-1,Num)==0
FX(i,j)=FX(i,j)-(f(i-1,j-1,7)+ftemp(i-1,j-1,5));
FY(i,j)=FY(i,j)-(f(i-1,j-1,7)+ftemp(i-1,j-1,5));
Balayage=Balayage+1;
end
if bin(i,j-1,Num)==0
FY(i,j)=FY(i,j)-(f(i,j-1,4)+ftemp(i,j-1,2));
Balayage=Balayage+1;
end
if bin(i+1,j-1,Num)==0
FX(i,j)=FX(i,j)+(f(i+1,j-1,8)+ftemp(i+1,j-1,6));
FY(i,j)=FY(i,j)-(f(i+1,j-1,8)+ftemp(i+1,j-1,6));
Balayage=Balayage+1;
end
end
end
end

alpha=0;

Fx=sum(sum(FX));
Fy=sum(sum(FY));
Fx=-Fx;
Fy=-Fy;

% Rotation of the referential
F=sqrt(FxFx+FyFy);
Ang1=acosd(Fx/F);
Ang2=asind(Fy/F);
Fx=Fcosd(Ang1+alpha);
Fy=Fsind(Ang2+alpha);

% --------- %
% - Cd/Cl - %
% --------- %

Cd=2abs(Fx)/uref/uref/R/rhoref;
Cl=2Fy/uref/uref/R/rhoref;

cd= drag coefficient
cl = lift coefficient

#8

@sofiane: Thanks i will test this code after modifications because i worked with fullway bounce back .

#9

I used this method bounce back.
is there possibility to modify this program.
that already sent to bounce back.

can you send me your email for many information

i have program in matlab in lattice boltzmann but is difucult

#10

Hi sofiane,
I have a code which uses the geometry exterior to generate Cd, Cl in matlab…and its pretty easy to understand with just 5 lines of code…If you are interested just reply.

#11

Hello sidneni,

I am interested in your code, could you please share it with me?

#12

I have a code which uses the geometry exterior to generate Cd, Cl in matlab…and its pretty easy to understand with just 5 lines of code…If you are interested just reply.

#13

Hi I read your explanation about drag force. it was very informative. i am trying to use Lallemand interpolation method in order to calculate drag force over cylinder Re=100. unfortunately my Cd is always higher than literature. to calculate drag force i sum up all f’s toward the boundary before streaming and add it to the inverse direction after streaming on fluid node next to the solid node. i am not sure where is my problem. i would be grateful if you help me thanks

Lallemand, Pierre, and Li-Shi Luo. “Lattice Boltzmann method for moving boundaries.” Journal of Computational Physics 184.2 (2003): 406-421.

#14

Dear sofiane,
Currently, I have some troubles in using the momentum exchange method to calculate the drag and lift force, could you
please send me some relevant codes as you said?
Thanks you in advance,
kind regards,
zhilin

#16

Is it possible to send MATLAB code related to drag coefficient to my email?