核反应堆热工分析课程设计剖析

if m==0 if q==0

[q,qmax,ql,qlmax,m]=T32(qMAX,DNBR,Fq,Dcs,L,Fa,Nt,n); else

[m,q,qmax,ql,qlmax]=T31(q,Fq,Dcs,Fa,L,n,Nt); end else

[q,qmax,ql,qlmax]=T3(m,n,Dcs,L,Fa,Nt,Fq); End

if qmax*DNBR>qMAX

disp('警报：qmax 过高，可能会发生快速烧毁') End

Def=T5(s,n0,delta,m);

Tf2=T6(W,casi,Tin,P,Nt,Fh,deltaB);

Vf2=T41(P,Tp,W,casi,m,n,s,Dcs,delta,n0,deltaB);%Vf2 为半高度处的冷却剂流速 h2l=T7(s,Dcs,P,Tp,Vf2); tcsmax=T8(P,qmax,h2l,Tf2); tcimax=T9(tcsmax);

tumax=T10(qlmax,Dci,Du,tcimax); tomax=T11(qlmax,tumax); if tomax>2200

disp('堆芯温度过高，最高温度为'); disp(tomax); else

disp('堆芯最高温度为：'); disp(tomax); End

函数 T1：

function [Tin,Tout]=T1(Tp,delta) Tin=Tp-delta/2; Tout=Tp+delta/2;

%已知平均温度和温度上升值，求入口出口值

函数 T2：

function Tout=T2(Tin,P,Fa,Nt,W,casi) Tout(10)=0; Tf(10)=0; Tout(1)=Tin+45; for i=1:10

Tf(i)=0.5*(Tout(i)+Tin);

Cp= WASPCHS('pt2Cp',P,Tf(i))*1000;

Tout(i+1)=Tin+(Fa*Nt*1000000)/((W/3.6)*(1-casi)*Cp); e=abs(Tout(i+1)-Tout(i))/Tout(i); if e>0.001 continue else

Tout=Tout(i+1); break end end

%已知入口温度，出口温度求解 函数 T3

function [q,qmax,ql,qlmax]=T3(m,n,Dcs,L,Fa,Nt,Fq) F=m*n*pi*Dcs*L; q=Fa*Nt*1000000/F;

qmax=q*Fq; ql=q*pi*Dcs; qlmax=ql*Fq;

%两个系数，求热流密度 函数 T31

function [m,q,qmax,ql,qlmax]=T31(q,Fq,Dcs,Fa,L,n,Nt) q=q*1000000; ql=q*pi*Dcs; qlmax=ql*Fq; qmax=q*Fq;

N=Nt*1000000*Fa/(pi*Dcs*L*q); m=ceil(N/n); m=ceil((m-1)/4)*4+1;

%在不知道燃料组数 m，却知道平均热流密度时候，求解 m,qlmax,qmax,等 函数 T32

function [q,qmax,ql,qlmax,m]=T32(qMAX,DNBR,Fq,Dcs,L,Fa,Nt,n) qmax=qMAX/DNBR; q=qmax/Fq; ql=q*pi*Dcs; qlmax=ql*Fq;

N=Nt*1000000*Fa/(pi*Dcs*L*q); m=ceil(N/n); m=ceil((m-1)/4)*4+1;

%不知道燃料组件数不知道热流密度情况下求解 m,qlmax,qmax;

函数 T42

function Vf2=T41(P,Tp,W,casi,m,n,s,Dcs,delta,n0,deltaB) rou= 1/WASPCHS('pt2V',P,Tp); Cp=WASPCHS('pt2Cp',P,Tp)*1000;

At=m*n*(s*s-(pi/4)*(Dcs*Dcs))+m*4*n0*s*delta/2; Vf2=(W/3.6)*(1-casi)*(1-deltaB)/(At*rou); %求热管半高度处冷却剂流速 函数 T5

function Def=T5(s,n0,delta,m) T=(n0*s+delta)^2; Def=sqrt((m*T/(pi/4))); %求等效直径 函数 T6

function Tf2=T6(W,casi,Tin,P,Nt,Fh,deltaB) V= WASPCHS('pt2V',P,Tin); H= WASPCHS('pt2H',P,Tin);

Wef=(1-casi)*(1-deltaB)*W/(3.6*V); H2=H+(Nt*1000000*Fh)/(2*Wef); Tf2=WASPCHS('pH2T',P,H2);

%求取热管半高处水的焓值，求出半高处温度。1.67 为因子 函数 T7

function h2l=T7(s,Dcs,P,Tp,Vf2) De=4*(s*s-pi*Dcs*Dcs/4)/(pi*Dcs); U= WASPCHS('pt2U',15,305); RAMD= WASPCHS('pt2RAMD',15,305);

PRN= WASPCHS('pt2PRN',15,305); rou=1/WASPCHS('pt2v',15,305); Re=Vf2*De/U;

Nu=0.023*(Re^0.8)*(PRN^0.4); h2l=Nu*RAMD/De; %求半高度处的比焓 函数 T8

function tcsmax=T8(P,qmax,h2l,Tf2) deltaf=qmax*1000000/h2l; ts= WASPCHS('p2T',P);

deltaj=ts+25*(qmax)^0.25*(exp(-P/6.2))-Tf2; if deltaf>deltaj tcsmax=Tf2+deltaj; else

tcsmax=Tf2+deltaf; end

%计算燃料元件表面最高温度 tcsmax 函数 T9

function tcimax=T9(tcsmax) tc1=tcsmax+25; tc=tcsmax+35;

while abs(tc-tc1)>=0.001 tc=tc1;

k=0.00547*(1.8*tc+32)+13.8;

tcimax=tcsmax+0.0168*10^6*2.32*1.03*0.57/(k*pi*(9.5+8.36)/2); tc1=(tcimax+tcsmax)/2; end