A simple error. In addition, how to calculate the absolutely permeability of the porous medium?

[xiangWu-WW]

[xiangWu-WW]

Another one,why there are some negative numbers in the velosity,every direction has negative numbers. Thanks for your code and reply!

[yjhp1016]

Hi Xiangwu, thanks for your help on identifying the bug.

To calculate permeability is not difficult, it is based on darcy's law, as we are working on a discrete space, there are several ways to do it. I can share with you what I did before, I calculate total flux Q_total (basically sum up all velocities), then divide them by volume V_total (q = Q_total/V_total), if you have pressure difference dp, then you need to use dp/Lx or dp/Ly or dp/Lz to calculate the pressure gradient, if you're using bodyforce to drive the flow, then you can simply use fx, or fy, or fz for pressure gradient.

You can calculate permeability K = q u (viscosity)/(fx) for body force case, or K = q u(viscosity)/(dp/Lx) for pressure drop case.

I don't understand your second question very well, what do you mean there are some negative numbers in velocity? Could you be more specific? Thanks.

On Wed, 5 Jul 2023 at 04:11, xiangWu-WW @.***> wrote:

Another one,why there are some negative numbers in the velosity,every direction has negative numbers. Thanks for your code and reply!

— Reply to this email directly, view it on GitHub https://github.com/yjhp1016/taichi_LBM3D/issues/15#issuecomment-1620953864, or unsubscribe https://github.com/notifications/unsubscribe-auth/AJEDKQBAPRGPYMBPA3ASMCDXOTLPVANCNFSM6AAAAAAZ6K3ZFM . You are receiving this because you are subscribed to this thread.Message ID: @.***>

[xiangWu-WW]

Thanks for your timely reply. As you told, i calculate the absolutely permeability of a Plane model which width is 10px. Theoretically，the result is 1010(pixel resolution)/12=8.33 um2.（Assuming resolution is 1um/px）.But the result is 0.002573, i wonder how to contact the two different results. Further more, as you can see, in the Python printed,the velocity has negative value,but in the simulating result it is not clear, what does the negative values mean?

![Uploading PlaneModel.png…]()

[yjhp1016]

Hi XiangWu,

I think there are several things you can check, first of all, v is an "array" of dimension of 4, so if you want to show the min of x direction , it should be something like np.min(vel[:, :, :, 0]).

And we calculate K in three directions, so we should use something like np.sum(vel[:, :, :, 0]) in the permeability calculation, and you need to use corresponding bodyforce/pressure gradient of same direction

Last thing is more about LBM boundary conditions, several personal comments:

• LBM pressure boundary condition is reliable when the boundary rho is around 1, I saw that you are using rho=0.1 for one of the boundaries, it will make the simulation not very stable, my personal experience is that I only use rho between 0.95-1.05.
• Big pressure gradient need longer time to reach equilibrium,
• I recommend you only use uniform bodyforce + periodic boundary condition (default) to drive the flow, this is the most accurate and stable solution

Thanks

Jianhui

On Thu, 6 Jul 2023 at 06:13, xiangWu-WW @.***> wrote:

Thanks for your timely reply. As you told, i calculate the absolutely permeability of a Plane model which width is 10px. Theoretically，the result is 1010(pixel resolution)/12=8.33 um2.（Assuming resolution is 1um/px）.But the result is 0.002573, i wonder how to contact the two different results. Further more, as you can see, in the Python printed,the velocity has negative value,but in the simulating result it is not clear, what does the negative values mean? [image: PlaneModel] https://user-images.githubusercontent.com/83948137/251349662-35ac64c7-adae-4def-920d-4680f3e44a2f.png [image: codes20230706124853] https://user-images.githubusercontent.com/83948137/251349611-23b512ad-2b0c-4a9f-b97c-bc4ff809368e.jpg [image: Uploading PlaneModel.png…] [image: result1] https://user-images.githubusercontent.com/83948137/251349709-5e1aa049-14ff-4b9d-aa7d-bf35cd3dbb14.png [image: result2] https://user-images.githubusercontent.com/83948137/251349747-eec7001e-90d8-486f-b143-b8b8635cc81f.png

— Reply to this email directly, view it on GitHub https://github.com/yjhp1016/taichi_LBM3D/issues/15#issuecomment-1623000896, or unsubscribe https://github.com/notifications/unsubscribe-auth/AJEDKQHNO2QZDTDQSSXVJPLXOZCOBANCNFSM6AAAAAAZ6K3ZFM . You are receiving this because you commented.Message ID: @.***>

[xiangWu-WW]

Hi Dr.Yang

Thank for your response！ The BC is: lb3d.set_bc_rho_z1(0.9) lb3d.set_bc_rho_z0(1.0)

Then i change the formula to calculate the K(K = q u(viscosity)/(dp/Lx)): vel = self.v.to_numpy() K = np.sum(vel[:,:,:,2])/(self.nxself.nyself.nz)self.niu/((3/8(1-0.9))/self.nz)#rho_in=1,rho_out=0.9,P=3/8rho The result converge at 1.122(Platemodel, width=10px）。What is the unit of the result, 1 Darcy? I wonder is there any problem in the formula.

Further more, where the rho has a berak in the begin and end: slice1 1 slice2 0.960676 … slice-2 0.93958 slice-1 0.9

Finally,as you mentioned,add body force to force the fluid,how to change the code,is it simplely change: self.fx,self.fy,self.fz = 0.0e-6,0.0,0.0 i tried this,but the result hardly convergent.

Thanks！

[yjhp1016]

Hi xiangWu,

1.0 and 0.9 looks too big difference for me. your mesh resolution can't handle it. Please try 1.0 and 0.99 if you want to use pressure BC

the unit for your case if your resolution is 1um/dx, then calculated K unit should be um^2, because you are calculating K in lattice unit, but the dimensionality of permeability is L^2 so to convert LBM K to real physics K is K_lbm/(resolution*resolution) Please refer to this documents for more information about unit conversion in LBM: https://www.biofm-research.com/wp-content/uploads/2021/07/Krueger_Edmonton_scaling.pdf

When you use body force, have you remove these two lines? lb3d.set_bc_rho_z1(0.9) lb3d.set_bc_rho_z0(1.0)

On Thu, 6 Jul 2023 at 10:58, xiangWu-WW @.***> wrote:

Hi Dr.Yang

Thank for your response！ The BC is: lb3d.set_bc_rho_z1(0.9) lb3d.set_bc_rho_z0(1.0)

Then i change the formula to calculate the K(K = q u(viscosity)/(dp/Lx)): vel = self.v.to_numpy() K = np.sum(vel[:,:,:,2])/(self.nxself.nyself.nz)self.niu/((3/8 (1-0.9))/self.nz)#rho_in=1,rho_out=0.9,P=3/8rho The result converge at 1.122(Platemodel, width=10px）。What is the unit of the result, 1 Darcy? I wonder is there any problem in the formula.

Further more, where the rho has a berak in the begin and end: slice1 1 slice2 0.960676 … slice-2 0.93958 slice-1 0.9

Finally,as you mentioned,add body force to force the fluid,how to change the code,is it simplely change: self.fx,self.fy,self.fz = 0.0e-6,0.0,0.0 i tried this,but the result hardly convergent.

Thanks！

— Reply to this email directly, view it on GitHub https://github.com/yjhp1016/taichi_LBM3D/issues/15#issuecomment-1623385674, or unsubscribe https://github.com/notifications/unsubscribe-auth/AJEDKQDTOEHWI3KPTJM6DRDXO2D3TANCNFSM6AAAAAAZ6K3ZFM . You are receiving this because you commented.Message ID: @.***>

[xiangWu-WW]

Hi Dr.Yang,

I set the body force as you told me.The model is a 128^3 cubic and there is a plate(Aperture = 20, assert the pixel resolution is 1um/px ) in the model.Theoretically，the absolutely permeability is 400/12=33.3um^2. The formula to calculate the K is qu/f resolution^2. But the result is 1.69um^2， i wonder the q=Qsum/Vtotal or q = Qsum(Pore) / Vtotal(Pore), If it is the latter，the result will be 1.69/(20/128)=10.8um^2. Additionaly, why the rho distribution is so disordered.

[yjhp1016]

Hi XiangWu,

It looks to me that the calculation of absolute premability of your model is independent of your model size? ... Which looks not quite right for me

I think if you want to check if the model is doing the right thing, the best thing is to check the velocity field, in your model, it is a poiseuille flow, you can easily find analytical solutions for this problem, then you can compare to LB simulations. If you can confirm the velocity calculation is correct, there is no way the permeability calculation is wrong...

You can refer to this book for more information about poiseuille flow: Lattice Boltzmann Modeling: An Introduction for Geoscientists and Engineers by Michael C. Sukop, Daniel T. Thorne. I also give you some snapshot for reference: [image: image.png] [image: image.png]

The last question about the rho, if we use body force, we normally don't change the rho too much, it means we would actually have a very uniform rho over the computing domain. If you check your visualization scale bar, you can find rho actually doesn't change too much, the difference of rho over space is more or less same as machine precision level, so you got very rough visualization, but actually it's a uniform field.

On Sat, 22 Jul 2023 at 08:51, xiangWu-WW @.***> wrote:

Hi Dr.Yang,

I set the body force as you told me.The model is a 128^3 cubic and there is a plate(Aperture = 20, assert the pixel resolution is 1um/px ) in the model.Theoretically，the absolutely permeability is 400/12=33.3um^2. The formula to calculate the K is qu/f resolution^2. But the result is 1.69um^2， i wonder the q=Qsum/Vtotal or q = Qsum(Pore) / Vtotal(Pore), If it is the latter，the result will be 1.69/(20/128)=10.8um^2. Additionaly, why the rho distribution is so disordered. [image: _020230722154949] https://user-images.githubusercontent.com/83948137/255318404-d5706b25-e204-4fa5-b31d-1f158ced7c90.png [image: 2023-07-22_153703] https://user-images.githubusercontent.com/83948137/255318409-26d1f174-7df3-4fc3-97bb-0e1bffe13d33.png [image: 20230722153629] https://user-images.githubusercontent.com/83948137/255318411-a06ae74f-b7dc-46b8-b469-e04712987938.png [image: 20230722153553] https://user-images.githubusercontent.com/83948137/255318425-54e54074-c858-4c0e-a439-dab6afa6d2a0.png

— Reply to this email directly, view it on GitHub https://github.com/yjhp1016/taichi_LBM3D/issues/15#issuecomment-1646519214, or unsubscribe https://github.com/notifications/unsubscribe-auth/AJEDKQCHD3BT4IO6BLGXZ3LXROA75ANCNFSM6AAAAAAZ6K3ZFM . You are receiving this because you commented.Message ID: @.***>

[hangqqq]

I would like to ask how I should calculate the pressure gradient?

[xiangWu-WW]

I would like to ask how I should calculate the pressure gradient?

ΔP = cs^2 * (Δρ). ΔP represents the pressure difference. c_s is the speed of sound, typically a constant in the LBM model. cs = 3^(1/2)/2 Δρ represents the density difference between the inlet and outlet. ∇P = (ΔP) / L ∇P represents the pressure gradient. ΔP is the previously calculated pressure difference.

[hangqqq]

The pressure gradient is only related to density and length. Doesn't it have anything to do with the stamina or initial speed I set?

[yjhp1016]

Like other methods, if you use pressure boundary condition for a boundary, you can't apply velocity, the velocity will be calculated by the programme. The same for velocity condition, if you use velocity condition for a surface, the pressure will be calculated. The value in the input won't do anything then.

That's why you need to select which BC you want (for all 6 out surfaces), 0 = periodic, 1=pressure, 2=velocity

On Tue, 5 Dec 2023 at 02:10, hangqqq @.***> wrote:

The pressure gradient is only related to density and length. Doesn't it have anything to do with the stamina or initial speed I set? LNHFU2DRW_K056AL.XLFD.png (view on web) https://github.com/yjhp1016/taichi_LBM3D/assets/140054215/28699c07-d1b8-4c45-94e3-991d921d5a9f HH6J.C2Y2CWJ.6A.1U.WI.png (view on web) https://github.com/yjhp1016/taichi_LBM3D/assets/140054215/d1151535-e4e1-4d65-9890-a227ae406791

— Reply to this email directly, view it on GitHub https://github.com/yjhp1016/taichi_LBM3D/issues/15#issuecomment-1839879829, or unsubscribe https://github.com/notifications/unsubscribe-auth/AJEDKQE6IGWOSSUUH2PFZV3YHZ7CHAVCNFSM6AAAAAAZ6K3ZFOVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMYTQMZZHA3TSOBSHE . You are receiving this because you commented.Message ID: @.***>

[hangqqq]

Why is there such a sudden change in value?

[yjhp1016]

It means the program didn't converge, the boundary condition and flow properties are set too tough for this solver.

LBM use fix dx and dt, so need.to be careful to keep the simulation stable.

hangqqq @.***> 于 2024年1月11日周四 上午2:03写道：

image.png (view on web) https://github.com/yjhp1016/taichi_LBM3D/assets/140054215/e31cca7e-ea66-4b06-851a-ba97f03c6f5f image.png (view on web) https://github.com/yjhp1016/taichi_LBM3D/assets/140054215/ec253709-ea0b-487d-b274-151e22d262d1 Why is there such a sudden change in value?

— Reply to this email directly, view it on GitHub https://github.com/yjhp1016/taichi_LBM3D/issues/15#issuecomment-1886079028, or unsubscribe https://github.com/notifications/unsubscribe-auth/AJEDKQANXAPDVPGZTEAXZ5DYN5CAXAVCNFSM6AAAAAAZ6K3ZFOVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMYTQOBWGA3TSMBSHA . You are receiving this because you commented.Message ID: @.***>

[hangqqq]

Hi Dr.Yang,"I am using the LBM2phase program and I would like to know how to calculate the relative permeability of two phases?"

[yjhp1016]

Personally I don't think calculating rel perm in 2D is a quite good idea, as in 3D there will be some corner structures in porous medium, the wetting phase will occupy the corners but non-wetting phase still in center, so you will see wetting phase perm is still connected by these small corner channels even at low wetting phase saturation. But in 2D, the wetting phase need to build a connecting flow path between inlet and outlet, to have a wetting phase perm bigger than 0. So seems not very realistic compared to real scenario. But you can still try, using similar approach we cal 3D relative permeability. (I explain the method in another post, please refer to that post)

On Fri, 17 May 2024 at 08:03, hangqqq @.***> wrote:

Hi Dr.Yang,"I am using the LBM2phase program and I would like to know how to calculate the relative permeability of two phases?"

— Reply to this email directly, view it on GitHub https://github.com/yjhp1016/taichi_LBM3D/issues/15#issuecomment-2116887802, or unsubscribe https://github.com/notifications/unsubscribe-auth/AJEDKQBGE3ACSJ2JQYUCPZ3ZCWTSHAVCNFSM6AAAAAAZ6K3ZFOVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDCMJWHA4DOOBQGI . You are receiving this because you commented.Message ID: @.***>

[hangqqq]

Hi Dr.Yang， "I am using the LBM-3d-2phase program and giving different initial phase distributions with one initial invading phase occupying 10% and another occupying 25%. However, by the time I output after 800000 iterations, the invading phase occupies almost 90% in both cases. Is this normal?"

[yjhp1016]

It's possible, maybe your porous medium is good at connectivity and is water wet, the invading phase is easy to enter... But you need to check visualizations to see if the whole invading process makes sense

On Thu, 23 May 2024 at 09:20, hangqqq @.***> wrote:

Hi Dr.Yang， "I am using the LBM-3d-2phase program and giving different initial phase distributions with one initial invading phase occupying 10% and another occupying 25%. However, by the time I output after 800000 iterations, the invading phase occupies almost 90% in both cases. Is this normal?"

— Reply to this email directly, view it on GitHub https://github.com/yjhp1016/taichi_LBM3D/issues/15#issuecomment-2126512737, or unsubscribe https://github.com/notifications/unsubscribe-auth/AJEDKQEZAUZ67VJ4523YX33ZDWRFVAVCNFSM6AAAAAAZ6K3ZFOVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDCMRWGUYTENZTG4 . You are receiving this because you commented.Message ID: @.***>