0 引言

电晕放电引起的无线电干扰、可听噪声、电晕损失等效应,已成为制约特高压输电线路结构设计和导线选型的重要因素^[1-2]。电晕笼的结构便于拆卸和更换导线,兼具适用性和经济性,是各机构研究电晕效应的重要试验装置,在笼内试验导线上施加一定电压后能在导线表面产生较高的电场强度,可用于模拟实际特高压线路电晕放电^[3-5]。但电晕笼内试验导线长度较短,其表面附近电场强度分布与实际输电线路存在差异,特别是导线端部存在端部效应会改变笼内空间电场分布,对试验结果造成影响。目前对交流电压下电晕笼内导线表面附近电场分布的计算主要集中在2维空间内,鲜有考虑到导线有限长特性的影响以及电晕条件下带电离子运动的3维方向性。笼内3维电场强度分布与试验导线的电晕放电分布直接相关,当电晕笼试验采用较高电压等级时,有必要对笼内有限长导线表面附近电场分布进行深入研究。

Comber M G等人通过对1 100 kV输电线路和电晕笼的不同分裂导线的最大电场强度进行研究,最早提出了笼内表面场强的计算理论,结果表明交流导线电晕放电产生的空间电荷被束缚在导线周围很小的区域内^[6]。谭震宇等人在研究导线起晕电压时对分裂导线表面电场进行了分析,并测量得到了雾霾条件下不同分裂导线的起晕场强^[7]。上述研究均忽略了空间电荷的影响,但随着电压等级的提高,电晕产生的带电粒子对空间电场的影响不断增大,有必要进行考虑空间电荷的交流电晕笼电场分析。

J.J. Clade^[8]等人采用Kaptzov假设^[9],建立了电晕笼单根导线交流离子流场计算模型。Salam等人采用模拟电荷法建立了单相导线线板结构的离子流场的计算模型^[10],并提出了起晕电荷^[11]的概念,考虑了非对称几何结构中空间电荷对电场方向的影响。以上对电晕笼内离子流场的研究均在2维空间内,模型中采用的无限长直电荷无法考虑电晕笼导线的端部效应,其作用范围和影响程度尚待研究,因此需要建立更符合实际的3维离子流场计算模型。

为了研究交流电晕笼内导线表面附近3维电场分布特性,本文采用模拟电荷法建立了考虑空间电荷运动的交流电晕笼3维离子流场计算模型,模型考虑了导线表面场强不均匀性,以及导线起晕、空间电荷发射、迁移和损耗的过程。按此模型计算并分析了500 kV交流电压下特高压电晕笼内场强沿线分布规律及其附近空间场强分布规律,最后与ANSYS仿真结果对比完成了误差分析。所得结论对电晕笼试验结果分析及电晕特性研究有参考意义。

1 计算方法

1.1 假设条件

1）Kaptzov假设

Kaptzov假设认为导线起晕后其表面场强维持在起晕场强值不变^[9]。Salam在此基础上提出临界起晕电荷,即导线起晕后其电荷维持在起晕电荷不变。本文建模时沿用了临界起晕电荷的概念,在每个时间步长下计算分裂子导线表面各点的临界起晕电荷。

2）忽略电晕层厚度

考虑到电子运动速度较快且不稳定,而导线表面电晕层内涉及一系列电离、吸附等过程难以定量计算,且电晕层厚度相对子导线半径可忽略不计,故本文所建模型中导线表面电晕直接产生并发射离子。

3）忽略弧垂影响

我国现用的户外特高压电晕笼长度<40 m,笼内导线长度较实际线路可忽略不计,因此不考虑弧垂产生的影响。

4）忽略离子扩散作用

在离子流场中,交流电压产生的电场是空间离子运动的最主要原因,空间电荷因浓度差引起的扩散较小,可忽略不计。

1.2 交流电晕笼3维电场计算模型

模拟电荷法^[11-13]主要用于求解静电场,其基本思想是基于电磁场的唯一性定理,用一组离散化的模拟电荷来等效电极表面连续分布的自由电荷或介质分界面连续分布的束缚电荷。针对模拟电荷设置相应的匹配点,建立矩阵方程,求解模拟电荷的值,并应用叠加原理求解场域中的电场。

1.2.1 模拟电荷及其镜像电荷

求解中电位系数的计算量决定了模拟电荷法的计算量。而电位系数计算与模拟电荷的类型紧密相关。3维空间中电荷和场点的数量比2维时扩大了很多倍,而模拟电荷的数量、位置和类型的选取直接决定了空间场强计算结果的准确性,因此不仅要符合实际情况,也要兼顾计算的方便性。与2维模型中的无限长线电荷不同,3维电晕笼电场计算模型中的模拟电荷采用有限长的线电荷,分布在分裂子导线内部,方向与分裂子导线相同,均垂直于xy平面,通过z方向上的分段和编号实现分裂导线3维离散化,如图1所示。

实际中电晕笼体保持良好接地,故从右侧开始按逆时针方向采用镜像法依次对四边笼壁作一次镜像电荷,笼内子导线按顺时针依次编号1—8号,镜像电荷不作编号。综合考虑计算精度和计算量,本文选择每根子导线内平均分布8个模拟电荷,均匀分布在子导线半径的中点上,并在各子导线表面对应位置上设置8个起晕校核点。任意一段导线的模拟电荷、镜像电荷及其校核点分布在xy平面上,如图2所示。

1.2.2 电位系数和场强系数计算

3维空间中,电晕笼笼内各空间点的电位和场强受所有模拟电荷和空间电荷影响。根据静电场原理,图3中有限长线电荷q在空间点(x₀,y₀,z₀)处的电位和电场强度可以分别采用电位系数p和场强系数f_x、f_y、f_z表示：

式中：l为线电荷长度,m;(x_qa, y_qa, z_qa)和(x_qb, y_qb, z_qb)分别为线电荷两端端点坐标,m;d_a和d_b为线电荷两端点与空间点间的距离,m;ε为空气介电常数。

1.2.3 导线起晕与电荷发射

当分裂导线表面各点的电场强度超过起晕场强时,导线表面有电晕产生。电晕笼中分裂子导线表面场强受各子导线模拟电荷产生的屏蔽效应和空间电荷的影响而具有不均匀性,导线表面不同位置发射的电荷量也有所区别,所以在每一个时间步长上分别计算导线表面各点的起晕电荷^[14]。

导线模拟电荷及其镜像电荷采用矩阵Q_cond表示,空间电荷及其镜像电荷采用矩阵Q_space表示。导线表面的电位矩阵U由Q_cond产生的电位和Q_space产生的电位叠加而成,属于标量的叠加,用电位方程式(5)表示;模拟电荷所在处的场强矩阵E由3维方向上的场强叠加而成,属于矢量的叠加。而导线模拟电荷Q_cond和空间电荷Q_space在笼内任意点处产生的场强可用场强方程式(6)表示：

式中：P_cond和P_space分别是电荷Q_cond和Q_space对导线表面各点的电位系数矩阵;F_cond和F_space则分别是电荷Q_cond和Q_space对导线表面各点的场强系数矩阵,由3维方向的场强系数组成。

图1 3维电晕笼导线分段示意图 Fig.1 Three-dimensional view of corona cage wire segment

图2 模拟电荷和镜像电荷示意图 Fig.2 Simulating charge and image charge

图3 有限长线电荷和空间点位置示意图 Fig.3 Finite line charge and spatial point position schematic

首个交流正半周期内,升压未起晕时导线表面场强为模拟电荷产生的标称场强,当电压U升高至第k个点场强正好超过起晕场强E_onset±时,记录此时的模拟电荷Q_cond,k为k点的起晕电荷Q_onset,联立式(5)和式(6)后得到式(7),可以分别计算不同表面点的起晕电荷。正、负半周期的起晕场强E_onset+和E_onset-由Peek公式计算得到^[15]

如果空间中已经存在空间电荷,根据空间电荷的分布情况,可求出空间电荷对导线表面各点的电位U_space和场强E_space,经叠加后再按上述方法求得表面各点的起晕电荷量,联立式(5)—(7)得到

计算时,对比模拟电荷量Q_cond与起晕电荷Q_onset±,如果某点对应的模拟电荷量超出起晕电荷,就判定为该点起晕并发射空间电荷。空间电荷量为模拟电荷与起晕电荷的差值,按式(9)计算

式中：下标k表示第k个点的相应物理量。

1.2.4 空间电荷的迁移与复合

当空间电荷与模拟电荷同极性时,在电场力作用下远离导线;当2者异极性时则靠近导线。若某个时间空间电荷到达表面,将其电荷量置为0。运动中空间电荷的迁移速度为

式中：μ为正负离子迁移率,其中正离子迁移率为1.5×10^-4 m²/(V·s),负离子迁移率为1.8×10^-4 m²/(V·s);Δd为空间电荷运动的距离;Δt为运动时间步长。

导线表面附近的异性电荷在运动中有小概率复合,复合过程中电荷损耗的速度与局部电荷密度相关。某空间电荷在多个时间步长后,若其空间电荷密度降低了2个数量级,就认为其不再产生影响并将其除去。t时刻某空间电荷Q_t经过一个时间步长的损耗后其电荷量Q_t+Δt变为

式中：ρ_±为正负空间电荷所在区域的电荷密度;复合系数γ=1.5×10^-12 m²/s。

每个时间步长结束时,空间电荷沿其电场方向运动到新的位置,形成新的空间电荷分布。统计本周期内的空间电荷产生量Q_newspace,并与前一个周期内的该值比较,2者相差<1%时认为计算稳定。交流电晕笼3维电场计算模型的基本流程见图4。

2 交流电晕笼3维电场计算模型实例

计算中采用的电晕笼为特高压电晕笼,笼体全

图4 3维电场计算模型基本流程 Fig.4 Process of three-dimensional electric field strength calculation model

长35 m,包括两端各5 m防护笼和中部25 m测量笼,其截面均为边长8 m的正方形。笼体内安装型号为8×LGJ720的钢芯铝绞线,分裂间距为400 mm,分裂导线半径为18.12 mm。导线轴向剖分为每段3.5 m,共10段。试验导线上施加的交流电压为500 kV。

减小计算时间步长可以更准确计算各变量的时域值,为将时间间隔引起的误差控制在3%范围内,将每个交流周期均分为40份,计算时间步长取0.5 ms。

3 空间电场计算结果及分析

采用本文所建模型计算特高压电晕笼实例,研究笼内导线表面附近空间电场的变化规律。3维空间内远离导线方向为电场强度的正方向。经7~8个交流周期后,计算达到稳定,考虑空间电荷影响时取第9个周期内的电场强度计算结果,不考虑空间电荷影响时取第1个周期内的电场强度计算结果。

3.1 工频周期不同时刻沿子导线表面最大场强分布

根据文献[16]的研究,导线某处各分裂子导线表面场强绕子导线一周的分布规律仅存在相位差异,其最大场强均出现在分裂子导线外径点上,且相差≤3%,因此本文选择1号子导线作为沿导线表面最大场强分布的研究对象,如图5所示。

选取正、负半周期测量笼内导线未起晕(瞬时电压200 kV)和电压峰值时刻1号子导线表面外径点上的沿线场强进行比较,结果如图6所示。电压峰值时刻导线已起晕。

由图6可见,由于笼体结构上的对称,中心点两边的导线表面校核点上的场强也是对称的,导线未起晕时刻防护笼内导线表面场强发生畸变^[16],峰值时刻畸变减弱,表面场强沿线分布均匀。

由于正、负半周期起晕场强不同,正半周期测量笼内导线未起晕和电压达到峰值时导线表面场强均比负半周期要大,电晕时发射的空间电荷抑制了表面场强继续增大,正极性未起晕时上一个交流周期内残留的空间负电荷更多。导线未起晕时,中部导线表面场强分布均匀,由于导线端部效应,端部5 m内导线表面场强分别高出测量笼内该值72.3%和75.8%,此时残留的空间电荷返回导线表面。电压达到正、负半周峰值时,沿线场强稍大于起晕场强,端部导线表面场强畸变减弱,其原因是由于导线起晕发射空间电荷后模拟电荷减少至起晕电荷,起晕产生的同极性空间电荷抑制了端部导线表面场强增大。

选取1号子导线中部15.75 m处的外径点,一个交流周期内的场强时域分布如图7所示。由于残留的异极性空间电荷增大了导线表面场强,子导线外径点场强在时域上领先于电压。正半周期内导线在3.5 ms起晕后模拟电荷减少至起晕电荷,导线表面场强维持在起晕场强附近,维持2 ms后停止起晕;由于负极性起晕场强比正极性时小,负半周期内子导线外径点在13 ms时开始起晕,维持3 ms后停止起晕。

3.2 电压峰值时中部空间电场分布

正半周期峰值时刻,远处残留的空间负电荷正在返回导线表面,此时导线全线起晕,并已发射空间正电荷。以点表示有限长空间线电荷,电晕笼内空间电荷分布如图8所示。

图8中空间电荷分布在分裂导线外侧,第2—8段导线表示的测量笼内电荷分布情况相似,而防护笼体内空间电荷数量更多,分布半径更大。选取第5段导线中点所在径向截面,此时残留的负极性空

图5 子导线外径点示意图 Fig.5 Outside point of sub-conductor surface

图6 子导线(1号)场强沿线分布 Fig.6 Electric field intensity distribution along 1th sub-conductor

图7 中部子导线外径点场强时域分布图 Fig.7 Temporal distribution of electric field intensity on middle sub-conductor surface

间电荷已返回至导线表面附近,而正极性空间 电荷正背离导线运动至远处,空间电荷分布如图9所示。

在当前3维电荷分布下,按前述方法计算图9所在截面上的空间场强分布如图10所示,鉴于空间电荷集中在子导线附近一定范围内,只展示分裂导

图8 正极性电压峰值时3维电荷分布 Fig.8 Three-dimensional charge distribution in positive polarity voltage peak

线中心2 m×2 m范围内的场强分布情况。

从图10可知,10 kV/cm以上的强电场区域集中在导线表面外侧半径0.1 m范围内,此区域内最接近模拟电荷且不受屏蔽作用影响。空间电荷主要分布于导线外半径0.1~0.5 m的环形区域内,该区域内的场强受空间电荷影响,由内至外从10 kV/cm减小至4 kV/cm;由于屏蔽效应,分裂子导线内侧区域内场强较小,子导线内表面往内0.1 m距离空间电场就已经<2 kV/cm,这与文献[17]中的规律

相同。

3.3 与ANSYS仿真结果对比分析

3.3.1 计算误差分析

由于对分裂子导线进行了电荷离散化处理,所以需要验证其计算误差。有限元仿真软件ANSYS以Maxwell方程为基础,在导线表面施加电压载荷,经过有限元剖分、求解和电位梯度提取就能准确计算静电场域内的空间电场。文献[18]在ANSYS中搭建了小电晕笼仿真模型,计算了不考虑电晕时导线表面的电场分布,参考该方法在ANSYS中搭建了特高压电晕笼仿真模型,在导线表面施加不同电压,求解导线表面附近的电场域,并取笼体中部1号子导线表面的最大场强,与本文中无残余空间电荷存在时的结果对比,如图11所示。

图11中,未起晕时本文计算结果与仿真结果一致性较高,与电压呈现出线性增长的关系,并在电压574 kV相差最大,此时误差来源主要为计算误差。笼内导线起晕后,子导线表面发生电晕处的模拟电荷减少至起晕电荷,新发射的空间正电荷对导线表面场强有一定的抑制,且随电压升高电晕范围

图9 正极性电压峰值时第5段导线中点径向电荷分布 Fig.9 Charge distribution of the 5th conductor center’s section in positive polarity voltage peak

图10 笼体中部横截面场强分布云图 Fig.10 Contours of electric field intensity distribution in the central cage section

图11 不同电压下模拟电荷法与ANSYS计算笼体中部导线表面场强比较 Fig.11 Comparison of field strength on central wires’ surface between charge simulation method and ANSYS at different voltages

扩大,模拟电荷法计算场强逐渐趋于平稳,略高于起晕场强;ANSYS仿真模型无法考虑导线起晕和电荷发射,计算所得场强与电压保持线性关系,故导线起晕后其值高于模拟电荷法。

3.3.2 沿分裂子导线中心径向场强对比分析

计算结果稳定后,当电压达到正极性峰值时,选取笼体15.75 m处截面上1号分裂子导线中心所在径向,计算此路径上的电场强度分布,并与施加同等电压下的ANSYS仿真结果对比,如图12所示。

由图12可见,分裂导线中心至子导线内表面点间,因屏蔽效应影响,此路径上的场强低于起晕场强,导线内表面点场强最大但未起晕,场强

>1 kV/cm的区域仅存在于子导线内表面点外0.2 m范围内,且越靠近子导线内表面场强随距离变化的梯度越大,曲线光滑说明此范围内没有空间电荷存在;距分裂导线中心0~2 m的路径上,最大场强出现在分裂子导线外表面点,此点外1 m范围内模拟电荷法计算所得场强受空间电荷影响出现波动,正极性空间电荷靠近导线侧的场强有所减小,远离导线侧的场强有所增大,负极性空间电荷则反之;而距分裂导线中心1.5 m以外的区间上,空间电荷影响较小,2种方法计算所得场强接近,相差8.02%。

3.3.3 与考虑均压环仿真结果对比分析

考虑到实际笼体两侧均压环的影响,另在ANSYS计算模型中防护笼外侧3 m处增加直径为1.5 m的均压环,重新计算无空间电荷且电压为400 kV时沿1号子导线表面最大场强,与本文所得结果进行比较,如图13所示。

图13中不带均压环情况下ANSYS计算所得沿线场强分布与本文结果相符合,测量笼内导线表面最大电场强度相差≤0.5 kV/cm,防护笼内相差≤1.5 kV/cm,其端部5 m往外方向上导线表面最大场强逐渐增大,这与带均压环的情况所有不同。考虑均压环影响后,防护笼内导线表面电场强度在端部5 m内先减小至10 kV/cm以下,仅在两侧0.5 m内增大至22 kV/cm,可见均压环改善了端部效应,使导线表面起晕点集中在测量笼内。电晕笼试验中均压环的存在改变了沿导线产生的电晕电流分布,后续计算中有必要考虑其影响。

4 结论

1）采用模拟电荷法建立了交流电晕笼电场强度3维计算模型。本模型将每根分裂子导线分段后

图12 模拟电荷法与ANSYS计算场强沿径向对比 Fig.12 Comparison of field strength between charge simulation method and ANSYS radially

图13 模拟电荷法与ANSYS带均压环、不带均压环计算沿线场强对比 Fig.13 Comparison of field strength along the line between charge simulation method, ANSYS with and without grading ring

替换为8根有限长模拟线电荷,按起晕电荷分别判定导线表面各点起晕,并在3维空间内计及空间电荷的发射、迁移和损耗,更接近物理实际。

2）按所建模型计算了特高压电晕笼8×LGJ720导线的离子流场。结果表明,沿线方向上端部效应增大了防护笼内导线表面场强,使端部2 m内导线比中部导线提前起晕;在径向方向上,场强>10 kV的区域主要集中在子导线表面外侧0.1 m范围内,受空间电荷影响,其外侧1 m内的空间场强有-5~7 kV/cm的波动。

3）与电晕笼仿真计算结果对比,误差分析表明本文模型与仿真计算结果相差<1.5 kV/cm,但有必要考虑均压环对电晕笼导线附近场强分布的影响。本文结果为下一步电晕笼导线电晕损失的计算打下了良好的基础。

参考文献

[1] 舒印彪,张文亮. 特高压输电若干关键技术研究[J]. 中国电机工程学报,2007,27(31):1-6. SHU Yinbiao, ZHANG Wenliang.Research of key technologies for UHV transmission[J]. Proceedings of the CSEE, 2007, 27(31): 1-6.

[2] 谢辉春,崔翔,刘华钢,等. 特高压同塔双回交流输电线路采用(8+2)分裂导线与常规导线电晕特性对比[J]. 高电压技术,2016,42(3):966-972. XIE Huichun, CUI Xiang, LIU Huagang, et al.Corona characteristic comparison of (8+2) bundled wire and normal wire applying in UHV AC double-circuit transmission line on the same tower[J]. High Voltage Engineering, 2016, 42(3): 966-972.

[3] 杨兰均,王维,林岑,等. 电晕放电离子风实验与理论研究进展及应用发展前景[J]. 高电压技术,2016,42(4):1100-1108. YANG Lanjun, WANG Wei, LIN Cen, et al.Experimental and theoretical research progress in ionic wind produced by corona discharge and its application[J]. High Voltage Engineering, 2016, 42(4): 1100-1108.

[4] 裴春明,陈豫朝,陈小月,等. 交流电晕笼无线电干扰试验校准系数研究[J]. 中国电机工程学报,2015,35(24):6566-6572. PEI Chunming, CHEN Yuchao, CHEN Xiaoyue, et al.Research on measurement calibration coefficient of radio interference with AC corona cages[J]. Proceedings of the CSEE, 2015, 35(24): 6566-6572.

[5] 刘元庆,郭剑,陆家榆. 特高压直流导线在单双极电晕笼中的可听噪声测量与分析[J]. 电网技术,2012,36(5):1-5. LIU Yuanqing, GUO Jian, LU Jiayu.Measurements and analysis on UHV DC conductor’s audible noise in unipolar and bipolar corona cages[J]. Power System Technology, 2012, 36(5): 1-5.

[6] COMBER M G, GARBERRY R E, CHARTIER V L, et al.A comparison of methods for calculating audible noise of high voltage transmission lines[J]. IEEE Transactions on Power Apparatus and System, 1982, 101(10): 4090-4099.

[7] 谭震宇,杨金洪,徐明铭,等. 雾霾对高压直流输电线路电晕离子流场的影响[J]. 高电压技术,2016,42(12):244-249. TAN Zhenyu, YANG Jinhong, XU Mingming, et al.Influence of fog-haze on corona ion flow field of HVDC transmission lines[J]. High Voltage Engineering, 2016, 42(12): 244-249.

[8] CLADE J J, GARY C H, LEFEVRE C A.Calculation of corona losses beyond the critical gradient in alternating voltage[J]. IEEE Transactions on Power Apparatus and Systems, 1969, 88(5): 695-703.

[9] CLADE J J, GARY C H.Predentermination of corona losses under rain: influence of rain intensity and utilization of a universal chart[J]. IEEE Transactions on Power Apparatus and Systems, 1970, 89(6): 1179-1185.

[10] ABDEL-SALAM M, SHAMLOUL D.Computation of ion-flow fields of AC coronating wires by charge simulation techniques[J]. IEEE Transactions on Electrical Insulation, 1992, 27(2): 352-361.

[11] ABDEL-SALAM M, ABDEL-AZIZ E Z. A charge-simulation-based method for calculating corona loss on AC power transmission lines[J]. Journal of Physics D: Applied Physics, 1994(27): 2570-2579.

[12] 张宇,魏远航,阮江军. 高压直流单极离子流场的有限元迭代计算[J]. 中国电机工程学报,2006,26(23):158-162. ZHANG Yu, WEI Yuanhang, RUAN Jiangjun.Finite element iterative computation of high voltage direct current unipolar ionized field[J]. Proceedings of the CSEE, 2006, 26(23): 158-162.

[13] 杨文翰,吕英华. 用模拟电荷法求解高压输电线路附近电磁场[J]. 电网技术,2008,32(2):47-50. YANG Wenhan, LÜ Yinghua.Application of emulation charge method in calculation of electromagnetic environment near to HV transmission lines[J]. Power System Technology, 2008, 32(2): 47-50.

[14] 黎振宇,蒋兴良,李立浧,等. 涂覆导线表面的淋雨电晕特性[J]. 高电压技术,2017,43(11):3740-3747. LI Zhenyu, JIANG Xingliang, LI Licheng, et al.Corona performance of coated conductors under raining condition[J]. High Voltage Engineering, 2017, 43(11): 3740-3747.

[15] SARMA M P.Corona performance of high-voltage transmission lines[R]. New York, USA: Research Studies Press Company Limited, 2000: 77-83.

[16] 李威. 特高压电晕笼交流离子流场与电晕损失数值仿真研究[D]. 武汉:武汉大学,2017:26-29. LI Wei.Numerical simulation research of AC ion flow field and corona loss in UHV corona cage[J]. Wuhan, China: Wuhan University, 2017: 26-29.

[17] 肖国洲,蓝磊,陈小月,等. 特高压电晕笼空间电场计算分析[J]. 水电能源科学,2015,33(10):182-185,140. XIAO Guozhou, LAN Lei, CHEN Xiaoyue, et al.Calculation and analysis of space electric field of UHV corona cage[J]. Water Resources and Power,2015, 33(10): 182-185, 140.

[18] 徐渝淇. 电晕笼内导线三维电场强度分布浅析[J]. 四川电力技术,2015,38(5):33-37. XU Yuqi.Calculation and analysis of three dimensional electrical field distribution in corona cage[J]. Sichuan Electrical Technology, 2015, 38(5): 33-37.