Learning Forum

[Haotian Wu]

I'm working on an evaporating law for my DPM simulation, but I'm currently stuck on how to call fuel mass fraction on droplet surface with macros. I did not find it in User Manual. Can somebody please help me? 

[Rob]

There isn't a droplet surface as such: it's all "droplet". The droplet mass fraction is what you're after, so look in the DPM macros. 

[Haotian Wu]

Hi Rob,

Thanks for replying. The droplet mass fraction is unity since the droplet has only one species. Let me explain the problem more clearly. The model I want to use included the Spalding mass transfer number which is calculated by the mass fraction of droplet species in the gas-vapor mixture surrounding the droplet. Is there a macro for that value?

Regards,

Haotian

[Rob]

Yes, it's the species mass fraction, https://ansyshelp.ansys.com/account/Secured?returnurl=/Views/Secured/corp/v231/en/flu_udf/flu_udf_DataAccessMacros.html%23x1-4610004.2.3   Click on Help in the solver & paste the link into that browser. 

The harder part is getting the cell that the particle is in to do that maths. That's covered in the DPM part, https://ansyshelp.ansys.com/account/Secured?returnurl=/Views/Secured/corp/v231/en/flu_udf/flu_udf_DPMDEFINE.html  

[Haotian Wu]

I have found the macro but the UDF code was causing a bad termination.

Could you please give some instructions on the code?

#include "udf.h"

DEFINE_DPM_LAW(Boiling_Law, p, coupled)
{
    Material* m = P_MATERIAL(p);
    Thread* t0 = P_CELL_THREAD(p);
    cell_t c = P_CELL(p);
    const cphase_state_t* c0 = &(p->cphase[0]);

    int i = P_EVAP_SPECIES_INDEX(p);

    real Tp = P_T(p);
   
    real P_sat = DPM_VAPOR_PRESSURE(p, m, Tp);

    real T_ambi = C_T(c, t0);
    real Tb = DPM_BOILING_TEMPERATURE(p, m);//SATURATED_TEMPERATURE
    real T_spe = Tb + (T_ambi - Tb) / 3;

    real rho = C_R(c, t0) + (C_R(c, t0) - DPM_RHO(p, m, Tb)) / 3;

    real Ap = M_PI * P_DIAM0(p) * P_DIAM0(p);//DPM_AREA(Dp)
    real LTb = DPM_LATENT_HEAT(p, m);
    real lamda = C_K_L(c, t0); //thermal conductivity of gas
    real cpp = DPM_SPECIFIC_HEAT(p, T_spe);
    real cpg = C_CP(c, t0); /*specific heat of gas*/
    real r0 = 0.5 * P_DIAM0(p);
    real h_infi = T_ambi * C_CP(c, t0);//enthaphy of gas
    real h_b = cpp * Tp;//enthaphy of particle
    real visg = DPM_Liquid_Property(p, m, T_spe); /*turbulent viscosity of gas*/

    real B_T = 2.0, B_M, F_M, F_P;//spalding_heat_transfer_number;
    real epsilon = DPM_EMISSIVITY(p, m);
    real sigma = 5.670e-8;//STEFAN-BOLTZMANN_CONSTANT

    real vel = sqrt(pow((P_VEL(p)[0] - C_U(c, t0)), 2.0) + pow((P_VEL(p)[1] - C_V(c, t0)), 2.0) + pow((P_VEL(p)[2] - C_W(c, t0)), 2.0));
    real Re = ((vel * P_DIAM0(p) * C_R(c, t0)) / visg);//Re num
    real Pr = (cpp * visg) / lamda;//Prandtl num Pr = (c0->sHeat)*(c0->mu)/(c0->tCond);
    real Sc = visg / (P_RHO(p) * DPM_DIFFUSION_COEFF(p, T_spe));//Sc = (c0->mu)/((c0->rho)*(diffusion_coeff));
    real phi_radi, alphaf;
    real Nu, Sh;
    real m_flash, m_heat, m_radi;
    real mp_dot;

    /*compute m_flash*/
    if (0 <= (P_T(p) - Tb) && (P_T(p) - Tb) <= 5) {
        alphaf = 760 * pow((P_T(p) - Tb), 0.26);
    }
    else if ((P_T(p) - Tb) > 5 && (P_T(p) - Tb) < 25) {
        alphaf = 27 * pow((P_T(p) - Tb), 2.33);
    }
    else if ((P_T(p) - Tb) >= 25) {
        alphaf = 13800 * pow((P_T(p) - Tb), 0.39);
    }

    m_flash = (alphaf * Ap * (P_T(p) - Tb)) / LTb;

    /* compute surface fuel mass fraction */
    real Y_fuel_surf = C_DPMS_SURF_YI(c, t0, i);// P_sat/C_P(c,t0)
    real Y_fuel_ambi = C_YI(c, t0, i);
    /* compute initial values of Sherwood and Nusselt numbers using Frossling correlations */
    Sh = 2.0 + 0.552 * pow(Re, 0.5) * pow(Sc, 1.0 / 3.0);
    Nu = 2.0 + 0.552 * pow(Re, 0.5) * pow(Pr, 1.0 / 3.0);

    /* calculate B_M and F_M */
    B_M = (Y_fuel_surf - Y_fuel_ambi) / (1.0 - Y_fuel_surf);
    F_M = pow(1.0 + B_M, 0.7) * log(1.0 + B_M) / B_M;

    real Le = (lamda / (rho * cpp * DPM_DIFFUSION_COEFF(p, T_spe)));

    /*spalding heat transfer number*/
    real F_T, B_T_prev = 0.0, phi;
    while (fabs(B_T - B_T_prev) > 1e-6) {
        /* update the previous value of B_T */
        B_T_prev = B_T;

        /* calculate the current values of Sh and Nu */
        F_T = pow(1.0 + B_T, 0.7) * log(1.0 + B_T) / B_T;
        Sh = 2.0 + (Sh - 2.0) / F_M;
        Nu = 2.0 + (Nu - 2.0) / F_T;

        /* calculate phi and update B_T */
        phi = (cpp / cpg) * (Sh / Nu) / Le;
        B_T = pow(1.0 + B_M, phi) - 1.0;
    }

    real m_heat_old = 0.0;
    real remaining_mass = 0.0;

    BoilingLaw(p);
    real dt = P_DT(p);
    while (fabs(m_heat - m_heat_old) > 1e-4) {
        m_heat_old = (P_MASS0(p) - P_MASS(p)) / dt;
        m_heat = 2 * M_PI * lamda / cpp * r0 * Nu / (1 + m_flash / m_heat_old)
            * log(1 + (1 + m_flash / m_heat_old) * (h_infi - h_b) / LTb);
        m_heat_old = m_heat;
    }

    phi_radi = epsilon * sigma * Ap * (pow(T_ambi, 4.0) - pow(Tb, 4.0));
    m_radi = phi_radi / LTb;

    mp_dot = m_heat + m_flash + m_radi;
    if (mp_dot < 0.)
        mp_dot = 0.0;

      UpdateTemperature(p, -mp_dot * LTb, P_MASS(p), 1.);
    remaining_mass = P_MASS(p) - mp_dot * dt;
    /* if the particle is very small, or the particle temperature
    exceeds the boiling temperature dump the remaining mass */
    if ((remaining_mass < TP_DPM_MINIMUM_MASS(p)) || (P_DIAM(p) <= DPM_LOWEST_DIAMETER))
    {
        P_DIAM(p) = 0.0;
        P_T(p) = Tp;
        P_MASS(p) = 0.;
    }
    else
    {
        P_MASS(p) -= mp_dot * dt;
        P_DIAM(p) = DPM_DIAM_FROM_VOL(P_MASS(p) / P_RHO(p));
        P_T(p) -= LTb / P_MASS(p) * DPM_SPECIFIC_HEAT(p, P_T(p));
    }
}

[Rob]

Other than recommending indents and comments to identify & record what each section is for, no. If you are seeing a compilation failure read the error carefully as they're fairly useful. 

[Haotian Wu]

The problem is not during compiling, it was after the calculation started. The only message was:

BAD TERMINATION OF ONE OF YOUR APPLICATION PROCESSES

RANK 5 PID 47124 
EXIT STATUS: -1

and: BAD TERMINATION OF ONE OF YOUR APPLICATION PROCESSES
RANK 5 PID 112408 
EXIT STATUS: -1

After these messages, the fluent was terminated. 

[Rob]

Debug on a smaller model. Then turn bits off until it works. Turn bits on until it doesn't. That's pretty much how I debug UDFs in the event I can't get a colleague to write them for me. 

[Haotian Wu]

Thanks for the tips.