CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   Fluent UDF and Scheme Programming (https://www.cfd-online.com/Forums/fluent-udf/)
-   -   Gidaspow udf (https://www.cfd-online.com/Forums/fluent-udf/116318-gidaspow-udf.html)

ravindra April 17, 2013 01:26

Gidaspow udf
 
hey hi friends,
i have a little problem but i am unable to figure it out
i am using this gidaspow drag model for my fluidised bed simulation but i am getting different results when i am incorporating this drag law through inbuilt gidaspow law and by writing a UDF for it.
i have seen from fluent manual that the equations which i have used are same in the udf as given in manual.

Can anybody please tell me what may be the reason behind this?

Thank you in advance.:)

ghost82 April 17, 2013 03:02

Maybe it could be of great help if you attach your udf.

Daniele

ravindra April 17, 2013 11:18

this is my udf for gidaspow model
#include "udf.h"

#define diam2 875.e-7


DEFINE_EXCHANGE_PROPERTY(custom_drag,cell,mix_thre ad,s_col,f_col)
{
Thread *thread_g, *thread_s;
real x_vel_g, x_vel_s, y_vel_g, y_vel_s, abs_v, slip_x, slip_y,rho_g, rho_s, mu_g, reyp;
real void_g,void_s,fdrgs,k_g_s;

/* find the threads for the gas (primary) */
/* and solids (secondary phases) */

thread_g = THREAD_SUB_THREAD(mix_thread, s_col);/* gas phase */
thread_s = THREAD_SUB_THREAD(mix_thread, f_col);/* solid phase*/

/* find phase velocities and properties*/

x_vel_g = C_U(cell, thread_g);
y_vel_g = C_V(cell, thread_g);

x_vel_s = C_U(cell, thread_s);
y_vel_s = C_V(cell, thread_s);

slip_x = x_vel_g - x_vel_s;
slip_y = y_vel_g - y_vel_s;

rho_g = C_R(cell, thread_g);
rho_s = C_R(cell, thread_s);

mu_g = C_MU_L(cell, thread_g);

/*compute slip*/
abs_v = sqrt(slip_x*slip_x + slip_y*slip_y);

/*compute Reynold's number*/
reyp = rho_g*abs_v*diam2/mu_g;

void_g = C_VOF(cell, thread_g);/* gas vol frac*/
void_s=1.-void_g;

if(void_g>0.8)
{
if(reyp<=1000)
{
fdrgs=(24/(void_g*reyp))*(1+0.15*(pow((void_g*reyp),0.687))) ;
}
else
{
fdrgs=0.44;
}

k_g_s=(3/4)*(fdrgs*void_s*void_g*rho_g*abs_v)*pow(void_g,-2.65)/diam2;
}
else
{
k_g_s=150*void_s*void_s*mu_g/(void_g*diam2*diam2)+1.75*void_s*rho_g*abs_v/diam2;
}
return k_g_s;
}

ghost82 April 18, 2013 04:31

Hi,
I read the fluent user guide; when you compute fdrgs you are discriminating For the particle Reynolds number, but as you can see here (https://www.sharcnet.ca/Software/Flu...mp-cd-gidaspow) the Gidaspow model calculates fdrgs with one formula:

fdrgs=(24/(void_g*reyp))*(1+0.15*(pow((void_g*reyp),0.687)))

Accordingly to the user guide your calculations of Cd refer to the "symmetric" or to the "Shiller and Naumann" models.
Maybe this is the "error" in your udf.

Daniele

ghost82 April 18, 2013 04:41

In addition I think it is better to add a fabs:

Code:

abs_v = sqrt(fabs(slip_x*slip_x + slip_y*slip_y));

ravindra April 19, 2013 01:26

hey Daniele,
i have also tried in which i didnot discriminate for reynolds number...but still the results are different for cases with udf and without udf.

i will do what you have suggested for using fabs() function and get back to you...

thanks

ravindra April 25, 2013 05:07

hey daniele,
i am using this udf but results are still different with those obtained without udf.
can you suggest what is wrong.

#include "udf.h"

#define diam2 875.e-7


DEFINE_EXCHANGE_PROPERTY(custom_drag,cell,mix_thre ad,s_col,f_col)
{
Thread *thread_g, *thread_s;
real x_vel_g, x_vel_s, y_vel_g, y_vel_s, abs_v, slip_x, slip_y,rho_g, rho_s, mu_g, reyp;
real void_g,void_s,fdrgs,k_g_s;

/* find the threads for the gas (primary) */
/* and solids (secondary phases) */

thread_g = THREAD_SUB_THREAD(mix_thread, s_col);/* gas phase */
thread_s = THREAD_SUB_THREAD(mix_thread, f_col);/* solid phase*/

/* find phase velocities and properties*/

x_vel_g = C_U(cell, thread_g);
y_vel_g = C_V(cell, thread_g);

x_vel_s = C_U(cell, thread_s);
y_vel_s = C_V(cell, thread_s);

slip_x = x_vel_g - x_vel_s;
slip_y = y_vel_g - y_vel_s;

rho_g = C_R(cell, thread_g);
rho_s = C_R(cell, thread_s);

mu_g = C_MU_L(cell, thread_g);

/*compute slip*/
abs_v = sqrt(fabs(slip_x*slip_x + slip_y*slip_y));

/*compute Reynold's number*/
reyp = rho_g*abs_v*diam2/mu_g;

void_g = C_VOF(cell, thread_g);/* gas vol frac*/
void_s=1.-void_g;

if(void_g>0.8)
{

fdrgs=(24/(void_g*reyp))*(1+0.15*(pow((void_g*reyp),0.687))) ;

k_g_s=(3/4)*(fdrgs*void_s*void_g*rho_g*abs_v)*pow(void_g,-2.65)/diam2;
}
else
{
k_g_s=150*void_s*void_s*mu_g/(void_g*diam2*diam2)+1.75*void_s*rho_g*abs_v/diam2;
}
return k_g_s;
}

mohsen0488 October 14, 2013 03:27

Hi
i want to write a udf about "gidaspow drag model with switch function" can you help me please


thanks a lot in advance.

Musa October 15, 2013 06:19

@ Mohsen it's very easy to do that. Bear in mind that there are several stitching functions have been proposed for the Gidaspow model with the most commonly used being that due to Hulin and Gidaspow. The model is actually available in ANSYS FLUENT 14. Alternatively you can browse the source codes for the drag models in MFIX code and get an idea how to code it for FLUENT. If you still have trouble or doubts then feel free to contact me and I will you the specific details.

mohsen0488 October 15, 2013 06:38

Thanks a lot Mr Musa for your reply.
I will try and if I had a problem I ask you.

Thanks again.

Musa October 15, 2013 06:53

@ Mohsen you are welcome. Look at this paper Leboreiro, J., Joseph, G. G., & Hrenya, C. M. Revisiting the standard drag law for bubbling, gas-fluidized beds. Powder Technol 2008, 183, 385-400. for the discussion regarding the different stitching functions to address the discontinuity of the Gidaspow model. For browsing the MFIX code you have to be a member to have access. It's easy and free to join. Membership will usually take a day or two to be activated uopn application. Though the source code is in FORTRAN it can be easily inferred for use in C.

mohsen0488 October 15, 2013 07:14

Ok.Thanks.

mohsen0488 October 17, 2013 11:48

Dear Musa
I read article that said you, and I wrote the code as follows:
But my results are different from other authors , and i don't know what's the problem can you help me please? Thanks

#include "udf.h"
#define pi 4.*atan(1.)
#define diam2 1.41e-3
DEFINE_EXCHANGE_PROPERTY(custom_drag, cell, mix_thread, s_col, f_col)
{
Thread *thread_g, *thread_s;
real x_vel_g, x_vel_s, y_vel_g, y_vel_s, abs_v, slip_x, slip_y,
rho_g, rho_s, mu_g, reyp, cd,
void_g,void_s, vfac, fi_gs, k_g_s, k_g_s_eg, k_g_s_wy;
/* find the threads for the gas (primary) */
/* and solids (secondary phases) */
thread_g = THREAD_SUB_THREAD(mix_thread, s_col);/* gas phase */
thread_s = THREAD_SUB_THREAD(mix_thread, f_col);/* solid phase*/
/* find phase velocities and properties*/
x_vel_g = C_U(cell, thread_g);
y_vel_g = C_V(cell, thread_g);
x_vel_s = C_U(cell, thread_s);
y_vel_s = C_V(cell, thread_s);
slip_x = x_vel_g - x_vel_s;
slip_y = y_vel_g - y_vel_s;
rho_g = C_R(cell, thread_g);
rho_s = C_R(cell, thread_s);
mu_g = C_MU_L(cell, thread_g);
/*compute slip*/
abs_v = sqrt(slip_x*slip_x + slip_y*slip_y);
void_g = C_VOF(cell, thread_g);/* gas vol frac*/
void_s = C_VOF(cell, thread_s);/* solid vol frac*/
/*compute reynolds number*/
reyp = rho_g*abs_v*diam2/mu_g;
/*compute Cd*/
if(reyp<=1000)
cd=24*(1+0.15*pow(void_g*reyp,0.687))/void_g*reyp;
else
cd=0.44;
fi_gs=atan(150*1.75*(void_g-0.8))/pi+0.5;
/*compute drag and return drag coeff, k_g_s*/
k_g_s_eg=(150.*(pow(void_s, 2)*mu_g)/(void_g*pow(diam2, 2)))+(1.75*void_s*rho_g*abs_v)/diam2;
k_g_s_wy=((3./4.)*cd*void_s*void_g*rho_g*abs_v*pow(void_g, -2.65))/diam2;
k_g_s=fi_gs*k_g_s_wy+(1-fi_gs)*k_g_s_eg;
return k_g_s;
}

Musa October 17, 2013 19:21

Hi Mohsen you forgot the void fraction condition that switches from Ergun to Wen and Yu. Here is my code below. Also notice how I define the particle diameter inside the code to avoid changing the particle diameter every time I have a new particle size.

#include "udf.h"
#define pi 4.*atan(1.)


DEFINE_EXCHANGE_PROPERTY(gidaspow_blend,cell,mix_t hread,s_col,f_col)
{
Thread *thread_g, *thread_s;
real x_vel_g, x_vel_s, y_vel_g, y_vel_s, abs_v, slip_x, slip_y,
rho_g, rho_s, mu_g, reyp, Res, afac, Cd,
vof_g, vof_s, k_g_s,Ergun,WenYu,phi_gs, diam2;

/* find the threads for the gas (primary) */
/* and solids (secondary phases) */

thread_g = THREAD_SUB_THREAD(mix_thread, s_col);/* gas phase */
thread_s = THREAD_SUB_THREAD(mix_thread, f_col);/* solid phase*/

/* find phase velocities and properties*/

x_vel_g = C_U(cell, thread_g);
y_vel_g = C_V(cell, thread_g);

x_vel_s = C_U(cell, thread_s);
y_vel_s = C_V(cell, thread_s);

slip_x = x_vel_g - x_vel_s;
slip_y = y_vel_g - y_vel_s;

rho_g = C_R(cell, thread_g);
rho_s = C_R(cell, thread_s);

mu_g = C_MU_L(cell, thread_g);

diam2 = C_PHASE_DIAMETER(cell, thread_s);

/*compute slip*/
abs_v = sqrt(slip_x*slip_x + slip_y*slip_y);

/*compute Reynold's number*/

reyp = rho_g*abs_v*diam2/mu_g;

/* compute particle relaxation time */

vof_g = C_VOF(cell, thread_g);/* gas vol frac*/
vof_s = C_VOF(cell, thread_s);/* particle vol frac*/

Res = vof_g * reyp;

if(Res>=1000)
Cd=0.44;
else
Cd=24/Res*(1+0.15*pow(Res,0.687));


if(vof_g <=0.8){
/* Ergun Drag Model */
Ergun = 150*vof_s*(1-vof_g)/vof_g*mu_g/diam2/diam2 + 1.75*vof_s*rho_g/diam2*abs_v;
}


else
{
/* Wen & Yu Model */
WenYu = 3/4*Cd*vof_g*vof_s*rho_g*abs_v/diam2*pow(vof_g,-2.65);
}

phi_gs = atan(150*1.75*(0.2-vof_s))/pi +0.5;

k_g_s = (1-phi_gs)*Ergun + phi_gs*WenYu;

return k_g_s;
}

Musa October 17, 2013 19:23

Try my code and let me know if you get the desired results. All the best.

mohsen0488 October 18, 2013 06:47

Dear Musa
Thank you so much
I tried with your code , but still the results are different from other authors
Do you think the code has still a problem?

Musa October 18, 2013 07:52

Hi Mohsen the code does not have a problem. I have coded it as given in several references. I also made reference to the MFIX source code of this drag model. Results might be different due to other reasons perhaps. You can go through the authors complete CFD code and see if matches with yours. Ensure that the boundary conditions and other parameters are the same. Check your Grid also. If you still don't get the same results then perhaps you can contact the authors of the work you are trying to mimic and ask them for the specifics. Good luck!

mohsen0488 October 18, 2013 10:52

hi musa
this paper " CFD modeling of the gas–particle flow behavior in spouted beds , Wu Zhonghua ⁎, Arun S. Mujumdar, " used this code as below: i tried by this code and my results are same.but i think it has a problem in solid reynolds definition do you think too?
code as follows:
#include "udf.h"
#define pi 4.*atan(1.)
#define diam2 1.41e-3
DEFINE_EXCHANGE_PROPERTY(custom_drag, cell, mix_thread, s_col, f_col)
{
Thread *thread_g, *thread_s;
real x_vel_g, x_vel_s, y_vel_g, y_vel_s, abs_v, slip_x, slip_y,
rho_g, rho_s, mu_g, reyp, cd,
void_g, vfac, fi_gs, k_g_s, k_g_s_eg, k_g_s_wy;
/* find the threads for the gas (primary) */
/* and solids (secondary phases) */
thread_g = THREAD_SUB_THREAD(mix_thread, s_col);/* gas phase */
thread_s = THREAD_SUB_THREAD(mix_thread, f_col);/* solid phase*/
/* find phase velocities and properties*/
x_vel_g = C_U(cell, thread_g);
y_vel_g = C_V(cell, thread_g);
x_vel_s = C_U(cell, thread_s);
y_vel_s = C_V(cell, thread_s);
slip_x = x_vel_g - x_vel_s;
slip_y = y_vel_g - y_vel_s;
rho_g = C_R(cell, thread_g);
rho_s = C_R(cell, thread_s);
mu_g = C_MU_L(cell, thread_g);
/*compute slip*/
abs_v = sqrt(slip_x*slip_x + slip_y*slip_y);
void_g = C_VOF(cell, thread_g);/* gas vol frac*/
/*compute reynolds number*/
reyp = rho_g*abs_v*diam2/mu_g;
/*compute Cd*/
if(reyp<=1000)
cd=24*(1+0.15*pow(reyp,0.687))/reyp;
else
cd=0.44;
fi_gs=atan(150*1.75*(void_g-0.8))/pi+0.5;
/*compute drag and return drag coeff, k_g_s*/
k_g_s_eg=(150*(1-void_g)+1.75*reyp)*(1-void_g)*mu_g/void_g/diam2/diam2;
k_g_s_wy=0.75*cd*(1-void_g)*mu_g*reyp*pow(void_g,-2.65)/diam2/diam2;
k_g_s=fi_gs*k_g_s_wy+(1-fi_gs)*k_g_s_eg;
return k_g_s;
}

Musa October 18, 2013 17:46

Mohsen you are absolutely right the reyp has to be multiplied by void_g in the equation containing Cd but this has not been done in the code although it has been shown in the paper . Secondly there is no condition for switching from one regime to the other i.e from Ergun for the dense bed (void_g > 0.8) to WenYu for ( void_g<= 0.8). The stitching function simply ensures smooth transition from one regime to the other meaning that the condition of regime change must still be given.

mohsen0488 October 18, 2013 18:20

Musa I tried with that code and the results were similar results of paper
so do you think results of this paper are incorrect?

Musa October 18, 2013 22:02

@Mohsen certainly the drag model has not been properly coded meaning the results obtained in the paper are questionable. Did the authors of the paper use the same code? If so bring it to their attention that the code has some errors! i.e. the definition of the particle reynolds number is not correct and the condition for switching to different regimes is missing.

mohsen0488 October 19, 2013 06:48

I read thesis of the paper and in appendix of the thesis that code was written for gidaspow with switch function.and i think that the code has some errors.


thanks a lot for your replies Dear Musa ,You've helped me a lot.

I wish the best for you.

Musa October 19, 2013 07:09

You are welcome buddy!

temo October 24, 2013 15:18

aoa, i have used this udf for gas solid fluidized bed and gives accurate results now i want to change from gas to liquid(water), tell me what are the changes required to do in code for this, waiting for ur reply
udf for customizing syamlal drag law in Fluent

#include "udf.h"
#define pi 4.*atan(1.)
#define diam2 3.e-4
DEFINE_EXCHANGE_PROPERTY(custom_drag, cell, mix_thread, s_col, f_col)
{
Thread *thread_g, *thread_s;
real x_vel_g, x_vel_s, y_vel_g, y_vel_s, abs_v, slip_x, slip_y,
rho_g, rho_s, mu_g, reyp, afac,
bfac, void_g, vfac, fdrgs, taup, k_g_s;
/* find the threads for the gas (primary) */
/* and solids (secondary phases) */
thread_g = THREAD_SUB_THREAD(mix_thread, s_col);/* gas phase */
thread_s = THREAD_SUB_THREAD(mix_thread, f_col);/* solid phase*/
/* find phase velocities and properties*/
x_vel_g = C_U(cell, thread_g);
y_vel_g = C_V(cell, thread_g);
x_vel_s = C_U(cell, thread_s);
y_vel_s = C_V(cell, thread_s);
slip_x = x_vel_g - x_vel_s;
slip_y = y_vel_g - y_vel_s;
rho_g = C_R(cell, thread_g);
rho_s = C_R(cell, thread_s);
mu_g = C_MU_L(cell, thread_g);
/*compute slip*/
abs_v = sqrt(slip_x*slip_x + slip_y*slip_y);
/*compute Reynold's number*/
reyp = rho_g*abs_v*diam2/mu_g;
/* compute particle relaxation time */
taup = rho_s*diam2*diam2/18./mu_g;
void_g = C_VOF(cell, thread_g);/* gas vol frac*/
/*compute drag and return drag coeff, k_g_s*/
afac = pow(void_g,4.14);
if(void_g<=0.85)
bfac = 0.281632*pow(void_g, 1.28);
else
bfac = pow(void_g, 9.076960);
vfac = 0.5*(afac-0.06*reyp+sqrt(0.0036*reyp*reyp+0.12*reyp*(2.*bfac-afac)+
afac*afac));
fdrgs = void_g*(pow((0.63*sqrt(reyp)/
vfac+4.8*sqrt(vfac)/vfac),2))/24.0;
k_g_s = (1.-void_g)*rho_s*fdrgs/taup;
return k_g_s;
}

temo October 24, 2013 15:20

i have used this udf for gas solid fluidized bed and gives accurate results now i want to change from gas to liquid(water), tell me what are the changes required to do in code for this, waiting for ur reply
udf for customizing syamlal drag law in Fluent

#include "udf.h"
#define pi 4.*atan(1.)
#define diam2 3.e-4
DEFINE_EXCHANGE_PROPERTY(custom_drag, cell, mix_thread, s_col, f_col)
{
Thread *thread_g, *thread_s;
real x_vel_g, x_vel_s, y_vel_g, y_vel_s, abs_v, slip_x, slip_y,
rho_g, rho_s, mu_g, reyp, afac,
bfac, void_g, vfac, fdrgs, taup, k_g_s;
/* find the threads for the gas (primary) */
/* and solids (secondary phases) */
thread_g = THREAD_SUB_THREAD(mix_thread, s_col);/* gas phase */
thread_s = THREAD_SUB_THREAD(mix_thread, f_col);/* solid phase*/
/* find phase velocities and properties*/
x_vel_g = C_U(cell, thread_g);
y_vel_g = C_V(cell, thread_g);
x_vel_s = C_U(cell, thread_s);
y_vel_s = C_V(cell, thread_s);
slip_x = x_vel_g - x_vel_s;
slip_y = y_vel_g - y_vel_s;
rho_g = C_R(cell, thread_g);
rho_s = C_R(cell, thread_s);
mu_g = C_MU_L(cell, thread_g);
/*compute slip*/
abs_v = sqrt(slip_x*slip_x + slip_y*slip_y);
/*compute Reynold's number*/
reyp = rho_g*abs_v*diam2/mu_g;
/* compute particle relaxation time */
taup = rho_s*diam2*diam2/18./mu_g;
void_g = C_VOF(cell, thread_g);/* gas vol frac*/
/*compute drag and return drag coeff, k_g_s*/
afac = pow(void_g,4.14);
if(void_g<=0.85)
bfac = 0.281632*pow(void_g, 1.28);
else
bfac = pow(void_g, 9.076960);
vfac = 0.5*(afac-0.06*reyp+sqrt(0.0036*reyp*reyp+0.12*reyp*(2.*bfac-afac)+
afac*afac));
fdrgs = void_g*(pow((0.63*sqrt(reyp)/
vfac+4.8*sqrt(vfac)/vfac),2))/24.0;
k_g_s = (1.-void_g)*rho_s*fdrgs/taup;
return k_g_s;
}

Musa October 24, 2013 20:48

@ Temour buddy sorry but I do not have any experience with gas-liquid fluidized beds!

eddy123 July 19, 2016 08:42

Quote:

Originally Posted by mohsen0488 (Post 457741)
Musa I tried with that code and the results were similar results of paper
so do you think results of this paper are incorrect?

How to write the UDF for DPM drag. I tried to write the UDF for DPM-Gidaspow drag model but it did not work properly. Any one can help me where I am wrong. Any help will be appreciated. Below is my UDF

#include "udf.h"
#include "dpm.h"
#define ZERO 0.0000001

DEFINE_DPM_DRAG(dpm_gidaspow_drag,Re,p)
{
Thread *thread_g, *thread_s, *mix_thread;
real x_vel_g, x_vel_s, y_vel_g, y_vel_s, abs_v, slip_x, slip_y,
rho_g, rho_s, mu_g, dp, reyp, reyp1,void_g, void_s, cd, D;
cell_t cell;
real left_term;
real dragcoefficient;

/*find the cell index and thread of the cell that the particle is currently in*/
cell=P_CELL_THREAD(p);
mix_thread=THREAD_SUB_THREAD(P_CELL_THREAD(p),0);

/* find the threads for the gas (primary) */

thread_g = THREAD_SUB_THREAD(mix_thread,0);/* gas phase thread */
thread_s = THREAD_SUB_THREAD(mix_thread,1);/* solid phase thread */
/* find phase velocities and properties*/

x_vel_g = C_U(cell, thread_g);
y_vel_g = C_V(cell, thread_g);

x_vel_s=P_VEL(p)[0];
y_vel_s=P_VEL(p)[1];

slip_x = x_vel_g - x_vel_s;
slip_y = y_vel_g - y_vel_s;

rho_g = C_R(cell, thread_g); /*density*/
rho_s =P_RHO(p);

mu_g = C_MU_L(cell, thread_g); /*laminar viscosity*/

dp =P_DIAM(p);

/*compute slip*/
abs_v = sqrt(slip_x*slip_x + slip_y*slip_y); /*absolute value of slip velocity*/

void_g = C_VOF(cell, thread_g);/* gas vol frac*/

Message ("Voidage of gas phase %f \n", void_g);

void_s = 1-void_g;/*particle vol frac*/

/*compute reynolds number*/
reyp =rho_g*abs_v*dp/mu_g; /*no volume fraction???*/
reyp1 = void_g*reyp;


if(reyp1<1000)
cd=24/reyp1*(1+0.15*pow(reyp1,0.687));
else
cd = 0.44;
if (void_g <=0.8)
{
/* Ergun Drag Model */
D = 150.0*void_s*void_s*mu_g/void_g/dp/dp+1.75*void_s*rho_g*abs_v/dp;
}
else
{
/* Wen & Yu Model */
D = 0.75*cd*void_g*void_s*rho_g*abs_v/dp/pow(void_g,2.65);
}

left_term=rho_s*pow(dp,2)/mu_g;
dragcoefficient=left_term*D/(abs_v+ZERO);
return dragcoefficient;

}

The_seeker December 13, 2021 04:43

Hi All,

I am trying to model 2D fluidized bed with laminar flow conditions and validate it with a paper. Following is the Gidaspow model that I wrote in UDF. Looking at the discussion in this thread, I feel that the code that I have written is correct but I am not able to get good convergence for continuity. I would appreciate if someone can guide me on this.


DEFINE_EXCHANGE_PROPERTY (custom_drag, cell, mix_thread,s_col,f_col)

{
Thread *thread_g, *thread_s;
real x_vel_g, x_vel_s, y_vel_g, y_vel_s, abs_v, slip_x, slip_y, rho_g, mu_g, reyp, void_g, k_g_s, cd;

/*find the threads for the gas (primary) */
/*and solid (secondary phases) */

thread_g =THREAD_SUB_THREAD (mix_thread, s_col); /*gas phase */
thread_s=THREAD_SUB_THREAD (mix_thread, f_col);/*solid phase*/

/*find phase velocities*/

x_vel_g = C_U(cell, thread_g);
y_vel_g = C_V(cell, thread_g);

x_vel_s = C_U(cell, thread_s);
y_vel_s = C_V(cell, thread_s);

slip_x = x_vel_g-x_vel_s;
slip_y = y_vel_g-y_vel_s;

rho_g = C_R (cell, thread_g);
mu_g = C_MU_L (cell, thread_g);
void_g = C_VOF (cell, thread_g); /*gas vol fraction*/

/*compute slip*/
abs_v=sqrt(slip_x*slip_x+slip_y*slip_y);


/*compute Reynolds number*/
reyp=rho_g*abs_v*diam2/mu_g;

/*compute drag coefficient*/
if (reyp <= 1000)
cd = (24/(void_g*reyp)*(1+0.15*pow((void_g*reyp), 0.687)));
else
cd = 0.44;
if (void_g <= 0.8)
k_g_s = ((150*(1-void_g)*(1-void_g)*mu_g)/(void_g*diam2*diam2))+((1.75*rho_g*abs_v)/(diam2));
else
k_g_s = (3/4)*cd*((void_g*(1-void_g)*rho_g*abs_v*pow(void_g,-2.65))/diam2);
return k_g_s;
}

AlexanderZ December 13, 2021 22:45

if UDF works without problems, your convergence issue is related to other factors (most likely):
1. discretization scheme
2. time step resolution
3. mesh resolution
4. boundary conditions

Amin Gh February 17, 2024 08:16

asking some questions
 
Quote:

Originally Posted by Musa (Post 457739)
Mohsen you are absolutely right the reyp has to be multiplied by void_g in the equation containing Cd but this has not been done in the code although it has been shown in the paper . Secondly there is no condition for switching from one regime to the other i.e from Ergun for the dense bed (void_g > 0.8) to WenYu for ( void_g<= 0.8). The stitching function simply ensures smooth transition from one regime to the other meaning that the condition of regime change must still be given.

Hi Musa,
I'm trying to write the UDF of gidaspow drag model then I want to write my own drag model UDF buy I don't know how to do it would you please help me?


All times are GMT -4. The time now is 00:53.