|
马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?注册
x
我做的是等离子体发生器的数值模拟,用UDF编写ar等离子体物性,通过双线性插值法编写物性,通过1atm,0.1atm,0.01atm,0.001atm压力下,温度范围在300K~30000k,间隔50K,通过插值得到粘性,密度,热导率,声速,但在初始化的时候出现问题。
#include "udf.h"
float gasprop(float t,float table[595]);
float rho,vis,ak,ah,cp,sig,speed; //密度 粘性系数 导热系数 比焓 定压比热 电导率 声速
float tbt0[595],tbrho0[595],tbvis0[595],tbak0[595],tbsig0[595],tbah0[595],tbcp0[595],tbspeed0[595]; //温度 密度 粘性系数 导热系数 电导率 比焓 定压比热 声速
float tbt1[595],tbrho1[595],tbvis1[595],tbak1[595],tbsig1[595],tbah1[595],tbcp1[595],tbspeed1[595];
float tbt2[595],tbrho2[595],tbvis2[595],tbak2[595],tbsig2[595],tbah2[595],tbcp2[595],tbspeed2[595];
float tbt3[595],tbrho3[595],tbvis3[595],tbak3[595],tbsig3[595],tbah3[595],tbcp3[595],tbspeed3[595]; //定义用于存储物性的数组;
float temperature; //同fluent联系
float pressure; //同fluent联系
float rhoh,rhol,vh,vl,akh,akl,cph,cpl,sigh,sigl,speedh,speedl;
float p0,p1,p2,p3;
//////////////////////////////////////////读取物性参数////////////////////////////////////////////////////////
DEFINE_ON_DEMAND(Property_data)
{
FILE *fp;
int i;
fp = fopen("1.0atm.DAT", "r");
for(i=0;i<595;i++)
fscanf(fp,"%e %e %e %e %e %e %e %e",&tbt0,&tbrho0,&tbvis0,&tbak0,&tbsig0,&tbah0,&tbcp0,&tbspeed0); //温度 密度 粘性系数 导热系数 电导率 比焓 定压比热 声速
fclose(fp);
fp = fopen("0.1atm.DAT", "r");
for(i=0;i<595;i++)
fscanf(fp,"%e %e %e %e %e %e %e %e",&tbt1,&tbrho1,&tbvis1,&tbak1,&tbsig1,&tbah1,&tbcp1,&tbspeed1);
fclose(fp);
fp = fopen("0.01atm.DAT", "r");
for(i=0;i<595;i++)
fscanf(fp,"%e %e %e %e %e %e %e %e",&tbt2,&tbrho2,&tbvis2,&tbak2,&tbsig2,&tbah2,&tbcp2,&tbspeed2);
fclose(fp);
fp = fopen("0.001atm.DAT", "r");
for(i=0;i<595;i++)
fscanf(fp,"%e %e %e %e %e %e %e %e",&tbt3,&tbrho3,&tbvis3,&tbak3,&tbsig3,&tbah3,&tbcp3,&tbspeed3);
fclose(fp);
}
//////////////////////////////////////////////////定义粘性系数/////////////////////////////////////
DEFINE_PROPERTY(viscosity,c,t)
{
real mu_lam;
p0 = 101325 * 1;
p1 = 101325 * 0.1;
p2 = 101325 * 0.01;
p3 = 101325 * 0.001;
temperature= C_T(c,t);
pressure=C_P(c,t);
if (pressure>=p0)
{
vh= gasprop(temperature,tbvis0);
vis=vh;
}
else if (pressure>=p1&&pressure<p0)
{
vh= gasprop(temperature,tbvis0);
vl=gasprop(temperature,tbvis1);
vis=vh+(vl-vh)/(p1-p0)*(pressure-p0);
}
else if (pressure>=p2&&pressure<p1)
{
vh= gasprop(temperature,tbvis1);
vl=gasprop(temperature,tbvis2);
vis=vh+(vl-vh)/(p2-p1)*(pressure-p1);
}
else if (pressure>=p3&&pressure<p2)
{
vh= gasprop(temperature,tbvis2);
vl=gasprop(temperature,tbvis3);
vis=vh+(vl-vh)/(p3-p2)*(pressure-p2);
}
else
{
vh= gasprop(temperature,tbvis3);
vis=vh;
}
mu_lam=vis;
return mu_lam;
}
//////////////////////////////////////////////////定义热导率/////////////////////////////////////
DEFINE_PROPERTY(thermal_conductivity,c,t)
{
real ther;
p0 = 101325 * 1;
p1 = 101325 * 0.1;
p2 = 101325 * 0.01;
p3 = 101325 * 0.001;
temperature= C_T(c,t);
pressure=C_P(c,t);
if (pressure>=p0)
{
akh=gasprop(temperature,tbak0);
ak=akh;
}
else if (pressure>=p1&&pressure<p0)
{
akh= gasprop(temperature,tbak0);
akl=gasprop(temperature,tbak1);
ak=akh+(akl-akh)/(p1-p0)*(pressure-p0);
}
else if (pressure>=p2&&pressure<p1)
{
akh= gasprop(temperature,tbak1);
akl=gasprop(temperature,tbak2);
ak=akh+(akl-akh)/(p2-p1)*(pressure-p1);
}
else if (pressure>=p3&&pressure<p2)
{
akh= gasprop(temperature,tbak2);
akl=gasprop(temperature,tbak3);
ak=akh+(akl-akh)/(p3-p2)*(pressure-p2);
}
else
{
akh=gasprop(temperature,tbak3);
ak=akh;
}
ther=ak;
return ther;
}
//////////////////////////////////////////////////定义密度/////////////////////////////////////
DEFINE_PROPERTY(density,c,t)
{
real dens;
p0 = 101325 * 1;
p1 = 101325 * 0.1;
p2 = 101325 * 0.01;
p3 = 101325 * 0.001;
temperature=C_T(c,t);
pressure=C_P(c,t);
if (pressure>=p0)
{
rhoh=gasprop(temperature,tbrho0);
rho=rhoh;
}
else if (pressure>=p1&&pressure<p0)
{
rhoh=gasprop(temperature,tbrho0);
rhol=gasprop(temperature,tbrho1);
rho=rhoh+(rhol-rhoh)/(p1-p0)*(pressure-p0);
}
else if (pressure>=p2&&pressure<p1)
{
rhoh= gasprop(temperature,tbrho1);
rhol=gasprop(temperature,tbrho2);
rho=rhoh+(rhol-rhoh)/(p2-p1)*(pressure-p1);
}
else if (pressure>=p3&&pressure<p2)
{
rhoh=gasprop(temperature,tbrho2);
rhol=gasprop(temperature,tbrho3);
rho=rhoh+(rhol-rhoh)/(p3-p2)*(pressure-p2);
}
else
{
rhoh=gasprop(temperature,tbrho3);
rho=rhoh;
}
dens=rho;
return dens;
}
//////////////////////////////////////////////////定义声速/////////////////////////////////////
DEFINE_PROPERTY(sonic_speed,c,t)
{
real sonic;
p0 = 101325 * 1;
p1 = 101325 * 0.1;
p2 = 101325 * 0.01;
p3 = 101325 * 0.001;
temperature= C_T(c,t);
pressure=C_P(c,t);
if (pressure>=p0)
{
speedh=gasprop(temperature,tbspeed0);
speed=speedh;
}
else if (pressure>=p1&&pressure<p0)
{
speedh=gasprop(temperature,tbspeed0);
speedl=gasprop(temperature,tbspeed1);
speed=speedh+(speedl-speedh)/(p1-p0)*(pressure-p0);
}
else if (pressure>=p2&&pressure<p1)
{
speedh=gasprop(temperature,tbspeed1);
speedl=gasprop(temperature,tbspeed2);
speed=speedh+(speedl-speedh)/(p2-p1)*(pressure-p1);
}
else if (pressure>=p3&&pressure<p2)
{
speedh=gasprop(temperature,tbspeed2);
speedl=gasprop(temperature,tbspeed3);
speed=speedh+(speedl-speedh)/(p3-p2)*(pressure-p2);
}
else
{
speedh=gasprop(temperature,tbspeed3);
speed=speedh;
}
sonic=speed;
return sonic;
}
float gasprop(float t,float table[595])
{
int k;
float delt,gasp;
if(t<300) t=300;
if(t>30000) t=30000;
k=(int)((t-300)/50);
delt=(t-k*50-300);
gasp=table[k]+delt/50*(table[k+1]-table[k]);
return (gasp);
}
初始化时,出现这个问题
Error:
FLUENT received fatal signal (ACCESS_VIOLATION)
1. Note exact events leading to error.
2. Save case/data under new name.
3. Exit program and restart to continue.
4. Report error to your distributor.
Error Object: ()
我曾尝试之导入粘性和热导率的UDF都没有问题,但是只要密度的UDF一加入就会出这个错误,但是我找不出原因!!
[ 本帖最后由 jony322 于 2013-5-18 13:02 编辑 ] |
|