FLUENT 用UDF写的PMV和PPD的问题

shuibixiansheng

我现在需要用fluent计算PMV,PPD,PD, 刚开始学习UDF,在网上找了别人编写的程序想试试正确与否。用interpreted udf. 编译后，云图显示的全是0，有人知道是啥问题吗？
下面是程序：
#include "udf.h"

DEFINE_ON_DEMAND(comfort)
{
Domain *d;
Thread *t;
cell_t c;
real CLO; // Clothing (clo)
real MET; // Metabolic rate (met)
real WME; // External work, normally around 0 (met)
real TA; // Air temperature (°C)
real TR; // Mean radiant temperature (°C)
real VEL; // Relative air velocity (m/s)
real RH; // Relative humidity (%)
real PA; // Water vapour pressure (Pa)
real TU; // Turbulence intensity (%)
real FNPS; // Saturated vapour pressure, kPa
real ICL; // Clothing insulation, m2K/W
real TCL; // Clothing surface temperature, °C
real M; // Metabolic rate, W/m2
real W; // External work, W/m2
real MW; // Internal heat production in the human body, W/m2
real FCL; // Clothing area factor
real HCF; // Forced convective heat transfer coefficient, W/m2K
real HCN; // Natural convective heat transfer coefficient, W/m2K
real HC; // Convective heat transfer coefficient, W/m2K
real TAA; // Air temperature, K
real TRA; // Mean radiant temperature, K
real HL1; // Heat loss difference through skin
real HL2; // Heat loss by sweating (comfort)
real HL3; // Latent respiration heat loss
real HL4; // Dry respiration heat loss
real HL5; // Heat loss by radiation
real HL6; // Heat loss by convection
real TS; // Thermal sensation transfer coefficient
real EPS; // Small value
real i; // Iterator
real TCLA; // First guess for surface temperature of clothing
real P1,P2,P3,P4,P5; // Calculation terms
real XN,XF; // Iterative solutions
real PMV,ePMV, PPD, DR; // Predicted mean vote, extended predicted mean vote, predicted percentage dissatisfied and draft
d=Get_Domain(1);
EPS=0.00015;
CLO=0.5;
MET=1.2;
WME=0.0;
RH=75.0;
ICL=0.155*CLO;
M=MET*58.15;
W=WME*58.15;
MW=M-W;
if (ICL<=0.078) FCL=1.0+1.29*ICL; else FCL=1.05+0.645*ICL;
thread_loop_c(t,d)
{
begin_c_loop(c,t)
{
//TA=C_T(c,t)-273;
TA=299-273;
TR=26;
VEL=sqrt(pow(C_U(c,t),2)+pow(C_V(c,t),2)+pow(C_W(c,t),2));
FNPS=exp(16.6536-4030.183/(TA+235));
PA=RH*10.0*FNPS;
HCF=12.1*sqrt(VEL);
TAA=TA+273.0;
TRA=TR+273.0;
TCLA=TAA+(35.5-TA)/(3.5*(6.45*ICL+0.1));
P1=ICL*FCL;
P2=P1*3.96;
P3=P1*100;
P4=P1*TAA;
P5=308.7-0.028*MW+P2*pow(TRA/100,4);
XN=TCLA/100.0;
XF=XN;
for (i=0;i<150;i++)
{
XF=(XF+XN)/2.0;
HCN=2.38*pow(fabs(100.0*XF-TAA),0.25);
if (HCF>HCN) HC=HCF; else HC=HCN;
XN=(P5+P4*HC-P2*pow(XF,4))/(100+P3*HC);
if (fabs(XN-XF)<EPS) break;
}
TCL=100.0*XN-273.0;
HL1=3.05*0.001*(5733.0-6.99*MW-PA);
if (MW>58.15) HL2=0.42*(MW-58.15); else HL2 = 0;
HL3=1.7*0.00001*M*(5867.0-PA);
HL4=0.0014*M*(34.0-TA);
HL5=3.96*FCL*(pow(XN,4)-pow(TRA/100.0,4));
HL6=FCL*HC*(TCL-TA);
TS=0.303*exp(-0.036*M)+0.028;
PMV=TS*(MW-HL1-HL2-HL3-HL4-HL5-HL6);
ePMV=0.9*PMV;
PPD=100.0-95.0*exp(-0.03353*pow(ePMV,4)-0.2179*pow(ePMV,2));
TU=pow(2*C_K(c,t)/3,0.5)/VEL;
if (VEL<0.05) DR=0; else DR=(34-TA)*pow(VEL-0.05,0.62)*(0.37*VEL*TU*100+3.14);
if (DR>100) DR=100;
C_UDMI(c,t,0)=ePMV;
C_UDMI(c,t,1)=PPD;
C_UDMI(c,t,2)=DR;
}
end_c_loop(c,t)
}
}