[Juan Torres de Lima]

I need help implementing a UDF to insert the gas phase momentum source term.
I have three phases:
-Primary phase - Alkaline water
-Secondary phase - Hydrogen
-Secondary phase - Oxygen

1) Is it possible to use this term as an input parameter through "Named expressions"?
2) Any alternative or reason why the code is not even initializing (it appears to be a memory error as the multiphase model is only valid in parallel after version 2018/2).


It is a dispersion force similar to the turbulent dispersion model of Burns et al. It is calculated by:



Fd_x = - [(Fraction of liquid or gas) * (density_liquid) * (Coef_Dispersion) * (Mag_vel_relative) * (Gradient_void fraction)] / (db_diameterbubble)








#include "udf.h"

DEFINE_SOURCE(disperse_force_x, c, t, dS, eqn)
{
real F_DB_x;
real epsilon_water, epsilon_h2, epsilon_o2, density_liq, K_h2, K_o2, Vr_mag, grad_epsilon;
real db = 0.000058; // Bubble diameter [m]
#if !RP_HOST
Thread* liq_thread, * h2_thread, * o2_thread;
int phase_index = -1;
#endif

// Coefficients and problem properties
epsilon_water = C_VOF(c, t); // Water fraction
liq_thread = THREAD_SUB_THREAD(t, 0); // Liquid phase thread

h2_thread = THREAD_SUB_THREAD(t, 1); // Hydrogen phase thread
epsilon_h2 = C_VOF(c, h2_thread); // Hydrogen fraction

o2_thread = THREAD_SUB_THREAD(t, 2); // Oxygen phase thread
epsilon_o2 = C_VOF(c, o2_thread); // Oxygen fraction

density_liq = C_R(c, liq_thread); // Liquid phase density
K_h2 = 9.0 * db; // Hydrogen dispersion coefficient (K_h2/db)
K_o2 = 4.0 * db; // Oxygen dispersion coefficient (K_o2/db)

// Obtaining the phase index of the cell using THREAD_STORAGE
phase_index = *((int*)THREAD_STORAGE(t, SV_CELL));

#if !RP_HOST
// Calculate the gradient according to the cell phase
if (phase_index == 1) { // If the cell contains H2
grad_epsilon = C_VOF_G(c, h2_thread)[0]; // Gradient for hydrogen
}
else if (phase_index == 2) { // If the cell contains O2
grad_epsilon = C_VOF_G(c, o2_thread)[0]; // Gradient for oxygen
}
else { // If the cell contains water or other phases
grad_epsilon = C_VOF_G(c, liq_thread)[0]; // Gradient for water
}

// Calculate the relative velocity between gas and liquid phases
if (phase_index == 1) { // If the cell contains H2
Vr_mag = sqrt(pow(C_U(c, t), 2) + pow(C_V(c, t), 2)); // Module of the relative velocity for hydrogen
F_DB_x = -epsilon_h2 * density_liq * K_h2 * Vr_mag * grad_epsilon / db;
}
else if (phase_index == 2) { // If the cell contains O2
Vr_mag = sqrt(pow(C_U(c, t), 2) + pow(C_V(c, t), 2)); // Module of the relative velocity for oxygen
F_DB_x = -epsilon_o2 * density_liq * K_o2 * Vr_mag * grad_epsilon / db;
}
else { // If the cell contains water or other phases
Vr_mag = 0.0; // Set the relative velocity to zero for the water phase
F_DB_x = -epsilon_water * density_liq * K_h2 * Vr_mag * grad_epsilon / db;
}
#endif

// Define the momentum source term in the x-direction momentum equation
dS[eqn] = 0.0; // There are no momentum sources on the surface
if (eqn == 0) // Momentum equation in the x-direction
return F_DB_x;
else
return 0.0;
}