ANSYS中Fourier级数的分解

1 问题提出：

用Ansys作电气设备场分析时，常需要对磁密等物理量做谐波场分析。此时，必须对这些数据要进行Fourier变换，包括Fourier正变换和反变换。

对给定的数据，或者计算得到的数据进行Fourier分解可采用Matlab的FFT指令和Ansys的*Mfourier指令。

下面分别就这两种软件的指令用法进行简述，并通过算例说明。

2.指令：

*MFOURI, Oper, COEFF, MODE, ISYM, THETA, CURVE（计算一个数据的傅立叶的系数或者根据傅里叶系数求傅里叶反变换）

(1) Oper:

Oper ：为FIT，表示根据 MODE, ISYM, THETA, CURVE 求傅里叶系数COEFF。

Oper ：为EVAL，表示根据COEFF, MODE, ISYM, THETA求傅里叶反变换数据CURVE。

(2) COEFF：用来存放计算得到傅里叶系数的数组参数名。

(3) MODE：预期要分解的傅里叶系数的次数的数组参数名。

(4) ISYM： 包含着相应次傅立叶系数的对称信息的数组参数名；

(5) THETA：原始数据对的角度数组参数名。

(6) CURVE：原始数据或者傅里叶反变换得到数据的数组参数名；

注意点

(1) 数组COEFF, MODE, ISYM的长度必须相同且为预期分析的傅里叶次数的两倍，因为每次傅里叶系数需要有两个系数组成(cos和sin前的系数)。

(2) 数组ISYM的第一项一般是0，表示分解出一个直流分量，其实本人试验过取1和-1对结果没有影响；从第二项开始设置ak和bk放置的顺序，0或1表示cos，说明那个系数是ak，-1表示用sin项，表示那个系数是bk。

(3) 数组Mode长度确定的一般原则是：一般由N个数据表示的波形，其最多可得到N/2(N为偶数)或者(N-1)/2(N为奇数)次傅里叶系数，所以数组Mode的长度应小于N/2或者(N-1)/2，否则会报错。

(4) 原创经验：COEFF的至少要大于Mode维数1，否则出错。

(5) 因为Mode要包含直流和ak和bk，所以Mode最好为奇数；对应地，Coeff数组即为偶数。

3 算例验证：

(1) 三角波程序

*DIM,COEFF,Aarry,24            !求24次傅里叶系数

*DIM,MODE,TABLE,24          !傅里叶系数次数的数组

*DIM,ISYM,TABLE,24           !对应次傅里叶级数项的对称性

*DIM,THETA,TABLE,121         !被分析数据点的横坐标，以角度表示，也可不是

*DIM,CURVEI,TABLE,121        !被分析数据点的纵坐标

*DIM,CURVEO,TABLE,121       !傅里叶反变换得到的数据点的纵坐标

*VFILL,MODE(2),RAMP,1,2       !指定预分析傅里叶级数的次数

*VFILL,ISYM(2),RAMP, -1,0       !三角波是关于y轴反对称的，故从第二项开始取-1

*VFILL,THETA(1),RAMP,0,3       !每隔三度选取一个数据点

*VFILL,CURVEI(1),RAMP,0,1/30    !0到90度内数据点

*VFILL,CURVEI(31),RAMP,1,-1/30  !90到270度内数据点

*VFILL,CURVEI(91),RAMP,-1,1/30  !270到360度内数据点

*MFOURI,FIT, COEFF(1),MODE(1),ISYM(1),THETA(1),CURVEI(1) !傅里叶变换

*MFOURI,EVAL,COEFF(1),MODE(1),ISYM(1),THETA(1),CURVEO(1) !傅里叶反变换

(2) 以f(t)=10+5cos(t)+2sin(t)+8cos(5t)-3sin(5t)来验证

IFIN=49                    !谐波次数,(49-1)/2=24

pi=4*atan(1)

N=122                      !采样点数

*DIM,COEFF,Array,IFIN

*DIM,dMODE,TABLE,IFIN

*DIM,ISYM,TABLE,IFIN

*DIM,THETA_F,TABLE,N

*DIM,CURVEI,TABLE,N                  !CURVE INPUT TO PROGRAM

*DIM,CURVEO,TABLE,N

*DO,I,1,N

       THETA_F(I)=(I-1)*360/(N-1)           !每个采样点的角度

       THETA=(I-1)*2*pi/(N-1)

       CURVEI(I)=10+5*cos(THETA)+2*sin(THETA)+8*cos(5*THETA)-3*sin(5*THETA)

*EndDo

MODE(1)=0                               !计算傅里叶系数 

ISYM(1)=1

ISTART=2

*DO,I,ISTART,IFIN,2

        MODE(I)=I/2                     !FILL EVEN INDICIES OF {MODE}

        ISYM(I)=1

*ENDDO

ISTART=3                                !FILL ODD INDICIES OF {MODE}

*DO,I,ISTART,IFIN,2

        MODE(I)=(I/2)-.5     

        ISYM(I)=-1

*ENDDO

*MFOURI,FIT, COEFF(1),MODE(1),ISYM(1),THETA_F(1),CURVEI(1)

*MFOURI,EVAL,COEFF(1),MODE(1),ISYM(1),THETA_F(1),CURVEO(1)

*dim,C,,N                                  !输入与反变换后曲线的误差

    *DO,I,1,N

        C(I)=CURVEO(I)-CURVEI(I)

    *ENDDO

计算结果：

直流分量幅值为10，1次谐波cos的幅值为5，sin的幅值为2，5次谐波cos的幅值为8，sin的幅值为-3，与所验证函数对应。