#include<iostream>
#include<math.h>
using namespace std;
double VPAD(long PNo, double pl);
double MQV(double to,int ng,int nf[],int FNo[],int k[],int PNo,double Ap,int ni,double Ai[],double mui[]);
int main()
{double pl=54.0;int PNo=400110;double to=33.2;int ng=2;int k[]={1,1};
int nf[]={1,2};long FNo[]={900057,900080};double Ap=154.0;
int ni=5;double Ai[]={1,2,3,4,5};double mui[]={0.5,0.5,0.5,0.5,0.5};
double Qf;
Qf=MQV(to,ng,nf,FNo,k,PNo,Ap,ni,Ai,mui);
cout<<"Qf="<<Qf<<endl;
return 0;
}
//通风机风量函数QFAN的变量说明：
//FNo—风机的代号，规定如下：
// 9FJ5.6 ——900056
// 9FJ6.0 ——900060
// 9FJ7.1 ——900071
// 9FJ9.0 ——900090
// 9FJ10.0——900100
// 9FJ12.5——900125
// 9FJ14.0——900140
// 如非以上型号风机,则以风机叶轮直径为代号,(单位cm,如有小数,取整)
// 但这样做实际上是近似按9FJ系列风机确定其性能，可能会有较大误差，
// 并且，这样不能确定直径560～1400mm范围以外风机的性能，
// 所以，为准确起见，应尽可能补充该种风机的性能资料
// 补充其他风机资料时，风机代号应为10000～900000之间的一整数，
//pl—通风阻力，Pa;
//通风机风量函数QFAN程序功能与使用说明：
//①函数QFAN执行结果返回单台风机的风量，m3/h;
//②目前仅有9FJ系列风机资料，函数仅针对该系列风机进行计算，
// 但今后可不断补充其他风机的资料，只需在程序中增加语句如下：
// if(FNo==风机代号) q=该风机风量的经验回归公式;
//③调用QFAN函数的实参依次为：
// 风机代号FNo；
// 风机的工作静压(通风阻力)pl，Pa；
//④如输入风机代号错误，函数QFAN执行结果给出提示，并返回值－300；
//⑤函数QFAN适用的静压为0～80Pa
// 如所给静压小于零，函数QFAN执行结果给出提示，并返回值－300；
// 如所给静压大于80Pa,函数QFAN仍给出计算结果，但给出提示；
// 如所给静压过大导致计算所得风机风量小于零，函数QFAN取风机风量为零。
//查算各型号通风机风量的函数
double QFAN(long FNo[], double pl)
{
double x,q,q1,q2;
if(pl<0.0)
{
printf("\n无该风机反向运行的风量资料，请检查输入静压值是否正确");
return (-300.0);
}
if(FNo<10000) FNo=FNo+900000;
x = pl;q=1000000.0;
if(FNo>=900056 && FNo<900060)
{
q1=-0.3198*x*x - 26.723*x + 10526;
q2=-0.2308*x*x - 29.875*x + 11970;
q=q1+(q2-q1)*(FNo-900056)/(900060-900056);
}
if(FNo>=900060 && FNo<900071)
{
q1=-0.2308*x*x - 29.875*x + 11970;
q2=-0.1699*x*x - 24.312*x + 13732;
q=q1+(q2-q1)*(FNo-900060)/(900071-900060);
}
if(FNo>=900071 && FNo<900090)
{
q1=-0.1699*x*x - 24.312*x + 13732;
q2=-0.1546*x*x - 82.548*x + 20081;
q=q1+(q2-q1)*(FNo-900071)/(900090-900071);
}
if(FNo>=900090 && FNo<900100)
{
q1=-0.1546*x*x - 82.548*x + 20081;
q2=-0.4826*x*x - 97.046*x + 26010;
q=q1+(q2-q1)*(FNo-900090)/(900100-900090);
}
if(FNo>=900100 && FNo<900125)
{
q1=-0.4826*x*x - 97.046*x + 26010;
q2=-3.5323*x*x - 17.305*x + 32765;
q=q1+(q2-q1)*(FNo-900100)/(900125-900100);
}
if(FNo>=900125 && FNo<=900140)
{
q1=-3.5323*x*x - 17.305*x + 32765;
q2=-3.0243*x*x - 27.947*x + 56682;
q=q1+(q2-q1)*(FNo-900125)/(900140-900125);
}
// 增补风机的位置及句式：if(FNo==风机代号) q=该风机风量的经验回归公式;
if(q==1000000.0)
{
printf("\n无该风机的资料，请检查输入的风机代号或叶轮直径是否正确");
return (-300.0);
}
if(pl>80.0) printf("\n静压值超过风机风量的经验公式范围，计算风量误差较大");
if(q<0.0) return (0.0);
return q;
}
//湿垫过垫风速函数VPAD的变量说明：
//PNo—湿垫的代号, 规定如下：
// 富通湿垫100mm厚——400100
// 富通湿垫125mm厚——400125
// 富通湿垫140mm厚——400140
// 富通湿垫165mm厚——400165
// 富通湿垫200mm厚——400200
// 如非以上型号湿垫,则以湿垫厚度为代号,(单位mm,如有小数,取整)
// 但这样做实际上是近似按富通湿垫确定其性能，可能会有较大误差，
// 并且，这样不能确定厚度100～200mm范围以外湿垫的性能，
// 所以，为准确起见，应尽可能补充该种湿垫的性能资料
// 补充其他湿垫资料时，湿垫代号应为10000～400000之间的一整数，
//pl—温室内外静压差(通风阻力)，Pa;
//vp—平均过垫风速，m/s;
//湿垫过垫风速函数VPAD的程序功能与使用说明：
//①函数VPAD执行结果返回过垫平均风速，m/s;
// 湿垫性能回归是以过垫风速1m/s为分界点，分二段给出计算式；
//②目前仅有富通湿垫资料，函数仅针对该湿垫进行计算，
// 但程序易于增补扩充，补充新的湿垫资料，只需在程序中增加语句如下：
// if(PNo==湿垫代号) vp=该湿垫过垫风速的经验回归公式;
//③如输入风机代号错误，函数QFAN执行结果给出提示，并返回值－300；
//④函数VPAD适用于过垫风速0～3m/s的范围；
// 温室内外静压差小于零，将给出出错提示，并返回值－300.
//查算各型号湿垫过垫风速的函数
double VPAD(long PNo, double pl)
{
double y,vp,v1,v2;
if (pl<0.0)
{
printf("\n温室内外气压差(湿垫阻力)不能为负值，请检查输入数据是否正确");
return (-300.0);
}
if(PNo<10000) PNo=PNo+400000;
y = pl; vp=1000000.0;
if(PNo>=400100 && PNo<400125)
{
if(y<=10.80) v1=y/10.80;
else v1=(-3.9+sqrt(15.21+28.0*(y+0.1)))/14.0;
if(y<=13.30) v2=y/13.30;
else v2=(-6.2+sqrt(38.44+32.0*(y+0.9)))/16.0;
vp=v1+(v2-v1)*(PNo-400100)/(400125-400100);
}
if(PNo>=400125 && PNo<400140)
{
if(y<=13.30) v1=y/13.30;
else v1=(-6.2+sqrt(38.44+32.0*(y+0.9)))/16.0;
if(y<=16.24) v2=y/16.24;
else v2=(-10.657+sqrt(113.57+33.143*(y+2.7)))/16.571;
vp=v1+(v2-v1)*(PNo-400125)/(400140-400125);
}
if(PNo>=400140 && PNo<400165)
{
if(y<=16.24) v1=y/16.24;
else v1=(-10.657+sqrt(113.57+33.143*(y+2.7)))/16.571;
if(y<=19.40) v2=y/19.40;
else v2=(-11.2+sqrt(125.44+40.0*(y+1.8)))/20.0;
vp=v1+(v2-v1)*(PNo-400140)/(400165-400140);
}
if(PNo>=400165 && PNo<=400200)
{
if(y<=19.40) v1=y/19.40;
else v1=(-11.2+sqrt(125.44+40.0*(y+1.8)))/20.0;
if(y<=23.17) v2=y/23.17;
else v2=(-1.8286+sqrt(3.3438+60.572*(y-6.2)))/30.286;
vp=v1+(v2-v1)*(PNo-400165)/(400200-400165);
}
// 增补湿垫的位置及句式：if(PNo==湿垫代号) vp=该湿垫过垫风速的经验回归公式;
if(vp==1000000.0)
{
printf("\n无该湿垫的资料，请检查输入的湿垫代号或湿垫厚度是否正确");
return (-300.0);
}
if(vp>5.0) printf("\n室内外压差过大，超过经验公式范围，通过湿垫的风速计算误差较大");
return vp;
}
//机械通风量计算函数MQV程序功能与使用说明：
//①函数MQV执行结果返回负压机械通风的温室通风量，单位为m3/s；
// 通风量统一按室外气温为准计算；
//②适用于温室具有湿垫与数个进风口、均可任意开闭的情况，风机可分组运行；
//③调用函数MQV的实参依次为：
// 室外空气温度to,℃;
// 风机的分组数量ng;
// 各组风机中的风机台数nf[];
// 各组风机的风机型号的代号FNo[];
// 各组风机开停状态k[];
// 湿垫型号的代号PNo;
// 湿垫的面积Ap，m2,如湿垫关闭，取为0；
// 进风口数量ni;
// 各进风口面积Ai，m2，注意取实际开启的面积，如某进风口关闭，取为0；
// 各进风口流量系数mui[];
//④函数MQV执行中需要调用通风机风量函数QFAN与湿垫过垫风速函数VPAD。
//⑤如因调用函数QFAN或函数VPAD发生错误（主要是因为输入实参有错误），
// 程序将终止执行，并返回值－300.0。
//机械通风量计算函数MQV
double MQV(double to,int ng,int nf[],int FNo[],int k[],int PNo,double Ap,int ni,double Ai[],double mui[])
{
double e1,e2,et,ro,pl1,pl2,plt,Qf,Qi,temp;
int e0,i,j;
e0=1; e1=1000.0;
ro=353.0/(to+273.15); pl2=25.0;
while (fabs(e1)>0.01)
{
Qf=0.0;
for(j=0;j<ng;j++)
{temp=QFAN(FNo[j],pl2);
if(temp==-300.0){printf("\n计算风机风量时输入数据错误，不能计算出机械通风量");return (-300.0);}
Qf=Qf+k[j]*nf[j]*temp/3600.0;
}
temp=VPAD(PNo,pl2);
if(temp==-300.0){printf("\n计算湿垫过垫风速时输入数据错误，不能计算出机械通风量");return (-300.0);}
Qi=Ap*temp;
for(i=0;i<ni;i++)
Qi=Qi+mui[i]*Ai[i]*sqrt(2*pl2/ro);
e2=Qf-Qi;
if(e0==1) {e1=e2;pl1=pl2;if(e2>0) pl2=pl2+10.0; else pl2=pl2-10.0;}
if(e0==2 && fabs(e1)<fabs(e2)) {et=e2;e2=e1;e1=et;plt=pl2;pl2=pl1;pl1=plt;}
if(e0>=2) {plt=pl2;pl2=pl1-e1*(pl1-pl2)/(e1-e2);e1=e2;pl1=plt;}
if(pl2<0.0) pl2=0.0;
e0=e0+1;
printf("\npl1=%9.4f pl2=%9.4f",pl1,pl2);
printf("\nQf=%9.4f Qi=%9.4f e1=%10.4f e2=%10.4f",Qf,Qi,e1,e2);
printf("\n\ni= "); scanf("%d",&i);
}
return Qf;
}
