马 捷,孔巧玲,刘雁集
(上海交通大学海洋工程国家重点实验室,上海 200030)
能量是物质的状态,物质的一定状态对应着一定的能量。物质储存了能量,就是储存了运动的能力。
物质储存能量的方式各有不同,通过相变的方式储存能量,称为相变储能。相变储能的机理,在于某些物质在物相变化过程中从环境吸收热量或是冷量,也可以向环境放出热量或是冷量,以实现能量储存和释放。借助这个过程,可以达到调节能量需求和实施能量供给的目的。相变是物质不同相之间的相互转变。物质在相的转变过程中伴随吸热和放热,就能发挥储能的功能。相变过程中,单位质量的物质吸收或放出的热量,就是相变潜热。
在相变储能过程中,物质温度基本保持不变,这时,传热的温差小,能够有效减少换热过程中的不可逆损失,有助于提高换热效率[1]。常见的相变方式有固-液相变、气-液相变和固-气相变。其中,固-液相变具有较大的承压能力。固-液相变储能问题即运动边界问题,亦称Stefan问题,因Stefan于1891年研究北极冰层厚度时提出该问题,因而得名。常温下,固-液相变过程的辐射效应可以忽略不计。固态的相变材料(phase change material,PCM)内,热传递的方式主要是热传导,在液态的相变材料内,热传递方式既有热传导,又有对流换热。
根据传热方式,固-液相变储能可分为两种。一种是导热控制的固-液相变,过程主要受固相区导热的影响,多存在于液相区处于凝固温度的场合。研究具有运动固-液相界面的单相区域的非稳态导热,可确定凝固层厚度与时间的关系。如果是简单几何形状下的凝固或融化问题,是可以用解析法进行求解的。另一种是固-液相变储能,过程是导热和对流耦合的固-液相变,液相区有明显的流动,固-液界面的形状和位置受固相区导热和液相区对流的共同影响。它与常规的对流问题有所不同,因为伴随相变的对流效应,受固相区域变化的影响,也受固-液界面边界形状和热流密度的影响,也就是说,存在导热和对流的耦合。该问题是近年来相变储能领域的研究热点。
20世纪70年代以来,固-液相变材料广泛应用于储能领域,其研究和应用涉及材料学、太阳能、工程热物理、空调采暖、工业废热利用、建筑物保温、航空航天等诸多学科及领域,研究内容也涉及多个方面[2]。
首先是相变材料的研究。相变材料是蓄热主体,寻找理想的相变材料一直是研究的重点。相变材料的选择和研究,包括材料的物性,诸如熔点、潜热、比热容、导热系数,还有复合材料、可变熔点相变材料,也涉及材料热物性的测定、热分析以及物理化学领域。常见的热物性测量方法,有热量法、差热扫描法、差热分析法。差热扫描法应用较多,Flaherty[3]用差热扫描法测量了碳氢化合物和石蜡的热物性,取得成果;Giavarini等[4]用差热扫描法成功分析了石油产品的热物性;阮德水[5]借助差热扫描法研究了多元醇、层状钙矿等几种相变材料的热物理性能。Baran[6]测出二元系统的成分比例并用差热扫描法测定热物性;Zhang等[7]发现D差热扫描法的缺陷,并提出了一种简便的时间历程法来分析相变材料的热物性;Marin等[8]提出测量焓和比热容随温度变化的有效方法,得到焓随温度变化的关系式;Manoo等[9]确定了部分石蜡类材料的导热系数、比热容以及焓随温度变化的关系式。Delaunay等[10]分析了一维圆柱导热,以此为基础测量了相变温度附近的相变材料导热系数;陈则韶等[11-12]在相变材料热物性测定和相变过程导热分析方面,也作了一系列工作。
其次,是储能的热物理问题研究。该研究包括相变储能过程机理、相变蓄热器设计、提高相变材料传热能力措施、蓄热装置传热过程强化及运行工况控制等。相变储能过程的机理研究方面,Cheng[13]将拉氏变换和控制容积法给合起来,求解了合金及有机物糊态区的二维相变问题;Giangi等[14]建立基于焓法的多孔介质模型,结合对流换热的影响,完成对蓄热器进行瞬态分析;Ismail等[15]采用控制容积法,对垂直管相变蓄热器进行瞬态数值模拟;Qarnia[16]对矩形通道的热能储存系统的瞬态热特性完成了数值模拟;Halawa等[17]在对流边界条件下,对于变壁温的储能单元的PCM相变过程进行了数值分析;Liu等[18]研究了第二类边界条件下的熔解过程,分析了储能系统的储热;Saitoh等[19]建立了瞬态熔化过程模型,研究了球形容器的自然对流和接触传热结合的熔解过程;Costa等[20]建立了焓法模型,采用完全隐式差分方法,研究了带肋片和不带肋片的潜热储能系统的热性能;Lamberg等[21]研究了内侧加肋片的凝固传热过程。管建春等[22]对圆球形相变蓄冷器的相变储能进行了数值模拟;贺友多等[23]对相变热源项进行处理,研究了钢液的凝固传热;柯秀芳等[24]建立了柱坐标下温度和流动微分方程组,以研究圆柱形金属储热体凝固过程的传热;李震等[25]给出了计算相变温度为一个温度区间的非理想相变材料蓄换热特性的适用条件;汤勇等[26]建立了蓄热材料相变储能的三维数学模型;施伟等[27]对流体在管内流动,相变材料在管外发生相变的储能器热性能进行了数值模拟;王志峰[28]对二元固液相变进行瞬态全三维的数值模拟;邢玉明[29]对空间站太阳能吸热储热器完成了数值计算;李海梅等[30]求解了相变过程的瞬态热传导问题。
在提高相变材料导热能力的措施方面,Sadasuke等[31]和Velraj[32]提出容器中加肋片的方法;Bugaje[33]在相变材料中加入金属颗粒;肖敏等[34]在复合相变材料中加入膨胀石墨,放热时间比纯石蜡缩短了 61%;王剑锋等[35]和 Gong[36]显著提高了系统传热效率;Fukai等[37]研究了布置碳纤维的方法,提高了蓄热器相变材料热传导性能;Cabeza等[38]研究了加入花岗石增强传热的效果,发现优于加入不锈钢片或铜块。
海洋资源开发与海洋环境探测的需求,促进了水下机器人及相关技术的发展。当前常见的水下机器人包括远程遥控机器人(remotely operated vehicle,ROV)和无人水下机器人(unmanned underwater vehicle,UUV),这两类机器人都不是完全自主的,极大地限制了机器人的活动领域及其工作效率。自治式水下机器人(autonomous underwater vehicle,AUV)具有水下活动范围大、机动性好、安全、结构简单等特点,海底地形地貌勘察、海洋资源及地质调查、海洋环境和水文参数测量、生物考察[39]都可采用。利用相变储能机理运行的水下机器人也称为水下滑翔机,它不需要自带能源就可以获得动力,与天空中运行的滑翔机很相像。
AUV的续航力主要由使用的能源决定,开发新能源是提高 AUV续航能力的必然趋势。近年来,低能耗、低成本和高续航能力的浮力驱动式水下滑翔机脱颖而出。这是一种依靠自身浮力驱动的水下机器人[40]。此类滑翔机通过改变机体外侧的皮囊体积,调节自身浮力。当水下机器人上浮或下沉时,作用在机翼上的升力在水平方向上的分力推动滑翔机向前航行。海洋温差能驱动的水下滑翔机,借助相变储能过程将海洋温差能转化为压力能,驱动滑翔机航行。温差能驱动水下滑翔机的动力能源来自海洋,不需要自带能源,具有续航能力强、体积小、机动灵活、制造成本低、噪音小的优点,在军事上和海洋探索研究上具有重要应用价值[41]。
水下机器人的结构外形见图1。图1(a)为Slocumm Thermal滑翔机,感温工质置于机体外侧的长圆柱形的储能管中[42];图1(b)为天津大学研制的热滑翔机,感温工质置于与机体制成一体的薄壁圆筒中[40]。
图2是Slocumm Thermal滑翔机的动力系统工作原理简图。其动力系统主要组成为:暴露在海水中的刚性储能管,是系统获取温差能的关键部件;蓄能器,存储和释放储能管中感温工质体积变化所传递的能量;体积可变的外胆,蓄能器的传递液体在其中存储和释放,以改变自身体积和滑翔机净浮力的大小;体积可变的内胆,存储和释放外胆中的传递液体;两个单向阀,单向阀 1控制储能管与蓄能器之间能量传递液体的流动;单向阀2控制储能管与内胆之间能量传递液体的流动;三通阀,控制蓄能器、内胆、外胆之间的液体流动。
图1 温差能驱动的水下滑翔机的结构外形图Fig.1 Geometry of underwater glider propelled by ocean thermal energy
图2 温差能动力系统原理图Fig.2 Principle diagram of the thermodynamic driving system
滑翔机在下潜和上浮中,调节内部质量块的位移以改变滑翔姿态角,按一定角度倾斜向下或倾斜向上运行。机翼上升力的水平方向分力,使滑翔机向水平方向运动,使滑翔机在水下以“之”字形在海洋温跃层间上下运行。获取温差能依靠相变装置,它是整个动力系统中的主要部件。在发生固-液或液-固相变时,内部的感温工质产生体积变化,这种体积变化会转化为滑翔机自身的体积变化,结果是上下的沉浮运动。感温工质的相变工作特性、动力装置的结构尺寸、外壁的边界条件等各项参数,会影响动力系统的传热速率,影响动力系统的输出功率。
依赖相变储能机理运行的水下运载器即水下滑翔机的动力系统,其储能材料的相变受制于多种因素。相变动力系统的结构尺寸及材料确定后,外界海水温度是影响相变工质体积变化及相变时间的主要因素。
水下滑翔机在水下运行,潜深变化,外界海水的温度也变化;在不同海域,海水的温度分布又不同。研究滑翔机在不同海洋温跃层的工作特性,就要研究不同海洋温跃层中储能材料的相变特性,得到滑翔机在冷、暖海水层间停留时间、储能材料体积膨胀率、储能材料相变时间,为动力系统的阀门控制和滑翔机的姿态调整提供确切依据。
初始时刻,滑翔机漂浮在水面,动力系统储能装置中,相变材料处于液态,初始温度等于海面上海水的温度;滑翔机下潜,滑翔机接触的外界海水温度逐渐降低,海水温度低于相变温度,感温工质开始凝固,液态逐步变为固态,体积收缩;滑翔机运动到最低水深,相变材料凝固,体积最小;调节浮力,滑翔机开始上浮,接触的外界海水温度高于相变温度时,感温工质开始熔化,固态变为液态,体积膨胀;滑翔机上升到水面,相变材料完全熔化,成为液态,完成一个循环。
海水层温度梯度大于深水跃层临界值,属于温度跃层。海水层的顶部水深称作跃层上界,海水层的厚度称作跃层厚度,海水层垂向温度梯度定义为跃层强度。
在不同的海洋温跃层中,温度分布也不同,其间运行的滑翔机,须设计不同的冷水层、暖水层停留时间,使相变材料完全相变。以稳定的深海跃层型温跃层[43]为研究对象,考察滑翔机运行时动力系统的工作,温跃层的数据采自赤道附近海域,温度分布如图3(a)所示,分段函数拟合后曲线如图3(b)所示。
图3 温跃层温度分布图[44]Fig.3 Temperature profile of the ocean[44]
温度随深度变化的函数为
式中,z为深度,m;T为温度,℃。
设定滑翔机处于中性浮力状态,最大体积变化350 mL。参照文献[40]滑翔机外型,用3根外置式长圆柱管为滑翔机动力装置储能管,直径0.03 m,长1.6 m,总传热面积为0.45 m2,总体积3393 mL。用正十六烷为相变材料,储能管沉浸海水中,忽略管与机体、管与管间传热影响,认为各管与外界海水传热性能相同。取一根管,进行相变传热分析,滑翔机速度0.25 m/s,稳定运行在赤道附近温跃层,海洋温跃层的温差不大,于对流换热系数影响小,取海水平均温度15℃,计算得对流换热系数297 W/m2·K。
由赤道附近海域温度分布,距离水面96 m处,海水温度等于正十六烷相变温度。离水面96 m处,海水层为冷暖水层分界面。滑翔机在从水面到水下96 m之间水域,相变为熔解过程;滑翔机在离水面96 m以下水域运行,相变为凝固过程。
初始时刻,滑翔机处于正浮力,浮于水面,储能装置中的十六烷初始温度等于外围水面温度(26 ℃),呈液态。
讨论两种工况。滑翔机运行到最大潜深,十六烷完全凝固,改变浮力大小,开始上浮行程,到达水面,完成一个循环。随后,改变浮力,开始新的下潜。
图4显示滑翔机工作过程,可见滑翔潜深、十六烷液相分数、中心温度和体积变化率随时间的变化曲线。图 4(a)显示滑翔机潜深随运行时间变化,滑翔机下潜到最大深度后返回水面为一循环。完成一循环所需时间193 min,最大潜深1215 m。
图4 滑翔机工作过程分析Fig.4 Analysis of gliding process for underwater glider
图4(b)为滑翔机运行过程中相变材料液相分数随时间变化曲线。第一个下潜行程,初始时刻,相变材料液相分数为 1。滑翔机运行到冷暖水层分界面以下,外界的海水稍低于相变温度,十六烷凝固。在这段行程内,十六烷保持液态,液相分数等于1,潜深增大,海水温度降低,液相分数减小,到达最大潜深,液相分数为 0,十六烷完全凝固。调节浮力大小,滑翔机上浮,滑翔机运行到冷暖水层分界面前,海水温度小于相变温度,十六烷保持固态,液相分数为 0;滑翔机继续上浮,海水温度稍高于相变温度,十六烷熔解。滑翔机到达水面,液相分数0.17,十六烷部分熔解。滑翔机又下潜,到达冷暖水层分界处之前的行程内,外界海水温度高于相变温度,十六烷熔解,液相分数达到0.22;滑翔机继续下潜,海水温度低于相变温度,十六烷凝固。同样行程,22%为液态十六烷,只需要22 min就完全凝固。在随后下潜行程,十六烷保持固态。在接下来的几个循环中,十六烷的液相分数都在0.22~0之间变化。即只有22%的材料参与相变。
图4(c)为十六烷体积变化率随时间变化曲线。可见:初始时刻,十六烷体积处于体积最大的膨胀状态,滑翔机下潜,液相分数保持为1时,体积变化率不变;海水温度降低到相变温度以下,体积收缩,十六烷完全凝固,体积达最小;上浮一段行程,海水温度仍低于十六烷相变温度,故很长一段时间十六烷呈固态,体积不变。滑翔机运行到水面,部分十六烷熔解,体积膨胀率3.9%。随后几个循环中,体积变化率在0~4.4%之间。
图4(d)给出容器中心温度随时间的变化。在初始时刻,十六烷的温度等于外界海水;随潜深增加,海水温度降低,导热使十六烷放出显热,温度随着降低;温度低于相变温度后,十六烷开始凝固。完全凝固前温度不变。
滑翔机运行到最大潜深,十六烷完成凝固,在上浮行程,十六烷温度高于海水温度,十六烷继续放热,温度迅速降低;随后,十六烷的温度随着海水温度而升高,海水温度超过相变温度,开始熔解,十六烷温度不变;到达水面,滑翔机开始新的下潜,十六烷处于部分熔解状态;外界海水温度高于相变温度,十六烷继续熔化,温度等于相变温度,保持不变;海水温度低于相变温度,十六烷放出热量,开始凝固。
从图4(d)可知,部分熔解的相变材料凝固后,十六烷的温度与海水温度相同,接下来,十六烷在相变温度与最低温度间变化。
由上分析,十六烷材料利用率低,22%十六烷发生相变,最大体积膨胀率只有4.4%,不能按照初始时设定的体积膨胀率进行循环,滑翔在水下不能稳定运行。
相变是物质系统不同相之间的相互转变。物质在相互转变时,伴有吸热、放热,就具备了储能效应,当同时伴有体积变化现象时,则可以应用于动力系统。相变储能过程中,物质温度基本保持不变,传热温差小,减少了换热过程中的不可逆损失,有助于提高换热效率。
固-液相变储能过程具有较大的承压能力。利用海洋温差能驱动的水下滑翔机工作在水下几百米深的海水中,承受巨大的外界海水压力,选用固-液相变材料作为感温工质-储能材料,可以保证相变储能动力系统稳定可靠地工作。
[1] Zhang Yanping(张寅平),Hu Hanping(胡汉平).Phase Change Energy Storage—Theory and Application(相变贮能——理论与应用)[M].Hefei:University of Science and Technology of China Press,1996:50-100 .
[2] Xie Chungang(谢春刚).Design and experiments on underwater gliders propelled by thermal engine[D].Tianjin:Tianjin University,2005.
[3] Flaherty B.Characterisation of waxes by differential scanning calorimetry[J].Appl.Chem.Biotechnol.1971,21(5):144-148.
[4] Giavarini C,Pochetti F.Characterization of petroleum products by DSC analysis[J].Thermal Anal.,1973(5):83-94.
[5] Ruan Deshui(阮德水).DSC research on phase change energy storage material[J].Acta Energiae Solaris Sinica(太阳能学报),1994,15(1):40-45.
[6] Baran G,Sari A.Phase change and heat transfer characteristics of a eutectic mixture of palmitic and stearic acids as PCM in a latent heat storage system[J].Energy Conversion and Management,2003,44:3227-3246.
[7] Zhang Yinping,Jiang Yi.A simple method, the T-history method, of determining the heat of fusion, specific heat and thermal conductivity of phase-change materials[J].Measurement Sci.Technol.,1999,10(3):201-205.
[8] Marin J M,Zalba B,Cabeza L F,Mehling H.Determination of enthalpy-temperature curves of phase change materials with the T-history method:Improvement to temperature dependent properties[J].Measurement Sci.Technol.,2003,14(2):184-189.
[9] Manoo A,Hensel E.One-dimensional two-phase moving boundary problem[C]//Phase Change Heat Transfer.ASME,1991(159):97-102.
[10] Delaunay D,Carre P.Dispositifs de mesure automatique de la conductivite thermique des materiaux a changement de phase[J].Rev.Phys.Appl,1982,17(4):209-215.
[11] Chen Zeshao(陈则韶).Thermophysical properties of phase change storage materials and the relationship between the enthalpy increment and relative volume increment for paraffin wax[J].Acta Energiae Solaris Sinica(太阳能学报),1983,4(1):9-15.
[12] Chen Zeshao(陈则韶).A simple heat-resistance method for the solution to heat conduction undergoing solidification[J].Journal of University of Science and Technology of China(中国科学技术大学学报),1991,21(3):69-77.
[13] Cheng T F.A numerical simulation for two-dimensional moving boundary problems with a mushy zone[J].Computational Mechanics,1999,23(5-6):440-447.
[14] Giangi M,Stella F,Kowalewski T A.Phase change problems with free convection:Fixed grid numerical simulation[J].Comput.Visual.Sci.,1999,2(2-3):123-130.
[15] Ismail K A R,Abugderah M M.Performance of a thermal storage system of the vertical tube type[J].Energy Conversion &Management,2000,41(11):1165-1190.
[16] Qarnia H E.Theoretical study of transient response of a rectangular latent heat thermal energy storage system with conjugate forced convection[J].Energy Conversion & Management,2004,45(9-10):1537-1551.
[17] Halawa E,Bruno F,Saman W.Numerical analysis of a PCM thermal storage system with varying wall temperature[J].Energy Conversion& Management,2005,46(15-16):2592-2604.
[18] Liu Zhongliang,Ma Chongfang .Numerical analysis of melting with constant heat flux heating in a thermal energy storage system[J].Energy Conversion & Management,2002,43(18):2521-2538.
[19] Saitoh T S,Hoshina H,Yamada K.Theoretical analysis and experiment on combined close-contact and natural convection melting in thermal energy storage spherical capsule[C]//Proceedings of the Intersociety Energy Conversion Engineering Conference,Energy Systems,Renewable Energy Resources, Environmental Impact and Policy Impacts on Energy,1997:1656-1661.
[20] Costa M,Buddhi D,Oliva A.Numerical simulation of latent heat thermal energy storage system with enhanced heat conduction[J].Energy Conversion & Management,1998,39(3-4):319-330.
[21] Lamberg P,Sirén K.Approximate analytical model for solidification in a finite PCM storage with internal fins[J].Applied Mathematical Modelling,2003,27(7):491-513.
[22] Guan Jianchun(管建春),Zhu Hua(朱华).Numerical analysis of heat transfer with phase change in cool storage hall[J].Journal of Zhejiang University(浙江大学学报),1999,33(1):85-89.
[23] Zhang Yin(张胤),He Youduo(贺友多),Li Shiqi(李士琦),Shen Yishen(沈颐身).A mathematical model for phase change heat transfer[J].Journal of Baotou University of Iron and Steel Technology(包头钢铁学院学报),1999,18(2):94-97.
[24] Ke Xiufang(柯秀芳),Zhang Renyuan(张仁元).Study on heat transfer of a cylindrical thermal storage unit with metal medium in the solidification process[J].New Energy(新能源),2000,22(3):23-27.
[25] Li Zhen(李震),Zhang Yinping(张寅平),Jiang Yi(江亿).Effect of specific heat of phase change material on heat charging or discharging performance[J].Acta Energiae Solaris Sinica(太阳能学报),2002,23(1):1-5.
[26] Tang Yong(汤勇),Zhang Juan(张娟),Zhou Deming(周德明),Wang Wei(王威),Wang Xiaowu(王小五).The analog and analysis on the solid/liquid interface of fiber-phase matrix composites forenergy storage[J].Journal of South China University of Technology(华南理工大学学报),2001,29(11):1-6.
[27] Shi Wei(施伟),Ge Xinshi(葛新石).Thermal performance of a latent heat thermal storage module with flow fluid inside tube and phase change material outside[J].Acta Energiae Solaris Sinica(太阳能学报),2004,25(4):497-502.
[28] Wang Zhifeng(王志峰).A numerical solution of solid-liquid binary phase change heat transfer[J].Acta Energiae Solaris Sinica(太阳能学报),2000,21(1):7-14.
[29] Xing Yuming(邢玉明).Numerical simulation of high temperature phase change heat storages system[J].Journal of Aerospace Power(航空动力学报),2002,17(2):231-235.
[30] Li Haimei(李海梅),Gu Yuanxian(顾元宪).Finite element solution of heat transfer with planar phasechange by equivalent heat capacity method[J].Journal of Dalian University of Technology(大连理工大学学报),2001,40(1):45-48.
[31] Sadasuke I,Naokatsu M.Heat transfer enhancement by fins in latentheat thermal energy storage devices[J].Solar Eng.,1991:223-228.
[32] Velraj R.Heat transfer enhancement in a latent heat storage system[J].Solar Energy,1999,65(3):171-180.
[33] Bugaje I M.Enhancing the thermal response of latent heat storage systems[J].Int.J.Energy Res.,1997,21(9):759-766.
[34] Xiao Min(肖敏),Gong Kecheng(龚克成).Preparation of a good thermal conductive, shape-stabilized phase change material and its performance study[J].Acta Energiae Solaris Sinica(太阳能学报),2001,22(4):427-431.
[35] Wang Jianfeng(王剑锋),Chen Guangming(陈光明).Study on charging and discharging rates of latent heat thermal energy storage systems employing multiple phase change materials[J].Acta Energiae Solaris Sinica(太阳能学报),1998,21(3):258-265.
[36] Gong Zhenxiang.Cyclic heat transfer in a novel storage unit of multiple phase change materials[J].Applied Thermal Engineering,1996,16(11):807-815.
[37] Fukai J,Kanou M.Thermal conductivity enhancement of energy storage media using carbon fibers[J].Energy Conversion &Management,2000,41(14):1543-1556.
[38] Cabeza L F,Mehling H.Heat transfer enhancement in water when used as PCM in thermal energy storage[J].Applied Thermal Engineering,2002,(22):1141-1151.
[39] Feng Xisheng(封锡盛),Liu Yongkuan(刘永宽).The research and development of autonomous underwater vehicle[J].High Technology Letters(高科技通讯),1999,(9):55-59.
[40] Wang Shuxin(王树新),Wang Yanhui(王延辉),Zhang Datao(张大涛).Design and trialon and underwater glider propelled by thermal engine[J].Ocean Technology(海洋技术),2006,25(1):1-5.
[41] Hui Shaotang(惠绍棠).Research on technology development of autonomus underwater vehicle[J].Ocean Technology(海洋技术),2001,20(4):13-14.
[42] Zhang Datao(张大涛).System design and experiments on profile propelled by thermal engine[D].Tianjin:Tianjin University,2005.
[43] Jenkins S A,Jolla L,Humphreys D E,et al.Alternatives for enhancement of transport economy in underwater gliders[C]//OCEANS,2003:948-950.
[44] Gwyn Griffiths.Technology and Applications of Autonomous Underwater Vehicles[M].Davis R E,Eriksen C E,Jones C P.London:Taylor and Francis,2002:324.