[Ding xiansheng]

According to the udf compiled by Johnson&Jackson wall boundary conditions, in the calculation of gas-solid two-phase flow in fluent, it always diverges. Can you show us the reason,thank you very much

#include "udf.h"
#include "sg_vof.h"
#include "sg_mphase.h"
#include "flow.h"
#include "mem.h"
#include "metric.h"
#define VOF_SF_MAX 0.63
#define SPE_COE 0.001
#define E_ss 0.9
#define E_w 0.8
#define mu_w 0.5

DEFINE_PROFILE(wall_shear_x, t, i)
{
cell_t c0;
face_t f;
Thread *t0 = THREAD_T0(t);
real wall_shear, goss;
real rho_s, vof_s, temp_s, p_sum, p_k, p_f;
real u_x, u_y, u_z, u_slip;
if (!Data_Valid_P()) return;
begin_f_loop(f, t)
{
c0 = F_C0(f, t);
u_x = C_U(c0, t0);
u_y = C_V(c0, t0);
u_z = C_W(c0, t0);
rho_s = C_R(c0, t0);
vof_s = C_VOF(c0, t0);
temp_s = C_GT(c0, t0);
p_sum = C_GP(c0, t0);
u_slip = sqrt(u_x*u_x + u_y*u_y + u_z*u_z);
goss = (1.0 - 7.0*vof_s / 16.0) / pow(1.0 - vof_s / VOF_SF_MAX, 2.0);
p_k = vof_s*rho_s*temp_s + 2.0*rho_s*(1.0 + E_ss)*vof_s*vof_s*goss*temp_s;
p_f = p_sum - p_k;
wall_shear = sqrt(3.0*temp_s)*M_PI*SPE_COE*rho_s*goss*vof_s*u_x / (6.0*VOF_SF_MAX) + mu_w*p_f*u_x / u_slip;
F_PROFILE(f, t, i) = wall_shear;
}
end_f_loop(f, t)
}

DEFINE_PROFILE(wall_shear_y, t, i)
{
cell_t c0;
face_t f;
Thread *t0 = THREAD_T0(t);
real wall_shear, goss;
real rho_s, vof_s, temp_s, p_sum, p_k, p_f;
real u_x, u_y, u_z, u_slip;
if (!Data_Valid_P()) return;
begin_f_loop(f, t)
{
c0 = F_C0(f, t);
u_x = C_U(c0, t0);
u_y = C_V(c0, t0);
u_z = C_W(c0, t0);
rho_s = C_R(c0, t0);
vof_s = C_VOF(c0, t0);
temp_s = C_GT(c0, t0);
p_sum = C_GP(c0, t0);
u_slip = sqrt(u_x*u_x + u_y*u_y + u_z*u_z);
goss = (1.0 - 7.0*vof_s / 16.0) / pow(1.0 - vof_s / VOF_SF_MAX, 2.0);
p_k = vof_s*rho_s*temp_s + 2.0*rho_s*(1.0 + E_ss)*vof_s*vof_s*goss*temp_s;
p_f = p_sum - p_k;
wall_shear = sqrt(3.0*temp_s)*M_PI*SPE_COE*rho_s*goss*vof_s*u_y / (6.0*VOF_SF_MAX) + mu_w*p_f*u_y / u_slip;
F_PROFILE(f, t, i) = wall_shear;
}
end_f_loop(f, t)
}

DEFINE_PROFILE(wall_shear_z, t, i)
{
cell_t c0;
face_t f;
Thread *t0 = THREAD_T0(t);
real wall_shear, goss;
real rho_s, vof_s, temp_s, p_sum, p_k, p_f;
real u_x, u_y, u_z, u_slip;
if (!Data_Valid_P()) return;
begin_f_loop(f, t)
{
c0 = F_C0(f, t);
u_x = C_U(c0, t0);
u_y = C_V(c0, t0);
u_z = C_W(c0, t0);
rho_s = C_R(c0, t0);
vof_s = C_VOF(c0, t0);
temp_s = C_GT(c0, t0);
p_sum = C_GP(c0, t0);
u_slip = sqrt(u_x*u_x + u_y*u_y + u_z*u_z);
goss = (1.0 - 7.0*vof_s / 16.0) / pow(1.0 - vof_s / VOF_SF_MAX, 2.0);
p_k = vof_s*rho_s*temp_s + 2.0*rho_s*(1.0 + E_ss)*vof_s*vof_s*goss*temp_s;
p_f = p_sum - p_k;
wall_shear = sqrt(3.0*temp_s)*M_PI*SPE_COE*rho_s*goss*vof_s*u_z / (6.0*VOF_SF_MAX) + mu_w*p_f*u_z / u_slip;
F_PROFILE(f, t, i) = wall_shear;
}
end_f_loop(f, t)
}

DEFINE_PROFILE(q_flux, t, i)
{
cell_t c0;
face_t f;
Thread *t0 = THREAD_T0(t);
real u_x, u_y, u_z;
real vof_s, temp_s, rho_s;
real goss, u_slip;
real q_w, q_w_1, q_w_2;
if (!Data_Valid_P()) return;
begin_f_loop(f, t)
{
c0 = F_C0(f, t);
vof_s = C_VOF(c0, t0);
rho_s = C_R(c0, t0);
temp_s = C_GT(c0, t0);
u_x = C_U(c0, t0);
u_y = C_V(c0, t0);
u_z = C_W(c0, t0);
u_slip = sqrt(u_x*u_x + u_y*u_y + u_z*u_z);
goss = (1.0 - 7.0*vof_s / 16.0) / pow(1.0 - vof_s / VOF_SF_MAX, 2.0);
q_w_1 = sqrt(3.0*temp_s)*M_PI*SPE_COE*rho_s*goss*vof_s*u_s lip*u_slip / (6.0*VOF_SF_MAX);
q_w_2 = sqrt(3.0*temp_s)*M_PI*rho_s*goss*vof_s*(1.0 - E_w*E_w)*temp_s / (4.0*VOF_SF_MAX);
q_w = q_w_1 - q_w_2;
F_PROFILE(f, t, i) = q_w;
}
end_f_loop(f, t)
}