永磁直流电机定子的振动分析
曹海翔 1 , 沈靖文 1 , 王善铭 1 , 牟佳男 2 , 杨占录 1 , 洪剑锋 1     
1. 清华大学 电机工程与应用电子系, 电力系统及发电设备控制和仿真国家重点实验室, 北京 100084;
2. 国家电网公司 国家电力调度控制中心, 北京 100031
摘要:永磁直流电机的定子振动具有鲜明的特性,对它的分析对于分析其他各类电机的振动具有借鉴作用。该文用解析法分析了电机内的磁场和径向电磁激振力,从物理概念上定性解释了直流电机振动的特征。特别地,指出永磁直流电机槽频电磁激振力的阶数与定子振动形态的空间阶数不同。为准确计算电机振动,利用电磁场有限元仿真了电机内的磁场与电磁力分布,并采用结构场有限元对定子振动进行了仿真。针对试验测试,指出为得到定子振动形态特征,需要在定子上测量足够多的测点。并在一台实际的直流电机上进行振动测试,试验结果与仿真结果相符。
关键词永磁直流电机     电磁激振力     振动     振动形态    
Analyses of permanent magnet DC motor stator vibrations
CAO Haixiang1, SHEN Jingwen1, WANG Shanming1, MU Jianan2, YANG Zhanlu1, HONG Jianfeng1     
1. State Key Laboratory of Control and Simulation of Power Systems and Generation Equipments, Department of Electrical Engineering, Tsinghua University, Beijing 100084, China;
2. National Electric Power Dispatching and Control Center, State Grid Corporation of China, Beijing 100031, China
Abstract:Stator vibrations in permanent magnet DC motors have distinct characteristics which are also useful for analyzing other types of motors. The magnetic fields and the electromagnetic forces in permanent magnet DC motors were analyzed to characterize the physical mechanisms driving the vibrations. The order of the slot-frequency electromagnetic exciting force is not equal to the spatial order of the stator slot-frequency vibration shape. Accurate analyses of the vibrations inputs the magnetic field and electromagnetic force distributions calculated by finite element simulations into a structural finite element simulation to predict the stator vibrations. Experiments must have enough sensors to obtain the vibration mode shape. The predictions agree well with data from DC motor tests.
Key words: permanent magnet DC motor     electromagnetic exciting force     vibration     mode of vibration    

随着技术的进步和人们对环境要求的不断提高,各种应用场合对电机的振动和噪声提出了越来越高的要求。电机气隙磁场在定转子上产生的电磁激振力是造成电磁振动和噪声的根源[1]。若要从源头上降低电机振动噪声,就要对电机磁场、电磁激振力和振动进行透彻分析,掌握其特征。

永磁直流电机由于其控制简单、转矩特性好的特点在很多场合中得到广泛应用。直流电机主磁场固定、相对定子不运动的特点使得直流电机定子受到的电磁激振力中最主要的分量为频率与槽数相关的槽频电磁激振力,而不存在频率与磁极数相关的极频激振力。这一特点使得永磁直流电机的振动分析更具针对性,符合从简单入手对问题进行分析的思路。因此本文对永磁直流电机的磁场、电磁激振力和定子振动进行分析,以得到有益的结论,进而指导各类电机的振动分析。

对永磁电机振动的已有研究主要集中于永磁同步电机以及无刷直流电机。文[2]对正弦以及变频器供电的永磁同步电机振动机理进行分析,并测量了振动频谱,与理论结果相对比。文[3]对分数槽永磁同步电机电磁振动进行分析,提出一种注入补偿电流减小电磁振动的抑制方法。文[4]首先采用基于二维有限元方法计算了永磁同步电机定子齿所受电磁力。然后将定子设为单环圆柱体,建立了预测定子齿径向位移的解析模型,并通过结构有限元分析和试验数据验证了解析模型。文[5]提出了分数槽无刷永磁电机的径向激振力的解析模型,利用有限元法进行验证。并利用解析模型研究了定子开槽、磁场切向分量、气隙半径等因素对径向力计算的影响。文[6]分析了无刷直流电机径向不平衡电磁力的产生原因,定量研究了电流换相时的径向不平衡电磁力。

对永磁直流电机振动的已有研究相对较少。文[7]对直流电机的振动进行的分析表明: 引起直流电机振动的主要原因是定子主极磁场与转子一阶齿谐波的相互作用。该文对作用在磁极上的交变径向力以及交变弯曲力矩进行了推导。文[8]建立了永磁直流电机的磁场—结构场—声场耦合的数值模型,对于直流电机的电磁力、振动形变以及噪声进行了计算。文[9]建立了永磁直流电机的有限元和边界元模型,通过数值计算来预测磁场、电磁力、模态参数、动态响应以及声压,并且研究了转子偏心以及不同胶水厚度对振动和噪声的影响。文[10]利用解析法建立了永磁直流电机的气隙磁场的模型,得到瞬时磁场表达式,并对电机的径向力以及切向力进行分析,得到直流电机的电磁力频率为转频与槽数乘积的整数倍的结论。文[11]分析了一台2极3槽直流电机的振动产生原因,并提出了通过在转子齿上开槽的方法减小不平衡力。对于永磁直流电机振动的研究,上述文献尚在以下2方面存在不足: 一是未从空间的角度对电磁激振力的阶数和电机振动形态加以阐述; 二是未从试验角度给出电机整体测试结果。以往的试验经常只测量一个或少数测点的振动加速度,不能体现完整的振动形态。

本文首先利用解析法分析了永磁直流电机的磁场与电磁力,从物理概念上阐述了直流电机振动的特征,指出其槽频电磁激振力的阶数与定子振动形态概念的差别。然后采用电磁场有限元计算了电机内的磁场与电磁力分布,利用结构场有限元对定子振动进行计算。最后对样机进行试验,对定子圆周不同位置处的振动同时进行测量,获取完整的振动形态的信息,验证了解析与仿真分析的结果。

1 永磁直流电机的磁场与电磁力 1.1 磁场分析

电枢反应会引起磁场变化,但是数值上比永磁体激励产生的磁场小很多[10],因此本文为简化分析将只考虑永磁体产生的磁场,而不考虑电枢电流的作用。图 1为永磁直流电机的截面图,以一对极的电机为例进行说明,在气隙中心的圆周上建立圆柱坐标系,以定子磁极中心线为原点,逆时针方向为正方向,磁动势以磁力线出定子进转子为正方向。永磁体产生的磁动势可表示为

$f\left( \theta ,t \right)=\left\{ \begin{array}{*{35}{l}} +{{H}_{\text{c}}}{{h}_{\text{m}}}, & \theta \in \left( -\frac{\alpha \pi }{2}, \right.\left. \frac{\alpha \pi }{2} \right]; \\ 0, & \theta \in \left( \frac{\alpha \pi }{2}, \right.\pi -\left. \frac{\alpha \pi }{2} \right]; \\ -+{{H}_{\text{c}}}{{h}_{\text{m}}}, & \theta \in \left( \pi -\frac{\alpha \pi }{2}, \right.\pi +\left. \frac{\alpha \pi }{2} \right]; \\ 0, & \theta \in \left( \pi +\frac{\alpha \pi }{2}, \right.2\pi -\left. \frac{\alpha \pi }{2} \right]. \\ \end{array} \right.$ (1)
图 1 永磁直流电机的截面图

其中: Hc为永磁体的矫顽力,hm为永磁体厚度,α为永磁体的极弧系数。该磁动势为一方波,每一对极下正负变化一次。

以转子齿中心与定子磁极中心线对齐时刻为计时起点,则在定子坐标系下的转子侧开槽时的气隙磁导可表示为

$\lambda \left( \theta ,t \right)=\frac{{{\mu }_{0}}}{\delta +{{h}_{\text{m}}}}\left( {{{\bar{\Lambda }}}_{0}}+\sum\limits_{k=1}^{\infty }{{{{\bar{\Lambda }}}_{k}}\cos kZ\left( \theta -{{\omega }_{\text{r}}}t \right)} \right).$ (2)

其中: μ0为真空磁导率; δ为气隙长度; ωr为转子旋转角速度; Z为转子槽数;${{{\bar{\Lambda }}}_{0}}$和${{{\bar{\Lambda }}}_{k}}$为平均磁导系数和第k次谐波磁导系数[8],${{{\bar{\Lambda }}}_{0}}\text{=}\frac{1}{{{K}_{\text{c}}}},{{{\bar{\Lambda }}}_{k}}=\left( {{K}_{\text{c}}}-1 \right)\left| \frac{\sin k\frac{{{K}_{\text{c}}}-1}{{{K}_{\text{c}}}}\pi }{k\frac{{{K}_{\text{c}}}-1}{{{K}_{\text{c}}}}\pi } \right|$,Kc为Carter系数。磁导系数中k=1时最大,后面各次系数较小,一般只考虑k=1。

永磁直流电机的空载气隙磁场由永磁磁动势与转子开槽的磁导作用产生。在忽略饱和时,气隙磁场表示为

$b\left( \theta ,t \right)=f\left( \theta ,t \right)\lambda \left( \theta ,t \right).$ (3)
$b\left( \theta ,t \right)=\left\{ \begin{array}{*{35}{l}} {{B}_{\delta }}{{{\bar{\Lambda }}}_{0}}+\sum\limits_{k=1}^{\infty }{{{B}_{\delta }}{{{\bar{\Lambda }}}_{k}}\cos \left( kZ\theta -kZ{{\omega }_{\text{r}}}t \right)} & \theta \in \left( -\frac{\alpha \pi }{2}, \right.\left. \frac{\alpha \pi }{2} \right]; \\ 0, & \theta \in \left( \frac{\alpha \pi }{2}, \right.\pi -\left. \frac{\alpha \pi }{2} \right]; \\ -{{B}_{\delta }}{{{\bar{\Lambda }}}_{0}}+\sum\limits_{k=1}^{\infty }{{{B}_{\delta }}{{{\bar{\Lambda }}}_{k}}\cos \left( kZ\theta -kZ{{\omega }_{\text{r}}}t \right)} & \theta \in \left( \pi -\frac{\alpha \pi }{2}, \right.\left. \pi +\frac{\alpha \pi }{2} \right]; \\ 0, & \theta \in \left( \pi +\frac{\alpha \pi }{2}, \right.\left. 2\pi -\frac{\alpha \pi }{2} \right]. \\ \end{array} \right.$ (4)

其中${{B}_{\delta }}=\frac{{{H}_{\text{c}}}{{h}_{\text{m}}}{{\mu }_{0}}}{\delta +{{h}_{\text{m}}}}$。 式(4)中不为零的2个表达式均由2项组成: 第1项表示转子开槽使平均磁阻增加时,由永磁体产生的主磁场,该项不随时间变化; 第2项表示由于转子开槽,气隙磁导发生周期性变化而导致主磁场变化的附加磁场分量。

1.2 电磁激振力分析

铁磁物质在电磁场中的受力可用Maxwell应力张量法进行计算。在一个磁极下的气隙磁场产生的电磁激振力密度为

$\begin{align} & {{p}_{n}}\left( \theta ,t \right)=\frac{{{b}^{2}}\left( \theta ,t \right)}{2{{\mu }_{0}}}= \\ & \frac{1}{2{{\mu }_{0}}}{{\left[ \pm {{B}_{\delta }}_{0}{{{\bar{\Lambda }}}_{0}}+\sum\limits_{k=1}^{\infty }{{{B}_{\delta }}{{{\bar{\Lambda }}}_{k}}\cos \left( kZ\theta -kZ{{\omega }_{\text{r}}}t \right)} \right]}^{2}}\text{=} \\ & \frac{B_{\delta }^{2}}{2{{\mu }_{0}}}\sum\limits_{i}{\sum\limits_{j}{\left[ {{{\bar{\Lambda }}}_{i}}{{{\bar{\Lambda }}}_{j}}\cos \left( iZ\theta -iZ{{\omega }_{r}}t \right)\cos \left( jZ\theta -jZ{{\omega }_{\text{r}}}t \right) \right]}}. \\ \end{align}$ (5)

式(5)中共有3项: 第1项表示永磁体产生的,不随时间变化的磁场自身作用得到的径向力密度,为恒定值; 第2项表示不随时间变化的磁场与周期性变化的附加磁场相互作用产生的,随时间变化的径向力密度; 第3项表示气隙磁导周期性变化产生的附加磁场之间作用产生的,随时间变化的径向力密度。由各项的系数可见,第3项的幅值由两次谐波磁导相乘得到,数值很小可以忽略。上述解析法定性分析了电磁激振力的主要成分以及来源,下面利用有限元进行定量计算。

1.3 电磁场有限元仿真结果

利用Maxwell 2D有限元软件计算电机内磁场。样机为2极24槽永磁直流电机,极弧系数为0.744。定子外径118 mm,内径108 mm,长度 195 mm; 永磁体内径79 mm,长度100 mm; 转子外径76 mm,转子齿根直径46 mm,转轴直径19 mm。

仿真工况的转子转速为3 000 r/min,计算得到的定子永磁体内表面上的径向磁通密度以及径向电磁力密度分布分别如图 23所示。在空间上,由式(5)中第2项所述,电磁力密度分布在360° 范围内有kZ个周期,k取1时即为槽数。由图 3中可知,若加上两磁极之间的区域,电磁力密度空间分布为24个周期。在时间上,转子转过一周的过程中,永磁体表面中心点的磁通密度波动24次(见图 4),即磁场和电磁激振力的频率均为电机转动频率与转子齿数Z的乘积的整数倍。

图 2 永磁体表面径向磁通密度分布

图 3 永磁体表面圆周上径向力密度分布

图 4 永磁体表面中心点磁通密度随时间变化

2 定子电磁振动的分析与仿真 2.1 振动频率分析

1) 极频振动。

由式(5)的第1项可见,磁导为平均磁导,数值比谐波磁导大,因此第1项幅值最大,但是该项不随时间变化,不引起振动。这完全不同于普通永磁电机,后者的振动频率为极数乘以转频,而直流电机中不存在该分量。

2) 槽频振动。

由式(5)的第2项可见,该力密度分量存在2个主要特征: 频率和阶数。该分量系数较大,对直流电机的振动影响最大,其频率为转子槽数与转子转速的乘积的整数倍,称为槽频:

$f=\frac{kZ{{\omega }_{\text{r}}}t}{2\pi },\text{ }k\in {{\mathbb{N}}_{+}}.$ (6)

该分量激振力的阶数,即为力密度的空间分布的周期数,为kZ阶。若直流电机每极下槽数较多,则该阶数较高,由振动理论可知,振动变形与激振力的阶数的4次方成反比,即阶数越高振动变形越小[12]。而实际试验中,槽频振动变形最为突出,该分析与试验结果存在较大矛盾。

2.2 槽频振动的振动形态分析

转子开槽会引起定子磁极表面磁场变化,而且该磁场随着转子转动而改变,引起槽频振动。该变化可采用凸极同步电机齿谐波电动势分析中的磁场纵振荡和横振荡进行定性分析[13]。由于磁通在空间分布不均匀地扫动而产生的齿谐波电势为横振荡; 由于磁通幅值变化而产生的齿谐波电势为纵振荡。

为了便于阐述清楚,以每个磁极下的转子槽数为整数和非整数这2种情况为代表进行说明,图 5为每极下为2.5个槽时的情况,在图 5a5b的2个转子位置时,磁极下的齿数分别为3和2,因此磁极下的磁通也随之发生变化,这种变化将使磁极表面的电磁激振力密度改变,引起振动。该磁场变化由磁场纵振荡引起,若每个极下情况相同,则该槽频振动的振动形态为椭圆形,固定节点数与磁极数相同。在图 5a5b的2个典型转子位置处,磁极受到的电磁激振力关于磁极中心线对称。

图 5 每个磁极下的转子槽数为2.5时的磁极受力

图 6为每极下槽数为2时的情况,转子在不同位置时,磁极下的齿数不变,不发生磁场的纵振荡。但是,由于齿相对于磁极位置会变化,导致磁极受到的电磁激振力一会偏在磁极的右侧,一会偏在磁极的左侧,对比图 3也可以看出每个磁极下的径向力分布情况,在齿尖受力较大而槽低受力较小。这种激振力位置相对磁极轴线的偏移将使磁极受到左右不平衡力的作用,再考虑到磁极之间的磁轭比磁极更薄,这种不平衡力引起的变形将更大。图 7a7b分别为两变形最大时刻的示意图。可见,该情况下振动的振动形态的固定节点数为2,与磁极数相同。需注意,与图 5不同的是,此时磁极受到的电磁激振力不再关于磁极中心线对称。

图 6 每个磁极下的转子槽数为2时的磁极受力

图 7 磁极受不平衡力产生的定子的变形

由上述2种情况可见,直流电机的振动形态与磁极数相关,而不是与式(5)中的槽数相关,这一点是直流电机振动分析中的重要特征,也是降低振动的指导方向。从频率上看,这种振动的频率依然为与槽数相关的槽频。当磁极下槽数不同、每极下情况不同时,也可以类似进行分析,得到振动形态特征。一般来说,电磁激振力波的空间阶数是与槽数相关,高阶数对应的振动模态频率非常高,不容易激发。而由于磁场与定子相对位置变化产生的纵振荡和横振荡易于产生低阶数的振荡,使得电磁激振力的阶数与定子振动形态的阶数不一致。

2.3 振动仿真

大量研究表明,电机的电磁振动主要是由定子振动引起的,转子的影响较小可以忽略。定子和永磁体的三维简化模型及网格划分模型如图 8所示,各部分材料的属性如表 1所示。

图 8 永磁直流电机定子三维简化模型及网格剖分

表 1 材料属性
材料密度/(kg·m-3)弹性模量/GPaPoisson比
结构钢7 8502000.3
永磁体4 6001800.28

利用Ansys Workbench有限元软件进行电磁场与结构场的耦合分析,将Maxwell 2D有限元软件中计算得到的转子处于不同位置的电磁力密度以谐波方式传递到结构场中,进行振动计算。利用谐响应分析求解振动。

在电机额定转速3 000 r/min下进行仿真,利用谐波加载方式,对Maxwell 2D有限元软件中得到的电磁力进行Fourier分解,在频域中分为实部力与虚部力,加载在永磁体上。仿真时设置定子圆环端面为弹性支承,模拟在试验过程中将电机吊装的试验条件。选择模态叠加法进行计算,得到图 9的槽频1 200 Hz的定子整体振动位移,图 9a9b分别为振动幅度最大的2个时刻,箭头的方向和长度分别表示位移的方向和大小。可以看出由于永磁直流电机磁钢刚度较大,电机振动形变主要发生在磁极间的机壳。定子的振动形态为椭圆形,2个磁极中心为振动位移最小的节点,节点数目与磁极数相同,与图 7相符合。

图 9 1 200 Hz的槽频定子振动

3 试验验证

对节1.3中的2极24槽的永磁直流样机进行了试验,采用悬挂安装的方式,利用Siemens LMS振动测试系统,测量电机定子机壳上的径向振动加速度。

在定子圆周360° 范围内均匀布置15个径向加速度传感器,随后在定子圆周180° 范围内均匀布置14个传感器,测量并比较振动加速度的幅值和相位。定子上的传感器位置和电机的悬吊方式如图 10所示。对于整个圆周,以磁极正中心点为起始点a1,按逆时针方向依次编号直至点a15。对于半个圆周,将点b1与点a12重合,按逆时针方向依次编号至点b14。测量得到机壳各个位置处的振动的频谱图。为便于显示,取具有代表性的点a7的振动加速度频谱图如图 11所示。在2 500 r/min转速下,对应槽频为

图 10 传感器位置分布

图 11 点a7的振动加速度频谱

$\begin{align} & f=k\cdot \frac{24\times 2\text{ }500}{60~s}=1\text{ }000k~Hz, \\ & k\in {{\mathbb{N}}_{+}}. \\ \end{align}$ (7)

图 11可知,槽频振动远大于其他频率下的振动,其次是2倍槽频下的振动,这与解析分析相一致。但是该频谱图不能显示不同位置处振动的相位信息。因此将各传感器相位与幅值综合,得到定子圆周各测点处的振动加速度变化情况。2个槽频变形最大时刻的各点处的振动加速度如图 12所示,其中虚线圆表示静止时的加速度传感器位置。可以看出各位置处的振动加速度具有相位差,并非同时达到幅值,此时槽频振动的振动形态为以两磁极中心点为不动节点的椭圆形,与仿真分析所得到的结果相一致。在180° 范围内布置的14个径向传感器的测量结果见图 13,结论得到进一步验证。

图 12 360° 范围内测量结果

图 13 180° 范围内测量结果

4 结 论

本文采用解析法定性分析了永磁直流电机的磁场以及电磁力,完整地建立了电机电磁场、结构场的计算模型,并得到了永磁直流电机的电磁力和振动的特征,阐明了振动分量中最为明显的槽频分量的来源并区分了电磁激振力的阶数和电机振动形态的阶数。之后通过样机试验进行了振动测试,采用了利用多个传感器同时测量的相位幅值信息得到定子振动形态的测试方法,使振动测试更为准确。得到的试验结果与解析分析以及仿真结果相符合,证明了分析的正确性,对直流电机进一步减振降噪具有指导意义。

参考文献
[1] 朱海峰. 异步电机电磁激振力分析[D]. 杭州:浙江大学, 2013. ZHU Haifeng. Analysis of Electromagnetic Exciting Force of Asynchronous Motors[D]. Hangzhou:Zhejiang University, 2013. (in Chinese)
[2] 陈秋明, 陈勇. 永磁同步电动机电磁振动噪声机理研究[J]. 微特电机, 2013, 41(8): 1–5. CHEN Qiuming, CHEN Yong. Study on theory of electro-magnetic vibration and noise of a permanent magnet synchronous motor[J]. Small & Special Electrical Machines, 2013, 41(8): 1–5. (in Chinese)
[3] 杨浩东, 陈阳生. 分数槽永磁同步电机电磁振动的分析与抑制[J]. 中国电机工程学报, 2011, 24: 83–89. YANG Haodong, CHEN Yangsheng. Electromagnetic vibration analysis and suppression of permanent magnet synchronous motor with fractional slot combination[J]. Proceedings of the CSEE, 2011, 24: 83–89. (in Chinese)
[4] Islam R, Husain I. Analytical model for predicting noise and vibration in permanent-magnet synchronous motors[J]. IEEE Transactions on Industry Applications, 2010, 46(6): 2346–1354. DOI:10.1109/TIA.2010.2070473
[5] Zhu Z Q, Xia Z P, Wu L J, et al. Analytical modeling and finite-element computation of radial vibration force in fractional-slot permanent-magnet brushless machines[J]. IEEE Transactions on Industry Applications, 2010, 46(5): 1908–1918. DOI:10.1109/TIA.2010.2058078
[6] 陈磊, 高宏伟, 柴凤, 等. 小型无刷直流电动机振动与噪声的研究[J]. 中国电机工程学报, 2006, 24: 148–152. CHEN Lei, GAO Hongwei, CHAI Feng, et al. Research on the vibration and noise of the minitype brushless DC motor[J]. Proceedings of the CSEE, 2006, 24: 148–152. (in Chinese)
[7] 陈永校, 诸自强, 应善成. 电机噪声的分析与控制[M]. 杭州: 浙江大学出版社, 1987. CHEN Yongxiao, ZHU Ziqiang, YING Shancheng. Analysis and Control of Motor Noise[M]. Hangzhou: Press of Zhejiang University, 1987. (in Chinese)
[8] Furlan M, Cernigoj A, Boltezar M. A coupled electromagnetic-mechanical-acoustic model of a DC electric motor[J]. The International Journal for Computation and Mathematics in Electrical and Electronic Engineering, 2003, 22(4): 1155–1165. DOI:10.1108/03321640310483075
[9] He G, Huang Z, Qin R, et al. Numerical prediction of electromagnetic vibration and noise of permanent-magnet direct current commutator motors with rotor eccentricities and glue effects[J]. IEEE Transactions on Magnetics, 2012, 48(5): 1924–1931. DOI:10.1109/TMAG.2011.2178100
[10] He G, Huang Z, Chen D. Two-dimensional field analysis on electromagnetic vibration-and-noise sources in permanent-magnet direct current commutator motors[J]. IEEE Transactions on Magnetics, 2011, 47(4): 787–794. DOI:10.1109/TMAG.2010.2103382
[11] Li Y B, Ho S L, Fu W N, et al. Analysis and solution on squeak noise of small permanent-magnet DC brush motors in variable speed applications[J]. IEEE Transactions on Magnetics, 2009, 45(10): 4752–4755. DOI:10.1109/TMAG.2009.2024882
[12] 舒波夫. 电机的噪声和振动[M]. 北京: 机械工业出版社, 1980. Shupov. Motor Noise and Vibration[M]. Beijing: China Machine Press, 1980. (in Chinese)
[13] 章名涛. 电机学[M]. 北京: 科学出版社, 1964. ZHANG Mingtao. Electric Machinery[M]. Beijing: Science Press, 1964. (in Chinese)