hi as an exercize, i tried to specify the drag closure (Morsi-Alexander) explicitly through a udf, for a gas-liquid two phase flow. I am using euler-euler model and k-epsilon for turbulence.

The problem is that the simulation runs for some time and then says Invalid No: Floating Point error after say 500 time steps

I am confused as to what could be the problem. this is the code.. someone please help me out #include "udf.h"

#define diam 0.001

DEFINE_EXCHANGE_PROPERTY(custom_drag,cell,mix_thre ad,s_col,f_col) {

Thread *thread_g, *thread_l

real x_vel_g, x_vel_l, y_vel_g, y_vel_l, abs_v, slip_x, slip_y,

rho_g, rho_l, mu_l, reyp, tau, void_g, void_l, a1, a2, a3,

cd, f_drag, k_drag;

thread_g = THREAD_SUB_THREAD(mix_thread, f_col);/*gas*/

thread_l = THREAD_SUB_THREAD(mix_thread, s_col);/*liq*/

/*geting the velocities and properties*/

x_vel_g = C_U(cell, thread_g);

y_vel_g = C_V(cell, thread_g);

x_vel_l = C_U(cell, thread_l);

y_Vel_l = C_V(Cell, thread_l);

slip_x = x_vel_g - x_vel_l;

slip_y = y_vel_g - y_vel_l;

rho_g = C_R(cell, thread_g);

rho_l = C_R(cell, thread_l);

mu_l = C_MU_L(cell, thread_l);

/*slip velocity*/

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

/*Relative Reynolds No:*/

reyp = rho_l*abs_v*diam/mu_l;

/*relaxation time*/

tau = rho_g*diam*diam/18/mu_l;

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

void_l = 1.0 - void_g; /*liq vol frac*/

if(reyp < 0.1)

{

a1 = 0.0;

a2 = 24;

a3 = 0.0;

}

else if(reyp < 1.0)

{

a1 = 3.69;

a2 = 22.73;

a3 = 0.0903;

}

else if(reyp < 10.0)

{

a1 = 1.222;

a2 = 29.1667;

a3 = -3.8889;

}

else if(reyp < 100.0)

{

a1 = 0.6167;

a2 = 46.5;

a3 = -116.67;

}

else if(reyp < 1000.0)

{

a1 = 0.3644;

a2 = 98.33;

a3 = -2778.0;

}

else if(reyp < 5000.0)

{

a1 = 0.357;

a2 = 148.62;

a3 = -47500.0;

}

else if(reyp < 10000.0)

{

a1 = 0.46;

a2 = -490.546;

a3 = 578700.0;

}

else

{

a1 = 0.5191;

a2 = -1662.5;

a3 = 5416700.0;

}

/* now calculating cd*/

cd = a1 + a2/reyp + a3/reyp/reyp;

/* now calculating f*/

f_drag = cd*reyp/24.0;

/*now calculate k*/

k_drag = void_l*void_g*rho_l*f_drag/tau;

return k_drag; /*return drag coeff*/ }

thanx in advance

hi it may be that you using a integer one place and float all the other places. it's here: if(reyp < 0.1)

{

a1 = 0.0;

a2 = 24;

a3 = 0.0; try using a2 = 24.; i had the same problem myself -lurk