When I run the GridRefinement3D example in the new version of the Palabos, several output files are generated. Seems they are sparse output option that is added to Palabos!

I could not open these separated files in Paraview for post processing!

Any helps for collecting these files as a one would be appreciated.

Cheers,

King]]>

I know the import export process is not perfected, but maybe there is a fix to this. when i load an LBM, it will work till i save and restart blender, when i restart blender, open the saved file, then try to render, i recieve this error

ERROR: Export aborted: 'bool' object has no attribute 'luxrender_texture'

its saying ""Failed to find texture "478ec31794c3cc473dfb6" in scene "" i dunno if that has anything to do with it

i should add that i delete that texture because it prohibits the rendering process and py lux from working.

Please help.

I didn't find the right solution from the Internet.

References:-

[www.luxrender.net]

Telematics Solution Video

Thanks!]]>

Thanks for your answer]]>

the general way to convert units in LBM to the physical units is using a dimensionless number like Re, Gr, Ra and etc.

the calculation you mentioned doesn't convert units "(u/Ucavity)_physical=(u/Ucavity)_LB)", it normalizes the velocity in both units.

for the temperature, you can use a dimensionless number like Ra or Gr. but if you know the high and the low value of temperature in the physical unit follow the below calculation:

T_phy = T_LBM * Diff_T_phy + T_low_phy.

for example:

T_low_phy = 20

T_high_phy = 100

Diff_T_phy = 100 - 20 = 80

T_low_LBM = 0

T_high_LBM = 1

T_LBM = 0.5

T_phy = 0.5 * 80 + 20 = 60]]>

For velocity, in a lid driven cavity problem for instance, we use (u/Ucavity)_physical=(u/Ucavity)_LB. Where Ucavity_LB is obtained by Re and nu_LB independent from its physical value. But no such ability is available for temperature to be computed.

I have not found an answer despite so many unit conversion discussions, except natural convection case in which we can use Ra number that has deltaT term in it.]]>

Even the dimensions do not match and that’s weird to me !

We have “length^2” vs. “second” (for dt) to be canceled and nu(Phys) dimension becomes equals to that of nu(LB).]]>

> Hi,

>

> umm, at least for the case of Rayleigh-Benard

>

> g*alpha*deltaT*L^3/nu*k=gr_lu*alpha_lu*deltaT_lu*(ly-1)^3/(nu_lu*k_lu)

>

> and then (maybe)

>

> deltaT_lu=

>

>

> where gr_lu will depend on delta_t and delta_x of

> alpha_lu=chosen constant.

> (ly-1)= number of lattice spaces in the vertical

> nu_lu= c2s*(tauNS-0.5) % kinematic viscosity in

> k_lu = c2s*(tauThermal-0.5) % thermal

>

> I hope it helps.. and that what I wrote make

>

> Andreay

Hello,

I’m a New user of LBM. my case is lid driven cavity flow (not affected by natural convection) with thermal boundary conditions.

Now that I can’t use Ra as my connection between LB and real parameters how can I input my temperatures into LB code?

I would be very thankful for any help .]]>

Additionally, the board delay is

The primary clock (80MHz/12.5ns) of FPGA is distributed from the internal clock divider and is used to generate the clock for the ADC (f_adc_clk=40MHz).

How do I specify setup and hold times for the data port of the FPGA with respect to ADC clock? I made several attempts without success. The spreadsheet/Timing Preferences view that i used to specify INPUT_SETUP is shown in the image below

However, the problem is that I can only use clk80 as a reference which is the output of the clock divider, and there is a considerable phase shift between this clock and the clock driving the AD converter as a result of delay from the clock divider to the output pin of the FPGA. How do I take this delay into account when specifying preferences?]]>

this is a really good question. The answer is not so trivial and longer than you think. I am assuming that T=time, M=mass, L=length. With this in mind, I end up with

[G] = L^3 * M^-1

which is the inverse of the units of the density. And this actually fits both equations: the interaction strength equation and the pressure equation.

First important thing is that the shan chen force is actually a volume force, thus, its units are Newton/meter^3 or in your notation M^1 * L^-2 * T^-2. Now, the force equation reads

F = -G \psi sum_i w_i psi(x + e_i) e_i

where the sum is a finite difference approximation for grad(psi). The problem now is that you can find two sets of weighting coefficients w_i in the literature. Take the following discretisation as an example:

e_i = (1,0) (0,1) (-1,0) (0,-1) (1,1) (-1,1) (-1,-1) (1,-1)

then usually the corresponding weights are

1/9, 1/9, 1/9, 1/9, 1/36, 1/36, 1/36, 1/36

However, these actually contain the squared lattice speed of sound! In the derivation of the finite difference stencils, the weights are usually normalized by the squared speed of sound, i.e.,

1/3, 1/3, 1/3, 1/3, 1/12, 1/12, 1/12, 1/12

(For the derivation, see e.g. Shan 2006 Phys Rev E 73, 047701 (2006) or Sbragaglia et al. Phys Rev E 75, 026702 (2007).

Depending on the weights you use, the weighting coefficients actually would have units.

Normalized -> w_i is dimensionless

Not normalized -> w_i has dimensions L^2*T^-2

Whats the consequence of this? If you use NORMALIZED weights, your finite difference formula will approximate

F = -G psi grad(psi)

On the other hand, if you use the non-normalized weights (1/9, 1/9, 1/9, 1/9, 1/36, 1/36, 1/36, 1/36), you will approximate

F = -G psi c_s^2 grad(psi)

This formula is actually used to derive the pressure equation in your post (p = c_s^2 rho + c_s^2*G*0.5*psi^2). If I now use the force equation with the squared speed of sound in it to obtain the units of G, I will end up with my result. Using the pressure equation, I will end up at the same result. The problem is that (depending on the weights) there will be a c_s^2 appearing in front of G in the pressure or not. This leads to the inconsistency.

I hope this is understandable.

Regards,

Knut]]>

I am in the process of converting an LBM code from simulation to physical units. I know the basics of the procedure, but I run into some difficulties when working with the Shan-Chen models for multiphase/multicomponent flows. In order to transform the

This reads: p = c_2^2*rho + c_s^2*G/2*psi(rho)^2. However, when I set G to achieve dimensional consistency, its units differ from those needed for dimensional consistency of the interaction force density. For simplicity, following Martys and Chen 1996, I set psi(rho) = rho. This gives units of G as L^5*M^{-1}*T^{-2}. Now, if I use the same effective mass (psi(rho)=rho) in the interaction strength equation, which reads F = -G*psi*sum_a(w_i*psi*e_i), G has units L^3*M^{-1}*T^{-1} ! In the force equation w_i is dimensionless and e_i has units of speed.

Has anyone had any similar problems? Any advice will be greatly appreciated.

Thank you,

Peter]]>

I am trying to understand static Smagorinsky model implemented in Palabos. As described in the documentation:

I could understand this from the relationship between the stress tensor and strain rate: PiQuoteIn the Palabos implementation, the strain-rate is computed from the stress tensor \Pi.

However, when I look at the file

.Language: C++static T computeOmega(T omega0, T preFactor, T rhoBar, Array<T,SymmetricTensor<T,Descriptor>::n> const& PiNeq) { T PiNeqNormSqr = SymmetricTensor<T,Descriptor>::tensorNormSqr(PiNeq); T PiNeqNorm = std::sqrt(PiNeqNormSqr); T alpha = preFactor * Descriptor<T>::invRho(rhoBar); T linearTerm = alpha*PiNeqNorm; T squareTerm = (T)2*alpha*alpha*PiNeqNormSqr; // In the following formula, the square-root appearing in the explicit form of // omega is developed to second-order. return omega0*(1-linearTerm+squareTerm); }

This calculation seems to be a different form compared to the one in the paper by Hou, et al. (1994)

Could anyone help explain this or provide any reference?

Thank you in advance.

With best,

Song

[1] Hou, et al. 1994, "A lattice Boltzmann subgrid model for high Reynolds number flows"

[2] Krafczyk 2003 "LARGE-EDDY SIMULATIONS WITH A MULTIPLE-RELAXATION-TIME LBE MODEL"]]>

Thank you very much. is this code based on the multi-block LBM?]]>

For

Re = PhyL * PhyU / PhyNu = Resolution * LatticeU / LatticeNu

Usually we use LatticeU= PhyU / M.

1. Must LatticeU be proportionable with PhyU? Can it be any value? And according to the comments of the code, it should be “proportional to Mach number”. I’ve no idea what it means.

2. In some of my cases of High Re with smagorinsky model, If I make my deltaX halved, I found that the original LatticeU can not be used in the new case for the simulating divergent. In this situation, I have to make my latticeU smaller. Why this happened?

3. I found that if the M is different in the same case (means LatticeU different), there are some differences in the results velocity. Was it right? If so, how should we determine the value of latticeU?

4. During the iteration, we use averaged rho to monitor the simulation. Usually it should be drifted around 1. But if the latticeU is not suitable, the ave rho will continually grow up, even from 1 to 4. Why will this happen? As we know that rho is the sum of the distribution functions, if it grow up, does it suit the mass conservation? And is this kind of result correct？

Maybe some of the questions are so primitive，thank you for your focus.

with best wishes,

steed188]]>

examples/codesByTopic/shanChenMultiPhase

Best regards

Orestis]]>

I have been recently trying to incorporate thermal effects to my multiphase simulations. I am still using the shan chen lbm as i am not able properly incorporate wall interaction when using other equations of state. in shan chen, i take T=1/G and use a separate mesh to model the transport of temperature as we usually do in thermal LBM. It works perfectly in single phase but doesnt work in multiphase. The simulation goes to NaN at the interfaces. I read in various papers about a heat source term which took phase change into account (involving rho, cv and delta t). But how do i implement this formula? I tried Cv=5.1 as given in one paper i saw, but it didnt work. I also tried converting Cv directly into lattice units but it still didnt solve my problem. Can anybody tell me what value to take for Cv and if it is fine to calculate delta T similar to how we calculate psi?

Thanks in Advance,

Regards,

Githin]]>

Please introduce me some research works on the area multi-block multi-phase LB (Shan Chen Method) simulation.

any help is appreciated

Thanks in advance]]>

I use same density for both fluids, which is 1. I want to get same results for two different system sizes. I mean same droplet size. based on the mentioned LB unit reference, I should change relaxation time (1 for 100x100 and 1.5 for 200x200) . However, for different system size results are completely different.

Please Help me

I am in rush :(]]>

I will be soon working on developing an unstructured LBM solver based on the paper given below:-

Finite volume TVD formulation of lattice Boltzmann simulation on unstructured mesh, Dhiraj V. Patil, K.N. Lakshmisha *

All the papers i found on unstructured LBM are on simple SCSP LBM. Are there any papers extending the unstructured LBM to multiphase simulations? If yes, it would be a great help if you could name a few. If no, is it possible to do this? Any ideas on how?

Thanks a lot in advance,

Githin]]>

I have a question about of U_ {LB} at palabos. I simulated a channel 2D (D2Q9) and 3D (D3Q27), where the initial velocity initial is zero, and I used a model with external force term. At the case 2D, my simulation woked with U_ {LB} 0.02, but at the case 3D this value does not work. I changed this value for 0.3.

The choice of these values was to see the examples of the website of Palabos.

But at this moment I still can not explain why it is necessary to make that change.

Can someone tell me the reason, or where can I find literature about this topic?

Thanks for your help]]>

My question is, when I modified resolution number's without changed my problem physical (physical velocity, Reynolds number and characteristic length) my velocity profile changed.

According to consistency of problem, this doesn't must change the solution because I would to change the resolution number, but that the effect of this change is increasing accuracy.

Please community someone can tell me why do it happen this, or where can I found literature about of this topic.

Thank you

I attache my code

******************************************************************

Language: C++[code="cpp" number] #include"palabos3D.h" #include"palabos3D.hh" #include <cstdlib> #include <vector> #include <cmath> #include <iostream> #include <sstream> #include <fstream> #include <iomanip> using namespace plb; using namespace std; using namespace plb::descriptors; typedef double T; #define DESCRIPTOR ForcedD3Q19Descriptor #define DYNAMICS GuoExternalForceBGKdynamics T channelForce(IncomprFlowParam<T> const& parameters) { return (T) 6e-6; } void channelSetup(MultiBlockLattice3D<T,DESCRIPTOR>& lattice, IncomprFlowParam<T> const& parameters, OnLatticeBoundaryCondition3D<T,DESCRIPTOR>& boundaryCondition) { const plint nx = parameters.getNx(); const plint ny = parameters.getNy(); const plint nz = parameters.getNz(); ============================ Box3D bottom(0, nx-1, 0, 0, 0, nz-1); Box3D top (0, nx-1, ny-1, ny-1, 0, nz-1); Box3D right (0, nx-1, 1, ny-2, nz-1, nz-1); Box3D left (0, nx-1, 1, ny-2, 0, 0); Box3D outlet(nx-1, nx-1, 1, ny-2, 0, nz-1); Box3D inlet (0, 0, 1, ny-2, 0, nz-1); ==================================== boundaryCondition.setVelocityConditionOnBlockBoundaries(lattice, inlet, boundary::outflow); initializeAtEquilibrium(lattice, inlet, (T)1., Array<T,3>((T)0.,(T)0.,(T)0.) ); boundaryCondition.setVelocityConditionOnBlockBoundaries(lattice, outlet, boundary::outflow); initializeAtEquilibrium(lattice, outlet, (T)1., Array<T,3>((T)0.,(T)0.,(T)0.) ); boundaryCondition.setVelocityConditionOnBlockBoundaries(lattice, bottom ); setBoundaryVelocity(lattice, bottom, Array<T,3>((T)0.,(T)0.,(T)0.)); initializeAtEquilibrium(lattice, bottom, (T)1., Array<T,3>((T)0.,(T)0.,(T)0.) ); boundaryCondition.setVelocityConditionOnBlockBoundaries(lattice, left, boundary::outflow); initializeAtEquilibrium(lattice, left, (T)1., Array<T,3>((T)0.,(T)0.,(T)0.) ); boundaryCondition.setVelocityConditionOnBlockBoundaries(lattice, right, boundary::outflow); initializeAtEquilibrium(lattice, right, (T)1., Array<T,3>((T)0.,(T)0.,(T)0.) ); boundaryCondition.setVelocityConditionOnBlockBoundaries(lattice, top ); setBoundaryVelocity(lattice, top, Array<T,3>((T)0.,(T)0., (T)0.) ); initializeAtEquilibrium(lattice, top, (T)1., Array<T,3>((T)0,(T)0., (T) 0.) ); Array<T,3> force(channelForce(parameters), (T)0.,(T)0); setExternalVector(lattice,lattice.getBoundingBox(),DESCRIPTOR<T>::ExternalField::forceBeginsAt,force); lattice.initialize(); } int main(int argc, char*argv[]){ plbInit(&argc,&argv); global::directories().setOutputDir("./tmp/"); IncomprFlowParam<T> parameters( (T) 3.350e-2, //Characteristic Velocity (S.I. units) (T) 0.095, //Lattice Boltzmann Velocity (T) 1000., //Reynolds (T) 0.3, //Characteristic Lenght (T) 48., //resolution number';s (T) 3., //Lx (T) 1., //Ly (T) 1. //Lz ); const T maxT = (T)4500; const T imSave = (T)maxT/(T)100; const T promedio = (T)1110; writeLogFile(parameters, "Parameters_Channel MRT"); plint nx = parameters.getNx(); plint ny = parameters.getNy(); plint nz = parameters.getNz(); MultiBlockLattice3D<T, DESCRIPTOR> lattice ( nx, ny, nz, new DYNAMICS<T, DESCRIPTOR>(parameters.getOmega())); lattice.periodicity().toggle(0,true); lattice.periodicity().toggle(2,true); OnLatticeBoundaryCondition3D<T,DESCRIPTOR>* boundaryCondition = createLocalBoundaryCondition3D<T,DESCRIPTOR>(); channelSetup(lattice, parameters, *boundaryCondition); Box3D Cavity(0, nx-1, 0, ny-1, 0, nz-1); plb_ofstream UPvsT("./tmp/VvsT_320.dat"); Box3D profileSectionHalf(0, nx-1, ny/2, ny/2,nz/2, nz/2); //axis X Box3D profileSectionMit(nx/2, nx/2, 0, ny - 1,nz/2, nz/2); //axis Y Box3D profileSectionFinal(nx/2, nx/2, ny/2, ny/2,0, nz - 1);//axis Z Box3D profileSectionUP1(nx/2, nx/2, ny/2, ny/2, nz/2, nz/2); Box3D profileSectionIniz(0, nx-1, ny - 1, ny - 1, 0, nz-1); Box3D profileSectionMitz(0, nx-1, ny - 1, ny - 1, 0, nz-1); Box3D profileSectionFinalz(0, nx-1, ny - 1, ny - 1, 0, nz-1); Box3D profileSectionDownz(0, nx-1, ny/4, ny/4, 0, nz-1); Box3D profileSectionHalfz(0, nx-1, ny/2, ny/2, 0, nz-1); Box3D profileSectionUpz(0, nx-1, 3*ny/4, 3*ny/4, 0, nz-1); Box3D StrengthTapUp(0,nx-1,ny-1,ny-1,0,nz-1); Box3D StrengthTapDown(0,nx-1,0,0,0,nz-1); plint counter = 0; MultiScalarField3D<T> Velocity_iter(nx,ny,nz); MultiScalarField3D<T> Velocidad_iter(nx,ny,nz); MultiScalarField3D<T> ShearUP_iter(nx,ny,nz); MultiScalarField3D<T> EsfuerzoUP_iter(nx,ny,nz); MultiScalarField3D<T> ShearDown_iter(nx,ny,nz); MultiScalarField3D<T> EsfuerzoDown_iter(nx,ny,nz); setToConstant(Velocidad_iter, Cavity, (T)0.); setToConstant(EsfuerzoUP_iter, Cavity, (T)0.); setToConstant(EsfuerzoDown_iter, Cavity, (T)0.); T unit_Nu = parameters.getLatticeNu() * parameters.getDeltaX() * parameters.getDeltaX() / parameters.getDeltaT(); for(plint iT=0; iT*parameters.getDeltaT()<1.01*maxT; ++iT){ counter += 1; double TimE = iT*parameters.getDeltaT(); Velocity_iter = *multiply (unit_velo, *computeVelocity(lattice)); Velocidad_iter = *add(Velocidad_iter, Velocity_iter); ShearUP_iter = *multiply (unit_Nu, *extractComponent( *computeStrainRateFromStress( lattice, StrengthTapUp), StrengthTapUp, 1)); EsfuerzoUP_iter = *add(ShearUP_iter,EsfuerzoUP_iter); ShearDown_iter = *multiply (unit_Nu, *extractComponent( *computeStrainRateFromStress ( lattice, StrengthTapDown), StrengthTapDown, 1)); EsfuerzoDown_iter = *add(ShearDown_iter,EsfuerzoDown_iter); if(counter%parameters.nStep(promedio)==0 && TimE > 2750){ T prom = counter; pcout << prom << endl; EsfuerzoUP_iter = *divide(EsfuerzoUP_iter, prom); std::string str=createFileName("./tmp/TapUp_Re_320_", iT, 6); const char * fileName = str.c_str(); plb_ofstream perfilUvsY1(fileName); perfilUvsY1 << setprecision(16) << EsfuerzoUP_iter << endl; //setToConstant(EsfuerzoUP_iter, Cavity, (T)0.); } if(counter%parameters.nStep(promedio)==0 && TimE > 2750){ T prom = counter; pcout << prom << endl; EsfuerzoDown_iter = *divide(EsfuerzoDown_iter,prom); std::string str=createFileName("./tmp/TapDown_Re_320_", iT, 6); const char * fileName = str.c_str(); plb_ofstream perfilUvsY2(fileName); perfilUvsY2 << setprecision(16) << EsfuerzoDown_iter << endl; //setToConstant(EsfuerzoDown_iter, Cavity, (T)0.); } if(counter%parameters.nStep(promedio)==0 && TimE > 2750){ T prom = counter; pcout << prom << endl; Velocidad_iter = *divide(*divide(Velocidad_iter,prom),parameters.getPhysicalU()); std::string str=createFileName("./tmp/ProfilesU_Re_320_", iT, 6); const char * fileName = str.c_str(); plb_ofstream perfilUvsY4(fileName); perfilUvsY4 << setprecision(16)<< Velocidad_iter << endl; //setToConstant(EsfuerzoRight_iter, Cavity, (T)0.); } UPvsT << setprecision(16) << TimE << " " << *computeVelocityNorm (lattice, profileSectionUP1) << endl; lattice.collideAndStream(); } successiveProfilesHalf << setprecision(16) << *computeVelocityNorm(lattice, profileSectionIni) << endl; successiveStrainUp << setprecision(16) << *multiply (unit_Nu, *extractComponent(*computeStrainRateFromStress(lattice, StrengthTapUp), StrengthTapUp, 1)) << endl; successiveStrainDown << setprecision(16) << *multiply (unit_Nu, *extractComponent(*computeStrainRateFromStress(lattice, StrengthTapDown), StrengthTapDown, 1)) << endl; delete boundaryCondition; }

[/code]]]>

My question is, when I modified resolution number's without changed my problem physical (physical velocity, Reynolds number and characteristic length) my velocity profile changed.

According to consistency of problem, this doesn't must change the solution because I would to change the resolution number, but that the effect of this change is increasing accuracy.

Please community someone can tell me why do it happen this, or where can I found literature about of this topic.

Thank you]]>

Did you write your own code or used opensource LBM solvers available?

Even I am trying to create an MCMP shan and chan model

Would like to discuss more how you achieved a density ratio of 100

Awaiting response

Regards

Sandeep]]>

I'm not sure if this answers you question. If not, let us know what do you want more specifically.]]>

I have started working on an LBM model with medium range repulsive force to simulate emulsions of heavy phase. I was inspired by the works of Chibbaro, 2008 (Phys Rev E 77, 036705) and Falcucci, 2010 (Soft Matter, 2010, 6, 4357–4365). For those of you who do not know what the repulsion model does, it is essentially incorporating a realistic component of repulsion between particles of the same phase, reducing the effective surface tension and promoting the formation of emulsions. This is compared to the original Shan-Chen model that only features short-range attractive forces.

My question is this. The repulsive model features two G parameters, uncoupling density ratio and surface tension. How are the values of these parameters determined? Does anyone have any experience? The expressions given in the two papers above as well as in other places differ significantly in both value and sign! Any help will be most appreciated.

Best Regards,

Peter]]>