第五章 結果討論
6.2 未來工作
現階段完成的成果，可以確認將汽電共生觀念應用於HTGR 循環與產氫循環 中，這個研究構想是可行的，亦作不少的分析，但在未來仍需要有更多的參數分 析才能完成這個研究工作，外來工作包括：
1.搜索更多產氫循環和新一代的核能電廠，因為產氫循環不只有碘硫循環一種，
而每種產氫循環皆有相對應最適合的核電廠。
2.考慮更多參數，較合乎實際情況，再藉由統計分析方法(例如：田口法)來評估 各參數的重要性。
3.除了利用 ORC 及布雷登循環進行廢熱回收外，可再增加更多廢熱回收策略。
參考文獻
[1] J. A. Lake, “The Fourth Generation of Nuclear Power,” Progress in Nuclear
Energy, Vol. 40 No.3~4, pp. 301~307, 2002.
[2] Y. Xu, and K. Zuo, “Overview of the 10 MW high temperature gas cooled reactor – test module project,” Nuclear Engineering and Design, Vol. 218, pp.
13-23, 2002.
[3] Z. Huang, J. Wang, and J. Li, “study on the thermodynamic cycle of HTR-10,”
2
ndInternational Topical Meeting on High Temperature Reactor Technology,
Beijing, China, September 22-24, 2004.[4] S. Hu, R. Wang, and Z. Gao, “Safety demonstration tests on HTR-10,” 2nd
International Topical Meeting on High Temperature Reactor Technology,
Beijing, China, September 22-24, 2004.[5] D.R. Nicholls, “Status of the Pebble Bed Modular Reactor,” Nuclear Energy, Vol. 39, No. 4, pp. 231, 2000.
[6] T. Takizuka, “Reactor Technology Development Under the HTTR Project,”
Progress in Nuclear Energy, Vol. 47, No. 1-4, pp. 283-291, 2005.
[7] X. Yan, K. Kunitomi, T. Nakata., “HTR300 design and development,” Nuclear
Engineering and Design, Vol. 222, pp. 247-262, 2003.
[8] C. H. Oh, C. B. Davis, S. R. Sherman., “Development of HyPEP, A Hydrogen Production Plant Efficiency Calculation Program,” INL, March 2006.
[9] A. Terada, J. Iwatsuki, S. Ishikural., “Development of Hydrogen Production Technology by Thermochemical Water Splitting IS Process Pilot Test Plan,”
Nuclear and science technology Vol.44, No.3, pp447-482, 2007.
[10] C. Huang, A. T-Raissi,”Analysis of sulfur-iodine thermochenical cycle for solar hydrogen production. PartI:decomposition of sulfuric acid,” solar
Energy,Vol.78,pp.632-646
[11] S. Kasahara, G. J.Hwang, H. Nakajima., ”Effects of Process Parameters of the IS Process on Total Thermal Efficiency to Produce Hydogen from Water,”
Chemical Engineering of Japan , Vol.36 No.7, pp. 887- 899,2003
[12] T. C. Hung, “Triple-Cycle: A Conceptual Arrangement of Multiple-Cycle toward Optimal Energy Conversion,” ASME Journal of Engineering for Gas
Turbines and Power, Vol. 124, pp. 429-436, 2002
.[13] T. C. Hung, T. Y. Shai, and S. K. Wang, “A review of organic Rankine cycles (ORCs) for the recovery of low-grade waste heat,” Energy, Vol. 22, No. 7, pp.
661-667, 1997.
[14] T. C. Hung, “Waste heat recovery of organic Rankine cycle using dry fluids,”
Energy Conversion and Management, Vol. 42, pp. 539-553, 2001.
[15] B. T. Liu, K. H. Chien, and C. C. Wang, “Effect of working fluids on organic Rankine cycle for waste heat recovery,” Energy, Vol. 29, pp.1207–1217, 2004.
[16] P. J. Mago, L. M. Chamra, and C. Somayaji, “Performance analysis of different working fluids for use in organic Rankine cycles,” Proceedings of the Institution
of Mechanical Engineers, part A: Journal of Power and Energy, Vol. 221, No. 3,
pp. 255~263, 2007.[17] D. Wei, X. Lu, Z. Lu., “Performance analysis and optimization of organic Rankine cycle (ORC) for waste heat recovery,” Energy Conversion and
Management, Vol. 48, pp.113–1119, 2007.
[18] R. E. Sonntag, C. Borgnakke and G. J. V. Wtlen , ”Fundamentals of
Thermodynamics 5
th,” John Wiley & Sons, NY, pp.359-387, 2003.[19] T. Schulenberg, H. Wider and M. A. Fütterer, “Electricity Production in Nuclear Power Plants - Rankine vs. Brayton Cycles,” Eisenhower's Global
Vision for Nuclear Energy, New Orleans, LA, United States, Nov 16-20,
2003.[20] A. Moisseytsev, J. J. Sienicki, “Transient accident a Analysis of a supercritical carbon dioxide Brayton cycle energy converter coupled to an autonomous
lead-cooled fast reactor,” Nuclear Engineering and Design, Vol. 238, pp.
2094–2105, 2008
[21] J. M. Smith, H. C. Van Ness, M. M. Abbott, ”Introduction to Chemical
Engineering Thermodynamics 5
th,” McGRAW-Hill, NY, pp.80-85,1996.[22] K. Tochigi, H. Futakuchi, K. Kojima. “Prediction of vapor–liquid equilibrium in polymer solutions using a Peng–Robinson group contribution model,” Fluid
Phase Equilibria, Vol.152, pp.209–217, 1998.
[23] K.A.M. Gasem, W. Gao, Z. Pan. “A modified temperature dependence for the Peng–Robinson equation of state,” Fluid Phase Equilibria, Vol.181,
pp.113–125, 2001.
[24] H. Lin, Y. Duan. “Empirical correction to the Peng–Robinson equation of state for the saturated region,” Fluid Phase Equilibria, Vol.233, pp.194–103, 2005.
[25] A.H. Farrokh-Niae, H. Moddarress and M. Mohsen-Nia. “A three-parameter cubicequation of state for prediction of thermodynamic properties of fluids,” J.
Chem. Thermodynamics, Vol.40, pp.84-95, 2008.
[26] P.K. Cheekatamarla, Caine M. Finnerty, Yanhai Du . “Advanced tubular solid oxide fuel cells with high efficiency for internal reforming of hydrocarbon fuels,”
Journal of Power Sources, Vol.188, pp.521-526, 2009.
[27] John A. Dean, “Lange’s Handbook of Chemistry 13th
, ” McGRAW-Hill,
pp.9.1-10.130, 1987.[28] Bruce E. Poling,John P.O’Connell, “The Properties of Gases and Liquids 5th
, ”
McGRAW-Hill, pp.3.1-8.198, 2001.附錄 A 計算流體熱力性值和各循環 Matlab 程式
clc;clear all;
Pin= input('please input pressure(MPa): ');
T= input('please input temperature(k): ');
R=8314.34;w=0.210;Pc=4.895*10^6;Vc=256*10^-3;Tc=562.05;
Hc=547.34;sc=1.1590;M=78.114;
P=Pin*10^5;
Tr=T/Tc;
Pvp=10^(3.98523-(1184.24/(T+217.572-273.15)))*10^5;
P0=Pvp;
Cp=(3.551-6.184*10^-3*T+14.365*10^-5*T^2-19.807*10^-8*T^3+8.234*10^-11*T^4)*R;
%%Peng-Robinson equation of state%%
B1=0.37464;B2=1.54230;B3=-0.26992;B4=0.5;B5=0.5;
b=0.007780*R*Tc/Pc;
a1=0.45724*(R^2)*(Tc^2)/Pc;
a2=(1+(B1+B2*w+B3*w^3)*(1-Tr^B4))^(1/B5);
a=a1*a2;
c1=P0;
c2=(P0*b-R*T);
c3=-3*P0*b^2-2*R*T*b+a;
c4=R*T*b^2-a*b-b^3*P0;
p=[c1 c2 c3 c4];
S=roots(p);
S1=max(S);
u=[1 -S1];
[q,r]=deconv(p,u);
S2=roots(q);
S4=max(S2);
S3=947.36*(Tr^6)-4073.9*(Tr^5)+7237.7*(Tr^4)-6796.6*(Tr^3)+3557.1*(Tr^2) -983.34*Tr+112.48;
Vl0=S3*Vc;
%%% get enthalpy and entropy B0=0.07780*P0*Tc/Pc/T;
Z0=P0*S1/R/T;
ZB0=(Z0+(1+2^0.5)*B0)/(Z0+(1-2^0.5)*B0);
bb=(1+k*(1-Tr^0.5))*log(ZB0);
aa=Tr*(Z0-1);
cc=2.078*(1+k);
deltaH1=R*Tc*(cc*bb-aa);
deltaH0=R*Tc*(2.078*k*((1+k)*(1+k*(1-(T/Tc)^0.5))*log(ZB0)-T/Tc*(Z0-1)))
;
Bc=0.07780;
Zc=Pc*Vc/R/Tc;
ZBc=(Zc+(1+(2)^0.5)*Bc)/(Zc+(1-(2)^0.5)*Bc);
%deltaH
deltaHc1=R*Tc*(2.078*(1+k)*log(ZBc)-Zc+1);
dd=(1+k*0)*log(ZBc);
ee=(Zc-1);
deltaHc=R*Tc*(dd*ff-ee);
S1;%vg Vl0;%vf
Pvp/10^6;%P(MPa)
CpTc=(3.551*Tc-6.184*(10^-3)/2*Tc^2+14.365*(10^-5)/3*Tc^3-19.807*(10^-8)/4*Tc^4+8.234*(10^-11)/5*Tc^5)*R;
CpT=(3.551*T-6.184*(10^-3)/2*T^2+14.365*(10^-5)/3*T^3-19.807*(10^-8)/4*
T^4+8.234*(10^-11)/5*T^5)*R;
Cpi=CpT-CpTc;
Dhv=deltaHc+Cpi-deltaH1;
dPvp=(3.98523*log(10)+log(10^5)-log(Pvp))*Pvp/(T+217.572-273.15);
dH=T*(S1-Vl0)*dPvp;
DHL=Dhv-dH;
Hv=Hc+Dhv/1000/M;
Hl=Hc+DHL/1000/M;
%deltas
gg=2.078*k*(((1+k)/Tr^0.5)-k)*log(ZB0);
deltas0=R*(gg-log(Z0-B0));
hh=2.078*k*(((1+k)/1^0.5)-k)*log(ZBc);
deltasc=R*(hh-log(Zc-Bc));
CpTTc=(3.551*log(Tc)-6.184*(10^-3)*Tc+14.365*(10^-5)/2*Tc^2-19.807*(10^
-8)/3*Tc^3+8.234*(10^-11)/4*Tc^4)*R;
CpTT=(3.551*log(T)-6.184*(10^-3)*T+14.365*(10^-5)/2*T^2-19.807*(10^-8)/
3*T^3+8.234*(10^-11)/4*T^4)*R;
CpTi=CpTT-CpTTc;
Dsl=Dsv-dH/T;
sv=sc+Dsv/1000/M sl=sc+Dsl/1000/M
%%%Organic Rankine cycle%%%
clc;clear all;
T3=input('請輸入 expander inlet 溫度: ');
T1=input('請輸入 pump inlet 溫度: ');
e=input('請輸入 expander efficiency: ');
e1=input('請輸入 heat exchanger efficiency: ');
Q=input('請輸入熱傳量:') M=78.114;
[Hv3 Hl3 sv3 sl3 Pvp3 Vv3 Vl03]=total(T3);
[Hv1 Hl1 sv1 sl1 Pvp1 Vv1 Vl01]=total(T1);
P2a=Pvp3;v2a=Vl01;h1a=Hl1;H3a=Hv3;P1a=Pvp1;P3a=Pvp3;s3a=sv3;
h2a=v2a*(P2a-P1a)/M/1000+h1a;
dx=1;ni=(T3-T1)/dx;
T(1:ni)=T1;
T(1)=T1;
for i=1:ni;
T=T1+(i-1)*dx;
[H4i s4i]=sh1(T1,T);
er=abs(s4i-s3a);
if er<=0.01;
break
end
H4a=H3a-e*(H3a-H4i);
eff=((H3a-H4a)+(h1a-h2a))/(H3a-h2a) m=Q*e1/(H3a-h2a)
Wnet=m*((H3a-H4a)+(h1a-h2a))
%%% hydrogen production%%%
%%%初始條件%%%
h3=0.0948097*1000;h4=0.0743523*1000;h5=0.0119523*1000;h6=0.0125493*1 000/2;h7=0.03189748*1000;
ef=0.9;mhe=27.5;cphe=5.1943;T0=863.8;
T=input('請輸入 SO3 裂解反應溫度( )=')℃ deltaT=input('請輸入熱交換溫差( )=')℃
%%%求出 SO3 裂解焓差%%%
t=(T+273.15)/1000;
%%A1~H1 SO3H 修正係數 A2~B2 SO2 修正係數 A3~B3 O2 修正係數 Z1~Z3 為標準生成焓%%
%%求出 SO3 裂解焓%%
Z1=-395.77;A1=24.02503;B1=119.4607;C1=-94.38686;D1=26.96237;E1=-0.11 7517;F1=-407.8526;G1=253.5186;H1=-395.7654;
Z2=-296.81;A2=21.43049;B2=74.35094;C2=-57.75217;D2=16.35534;E2=0.086 731;F2=-305.7688;G2=254.8872;H2=-296.8422;
Z3=0;A3=29.659;B3=6.137261;C3=-1.186521;D3=0.09578;E3=-0.219663;F3=-9.861391;G3=237.948;H3=0;
hr1=A1*t+B1*t^2/2+C1*t^3/3+D1*t^4/4-E1/t+F1-H1+Z1;
hr2=A2*t+B2*t^2/2+C2*t^3/3+D2*t^4/4-E2/t+F2-H2+Z2;
hr3=A3*t+B3*t^2/2+C3*t^3/3+D3*t^4/4-E3/t+F3-H3+Z3;
h1=abs(hr1-hr2-0.5*hr3);
%%%求出 SO3 升溫時所需的焓%%%定義 t1=(400+273.15)/1000%%%
t1=(400+273.15)/1000;
hr4=A1*t1+B1*t1^2/2+C1*t1^3/3+D1*t1^4/4-E1/t1+F1-H1+Z1;
%%%計算 SO3 質流率%%%
T1=T+deltaT;
mso3=ef*mhe*cphe*(T0-T1)/h1;
%%%求出各個物質質流率%%%
mh2so4=mso3;mso2=mso3;mi2=mso3;mhi=2*mso3;
%%%求出氦氣在 SO3 熱交換器出口溫度( ),℃ 求出熱交換至 SO3 之熱量
%%%
T2=T1-mso3*h2/mhe/cphe/ef Q2=mso3*h2;
%%%求出氦氣在 H2SO4 裂解反應出口溫度( ),℃ 求出熱傳量%%%
T3=T2-mh2so4*h3/mhe/cphe/ef Q3=mh2so4*h3;
%%%氦氣在 H2SO4 熱交換器出口溫度( )℃ ,求出熱傳量%%%
T4=T3-mh2so4*h4/mhe/cphe/ef Q4=mh2so4*h4;
%%%求出氦氣在 HI 裂解反應出口溫度( ),℃ 求出熱傳量%%%
T5=T4-mhi*h5/mhe/cphe/ef Q5=mhi*h5;
%%%氦氣在 HI 熱交換器出口溫度( ),℃ 求出熱傳量%%%
T6=T5-mhi*h6/mhe/cphe/ef Q6=mhi*h6;
%%%求出氦氣在 bunsen reaction 出口溫度( ),℃ 求出熱傳量%%%
T7=T6-mhi*h6/mhe/cphe/ef Q7=mhi*h7;
%%%helium Brayton cycle%%%
%%%初始條件%%%
clc;clear all;
T1=input('請輸入產氫反應出口溫度( )=');℃ e1=input('請輸入氣渦輪機效率(1<=e<=0)=');
e2=input('請輸入低壓段壓縮機效率(1<=e<=0)=');
e3=input('請輸入高壓段壓縮機效率(1<=e<=0)=');
T8=479.9;T4=25;T6=25;m=27.5;k=1.66618;cp=5.1943;
rtmax=4;dx=0.0001;rti=1;ni=(rtmax-rti)/dx+1;
Rt(1:ni)=rtmax;
Rt(1)=rti;
Rt(ni)=rti;
T7(1:ni)=T8;
T7(1)=T6;
T7(ni)=T6;
%%%計算輸入總熱量%%%
Qin=m*cp*(T1-T8);
for i=1:ni;
T2=T1-(T1+273.15)*(1-Rt(i)^((1-k)/k))*e1;
%%%計算出高壓段壓縮機出口溫度及作功及回熱器出口溫度%%%
T7=T6+(T6+273.15)*(r2(i)^((k-1)/k)-1)/e3;
T3A=T2+T7-T8;
Tef=abs(T3A-T3B);
if Tef<=0.01;
break end
end T3=T3A;
Whca=m*cp*(T6-T7) Whci=Whca*e3;
Wta=m*cp*(T1-T2) Wti=Wta/e1;
%%%計算低壓段壓縮機出口溫度及作功%%%
r1=r2(i);
T5=T4+(T4+273.15)*(r1^((k-1)/k)-1)/e3;
Wlca=m*cp*(T4-T5) Wlci=Wlca*e3;
%%%計算 precooler 及 intercooler 熱量%%%
Qp=m*cp*(T4-T3);
Qi=m*cp*(T6-T5);
%%%計算熱效率及靜功%%
Wnet=(Wta+Wlca+Whca) eff=(Wta+Wlca+Whca)/Qin
%%%Case 2%%%
%%%helium Brayton cycle%%%
clc;clear all;
T=input('請輸入 SO2 溫度=');
Q=input('請輸入 SO2 熱傳量=');
e1=input('請輸入氣渦輪機效率(1<=e<=0)=');
e2=input('請輸入高壓壓縮機效率(1<=e<=0)=');
e2=input('請輸入低壓機效率(1<=e<=0)=');
rtmax=input('請輸入最大壓縮比=');
;k=1.66618;cp=5.1943;rp0=1;R=8.314/4;
rtmax=20;dx=0.01;rti=2;ni=(rtmax-rti)/dx+1;
T4(1:ni)=25;
Rt(1:ni)=rtmax;
T1(1:ni)=105;
T2(1:ni)=T-15;
T6(1:ni)=25;
Wt(1:ni)=0;
Rt(1)=rti;
Rt(ni)=rti;
eff1(1:ni)=0;
eff(1)=1;
eff(ni)=1;
T11=max(T1);
T22=max(T2);
Qb(1:ni)=Q;
m=Q*0.9/cp/(T22-T11);
%%%記算輸入總熱量%%%
Rt(i)=rti+(i-1)*dx;
R1(i)=Rt(i)^0.5;
end
for i=1:ni
%%%計算出 turbine 作功及出口溫度%%%
T3(i)=T2(i)-(T2(i)+273.15)*(1-Rt(i)^((1-k)/k))*e1;
%%%計算出壓縮機出口溫度及作功%%%
T5(i)=T4(i)+(T4(i)+273.15)*(R1(i)^((k-1)/k)-1)/e2;
T7(i)=T6(i)+(T6(i)+273.15)*(R1(i)^((k-1)/k)-1)/e2;
Whca=m*cp*(T4(i)-T5(i));
Wlca=m*cp*(T6(i)-T7(i));
Wta=m*cp*(T2(i)-T3(i));
%%%計算低壓段壓縮機出口溫度及作功%%%
%%%計算 precooler 及 intercooler 熱量%%%
%Qp=m*cp*(T4-T3);
%Qi=m*cp*(T6-T5);
%%%計算熱效率及靜功%%
Wnet=(Wta+Whca);
eff=Wnet/Q;
eff1(i)=eff;
ds1(i)=cp*log((T2(i)+273.15)/(T1(i)+273.15))-R*log(rp0);
ds2(i)=cp*log((T3(i)+273.15)/(T2(i)+273.15))-R*log(1/Rt(i));
ds3(i)=cp*log((T4(i)+273.15)/(T3(i)+273.15))-R*log(rp0);
ds4(i)=cp*log((T5(i)+273.15)/(T4(i)+273.15))-R*log(R1(i));
ds5(i)=cp*log((T6(i)+273.15)/(T5(i)+273.15))-R*log(rp0);
I1(i)=m*(T4(i)+273.15)*(ds1(i)+ds2(i)+ds3(i)+ds4(i)+ds5(i)+ds6(i));
A(i)=(Qb(i)-I1(i))/Qb(i);
end
effa=max(eff1);
[AX,H1,H2] = plotyy(Rt,eff1,Rt, A,'plot');
%%%Case 1%%%
%%%helium Brayton cycle%%%
%%%初始條件%%%
clc;clear all;
T=input('請輸入 SO2 溫度=');
Q=input('請輸入 SO2 熱傳量=');
e1=input('請輸入氣渦輪機效率(1<=e<=0)=');
e2=input('請輸入壓縮機效率(1<=e<=0)=');
;k=1.66618;cp=5.1943;rp0=1;R=8.314/4;
rtmax=4;dx=0.01;rti=2;ni=(rtmax-rti)/dx+1;
T4(1:ni)=25;
Rt(1:ni)=rtmax;
T1(1:ni)=105;
T2(1:ni)=T-15;
T41(1:ni)=T4;
Wt(1:ni)=0;
Rt(1)=rti;
Rt(ni)=rti;
eff1(1:ni)=0;
eff(ni)=1;
T11=max(T1);
T22=max(T2);
Qb(1:ni)=Q;
m=Q*0.9/cp/(T22-T11)
%%%記算輸入總熱量%%%
T3(i)=T2(i)-(T2(i)+273.15)*(1-Rt(i)^((1-k)/k))*e1;
%%%計算出壓縮機出口溫度及作功%%%
T5(i)=T4(i)+(T4(i)+273.15)*(Rt(i)^((k-1)/k)-1)/e2;
Whca=m*cp*(T4(i)-T5(i));
Wta=m*cp*(T2(i)-T3(i));
%%%計算低壓段壓縮機出口溫度及作功%%%
%%%計算 precooler 及 intercooler 熱量%%%
Qp=m*cp*(T4-T3);
Qi=m*cp*(T6-T5);
%%%計算熱效率及靜功%%
Wnet=(Wta+Whca);
eff=Wnet/Q;
eff1(i)=eff;
ds2(i)=cp*log((T3(i)+273.15)/(T2(i)+273.15))-R*log(1/Rt(i));
ds3(i)=cp*log((T4(i)+273.15)/(T3(i)+273.15))-R*log(rp0);
ds4(i)=cp*log((T5(i)+273.15)/(T4(i)+273.15))-R*log(Rt(i));
I1(i)=m*(T41(i)+273.15)*(ds1(i)+ds2(i)+ds3(i)+ds4(i));
A(i)=(Qb(i)-I1(i))/Qb(i);
end
[effa n]=max(eff1);
Rp=rti+n*dx;
[AX,H1,H2] = plotyy(Rt,eff1,Rt, A,'plot');
xlabel('Rp(壓縮比)');
set(get(AX(1),'Ylabel'),'String','熱效率') set(get(AX(2),'Ylabel'),'String','有效率') Rtm=input('請輸入最佳壓縮比:');
%%%計算出 turbine 作功及出口溫度%%%
T3=T2-(T2+273.15)*(1-Rtm^((1-k)/k))*e1;
%%%計算出壓縮機出口溫度及作功%%%
T5=T4+(T4+273.15)*(Rtm^((k-1)/k)-1)/e2;
Whca=m*cp*(T4-T5);
Wta=m*cp*(T2-T3);
Wnet=(Wta+Whca);
Wnet1=max(Wnet) eff=Wnet/Q;
eff =max(eff);
T51=T5(n)