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/)
-   -   not working DEFINE_DPM_SOURCE (https://www.cfd-online.com/Forums/fluent-udf/231464-not-working-define_dpm_source.html)

loren November 6, 2020 06:53

not working DEFINE_DPM_SOURCE
 
Hi! I need your help, please:)
I want to model CO2 gas bubble dissolution in water. The VOF model is on to study the water-air inface; bubbles are modeled with dpm inert type particle; water is a mixture of water and dissolved CO2. To compute gas dissolution from the bubble and to add the dissolved gas to the mixture, I'm using DEEFINE_DPM_SOURCE macro. The source term so defined: S->species [co2_index] +=mp_dot;
I can compile and load the UDF but when I run calculations the following message appears:

DPM Iteration ....

================================================== ============================

Node 1: Process 3188: Received signal SIGSEGV.

================================================== ============================
The fl process could not be started.

Here is my UDF:

#include "dpm.h"
#include "dpm_types.h"
#include "dpm_laws.h"
#include "udf.h"

DEFINE_DPM_SOURCE(dissolution, c, t, S, strenght, tp)
{

real temp, P, Dco2, Sc, d, ki, kh, Ci, Cinf, massrate,mp_dot, Re, rel_vel;
int co2_index = 0;
temp = 288;
d = TP_DIAM(tp);
P = C_P(c, t);
rel_vel = sqrt(pow((TP_VEL(tp)[0] - C_U(c,t)),2.0) + pow((TP_VEL(tp)[1] - C_V(c,t)), 2.0) + pow((TP_VEL(tp)[2] - C_W(c,t)), 2.0));
Re = TP_RHO(tp) * rel_vel * d/0.001003;
Dco2 = (1.38*pow(10.0,-23.0)*temp)/(6.0*M_PI*0.001003*pow(10.0,-10.0));
Sc=0.001003/(TP_RHO(tp)*Dco2);
if (d<0.0035)
{
ki=(2+0.6*pow(Re,1./2.)*pow(Sc,1./3.))*Dco2/d;
}
if(d<=0.0045 && d>0.0035)
{
ki=0.2838*d-0.0009121;
}
if (d>=0.0035){
ki=2/pow(M_PI, 1./2.)*pow(Re*Sc,1./2.)*Dco2/d;
}
kh=44*1000/101325*pow(2.7183,2119.05/temp-17.4141);
Ci=P*kh;
Cinf=C_YI(c, t, co2_index);
massrate=M_PI*d*d*ki*(Ci-Cinf);
mp_dot=massrate*TP_DT(tp)*strenght;
S->species[co2_index] += mp_dot;
}

loren November 8, 2020 04:17

Help me please, it's for my thesis and I'm stuck on this UDF...

pakk November 8, 2020 07:57

Simplify your UDF in small steps until the error disappears. Find which part is the problem.

1.38*pow(10.0,-23.0) = 1.38e-23
0.001003*pow(10.0,-10.0) = 0.001003e-10
pow(2.7183,2119.05/temp-17.4141) = exp(2119.05/temp-17.4141)

1000/101325 will become zero, be careful!

loren November 8, 2020 09:02

Thankyou Pakk, I will try

loren November 9, 2020 09:23

Ok I think that I found the problem. If I don't use "Cinf=C_YI(c, t, co2_index);" and I use "Cinf=0;" Fluent doesn't crash, and gives some reasonable results.

Anyway the reason why I'm using TP_DT(tp) is because the mass rate is in kg/s and at each time step I'm computing the mass added to the continuous phase. I assume that the UDF is solved at each time step, so the mass rate computed in the time step is multiplied by the time step extension, is it correct?

Despite I get some results I see a particle diameter and mass reduction lower than expected and I don't understand why.

Here is the code:



#include "dpm.h"

#include "dpm_types.h"

#include "dpm_laws.h"

#include "udf.h"



DEFINE_DPM_SOURCE(dissolution, c, t, S, strenght, tp)

{



real temp, P, Dco2, Sc, d, ki, kh, Csol, Cinf, massrate,mp_dot, Re, rel_vel;

temp = 288;

d = TP_DIAM(tp);

P = C_P(c, t);

rel_vel = sqrt(pow((TP_VEL(tp)[0] - C_U(c,t)),2.0) + pow((TP_VEL(tp)[1] - C_V(c,t)), 2.0) + pow((TP_VEL(tp)[2] - C_W(c,t)), 2.0));

Re = (TP_RHO(tp) * rel_vel * d)/0.001003;

Dco2 = (1.38e-23*temp)/(1.8868e-12);

Sc=0.001003/(TP_RHO(tp)*Dco2);

if (d<0.0035)

{

ki=(2+0.6*pow(Re,1./2.)*pow(Sc,1./3.))*Dco2/d;

}

if(d<=0.0045 && d>0.0035)

{

ki=0.2838*d-0.0009121;

}

if (d>=0.0035){

ki=2/pow(3.14, 1./2.)*pow(Re*Sc,1./2.)*Dco2/d;

}

kh=(0.4342)*exp(2119.05/temp-17.4141);

Csol=P*kh;

Cinf=0;

massrate=3.14*d*d*ki*(Csol-Cinf);

TP_MASS(tp) = TP_MASS(tp) - massrate * TP_DT(tp);

mp_dot=massrate*TP_DT(tp)*strenght;

S->species[0] += mp_dot;

}


All times are GMT -4. The time now is 12:51.