CN1163435A - 利用矩量法模拟电磁场强度的设备和方法 - Google Patents

利用矩量法模拟电磁场强度的设备和方法 Download PDF

Info

Publication number
CN1163435A
CN1163435A CN97103317A CN97103317A CN1163435A CN 1163435 A CN1163435 A CN 1163435A CN 97103317 A CN97103317 A CN 97103317A CN 97103317 A CN97103317 A CN 97103317A CN 1163435 A CN1163435 A CN 1163435A
Authority
CN
China
Prior art keywords
frequency
transimpedance
moment
frequency spectrum
calculating
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN97103317A
Other languages
English (en)
Other versions
CN1118766C (zh
Inventor
大津信一
向井诚
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Fujitsu Ltd
Original Assignee
Fujitsu Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Fujitsu Ltd filed Critical Fujitsu Ltd
Publication of CN1163435A publication Critical patent/CN1163435A/zh
Application granted granted Critical
Publication of CN1118766C publication Critical patent/CN1118766C/zh
Anticipated expiration legal-status Critical
Expired - Fee Related legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G05CONTROLLING; REGULATING
    • G05BCONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B17/00Systems involving the use of models or simulators of said systems
    • G05B17/02Systems involving the use of models or simulators of said systems electric
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01RMEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R31/00Arrangements for testing electric properties; Arrangements for locating electric faults; Arrangements for electrical testing characterised by what is being tested not provided for elsewhere
    • G01R31/001Measuring interference from external sources to, or emission from, the device under test, e.g. EMC, EMI, EMP or ESD testing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/17Function evaluation by approximation methods, e.g. inter- or extrapolation, smoothing, least mean square method
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • G06F30/23Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00Details relating to CAD techniques
    • G06F2111/10Numerical modelling

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Engineering & Computer Science (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Physics (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Algebra (AREA)
  • Electromagnetism (AREA)
  • Automation & Control Theory (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明的设备包括:通过将波源的时间序列数据分段再加以付立叶变换得到频谱的变换装置;计算在采样频率的互阻抗、从算得的互阻抗和采样频率产生互阻抗的近似表达式、再用产生的近似表达式计算在变换装置得到的各频率的互阻抗的计算装置;按照矩量法从计算装置算得的互阻抗和变换装置得到的频谱求出在每个元流动的电源的频谱、从中求出电压、电场和磁场的频谱的模拟装置以及反变换装置。

Description

利用矩量法模拟电磁场 强度的设备和方法
本发明涉及利用矩量法模拟从电子设备辐射的电磁场的强度之类的模拟设备和模拟方法,具体地说,涉及利用矩量法高速模拟时域电磁场强度之类的模拟设备和模拟方法。
社会限制电子设备,不允许电子设备所产生的有害无线电波和噪声超过一定电平,一些国家为此还制定了严格的条例。
为了不违背这类无线电波条例,采取了一系列诸如屏蔽、滤波之类技术的抗辐射措施。在采取这些抗辐射措施时,必需开发能以定量方式模拟这些措施使无线电波或噪声降低的程度的模拟技术。
社会还要求电子设备不受其他电子设备所辐射的超过一定强度的无线电波或噪声的影响,一些国家也为此制字了严格的条例。
为了不违背这些无线电波条例,必需开发能分析电子设备所以辐射有害无线电波或噪声和所以受到无线电波或噪声干扰而不能工作的原因的模拟技术。
为了开发这种模拟技术,就需要能模拟随时间变化的电磁场强度的电磁场强度计算设备。然而,这种模拟随时间变化的电磁场强度的电磁场强度计算设备至今尚未实际投入使用,其原因将在下面加以说明。
任何形状的物体所辐射的电磁场的强度,如果知道在物体每个部位流动的电流和磁流的话,就可以很容易按照经典的理论方程加以计算。这些电流和磁流在理论上可以通过在一定的边界条件下解Maxwell电磁方程得到。
作为解这方程的一种方法,有所谓的矩量法。矩量法是一种解从Maxwell电磁方程导出的积分方程的方法,这种方法是通过将物体细分成一系列很小的元素计算电流和磁流,因此可以处理任何形状的三维物体。就这种矩量法,可参H.N.Wang,J.H.Richmond和M.C.Gilreath所撰写的“导电表面辐射和散射的正弦电抗表示法”(Sinusoidal reactionformulation for radiation and scattering from condrcting surface”,IEEE TRANSACTION  ANTENNAS PROPAGATION.Nol.ap-23,1975)。
但是,为了找出电子设备所以辐射有害无线电波或噪声和所以受到无线电波或噪声干扰而不能工作的原因,必需在时间域进行分析。这是因为使电子设备不能工作的大多数原因是由于脉冲式噪声的干扰。此外,电子设备不能工作往往是由于设备内的器件(如集成电路)工作不正常而造成的。因此,必需长时间连续监视这些器件,确认是否这些器件工作不正常。
两种有用的时域分析法是有限元法和有限差法。虽然能用有限元法和有限差法来进行时域分析,但是要处理含有诸如传输线、电缆、屏蔽罩之类各种构件的电子设备却是困难的。
这是因为有限元法和有限差法中必需对所要分析的物体及其周围的三维空间进行分割。如果对于诸如电缆终端部分之类的小部件进行细分的话,那么由于屏蔽罩和电缆周围的空间相当大,分割出来的元的数量就会过多而使计算机内存溢出。相反,如果对于电缆、屏蔽罩和其他构件进行粗分的话,那么就不能对在机理上起着重要作用、也是有害无线电波或噪声的大辐射源的电缆终端部分的影响进行分析。
此外,在有限元法和有限差法中,分割所用的座标系通常是正交座标系。然而,虽然设备的屏蔽罩可以是任何形状的,但机理上起着重要作用的电缆和电缆终端部分包括柱形件。因此,怎样分割所要分析的物体就是一个棘手的问题。
在这方面,矩量法就没有所有的这类问题,很适合处理包括诸如传输线、电缆、屏蔽罩之类各种构件的电子设备。
这是因为矩量法是一种边界元法,只要求对边界表面进行二维分割。此外,可以相当自由地确定分割节距,能对小的部分进行细分而对电缆和屏蔽罩进行粗分,从而分割出来的元在数量上要大大少于有限元法和有限差法所要求的。而且,由于可以采用任何形状的分割,因此也就不存在如何进行分割的问题。
所以,可以考虑用矩量法来模拟随时间变化的噪声电流、噪声电压和电磁场辐射强度。也就是说,如果给定了一个随时间变化的波源,就可以考虑采用这样一种方法:将这个波源变换到频率域,在所变换的频率域内模拟电磁场辐射强度,然后再将模拟值反变换回时间域。
然而,这种方法用有关技术也不能实现。这是因为矩量法需要对给定的频率计算各格形元之间的互阻抗、互导纳和互电抗,利用计算得的值解联立方程,但计算互阻抗等本身就要很长时间,而从时间域变换到频率域得出的频率的数量又相当大。由于必需分别对各个频率都计算互阻抗等,因此需要大量的处理时间。
也就是说,如果用矩量法计算电磁场辐射强度,正如下面将结合附图所作的说明那样,需要用非常多的时间去计算每一个互阻抗、互导纳和互电抗。此外,由于必需对于频率域内的每个频率都要执行这样的计算,因此在切合实际应用的时间内不可能完成。
具体些说,解矩量法联立方程,对于单一频率来说是几分钟,但计算互阻抗、互导纳和互电抗却要几个小时。由于必需分别对频率域的每个频率都要计算互阻抗等,因此在切合实际应用的时间内不可能完成这样的运算。
需要说明的是,互阻抗表示的是由在一个元中流动的电流感应的电场与在另一个元中流动的电流之间的关系,互导纳表示的是由流过一个元的磁流感应的磁场与流过另一个元的磁流之间的关系,而互电抗表示的是加在一个元上的电流(磁流)感应的磁场(电场)与加在另一个元上的磁流(电流)之间的关系。电流流过金属,而电流和磁流在介电体的表面上流动。
由上面的说明可见,目前还没有开发出能模拟随时间变化的磁场强度的实用电磁场强度计算设备可以用来分析含有印刷板、电缆和屏蔽罩的设备。
因此,考虑到这种情况而作出的本发明的目的是提出一种采用新的矩量法能高速模拟随时间变化的电磁场强度的模拟设备和方法。
为了达到这个目的,本发明所提出的利用新的矩量法能高速模拟时间域电磁场强度的模拟设备包括:一个变换装置,该装置通过分割波源的时间序列数据、进行傅立叶变换得到频谱;一个计算装置,该装置计算在采样频率的互阻抗、从计算出的互阻抗和采样频率产生互阻抗的近似表达式,再分别按变换装置得到的各个频率用所产生的近似表达式计算相应的互阻抗;一个模拟装置,该装置根据矩量法从计算装置计算出的互阻抗和变换装置得到的频谱求出在每个元流动的电流的频谱,再从求得的电流频谱求出相应的电压、电场和磁场的频谱;以及一个反变换装置,该装置对模拟装置求得的频谱进行付立叶反变换,输出反变换结果。
本发明的以上目的和优点从以下结合附图对优选实施例的说明中可以看得更加清楚。在这些附图中:
图1为本发明的配置的原理方框图;
图2为本发明的另一个配置的原理方框图;
图3为产生近似表达式的第一说明图;
图4A和4B为产生近似表达式的第二说明图;
图5A和5B为产生近似表达式的第三说明图;
图6为产生近似表达式的第四说明图;
图7为产生近似表达式的第五说明图;
图8为产生近似表达式的第六说明图;
图9为产生近似表达式的第七说明图;
图10为产生近似表达式的第八说明图;
图11A、11B和11C为产生近似表达式的第九说明图;
图12为产生近似表达式的第十说明图;
图13A和13B为产生近似表达式的第十一说明图;
图14为产生近似表达式的第十二说明图;
图15为产生近似表达式的第十三说明图;
图16为产生近似表达式的第十四说明图;
图17为产生近似表达式的第十五说明图;
图18为产生近似表达式的第十六说明图;
图19为本发明的实施例的流程图;
图20示出了输入波源的处理流程的实施例;
图21为输入波源的处理的说明图;
图22A和22B为脉冲测试波源的说明图;
图23A和23B为脉冲测试波源的说明图;
图24示出了计算近似表达式的处理流程的实施例;
图25示出了计算谐波电流的处理流程的实施例;
图26为矩量法联立方程的说明图;
图27示出了计算电流的处理流程的实施例;
图28示出了计算电压的处理流程的实施例;
图29示出了计算电磁场的处理流程的实施例;
图30为推导导体间电压的第一说明图;
图31A、31B和31C为推导时间电压的第二说明图;
图32为计算传输线耦合的处理的第一说明图;
图33为计算传输线耦合的处理的第二说明图;
图34为计算传输线耦合的处理的第三说明图;
图35为计算传输线耦合的处理的第四说明图;
图36A和36B为计算传输线耦合的处理的第五说明图;
图37为计算传输线耦合的处理的第六说明图;
图38为计算传输线耦合的处理的第七说明图;
图39A和39B为模拟的第一说明图;
图40为模拟的第二说明图;
图41为矩量法联立方程的说明图;
图42为天线辐射分析的第一说明图;
图43为天线辐射分析的第二说明图;
图44A和44B为天线辐射分析的第三说明图;以及
图45为矩量法的说明图;
在说明本发明的实施例之前,首先结合相应附应对有关技术及其缺点进行说明。
图45示出了已知的矩量法的流程图。
在用如图45所示流程的这种矩量法计算电磁场辐射强度时,首先读出需模拟的电子设备的划分成格形的外形信息,然后从需要加以计算的各频率中选择一个尚未处理的频率,用预定的计算方法对这个频率求出各格形元之间的互阻抗、互导纳和互电抗。求得的互阻抗等和由外形信息确定的波源一起代入这种矩量法的联立方程,解出在每个小区流动的电流和磁流。再用解出的这些值计算出各个观察点的电磁场辐射强度。然而,计算每个互阻抗、互导纳和互电抗需要很长的时间,并且要对频率域内的每一个频率进行这样的计算。如前面所提到的那样,这种运算在切实可行的短时间内是不能完成的。
也就是说,如前面所提到的那样,解矩量法的联立方程对于单一频率来说需要几分钟,但计算互阻抗、互导纳和互电抗则需要几个小时。由于这些必需对频率域内的每个频率都要进行计算,因此在切实可行的短时间内不可能完成。
本发明可以实现对随时间变化的电磁场强度之类的高速模拟。下面将对本发明进行详细说明。要注意的是,在象本发明那样用矩量法计算电磁场辐射强度的情况下,最好使用本发明者先前提出的矩量法高速处理技术(见日本专利申请No.7-298062)。因此,首先将对本发明者所提出的矩量法高速处理技术进行说明。
在矩量法高速处理中,需要注意的是,频率用f表示,波数用k(k=2πf/c,c为光速)表示,元间原始距离用r0表示,按矩量法分割成的元i和元j之间的互阻抗Zij可以近似表示为如下频率的指数多项式:
Zij=e-jkr0〔(a0+a1f2+a2f4+a3f6+a4f8+…)+j(b0f1+b1f+b2f3++b3f5+b4f7+)〕
互导纳Yij可以近似表示为如下频率的指数多项式:
Yij=e-jkr0〔(a0+a1f2+a2f4+a3f6+a4f8+…)+j(b0f1+b1f+b2f3++b3f5+b4f7+)〕
而互电抗Bij可以近似表示为如下频率的指数多项式:
Bij=e-jkr0〔(a0+a1f2+a2f4+a3f6+a4f8+…)+j(b0f+b1f3+b2f5+b3f7+b4f9+)〕
然后,首先在几个采样频率上按精确计算方法计算出互阻抗Zij、互导纳Gij和互电抗Bij。将所得到的值代入以上近似表达式,产生了求对以后求每个其他的互阻抗等值有用的系数ap和bp的联立方程,从这些方程就可解得系数ap和bp。此后,当给定了一个在这些采样频率以外的频率时,可用由系数ap和bp确定的以上近似表达式分别计算互阻抗Zij、互导纳Gij和互电抗Bij。这样就能实现矩量法高速处理。
注意,采样频率数由需要求得前面多少个系数ap和bp确定。例如,当采样频率数为5时,那就求得系数a0至a4和b0至b4
如即将要说明的那样,本发明采用这样一种方法:在给定了一个随时间变化的波源时,将这个波源变换到频率域,用矩量法在频率域内模拟各元的电流或电压和电磁场辐射强度,然后将模拟值反变换回时间域。此时,如果采用以上矩量法高速处理,则能以极高的速度计算出互阻抗Zij、互导纳Gij和互电抗Bij,因此就能以更高的速度执行矩量法。
这样,分析一个电子设备所以辐射有害无线电波或噪声和所以受到无线电波或噪声的干扰而不能工作的原因便成为可能。
本发明的原理和配置例示于图1和图2。
在这两个图中,标号1标的是本发明所提出的模拟设备。这种模拟设备将一个电子设备划分为一系列元,计算给定频率时这些元之间的互阻抗、互导纳和互电抗,按照矩量法从与各元分别相应的这些计算得到的值和波源模拟在每个元流动的电流和磁流。
图1所示的本发明的模拟设备1配备了管理装置10、变换装置11、计算装置12、模拟装置13和反变换装置14。
管理装置10管理需模拟的已划分为一系列元的电子设备的外形信息。变换装置11将波源的时间序列数据分段,执行付立叶变换,得出频谱。全部或部分频谱及其中每条谱线的频率由变换装置11设定为处理的对象。
计算装置12包括选择装置15、第一计算装置16、产生装置17和第二计算装置18。选择装置15根据变换装置11设定的频率选择采样频率。第一计算装置16计算在选择装置15选定的采样频率上的互阻抗、互导纳和互电抗。产生装置17从采样频率和第一计算装置16算得的互阻抗产生互阻抗的近似表达式,从采样频率和第一计算装置16算得的互导纳产生互导纳的近似表达式,以及从采样频率和第一计算装置16算得的互电抗产生互电抗的近似表达式。第二计算装置18利用产生装置17产生的这些近似表达式计算在变换装置11设定的各频率上的互阻抗、互导纳和互电抗。
模拟装置13按照矩量法从计算装置12算得的互阻抗、互导纳和互电抗以及变换装置11设定的频谱求出流过各元的电流频谱和磁流频谱,同时,如果需要的话再从电流频谱和磁流频谱求出电压频谱、电场频谱和电磁场频谱,将所求出的频谱中所需要的频率作为输出对象。反变换装置14用付立叶反变换将模拟装置13设定的需输出的频谱变换回时间域,输出变换结果。
模拟装置13按照操作员指示有时只进行求电流频谱和磁流频谱的处理。此时,图1的J/M装置19将从反变换装置14输出的时间域的电流和磁流直接求出时间域的电压、电场和磁场。
在本发明的模拟设备1中,当变换装置11通过将波源的时间序列数据分段进行付立叶变换而得到频谱时,全部或部分频谱及其中每条谱线的频率由变换装置11设定为处理的对象。选择装置15根据变换装置11设定的频率选择采样频率,第一计算装置16按照精确算法计算在选择装置15选定的采样频率上的互阻抗、互导纳和互电抗。
然后,产生装置17通过将互阻抗Zij近似为例如:
Zij=e-jkr0〔(a0+a1f2+a2f4+a3f6+a4f8+…)+j(b0f1+b1f+b2f3+b3f5+b4f7+…)〕
(其中:f为频率,k为波数,r0为元i和j之间的距离)再将采样频率和第一计算装置16算得的互阻抗代入这个近似表达式,列出以系数ap和bp为未知数的联立方程,从中解得互阻抗近似表达式的这些系数。
此外,产生装置17通过将互导纳Yij近似为例如:
Yij=e-jkr0〔(a0+a1f2+a2f4+a3f6+a4f8+…)+j(b0f1+b1f+b2f3+b3f5+b4f7+…)〕
(其中:f为频率,k为波数,r0为元i和j之间的距离)再将采样频率和第一计算装置16算得的互导纳代入这个近似表达式,列出以系数ap和bp为未知数的联立方程,从中解得互导纳近似表达式的这些系数。
此外,产生装置17通过将互电抗Bij近似为例如:
Bij=e-jkr0〔(a0+a1f2+a2f4+a3f6+a4f8+…)+j(b0f+b1f3+b2f3+b3f7+b4f9+…)〕
(其中:f为频率,k为波数,r0为元i和j之间的距离)再将采样频率和第一计算装置16算得的互电抗代入这个近似表达式,列出以系数ap和bp为未知数的联立方程,从中解得互电抗近似表达式的这些系数。
在以这种方式导出互阻抗、互导纳和互电抗的近似表达式后,第二计算装置18用这些近似表达式计算在变换装置11设定的每个频率(经L1)上的互阻抗、互导纳和互电抗。此后,模拟装置13按照矩量法从这些算得的互阻抗、互导纳和互电抗以及变换装置11设定的频谱(经L2)求出在每个元流动的电流频谱和磁流频谱。
然后,模拟装置13从以上求出的电流频谱和磁流频谱中规定希望模拟的电流频谱,求出希望模拟的电压频谱、电场频谱和磁场频谱,将这些频谱设定为输出的对象。反变换装置14通过进行付立叶反变换将设定的这些频谱变换回时间域,输出反变换结果。
这样,在图1所示的本发明的模拟设备1中,采用了一种用所似表达式的极高的速度计算出互阻抗、互导纳和互电抗的配置,因此,就能采用这样一种方法:在给定了随时间变化的波源后,将这个随时间变化的波源变换到频率域,用矩量法在频率域模拟电磁场强度等,再将模拟值反变换回时间域。
因此,使用本发明的模拟设备1就能在时间域模拟一个波源对电子设备的影响,从而能分析电子设备所以辐射有害无线电波或噪声和所以受到无线电波或噪声干扰而不能工作的原因。
另一方面,如图2所示的本发明的模拟设备1却配备有管理装置20、分解装置21、模拟装置22、计算装置23和计算装置24。
管理装置20管理划分为一系列元的电子设备的外形信息。在使用一个输出受调载波信号的波源的情况下,分解装置21将这个波源分解为三个具有不同的频率的波源,这几个不同的频率由载波信号的频率和调制信号的频率确定。模拟装置22模拟在每个元流动的频率域电流和频率域磁流。应注意,按照CISPR(国际无线电干扰特别委员会)规定,载波信号的频率为30MHz至1Ghz。
计算装置23计算时间域电流和时间域磁流。计算装置24从计算装置23求出的时间域电流和时间域磁流计算所希望得到的时间域电流、电压、电场和磁场。
在以一个对载波信号进行调制的源作为波源的图2所示本发明的模拟设备1中,分解装置21装这个波源分解为三个具有由载波信号频率和调制信号频率确定的不同频率的波源。接到这个信息后,模拟装置22对于分解装置21分解出的各波源应用矩量法分别模拟在每个元流动的相应的频率域电流和频率域磁流。
接收到这模拟后,计算装置23从模拟装置22得出的频率域电流计算出时间域电流,从模拟装置22得出的频率域磁流计算出时间域磁流。然后,计算装置24通过分别合并计算装置23求得的时间域电流和时间域磁流计算出时间域电流、电压、电场和磁场,通过使用合并后的时间域电流和时间域磁流求出时间域电压等,或者通过使用合并前的时间域电流以和时间域磁流求出时间域电压等,再将它们合并。
在一种配置中,如果波源是一个输出受调载波信号的波源,按照图1所示的本发明模拟设备那样通过将波源时间序列数据分段再进行付立叶变换得到频谱,那么分段数就会非常大。因此,这种方法对于一个输出受调载波信号并不适合。然而,使用图2所示的本发明的模拟系统1就能很方便地按照矩量法求出在每个元流动的时间域电流和时间域磁流。
下面将结合实施例对本发明进行更详细的说明。
本发明提供了一种如上所述能以极高的速度按照矩量法计算出电磁场强度的方法。
在这种方法中采用一种机理:将元i和元j之间的互阻抗(电流偶极子对电流偶极子)Zij近似为
Zij=e-jkro〔(a0+a1f2+a2f4+a3f6+a4f8+…)+j(b0f1+b1f+b2f3+b3f5b4f7+…)〕
将元i和元j之间的互导纳(磁流偶极子对磁流偶极子)Yij近拟为
Yij=e-jkr0〔(a0+a1f2+a2f4+a3f6+a4f8+…)+j(b0f1+b1f+b2f3+b3f5+b4f7+…)〕
将元i和元j之间的互电抗(磁流偶极子对磁流偶极子)Bij近拟为
Bij=e-jkr0〔(a0+a1f2+a2f4+a3f6+a4f8+…)+j(b0f1+b1f+b2f3+b3f5+b4f7+…)〕
其中f为频率,k为波数,(k=2πf/c,c为光速),而R0为两元之间的原始距离。
下面,将对互阻抗Zij可以用这种频率的指数多项式近似进行说明。
为了说明互阻抗Zij的近似表达式,首先考察图3所示的单极子情况。在图中,粗线表示单极子,而虚线表示下面还将提及的扩展函数J1和J2的形状。这里,单极子<1>与单极子<3>之间的倾斜角定义为φ1,而单极子<1>与单极子<4>之间的倾斜角定义为φ2。在图中,φ1和φ2由φ表示。
互阻抗Zij的通用方程示于图4A,其中ω为角频率,r为距离,ρ1=-1/jω×αJ1/αt,ρ2=-1/jω×αJ2/αt。此外,所以对S进行积分是因为考虑的不仅是对于单极子是线形的情况,而且也是对于单极子是面形的情况。
J1和J2是矩量法中的扩展函数。扩展函数示出了电流在单极子上分布的形状。对于不同的矩量法,扩展函数也不同,但下面的证明可以适用于所有的矩量法。也就是说,可以用正弦电流、三角形电流或脉冲函数电流中任何一种分布来作为扩展函数。此外,以下证明对于线形或面形单极子都适用。这里是用“分段”正弦矩量法进行证明的。注意:在图3中为了方便起见虚线所示的扩展函数是三角形的,但在以下说明中的扩展函数是正弦电流。
在分段正弦矩量法中,图3的单极子<1>至<4>的扩展函数可以表示为
电流单极子<1>J1=sink(z-z0)/sinkd1
电流单极子<2>J2=sink(z2-z)/sinkd2
电流单极子<3>J3=sink(t-t0)/sinkd3
电流单极子<4>J4=sink(t2-t)/sinkd4
其中,d1为单极子<1>的长度,d2为单极子<2>的长度,d3为单极子<3>的长度,d4为单极子<4>的长度。
使用这些扩展函数,首先求出单极子<1>和单极子<3>之间的互阻抗Z13以及单极子<1>和单极子<4>之间的互阻抗Z14。互阻抗Z13和Z14示于图4B。
考虑到单极子之间的距离r可以表示为“r(z2+t2-2ztcosφ+h2)1/2,因此互阻抗Z13和Z14的系数如图5所示。注意:α=cμ/4π(α为常数)。
此外,如果一个单极子的端点和另一个单极子的端点之间的距离表示为r0,则这两个单极子之间的距离r可以近似为
r=[r0 2+(r2-r0 2)]1/2=(r0 2+Δ)1/2=r0(1+Δ/r0 2)1/2
 r0(1+Δ/2r0 22/8r0 4+..)
 =r0+Δ/2r02/8r0 3+…=r0+d
(其中,Δ=(r2-r0 2),d=Δ/2r02/8r0 3+…,而r0=(z1 2+t1 2-2z1t1cosφ+h2)1/2)
由此可得出互阻抗Z13和Z14,如图5B中方程所示,其中A1=sinkd1·sinkd3
为了以简单的格式表示方程,令z-z0=u,t-t0=v,W=-t+t2。再在假设单极子很短的条件下利用多项式近似
sinkuku-(ku)3/6,cosku1-(ku)2/2
则互阻抗Z13的实部R1变换为图6所示方程,最终为图7所示方程。注意:式中略去了e-jkr0标记。
这样,在将k4等项的系数改写为P1至P8后,这方程可表示为
R1=(α/A1)[P1k4-P2k6+P3k8-P4k10)cosφ1
                          -P5k2+P6k4-P7k6+P8k8]
此外,如果采用诸如“A1=sinkd1·sinkd3=d1d3k2”的近似式,再将前面略去的e-jkr0补上,则互阻抗Z13的实部R1可表示为
R1=[αe-jkr0/(d1d3k2)][(P1k4-P2k6+P3k8-P4k10)cosφ1
                           -P5k2+P6k4-P7k6+P8k8]
  =[αe-jkr0/(d1d3)][(P1k2-P2k4+P3k6-P4k8)cosφ1
                           -P5+P6k2-P7k4+P8k6]
另一方面,互阻抗Z13的虚部I1表示为图8中的方程所示。由此,最终表示为图9中的方程所示。注意已略去了e-jkr0标记。
这样,在将k3等项的系数改写为Q1至Q8后,这方程可表示为
I1=(α/A1)[Q1k3-Q2k5+Q3k7-Q4k9)cosφ1
                    -Q5k+Q6k3-Q7k5+Q8k7]
此外,如果采用诸如“A1=sinkd1·sinkd3=d1d3k2”,再将前面略去的e-jkr0补上,则互阻抗Z13的虚部I1可表示为
I1=[αe-jkr0/(d1d2k2)][(Q1k3-Q2k5+Q3k7-Q4k9)cosφ1
                          -Q5k+Q6k3-Q7k5+Q8k7]
  =[αe-jkr0/(d1d3)][(Q1k-Q2k3+Q3k5-Q4k7)cosφ1
                          -Q5/k+Q6k-Q7k3+Q8k5]
由此,互阻抗Z13所表示为
Z13=R1+jI1
   =[αe-jkr0/(d1d3)][(P1k2-P2k4+P3k6-P4k8)cosφ1
                            -P5+P6k2-P7k4+P8k6]
      +j[αe-jkr0/(d1d3)][(Q1k-Q2k3+Q3k5-Q4k7)cosφ1
                            -Q5/k+Q6k-Q7k3+Q8k5]
类似,也可求得互阻抗Z14。在将K3等项的系数改写为R1至R8后,互阻抗Z14可表示为
Z14=[αe-jkr0/(d1d4)][(R1k2-R2k4+R3k6-R4k8)cosφ2
                            +R5-R6k2+R7k4-R8k6]
      +j[αe-jkr0/(d1d4)][(S1k-S2k3+S3k5-S4k7)cosφ2
                            +S5/k-S6k+S7k3-S8k5]
这样,可以如下面所示那样将互阻抗(Z13+Z14)用波数k的多项式来表示。在将k3等项的系数改写为C1至C9后,可得
Z13+Z14=e-jkr0[(C0+c1k2+c2k4+C3k6+C4k8+..)
             +j(C5k-1+C6k+c7k3+C8k5+C9k7+..)]
如果单极子<2>与单极子<3>之间的互阻抗用Z13表示,而单极子<2>与单极子<4>之间的互阻抗用Z24表示,则互阻抗(Z23+Z24)同样可以求得。因此,互阻抗(Z13+Z14+Z23+Z24)也可以用一个与以上方程类似的波数的多式表示。
也就是说,互阻抗Zij可以近似为
Zij=e-jkr0〔(a0+a1f2+a2f4+a3f6+a4f8+…)+j(b0f1+b1f+b2f3+b3f5b4f7+…)〕
可以证明,互导纳Yij可以用一个与上述互阻抗Zij类似的频率的多项式近似。而且,互导纳Yij的近似表达式与互阻抗Zij的近似表达式完全相同,因此不再重复予以证明。
也就是说,互导纳Yij可以近似为
Yij=e-jkr0〔(a0+a1k2+a2k4+a3k6+a4k8+…)+j(b0k-1+b1k+b2k3+b3k5+b4k7+…)〕
下面,将对互电抗Bij可以用频率的多项式近似这个问题进行说明。
为了说明互电抗Bij的近似表达式,考察如图10所示的单极子情况。在图中,粗线表示单极子,而虚线表示扩展函数的形状,这与图3所示的相同。
考虑电流源引起的一个磁场。这里,单位向量可以如图11A中所示。
在分段正弦矩量法中,图10的单极子<1>至<4>的扩展函数可以表示为
电流单极子<1>J1=sink(z-z0)/sinkd1
电流单极子<2>J2=sink(z2-z)/sinkd2
电流单极子<3>M3=sink(t-t0)/sinkd3
电流单极子<4>M4=sink(t2-t)/sinkd4
其中,d1为单极子<1>的长度,d2为单极子<2>的长度,d3为单极子<3>的长度,d4为单极子<4>的长度。
如果只在乙座标上有电流源,则电磁场呈现为圆对称(与φ座标无关),因此只存在Hφ。Hφ由图11B所示方程表示。
在图10中,单极子<3>的切向磁场分量Ht为-(h/p)×Hφsinφ,因此互电抗Bij由图11c中所示方程表示。
将图12中所示条件代入这方程后,即得出由图13A中所示方程表示的互电抗Bij
此外,单极子之间的距离r可以如前所说明的那样近似为
r=(z2+t2-2ztcosφ+h2)1/2r0+d,
因此互电抗Bij可以表示为图13B中所示的方程。
由此,单极子<1>和单极子<3>之间的互电抗B13以及单极子<1>和单极子<4>之间的互电抗B14可以表示为图14中所示的方程。
为了将这两个方程变换成比较简单的形式,令z-z0=u,t-t0=v,w=-t+t2。再在假设单极子的条件下利用多项式近似
sinkuku-(ku)3/6,cosku1-(ku)2/2
则互电抗B13的右侧第二项变换为图16所示方程。
这样,在将k2等项的系数改写为P1至P8后,该式可表示为
P1k2+P2k4+P3k6+P4k8+j(P5k3+P6k5+P7k7+p8k9)
同样,互电抗R13的右侧第一项变换为图17所示方程,最终表示为图18所示方程。
这样,在将k2等项的系数改写为Q1至Q8后,该式可表示为
Q1k4+Q2k6+Q3k8+Q4k10+j(Q5k3+Q6k5+Q7k7+Q8k9)
由此,互电抗B13可表示为
B13=[h1sinφ1/(4πsinkd1·sinkd3)]e-jkr0×[P1k2+P2k4
     +P3k6+P4k8+j(P5k3+P6k5+P7k7+P8k9)+Q1k4
     +Q2k6+Q3k8+Q4k10+j(Q5k3+Q6k5+Q7k7+Q8k9)]
此外,在采用了诸如“sinkd1·sinkd3d1d3k2”近似式后,则互电抗B13可表示为
B13=[h1sinφ1/(4πd1d3)]e-jkr0×{P1+(P2+Q1)k2
     +(P3+Q2)k4+(P4+Q3)k6+Q4k8+j[(P5+Q5)k
     +(P6+Q6)k3+(P7+Q7)k5+(P8+Q8)k7]}
   =e-jkr0[R1+R2k2+R3k4+R4k6+R5k8
     +j(R6k+R7k3+R8k5+R9k7)]
同样也可以求得互电抗B14如下
B14=e-jkr0[S1+S2k2+S3k4+S4k6+S5k8
                   +j(S6k+S7k3+S8k5+S9k7)]
这样,互电抗(B13+B14)就能用如下的波数k的多项式表示:
B13+B14=e-jkr0[C1+C2k2+C3k4+C4k6+C5k8
                   +j(C6k+C7k3+C8k5+C9k7)]
如果单极子<2>和单极子<3>之间的互电抗用B23表示,单极子<2>和单极子<4>之间的互是抗用B24表示,同样可求得互电抗(B23+B24),因此互电抗(B13+B14+B23+B24)也可以用与以上方程类似的波数k的多项式表示。
即互电抗Bij可以近似为
Bij=e-jkr0[(a0+a1f2+a2f4+a3f6+a4f8+..)
     +j(b0f+b1f3+b2f5+b3f7+b4f9+..)]
下面,将按照本发明的模拟设备1所执行的处理流程对本发明进行详细说明。注意。在下面说明的实施例中,为了方便起见假设模拟只考虑互阻抗。
图19示出了本发明的模拟设备1执行的处理流程的一个实施例。图中,管理文件100对需模拟的划分成一系列元的电子设备的外形信息进行管理。
如这个处理流程所示,本发明的模拟设备1首先执行步骤1(ST1)的“输入波源处理”,然后执行步骤2(ST2)的“计算近似表达式处理“,执行步骤3(ST3)的“计算谐波电波处理”,执行步骤4(ST4)的“计算电流处理”,执行步骤5(ST5)的“计算电压处理”,执行步骤6(ST6)的“计算电磁场处理”,执行步骤7的“计算导体间电压处理”,以及执行步骤8(ST8)的“计算传输线耦合处理”。
图20示出了图19的“输入波源处理”(ST1)流程的一个实施例。
如这个处理流程所示,在进行“输入波源处理”中,首先从管理文件100读出时间域的波源数据,然后将所读出的数据分段,进行付立叶变换,得出波源的频谱。
也就是如图21中所示,波源的时间序列数据分段后,通过付立叶变换,得出波源的频谱。此时,在这频谱的直流分量是必需时,按照欧姆定律,就会得到它的直流分量。
作为对电子设备的测试,可以是加一个高电压的脉冲或加一个大电流的脉冲,然后考察因此而产生的影响。为了模拟这种情况,必需将一个充有初始电压的电容假设为波源,或将一个加有初始电流的电感假设为波源。
如果波源是一个充有初始电压的电容,则这个波源可以变换成如图22A所示的等效电路。此外,由于进行付立叶变换要求有周期性,因此为了方便起见,这个等效电路所产生的阶越电压将处理成时间间隔t1和t2相当长的周期性脉冲电压,如图22B所示。电容初始电压所造成的影响通常呈现为一个衰减振荡现象,因此可以事先就预测出这个衰减振荡现象的收敛时间,从而令t1和t2的值大于这预测的收敛时间。
如果波源是一个加有初始电流的电感,则这个波源可以变换成如图23A所示的等效电路。此外,由于进行付立叶变换要求有周期性,因此为了方便起见,这个等效电路所产生的阶越电流将处理成时间间隔t1和t2相当长的周期性脉冲电流,如图23B所示,电感初始电流所造成的影响通常呈现为一个衰减振荡现象,因此可以事先就预测出这个衰减现象有收敛时间,从而令t1和t2的值大于这预测的收敛时间。
这样,在“输入波源处理”中,在从管理文件100读出时间域的波源数据后,将所读出的数据分段和付立叶变换,便得到波源的频谱,于是程序进至“计算近似表达式处理”,如前面所述(图19的ST2)。
图24示出了“计算近似表达式处理”流程的一个实施例。
如这个处理流程所示,在“计算近似表达式处理”中,首先在步骤1(ST1),从在“输入波源处理”中得到的波源频谱的那些正侧频率中规定需模拟的频率范围。如果要求精度高,则可将这频率范围规定为整个频谱的全部频率。如果首先要考虑的是减少处理时间,则可除掉一些谐波分量,而将模拟的频率限制在一个比较窄的范围内。
然后在步骤2(ST2),从所规定的范围内的这些频率中规定一些采样频率。例如,选定5个采样频率fs1、fs2、fs3、fs4和fs5。可以通过例如大致均分步骤1(ST1)中可规定的频率范围的方式从这个频率范围内的频率中选出这些采样频率。这里,如果要求精度高,则可增加采样数,而如果首先要考虑的是减少处理时间,则可减少采样数。
然后,在步骤3(ST3),从未处理的采样频率中选择一个采样频率;在步骤4(ST4),按照精确计算方法计算元i和元j(i,j=1至M)之间的对于所选采样频率的互阻抗Zij,将算得的值存入工作文件200后,处理程序返回步骤3(ST3)的开始处。
在步骤3(CT3),在已逐个选择了所有的采样频率时,确定已经完成了对于所有采样频率的互阻抗Zij(i,j=1至M)的计算,则处理程序进至步骤5(ST5),将互阻抗Zij近似为
Zij=e-jkr0[(a0+a1f2+a2f4+a3f6+a4f8+..)
     +j(b0f-1+b1f+b2f3+b3f5+b4f7+..)]
(其中:f为频率,k为波数,r0为元i与元j之间的距离)。将存储在工作文件200中的采样频率和互阻抗Zij代入这个近拟表达式,列出以系数a0至aL和b0至bL作为未知数的联立方程。从联立方程解出这个近似表达式的系数a0至aL和b0至bL。这些系数存入工作文件200。
例如,如果已求得在元i和元j之间对于5个采样频率fs1至fs5的互阻抗Zijs1至Zijs5
Zijs1=αijs1+jβijs1
Zijs2=αijs2+jβijs2
Zijs3=αijs3+jβijs3
Zijs4=αijs4+jβijs4
Zijs5=αijs5+jβijs5
则可列下如下联立方程αijs1=exp[-j(2πfs1/c)r0][a0+a1fs1 2+a2fs1 4+a3fs1 6+a4fs1 8]αijs2=exp[-j(2πfs2/c)r0][a0+a1fs2 2+a2fs2 4+a3fs2 6+a4fs2 8]αijs3=exp[-j(2πfs3/c)r0][a0+a1fs3 2+a2fs3 4+a3fs3 6+a4fs3 8]αijs4=exp[-j(2πfs4/c)r0][a0+a1fs4 2+a2fs4 4+a3fs4 6+a4fs4 8]αijs5=exp[-j(2πfs5/c)r0][a0+a1fs5 2+d2fs5 4+a3fs5 6+a4fs5 8]以及βijs1=exp[-j(2πfs1/c)r0][b0+b1fs1 2+b2fs1 4+b3fs1 6+b4fs1 8]βijs2=exp[-j(2πfs2/c)r0][b0+b1fs2 2+b2fs2 4+b3fs2 6+b4fs2 8]βijs3=exp[-j(2πfs3/c)r0][b0+b1fs3 2+b2fs3 4+b3fs3 6+b4fs3 8]βijs4=exp[-j(2πfs4/c)r0][b0+b1fs4 2+b2fs4 4+b3fs4 6+b4fs4 8]βijs5=exp[-j(2πfs5/c)r0][b0+b1fs5 2+b2fs5 4+b3fs5 6+b4fs5 8]
解以上联立方程,即可求得系数a0至a4和b0至b4,从而产生元i和元j之间的互阻抗Zij的近似表达式如下
 Zij=e-jkr0[(a0+a1f2+a2f4+a3f6+a4f8)
                +j(b0f-1+b1f+b3f5+b4f7)]
注意到执行这个“计算近似表达式程序”对于每个采样频率是独立的,因此如果使用一个并行计算机的话,就能执行更高速度的处理。
在通过“计算近似表达式处理”以这种方式求得表示元i和元j(i,j=1至M)之间的互阻抗Zij的近似表达式后,程序进至“计算谐波电流处理”(图19的ST3)。
图25示出了“计算谐波电流处理”流程的一个实施例。
如这个处理流程所示,在“计算谐波电流处理”中,如果输入了在“输入波源处理”中得到的具有正侧频率的波源的频谱,则首先在步骤1(ST1),从这些正侧频率中选择一个未处理的频率。
然后在步骤2(ST2),将所选频率代入通过“计算近似表达式处理”求出的互阻抗Zij的近似表达式,求出元j和元j(i,j=1至M)之间的对于所选频率的互阻抗Zij。这个互阻抗Zij的计算是一种简单的代入计算,因此可以极高的速度执行。
然后在步骤3(ST3),利用在步骤2(ST2)算得的互阻抗Zij(i,j=1至M)和由在步骤1(ST1)所选的fq规定的波源的频谱V(ω)(ω=2πfq)解图26所示矩量法联立方程,求出在每个元流动的频率为在步骤1(ST1)选定的频率fq的电流I1(fg)至IM(fg)。这种处理一直重复到在步骤1完成频率选择为止。因此,求出了在每个元流动的所有输入频率fq(q=1至N)的电流I1(fg)至IM(fg)。
这里,在图26所示的矩量法联立方程中,假设在元m存在波源v(ω)。注意,在解矩量法联立方程时候,如果采用例如本发明的发明者提出LU分解法(见日本专利申请No.7-342695),那就可以达到更高的处理速度。
通过“输入波源处理”得到的波源频谱包括具有负侧频率的部分,如图21所示。如上面所述,在“计算近似表达式处理”和“计算谐波电流处理”中,已经计算出的只是波源引起的具有正侧频率的电流。因此,必需计算波源引起的具有负侧频率的电流。
这个计算是利用在每个元m流动的电流Im(ω)和Im(-ω)之间存在的共轭关系执行的,这是因为Zij(ω)和Zij(-ω)之间有着如下那样的共轭关系
Zij(ω)=e-jkr0[(a0+a1f2+a2f4+a3f6+a4f8)
               +j(b0f-1+b1f+b2f3+b3f5+b4f7)]
Zij(-ω)=e-jkr0[(a0+a1f2+a2f4+a3f6+a4f8)
               -j(b0f-1+b1f+b2f3+b3f5+b4f7)]以及波源的v(ω)和v(-ω)的频谱v(ω)和v(-ω)之间有着共轭关系的缘故。也就是说,波源引起的具有负侧频率的电流Im(-ω)可以通过将波源引起的具有正侧频率的电流Im(ω)的虚部的符号颠倒求得。
注意到执行这个“计算谐波电流处理”对于各频率是独立的,因此如果使用一个并行计算机,那就能以更高的速度进行这个处理。
这样通过“计算谐波电流处理”求得在每个元流动的各频率的电流后,程序进至“计算电流处理”(图19的ST4)。
图27示出了“计算电流处理”流程的一个实施例。
如这个处理流程所示,在这个“计算电流处理中”,首先通过与使用者交互处理(ST1)规定需计算电流的点P(即元P),提取流过这个规定点p的各频率的电流Ip(f1)至Ip(fN)、IP(0)和Ip(-f1)至Ip(-fN),加以付立叶反变换(ST2),从而计算出流过这个规定点P的时间域电流。
这样通过“计算电流处理”求得流过规定点的时间域电流后,程序进至“计算电压处理”(图19的ST5)。
图28示出了“计算电压处理”流程的一个实施例。
如这个处理流程所示,在这个“计算电压处理”中,首先通过与使用者交互处理(ST1)规定需计算电压的点P(即元P),提取流过这个规定点p的各频率的电流Ip(fl)至Ip(fN)、IP(0)和Ip(-f1)至Ip(-fN)(ST2),然后将各提取的值分别乘以这个规定点P的相应频率的阻抗Zp(d1)至Zp(fn)、Zp(0)和Zp(-fl)至Zp(-fn),求得在这个规定点p产生的各频率的电压Vp(f1)至Vp(fN)、Vp(0)和Vp(-f1)至Vp(-fn)(ST3)。
然后,对因此求得的各频率的电压Vp(f1)至Vp(fN)、Vp(0)和Vp(-f1)至Vp(-fn)进行付立叶反变换(ST4),计算出在规定点P所产生的时间域电压。
这样通过“电压计算处理”求得在规定点所产生的时间域电压后,程序进至“计算电磁场处理”(图19的ST6)。
图29示出了“计算电磁场处理”流程的一个实施例。
如这个处理流程所示,在“计算电磁场处理中”,以计算电场强度情况作为例子加以说明。首先通过与使用者交互处理(ST1)规定观察点p,按照已知的电磁场理论公式求出在每个元流动的各频率的电流在观察点p产生的电场,从而求得在观察点p产生的各频率的电场Ep(fl)至Ep(fN)至Ep(fO)和Ep(-f1)至Ep(-fN)。
然后,对求得的各频率的电场Ep(f1)至Ep(fN)Ep(f0)和Ep(-f1)至Ep(-fN)进行付立叶反变换(ST4),计算出在观察点p产生的时间域电场。磁场计算亦可按此进行。
这样通过“计算电磁场处理”求得观察点的电磁场后,程序进至“计算导体间电压处理”(图19的ST7)。
在“计算导体间电压处理”中,通过与使用者交互处理(ST1)规定需处理的导体间位置p后,如果在元n流动的角频率为ω的电流表示为In(ω),而导体间位置p和元n之间角频率为ω的互阻抗表示为Zpn(ω),则导体间角频率为ω的电压Vp(ω)可以表示为公式1
Figure A9710331700281
由此,可求得导体间各频率的电压Vp(f1)至Vp(fN)、Vp(0)和Vp(-f1)至Vp(-fN)。
然后,对因此求得的导体间各频率的电压Vp(f1)至Vp(fN)、Vp(0)和Vp(-f1)至Vp(-fN)进行付立叶反变换,计算出在导体间位置p处所产生的导体间时间域电压。
所以得出公式1的逻辑说明如图30所示。如果在以P1和P2表示的导体之间插入一个电阻R,由于有着围绕每个导体的电场为零的边界条件,因此图31A所示方程成立。由此可得导体间电流IP为图31B中方程所示,从而导体间电压VP为图31c中方程所示。实际上导体之间并无电流流动,因此将R→∞,IP1、IP2→0代入图31中的方程,即可得到公式1。
也就是说,公式1是从由一个假想插在导体之间的电阻产生的电压频谱在这个电阻趋向无穷大时导出的。
这样通过“计算导体间电压处理”求得导体间规定位置处的电压后,程序进至“计算传输线耦合处理”(图19的ST8)。
执行“计算传输线耦合处理”是为了模拟由于离宽接地平面一定高度的传输线的外界电场而加到与传输线连接发接收电路之类上的电流或电压的情况。
假设传输线离宽接地平面的架设情况如图32所示。其中:发送端阻抗用Zd表示,接收端阻抗用zr表示,传输线长度用L表示,传输线高度用h表示,传输线传播常数用β表示,传输线特征阻抗用Z0表示,由于其他导体上的电流分布而投射到传输线的位置(x,z)上的电场的x分量用Ex i(x,z)表示,由于其他导体上的电流分布而投射到传输线的位置(x,z)上的电场的Z分量用Ez i(x,z)表示,而k(x)定义为
公式2
d(x)=Ex i(x,h)-Ex i(x,0)D定义为
公式3
D=(Z0Zd+Z0Zr)coSβL+j(Z0 2+ZdZr)sinβL
如图33所示,在这传输线以与矩量法相同的正弦扩展函数之类分段时,Ez i和k(x)可通过用矩量法计算出的其他导体的电流分布和互阻抗求得,如图34中的方程所示。其中,Znm的传输线的上导线和其他导体之间和的互阻抗,而Z′nm为传输线的下导线和其他导体之间的互阻抗。
另一方面,在存在这类外界电场Ex i和Ez i的情况下,传输线发送端处的电流Id(ω)已为C·D·Taylor求得,如图35中的方程所示。该方程中的A项经如图36A所示改写后,得到的在发送端处的电流Id(ω)如图36B中的方程所示。
此外,在存在外界电场的情况下,传输线接收端处的电流Ir(ω)也已为C.D.Taylor求得,如图37中的方程所示。该方程经改写后成为图38所示方程。
以上C.D.Taylor方程可参见C.D.Taylor,R.S.Satterwite和C.W.Harrison,Jr.撰写的“非均匀电磁场所引起的有载双导线传输线的响应”(“ The response of a terminated two-wire transmisson line excitedby a nonuniform electromagnetic field”,IEEE Trns.AntennasPropagation.,Ap-13,No.6,pp.987-989,Nov.1967”)。
由此,在得到对于传输线的外界电场Ex i和Ez i(亦即Ez i和k(x))时,就可按照C.D.Taylor得出的方程求出发送端处的电流Id(ω)和接收端处的电流Ir(ω),将所得到的值乘以相应的阻抗,即可求得发送端处的电压Vd(ω)和接收端处的电压Vr(ω)。
此后,在“计算传输线耦合处理”中,按照C.D.Taylor得出的方程求得流过发送端和接收端的各频率的电流IP(f1)至IP(fN)、IP(0)和IP(-f1)至Ip(-fN),再加以付立叶反变换后,计算出流过发送端和接收端的时间域电流。此外,按照C.D Taylor得出的方程求得在发送端和接收端产生的各频率的电压VP(f1)至Vp(fN)、VP(0)和VP(-fl)至VP(-fN),再加以付立叶反变换后,计算出在发送端和接收端产生的时间域电压。
通过这种“计算传输线耦合处理”,因此就能模拟外界电场加到传输线接收电路之类上的电流和电压的情况。
图39A和39B示出了采用“计算传输线耦合处理”的模拟结果。这个模拟是在图40所示的条件下进行的。
图39A示出了所假设的传输天线的传输电流,而图39B示出了由于这传输电流而引起的传输线的电流。在这个模拟中,模拟了电流以相同的模式流过离传输天线1米的传输线,时间滞后与无线电波传播时间相当,因此“计算传输线耦合处理”的合理性就可以得到了验证。
在以这种方式给定了一个随时间变化的波源时,本发明的模拟设备1将这个波源变换到频率域,利用近似表达式以极高的速度求出互阻抗,用矩量法在可变换的频率域内模拟电磁场强度等,然后将模拟结果反变换回时间域。
这个实施例将矩量法联立方程的解假设为如图26所示只考虑互阻抗的情况,但在使用如图41所示考虑互导纳和互电抗的矩量法联立方程的情况下,也可以利用前面所提到的近似表达式以极高的速度求出互导纳和互电抗。
在图41所示的方程中,Ic,n为流过金属物体的电流,Id,n为在介电物体表面上流动的电流,而Mn为在介电物体表面上流动的电流。右上角的上标0表示在空气中的值,右上角的上标d表示在介电物体中的值,下标c表示金属物体,而下标d表示介电物体。
此外,在这个实施例中,解出矩量法联立方程后,通过“计算电压处理”、“计算电磁场处理”、“计算导体间电压处理”和“计算传输线耦合处理”求得电压频谱等,再对所求得的电压频谱等进行付立叶反变换。但本发明并不局限于此,本发明能直接对解矩量法联立方程可得到的谐波电流进行付立叶反变换,将它变换成时间域电流,再在这时间域内执行以上在“计算电压处理”、“计算电磁场处理”、“计算导体间电压处理”和“计算传输线耦合处理”所进行的这些处理。
对于电子设备所进行的一种测试是测试电子设备是否受到天线所辐的电场的影响(即抗扰性测试)。
这种测试例如通过向电子设备辐射一个用频率为1KHz的调节信号对频率为30MHz至1GHz的载波信号进行80%的调制而得到的受调波来实现,如图42所示。
如果要用上面这种方法来进行这种测试,那么这个波源的时间序列数据就会非常多。这是因为,如果1GHz的一个周期分出20个点,那么1KHz就需要有多达20×106个点的数据。
因此,在这种情况下,如果载波信号的角频率表示为ωc,载波信号的振幅表示为Vc,调制信号的角频率表示为ωm,调制信号的调制度为m,注意到按照
V=Vc(1+msinωm)sinωct
 =Vcsinωct+(mVc/2)cos(ωcm)t
    -(mVc/2)cos(ωcm)t
天线波源可以分解为三个具有不同频率的波源(这里用了cost=jsint关系),如图43所示,于是可依次对这三个波源执行矩量法。
这样,在用这种矩量法求得在元m流动的频率域电流(如图44A所示)后,从求得的这些电流可求得时间域电流(如图44B所示),相加后所得出在元m流动的电流。
因此,不需要对一个输出一个通过对载波信号进行调制而得到的波的波源进行付立叶变换就能求得在每个元流动的时间域电流和时间域磁流。然后,就能从以上面这种方式得到的相加前和相加后的时间域电流和时间域磁流求得在规定元处产生的时间域电压和在规定点产生的时间域电场和磁场,因此就能模拟这个波源对电子设备的影响。
例如,通过合并求得的时间域电流求出在规定点产生的时间域电流,按照从由一个假想插在导体之间的电阻产生的电压在这电阻趋向无穷大时导出的计算方程、利用合并前的时域电流求出导体间的时间域电压,按照C.D.Taylor求得的方程计算出传输线的接收电路等接收到的频率域电流和电压后加以合并,从而计算出由于波源而影响传输线的接收电路等的时间域电流和电压。
采用这种方式就能快速模拟一个电子设备是否会受到图42中所示天线的电场辐射影响的试验。
如上所述,按照本发明,在执行矩量法时,利用近似表达式以极高的速度计算出互阻抗、互导纳和互电抗。然后,在给定一个随时间变化的波源时,将计算结果变换到频率域。这样就能采用一种方法,用矩量法模拟频率域电磁场强度等,再将这模拟结果反变换回时间域。
由此可见,采用本发明能模拟波源在时间域对电子设备的影响,从而就能分析一个电子设备可以辐射有害无线电波或噪声和所以受到无线电波或噪声的干扰而不能工作的原因。
所以,按照本发明,利用矩量法,能很容易地模拟一个输出一个通过对载波信号进行调制而得到的波的波源在时间域对电子设备的影响。
因此,采用本发明就能在桌上检验一个电子设备是否将会受到天线所投射的电场的影响而不能工作。

Claims (16)

1.一种采用矩量法的模拟设备,所述设备将一个电子设备划分成一系列元,计算在给定各频率时这些元之间的相应互阻抗,以及利用与所述各元分别相应的互阻抗和波源按照矩量法模拟在每个所述元流动的电流,所述设备的特征是它包括:
一个变换装置,其作用是通过将一个波源的时间序列数据分段、加以付立叶变换得到一个频谱,将这个频谱的全部或部分和其中所含的频率设定为处理对象;
一个计算装置,其作用是根据所述变换装置设定的频率采样频率,计算在采样频率的互阻抗,从算得的互阻抗和采样频率产生一个互阻抗的近似表达式,用这个近似表达式计算在所述变换装置设定的各频率的互阻抗;
一个模拟装置,其作用是按照矩量法从所述计算装置算得的互阻抗和所述变换装置设定的频谱求出在每个所述元流动的电流的频谱;以及
一个反变换装置,其作用是对所述模拟装置求得的电压频谱进行付立叶反变换,输出反变换结果。
2.一种采用矩量法的模拟设备,所述设备将一个电子划分成一系列元,计算在给定各频率时这些元之间的相应互阻抗,以及利用与所述各元分别相应的互阻抗和波源按照矩量法模拟在每个所述元流动的电流,所述设备的特征是它包括:
一个变换装置,其作用是通过将一个波源的时间序列数据分段、加工付立叶变换得到一个频谱,将这个频谱的全部或部分和其中所含的频率设定为处理对象;
一个计算装置,其作用是根据所述变换装置设定的频率选择采样频率,计算在采样频率的互阻抗,从算得的互阻抗和采样频率产生一个互阻抗的近似表达式,用这个近似表达式计算在所述变换装置设定的各频率的互阻抗;
一个模拟装置,其作用是按照矩量法从所述计算装置算得的互阻抗和所述变换装置设定的频谱求出在每个所述元流动的电流的频谱,同时,如果需要的话,从这电流频谱求出电压频谱、电场频谱和磁场频谱,将求得的这些频谱中所需要的频谱设定为输出对象;以及
一个反变换装置,其作用是对所述模拟装置设定的需输出的频谱进行付立叶反变换,输出反变换结果。
3.一种按权利要求2所述的采用矩量法的模拟设备,其中
所述模拟装置,在需要求出导体间电压频谱的情况下,按照从由一个假想插在导体之间的电阻产生的电压频谱在这个电阻趋于无穷大时导出的计算式计算导体间电压的频谱。
4.一种按权利要求2所述的采用矩量法的模拟设备,其中
所述模拟装置,在需要求出在传输线的终端的电流频谱或电压频谱的情况下,执行移除传输线的矩量法,按照执行结果求出由求得的各元的电流频谱加到终端的电场,再按照上述计算方程用求得的电场计算在传输线的终端的电流频谱和电压频谱。
5.一种按权利要求1至4中任何一个权利要求所述的采用矩量法的模拟设备,其中
所述变换装置,在需要将充在电容上的初始电压作为一个波源处理的情况下,将初始电压的周期信号分段,加以付立叶变换;而在需要将流过电感的初始电压作为一个波源处理的情况下,将初始电流的周期信号分段,加以付立叶变换;
6.一种按权利要求1至5中任何一个权利要求所述的采用矩量法的模拟设备,其中
所述计算装置除了计算互阻抗外还用与计算互阻抗相同的程序计算互导纳和/或互电抗;
所述模拟装置还结合所述计算装置算得的互导纳和互电抗考虑磁流频谱执行处理;以及
所述反变换装置执行与所述模拟装置配合的付立叶反变换。
7.一种按权利要求1至6中任何一个权利要求所述的采用矩量法的模拟设备,其中
所述计算装置将元i和元j之间的互阻抗Zij近似为zij=e-jkr0[(a0+a1f2+a2f4+a3f6+a4f8+..)
      +j(b0f-1+b1f+b2f3+b3f5+b4f7+..)]式中:f为频率,k为波数,r0为元i与元j之间的距离,而j为虚部符号
Figure A9710331700041
8.一种按权利要求6所述的采用矩量法的模拟设备,其中所述计算装置将元i和元j之间的互导纳Vij近拟为Yij=e-jkr0[(a0+a1f2+a2f4+a3f6+a4f8+..)
      +j(b0f-1+b1f+b2f3+b3f5+b4f7+..)]式中:f为频率,k为波数,r0为元i与元j之间的距离,而j为虚部符号
9.一种按权利要求6所述的采用矩量法的模拟设备,其中所述计算装置将元i和元j之间的互电抗Bij近似为Bij=e-jkr0[(a0+a1f2+a2f4+a3f6+a4f8+..)
      +j(b0f+b1f3+b2f5+b3f7+b4f9+..)]式中:f为频率,k为波数,r0为元i与元j之间的距离,而j为虚部符号
10.一种采用矩量法的模拟设备,所述设备将一个电子设备划分成一系列元,计算在给定各频率时这些元之间的相应互阻抗,以及利用与所述各元分别相应的互阻抗和波源按照矩量法模拟在每个所述元流动的电流,所述设备的特征是它包括:
一个分解装置,其作用是在将一个通过调制一个载波信号而得到的波用作一个波源的情况下,将这个波源分解为三个具有由载波信号频率和调制信号频率决定的不同频率的波源;
一个模拟装置,其作用是对所述分解装置分解得的这些波源采用矩量法求出在每个所述元流动的频率域电流;以及
一个计算装置,其作用是从所述模拟装置求得的频率域电流计算出时间域电流。
11.一种按权利要求10所述的采用矩量法的模拟设备,其中
所述模拟装置除了互阻抗外还考虑互导纳和/或互电抗采用矩量法求出频率域磁流;以及
所述计算装置还从所述模拟装置求得的频率域磁流计算出时间域磁流。
12.一种采用矩量法的模拟方法,所述方法将一个电子设备划分成一系列,计算在给定各频率时这些元之间的相应互阻抗,以及按照矩量法从与所述各元分别相应的互阻抗和波源模拟在每个所述元流动的电流,所述方法的特征是它包括下列处理过程:
一个第一过程,所述第一过程通过将一个波源的时间序列数据分段、加以付立叶变换得到一个频谱,将这个频谱的全部或部分和其中所含的频率设定为处理对象;
一个第二过程,所述第二过程根据所述第一过程设定的频率选择采样频率,计算在采样频率的互阻抗,从算得的互阻抗和采样频率产生一个互阻抗的近似表达式,用这个近似表达式计算在所述第一过程设定的各频率的互阻抗;
一个第三过程,所述第三过程按照矩量法从在所述第二过程算得的互阻抗和在所述第一过程设定的频谱求出在每个所述元流动的电流频谱;以及
一个第四过程,所述第四过程对所述第三过程求得的电流频谱进行付立叶反变换,输出反变换结果。
13.一种采用矩量法的模拟方法,所述方法将一个电子设备划分成一系列元,计算在给定各频率时这些元之间的相应互阻抗,以及按照矩量法从与所述各元分别相应的互阻抗和波源模拟在每个所述元流动的电流,所述方法的特征是它包括下列处理过程:
一个第一过程,所述第一过程通过将一个波源的时间序列数据分段、加以付立叶变换得到一个频谱,将这个频谱的全部或部分和其中所含的频率设定为处理对象;
一个第二过程,所述第二过程根据所述第一过程设定的频率选择采样频率,计算在采样频率的互阻抗,从算得的互阻抗和采样频率产生一个互阻抗的近似表达式,用这个近似表达式计算在所述第一过程设定的各频率的互阻抗;
一个第三过程,所述第三过程按照矩量法从在所述第二过程算得的互阻抗和在第一过程设定的频谱求出在每个所述元流动的电流的频谱,同时,如果需要的话,从这电流频谱求出电压频谱、电场频谱和磁场频谱,将求得的这些频谱中所需要的频谱设定为输出对象;以及
一个第四过程,所述第四过程对在所述第三过程设定的需输出的频谱进行付立叶反变换,输出反变换结果。
14.一种按权利要求12或13所述的采用矩量法的模拟方法,其中
在所述第二过程,除了互阻抗外,还用与计算互阻抗相同的程序计算互导纳和/或互电抗;
在所述第三过程,还结合在所述第二过程算得的互导纳和互电抗考虑磁流频谱执行处理;以及
在所述第四过程,执行与所述第三过程配合的付立叶反变换。
15.一种采用矩量法的模拟方法,所述方法将一个电子设备划分成一系列元,计算在给定各频率时这些元之间的相应互阻抗,以及按照矩量法从与所述各元分别相应的互阻抗和波源模拟在每个所述元流动的电流,所述方法的特征是它包括下列处理过程:
一个第一过程,所述第一过程在将一个通过调制一个载波信号而得到的波用作一个波源的情况下,将这个波源分解为三个具有由载波信号频率和调制信号频率决定的不同频率的波源;
一个第二过程,所述第二过程对在所述第一过程分解得的这些波源采用矩量法求出在每个所述元流动的频率域电流;以及
一个第三过程,所述第三过程从在所述第二过程求得的频率域电流计算出时间域电流。
16.一种按权利要求15所述的采用矩量法的模拟方法,其中
在所述第二过程,除了互阻互阻抗外还考虑互导纳和/或互电抗,采用矩量法求出频率域磁流;以及
在所述第三过程,还从在所述第二过程求得的频率域磁流计算出时间域磁流。
CN 97103317 1996-04-10 1997-03-17 利用矩量法模拟电磁场强度的设备和方法 Expired - Fee Related CN1118766C (zh)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
JP088213/1996 1996-04-10
JP8821396 1996-04-10
JP088213/96 1996-04-10

Publications (2)

Publication Number Publication Date
CN1163435A true CN1163435A (zh) 1997-10-29
CN1118766C CN1118766C (zh) 2003-08-20

Family

ID=13936634

Family Applications (1)

Application Number Title Priority Date Filing Date
CN 97103317 Expired - Fee Related CN1118766C (zh) 1996-04-10 1997-03-17 利用矩量法模拟电磁场强度的设备和方法

Country Status (3)

Country Link
US (1) US5903477A (zh)
CN (1) CN1118766C (zh)
DE (1) DE19710687B4 (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN100487490C (zh) * 2001-06-16 2009-05-13 维斯特恩格科地震控股有限公司 处理数据的方法
CN112730994A (zh) * 2020-12-22 2021-04-30 国网天津市电力公司电力科学研究院 基于matlab获取高压交流线路电场的方法及系统
CN113536229A (zh) * 2020-04-13 2021-10-22 富士通株式会社 采样装置、采样方法以及用于存储采样程序的存储介质

Families Citing this family (19)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH1115814A (ja) * 1997-06-26 1999-01-22 Fujitsu Ltd モーメント法を用いたシミュレーション装置及び方法並びにプログラム記憶媒体
JP3405905B2 (ja) * 1997-08-12 2003-05-12 富士通株式会社 電磁界強度算出装置及び方法並びにプログラム記憶媒体
JPH1194889A (ja) * 1997-09-19 1999-04-09 Fujitsu Ltd 多層基板からの放射電磁波解析装置
JP3488629B2 (ja) * 1998-04-07 2004-01-19 富士通株式会社 放射電磁界耐性算出装置及び方法並びにプログラム記録媒体
DE19964239B4 (de) * 1998-04-07 2008-06-26 Fujitsu Ltd., Kawasaki Vorrichtung zum Simulieren eines elektrischen Stroms, der aufgrund einer Funkwelle, die durch eine Antenne abgestrahlt wird, durch eine elektronische Vorrichtung fließt, und Programmspeichermedium
JPH11296499A (ja) * 1998-04-07 1999-10-29 Fujitsu Ltd モーメント法を用いたシミュレーション装置及び方法並びにプログラム記録媒体
JP2000214191A (ja) * 1999-01-27 2000-08-04 Fujitsu Ltd シミュレ―ション装置及び方法並びにプログラム記録媒体
US6289298B1 (en) * 1999-04-01 2001-09-11 Agere Systems Guardian Corp. Method and apparatus for quasi full-wave modeling of interactions in circuits
US6397171B1 (en) * 1999-04-01 2002-05-28 Agere Systems Guardian Corp. Method and apparatus for modeling electromagnetic interactions in electrical circuit metalizations to simulate their electrical characteristics
US6694289B1 (en) * 1999-07-01 2004-02-17 International Business Machines Corporation Fast simulation method for single and coupled lossy lines with frequency-dependent parameters based on triangle impulse responses
JP4001449B2 (ja) * 2000-03-08 2007-10-31 松下電器産業株式会社 不要輻射解析方法
US20020082812A1 (en) * 2000-12-26 2002-06-27 Fujitsu Limited Apparatus and method for simulating the receiving characteristic of radio waves
US7860496B2 (en) * 2003-06-27 2010-12-28 Robert Bosch Gmbh Device and method for simulating a data transmission via a data transmission network
US7222033B1 (en) 2003-08-18 2007-05-22 Steven Lynn Newson Electromagnetic emissions and susceptibility calculating method and apparatus
JP4499781B2 (ja) * 2005-03-28 2010-07-07 富士通株式会社 電磁界強度算出方法、電磁界強度算出装置、制御プログラム
JP2007263789A (ja) * 2006-03-29 2007-10-11 Fujitsu Ltd 電磁波解析プログラム、電磁波解析装置および電磁波解析方法
CN102752060B (zh) * 2012-05-23 2014-04-16 西北工业大学 一种基于微波暗室的动态干扰模拟方法
CN105388368B (zh) * 2015-12-02 2019-02-15 中国电力科学研究院 一种高压架空输电线路电磁散射的阻抗加载点选取方法
CN109460598B (zh) * 2018-10-30 2022-05-03 电子科技大学 电磁场在波导中传播的离散傅里叶逆变换解析验证方法

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US4855677A (en) * 1988-03-11 1989-08-08 Westinghouse Electric Corp. Multiple coil eddy current probe and method of flaw detection
US5270647A (en) * 1992-01-08 1993-12-14 Osaka Gas Company, Ltd. Pipe electromagnetic field simulation apparatus using Born's approximation rule
JP3091815B2 (ja) * 1994-02-25 2000-09-25 富士通株式会社 電磁界強度算出装置
JP3443162B2 (ja) * 1994-05-10 2003-09-02 富士通株式会社 電磁界強度算出装置
JP2768900B2 (ja) * 1994-05-10 1998-06-25 富士通株式会社 電磁界強度算出装置

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN100487490C (zh) * 2001-06-16 2009-05-13 维斯特恩格科地震控股有限公司 处理数据的方法
CN113536229A (zh) * 2020-04-13 2021-10-22 富士通株式会社 采样装置、采样方法以及用于存储采样程序的存储介质
CN112730994A (zh) * 2020-12-22 2021-04-30 国网天津市电力公司电力科学研究院 基于matlab获取高压交流线路电场的方法及系统

Also Published As

Publication number Publication date
DE19710687B4 (de) 2006-02-16
CN1118766C (zh) 2003-08-20
DE19710687A1 (de) 1997-10-16
US5903477A (en) 1999-05-11

Similar Documents

Publication Publication Date Title
CN1118766C (zh) 利用矩量法模拟电磁场强度的设备和方法
CN1323789C (zh) 焊接评估
CN1300573C (zh) 用于无线装置的相对吸收率测定装置
CN1188945C (zh) 放大器的校准装置及校准方法
CN100336066C (zh) 电阻值计算方法
CN1136459C (zh) 检测电力系统失调的方法和装置
CN1902506A (zh) 接收时刻计测装置以及使用该装置的距离计测装置
CN1645360A (zh) 信号处理方法、信号处理程序、记录介质及信号处理装置
CN1078779C (zh) 分集接收机
CN1190963C (zh) 数据处理装置和方法,学习装置和方法
CN1168198C (zh) 倍频电路
CN1488209A (zh) 多径干扰消除设备和多径干扰消除方法
CN1180271C (zh) 并行双输电线故障定位装置和方法
CN1265355C (zh) 音源矢量生成装置及语音编码/解码装置
CN1879103A (zh) 数字群时延补偿器
CN1890683A (zh) 信息处理装置、信息提供装置、信息处理方法及信息提供方法
CN1204187A (zh) 无线电接收系统
CN1801183A (zh) 信息处理装置和方法以及程序
CN1910843A (zh) 光信号质量监控电路以及光信号质量监控方法
CN1723453A (zh) 用于处理声场表现的方法和系统
CN1435019A (zh) 模拟调制法和应用该方法的光学发射器
CN1460227A (zh) 信号处理设备
CN1765099A (zh) 能够容易控制协议消息的通信系统的模拟装置与模拟方法
CN1695408A (zh) 多层印刷电路板、电子设备、安装方法
CN101074987A (zh) 高压电能计量装置综合误差实时在线监测方法及监测设备

Legal Events

Date Code Title Description
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C06 Publication
PB01 Publication
C14 Grant of patent or utility model
GR01 Patent grant
C17 Cessation of patent right
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20030820

Termination date: 20140317