这是我的udf
#include "udf.h" /*head files*/
#include "mem.h"
#include "flow.h"
#define PI 3.141592653 /*constant*/
#define Mv 18.01534 /*molecular weight of vapor(kg/kmol)*/
#define Ma 28.966 /*molecular weight of air(kg/kmol) */
int *a;
real *b;
real *d;
real *q;
real *g;
real *l;
real s1;
real s2;
real s3;
real s0;
int counter=0; /*counter*/
int j; /*index of interface wall cell*/
int k;
real latent_heat(real T)
{
real t;
real Hlg;
t=T-273.15;
if(t<10.0)
Hlg=2477208.8;
if((t>=10.0)&&(t<20.0))
Hlg=2477208.8+(t-10.0)*(2453549.6-2477208.8)/10.0;
if((t>=20.0)&&(t<30.0))
Hlg=2453549.6+(t-20.0)*(2429882.8-2453549.6)/10.0;
if((t>=30.0)&&(t<40.0))
Hlg=2429882.8+(t-30.0)*(2406058.7-2429882.8)/10.0;
if((t>=40.0)&&(t<50.0))
Hlg=2406058.7+(t-40.0)*(2382038.2-2406058.7)/10.0;
if((t>=50.0)&&(t<60.0))
Hlg=2382038.2+(t-50.0)*(2357691.0-2382038.2)/10.0;
if((t>=60.0)&&(t<70.0))
Hlg=2382038.2+(t-60.0)*(2333080.9-2382038.2)/10.0;
if((t>=70.0)&&(t<80.0))
Hlg=2333080.9+(t-70.0)*(2308065.7-2333080.9)/10.0;
if((t>=80.0)&&(t<90.0))
Hlg=2308065.7+(t-80.0)*(2282560.3-2308065.7)/10.0;
if((t>=90.0)&&(t<100.0))
Hlg=2282560.3+(t-90.0)*(2256472.9-2282560.3)/10.0;
if((t>=100.0)&&(t<110.0))
Hlg=2256472.9+(t-100.0)*(2229704.3-2256472.9)/10.0;
if((t>=110.0)&&(t<120.0))
Hlg=2229704.3+(t-110.0)*(2202149.7-2229704.3)/10.0;
if((t>=120.0)&&(t<130.0))
Hlg=2202149.7+(t-120.0)*(2173700.0-2202149.7)/10.0;
if((t>=130.0)&&(t<140.0))
Hlg=2173700.0+(t-130.0)*(2144243.7-2173700.0)/10.0;
if((t>=140.0)&&(t<150.0))
Hlg=2144243.7+(t-140.0)*(2113667.6-2144243.7)/10.0;
if((t>=150.0)&&(t<160.0))
Hlg=2113667.6+(t-150.0)*(2081855.9-2113667.6)/10.0;
if((t>=160.0)&&(t<170.0))
Hlg=2081855.9+(t-160.0)*(2048687.3-2081855.9)/10.0;
if((t>=170.0)&&(t<180.0))
Hlg=2048687.3+(t-170.0)*(2014031.4-2048687.3)/10.0;
if((t>=180.0)&&(t<190.0))
Hlg=2014031.4+(t-180.0)*(1977745.0-2014031.4)/10.0;
if((t>=190.0)&&(t<200.0))
Hlg=1977745.0+(t-190.0)*(1939668.5-1977745.0)/10.0;
if((t>=200.0)&&(t<210.0))
Hlg=1939668.0+(t-200.0)*(1899623.0-1939668.0)/10.0;
if((t>=210.0)&&(t<220.0))
Hlg=1899623.0+(t-210.0)*(1857409.0-1899623.0)/10.0;
if((t>=220.0)&&(t<230.0))
Hlg=1857409.0+(t-220.0)*(1812799.75-1857409.0)/10.0;
if((t>=230.0)&&(t<240.0))
Hlg=1812799.75+(t-230.0)*(1765537.21-1812799.75)/10.0;
if((t>=240.0)&&(t<250.0))
Hlg=1765537.21+(t-230.0)*(1715325.26-1765537.21)/10.0;
if(t>=250.0)
Hlg=1715325.26;
return Hlg;
}
real s_enthalpy(real T)
{
real t;
real H;
t=T-273.15;
if(t<10.0)
H=2519000.0;
if((t>=10.0)&&(t<20.0))
H=2519000.0+(t-10.0)*(2538000.0-2519000.0)/10.0;
if((t>=20.0)&&(t<30.0))
H=2538000.0+(t-20.0)*(2556000.0-2538000.0)/10.0;
if((t>=30.0)&&(t<40.0))
H=2556000.0+(t-30.0)*(2574000.0-2556000.0)/10.0;
if((t>=40.0)&&(t<50.0))
H=2574000.0+(t-40.0)*(2596000.0-2574000.0)/10.0;
if((t>=50.0)&&(t<60.0))
H=2596000.0+(t-50.0)*(2609000.0-2596000.0)/10.0;
if((t>=60.0)&&(t<70.0))
H=2609000.0+(t-60.0)*(2609000.0-2609000.0)/10.0;
if((t>=70.0)&&(t<80.0))
H=2626000.0+(t-70.0)*(2643000.0-2626000.0)/10.0;
if((t>=80.0)&&(t<90.0))
H=2643000.0+(t-80.0)*(2660000.0-2643000.0)/10.0;
if((t>=90.0)&&(t<100.0))
H=2660000.0+(t-90.0)*(2676000.0-2660000.0)/10.0;
if((t>=100.0)&&(t<110.0))
H=2676000.0+(t-100.0)*(2691000.0-2676000.0)/10.0;
if((t>=110.0)&&(t<120.0))
H=2691000.0+(t-110.0)*(2706000.0-2691000.0)/10.0;
if((t>=120.0)&&(t<130.0))
H=2706000.0+(t-120.0)*(2720000.0-2706000.0)/10.0;
if((t>=130.0)&&(t<140.0))
H=2720000.0+(t-130.0)*(2734000.0-2720000.0)/10.0;
if((t>=140.0)&&(t<150.0))
H=2734000.0+(t-140.0)*(2747000.0-2734000.0)/10.0;
if((t>=150.0)&&(t<160.0))
H=2747000.0+(t-150.0)*(2758000.0-2747000.0)/10.0;
if((t>=160.0)&&(t<170.0))
H=2758000.0+(t-160.0)*(2769000.0-2758000.0)/10.0;
if((t>=170.0)&&(t<180.0))
H=2769000.0+(t-170.0)*(2778000.0-2769000.0)/10.0;
if((t>=180.0)&&(t<190.0))
H=2778000.0+(t-180.0)*(2786000.0-2778000.0)/10.0;
if((t>=190.0)&&(t<200.0))
H=2786000.0+(t-190.0)*(2793000.0-2786000.0)/10.0;
if((t>=200.0)&&(t<=210.0))
H=2793000.0+(t-200.0)*(2798000.0-2793000.0)/10.0;
if((t>=210.0)&&(t<=220.0))
H=2798000.0+(t-210.0)*(2802000.0-2798000.0)/10.0;
if(t>=220.0)
H=2802000.0;
return H;
}
DEFINE_DELTAT(my_timestep,d)
{
real time_step;
real flow_time=RP_Get_Real("flow-time" ;
if(flow_time<4.0)
time_step=0.01;
else
time_step=0.05;
return time_step;
}
DEFINE_PROFILE(cwalltemp,t,i)
{
real x[ND_ND];
real s,ma,mb,mc,md,me;
face_t f;
ma=8.284;
mb=-25.585;
mc=22.699;
md=-5.3462;
me=350.87;
begin_f_loop(f,t)
{
F_CENTROID(x,f,t);
s=x[1];
F_PROFILE(f,t,i) = ma*s*s*s*s+mb*s*s*s+mc*s*s+md*s+me;
}
end_f_loop(f,t)
}
DEFINE_PROFILE(mass_flow,t,i)
{
face_t f;
begin_f_loop(f,t)
{
F_PROFILE(f,t,i) =-s0;
}
end_f_loop(f,t)
}
DEFINE_ON_DEMAND(INIT)
{
int thread_ID;
Domain *domain;
Thread *tt;
int i;
face_t f; /*face index on the wall*/
cell_t c0;
real A[ND_ND];
i=0;
thread_ID=6;
domain=Get_Domain(1);
tt=Lookup_Thread(domain,thread_ID);
/*calculate the total number of faces on the wall*/
begin_f_loop(f,tt)
{
c0=F_C0(f,tt);
counter++;
}
end_f_loop(f,tt)
a=(int*)malloc(counter*sizeof(int));
b=(real*)malloc(counter*sizeof(real));
d=(real*)malloc(counter*sizeof(real)); /* temp info*/
q=(real*)malloc(counter*sizeof(real)); /* temp info at previous time*/
g=(real*)malloc(counter*sizeof(real)); /* mass fraction info at previous time*/
l=(real*)malloc(counter*sizeof(real)); /* pressure info at previous time*/
begin_f_loop(f,tt)
{
c0=F_C0(f,tt);
a=c0;
F_AREA(A,f,tt);
b=NV_MAG(A);
i++;
}
end_f_loop(f,tt)
}
DEFINE_ADJUST(k_init,domain)
{
k=0;
}
DEFINE_ADJUST(temp_at_interface,domain)
{
Thread *f_thread;
face_t face;
int thread_ID;
int i=0;
thread_ID=6; /*ID of g_s_wall is 6*/
f_thread=Lookup_Thread(domain,thread_ID);
begin_f_loop(face,f_thread)
{
d=F_T(face,f_thread);
i++;
}
end_f_loop(face,f_thread);
}
DEFINE_ADJUST(avg,domain)
{
int i=0;
real tavg = 0;
real wavg = 0;
real Pavg = 0;
real temp,volume,vol_tot=0,w_h2o,Pa,condmass;
real sum=0;
Domain *d;
Thread *t;
cell_t c;
real x[ND_ND];
d = Get_Domain(1);
thread_loop_c(t,d)
{
begin_c_loop(c,t)
{
C_CENTROID(x,c,t);
volume = C_VOLUME(c,t); /* get cell volume */
temp = C_T(c,t); /* get cell temperature */
w_h2o= C_YI(c,t,i); /* get cell mass fraction of h2o */
Pa=C_P(c,t)+402145.725; /* get cell Absolute Pressure */
condmass=C_UDMI(c,t,0);
if ((-0.001<x[1])&(x[1]<2.001))
{
volume = C_VOLUME(c,t);
temp = C_T(c,t);
w_h2o= C_YI(c,t,i);
Pa=C_P(c,t)+402145.725;
condmass=C_UDMI(c,t,0);
}
else
{
temp =0;
volume = 0;
w_h2o=0;
Pa=0;
condmass=0;
}
sum+=condmass;
vol_tot+=volume;
tavg+=temp*volume;
wavg+=w_h2o*volume;
Pavg+=Pa*volume;
}
end_c_loop(c,t)
tavg /= vol_tot;
wavg /= vol_tot;
Pavg /= vol_tot;
}
s0=sum;
s1=tavg;
s2=wavg;
s3=Pavg;
}
DEFINE_EXECUTE_AT_END(temp_info)
{
Domain *domain;
Thread *f_thread;
Thread *c_thread;
face_t face;
cell_t cell;
int thread_ID;
int i=0;
int mf=0;
thread_ID=6; /*ID OF g_s_wall is 6*/
domain=Get_Domain(1);
f_thread=Lookup_Thread(domain,thread_ID);
begin_f_loop(face,f_thread)
{
q=F_T(face,f_thread);
i++;
}
end_f_loop(face,f_thread)
i=0;
thread_ID=2;/*ID of gas is 2*/
c_thread=Lookup_Thread(domain,thread_ID);
begin_c_loop(cell,c_thread)
{
if(cellTest(cell))
{
g=C_YI(cell,c_thread,mf);
l=C_P(cell,c_thread);
i++;
}
}
end_c_loop(cell,c_thread);
}
DEFINE_SOURCE(vapor_m_source,c,t,dS,eqn)
{
real source=0;
int i; /*index of vapor*/
real w; /*vapor mass fraction*/
real w_delt; /*vapor mass fraction at previous time step*/
real area; /*area of wall boundary*/
real V; /*volume of neighbouring cell*/
real Tw; /*temperature of wall, K*/
real Tw_delt; /*temperature of fluid at previous time step, K*/
real T; /*temperature of fluid, K*/
real Hlg; /*latent heat of vapor(J/Kg)*/
real Cu; /*adjustable coefficient in Uchida correlation*/
real Tref; /*reference temperature*/
real Cp_mix=1500; /*spefic heat of vapor at constant pressure, J/(Kg.K)*/
real Cp_air; /*spefic heat of air at constant pressure and referenct temperature, J/(Kg.K)*/
real P; /*pressure of neighbouring cells*/
real P_delt; /*pressure of neighbouring cells at previous time step, Pa*/
real Pv; /*pressure of vapor pressure in the neighbouring cells*/
real Ps; /*saturation pressure of vapor*/
real hv; /*vapor specific enthalpy, J/Kg*/
real Vx; /*axial velocity*/
real Vy; /*radial velocity*/
Cu=0.98; /*adjustable coefficient*/
Tref=273.15;
Cp_air=1005.0;
i=0;
T=s1;
w=s2;
P=s3;
Vx=C_U(c,t);
Vy=C_V(c,t);
V=C_VOLUME(c,t);
if(cellTest(c))
{
area=b[j];
Tw=d[j];
Tw_delt=q[j];
w_delt=g[j];
w=C_YI(c,t,i);
Cp_mix=C_CP(c,t);
P_delt=l[j]+402145.725;
/*Pv=Ma*w*P/(Ma*w+Mv*(1-w));*/
Pv=Ma*w_delt*P_delt/(Ma*w_delt+Mv*(1-w_delt));
Ps=2000*exp(18.5916-3991.11/(Tw_delt-39.31))/15;
hv=s_enthalpy(Tw_delt);
if(Pv>=Ps)
{
Hlg=latent_heat(Tw);
source=-Cu*(10189.3+0.0904164*P-(4314.4+0.046537*P)*log10(100*(1-w)))*area*pow((T-Tw),0.4)/(Hlg*V);/*source=-Cu*pow(x/(1-x),0.8)*area*(T-Tw)/(Hlg*V)*/ /*Uchida correlation*/
dS[eqn]=0.;
}
else
{
source=0.;
}
}
else
{
source=0.0;
}
C_UDMI(c,t,0)=source*V;
C_UDMI(c,t,1)=source*(T*Cp_mix-Tref*Cp_air)*V; /*C_UDMI(c,t,1)=source*hv;*/
C_UDMI(c,t,2)=source*Vx*V; /*moment source in axis direction*/
C_UDMI(c,t,3)=source*Vy*V; /*moment source in xradius direction*/
return source;
DEFINE_SOURCE(vapor_species_source,c,t,dS,eqn)
{
real source=0;
real V;
V=C_VOLUME(c,t);
source=C_UDMI(c,t,0)/V;
dS[eqn]=0.0;
return source;
}
DEFINE_SOURCE(vapor_h_source,c,t,dS,eqn)
{
real source=0;
real V;
V=C_VOLUME(c,t);
source=C_UDMI(c,t,1)/V;
dS[eqn]=0.0;
return source;
}
DEFINE_SOURCE(Yaxi_mom_source,c,t,dS,eqn)
{
real source=0;
real V;
V=C_VOLUME(c,t);
source=C_UDMI(c,t,2);
dS[eqn]=0.0;
return source;
}
DEFINE_SOURCE(Xrad_mom_source,c,t,dS,eqn)
{
real source=0;
real V;
V=C_VOLUME(c,t);
source=C_UDMI(c,t,3);
dS[eqn]=0.0;
return source;
}
DEFINE_DIFFUSIVITY(vapor_diffusivity,c,t,i) /*vapor diffusivity in the air*/
{
real diff=0;
real T;
real P;
T=C_T(c,t);
P=C_P(c,t)+402145.725;
diff=4.7931*pow(10.0,-5.0)*pow(T,1.81)/P;
return diff;
}
int cellTest(int dd) /*judge if the given cell belongs to the neighbouring cells*/
{
int aa=0;
int i;
for(i=0;i<counter;i++)
{
if(dd==a)
{
aa=1;
j=i;
k++;
break;
}
}
return aa;
}
![udf编译出问题求大神帮助]()
K3@AGLHRO0FETY`GI2B10%O.png |