电磁场数值计算上机作业.docx

电磁场数值计算上机作业.docx

《电磁场数值计算上机作业.docx》由会员分享，可在线阅读，更多相关《电磁场数值计算上机作业.docx（16页珍藏版）》请在冰豆网上搜索。

电磁场数值计算上机作业

电磁场数值计算上机作业报告

一、有限差分法及原理

有限差分法基本思想是把连续的定解区域用有限个离散点构成的网格来代替,这些离散点称作网格的节点；把连续定解区域上的连续变量的函数用在网格上定义的离散变量函数来近似；把原方程和定解条件中的微商用差商来近似,积分用积分和来近似,于是原微分方程和定解条件就近似地代之以代数方程组,即有限差分方程组,解此方程组就可以得到原问题在离散点上的近似解.然后再利用插值方法便可以从离散解得到定解问题在整个区域上的近似解.在采用数值计算方法求解偏微分方程时,若将每一处导数由有限差分近似公式替代,从而把求解偏微分方程的问题转换成求解代数方程的问题。

求解拉普拉斯方程：

为简单起见，将场域分成足够小的正方形网格，网格线之间的距离为h，

。

节点0、1、2、3、4上的电位分别用

、

、

、

和

表示。

点1、点3在x0处可微，沿x方向在x0处的泰勒级数展开式为

点2、点4在y0处可微，沿y方向在y0处的泰勒级数展开式为

忽略高次项

稍作变化得到拉普拉斯方程的五点差分格式：

利用超松弛迭代法求解以上差分方程，二维场拉普拉斯方程等距剖分差分格式公式为：

其中

为超松驰因子，

.

计算流程如图1所示：

图1超松驰迭代法计算流程

利用上述方法求解两题。

二、第一题

计算长直接地金属槽中的电场分布。

金属槽横截面如图所示，其侧壁与底面电位均为零，顶盖电位相对值为10。

槽内电位函数满足拉普拉斯方程。

计算槽内电位分布。

图2题1示意图

要求：

（1）先用正方形网格粗分，每边取4个网格计算，取不同的松弛因子，比较其收敛速度。

取计算精度为千分之一。

（2）划分网格加倍，计算电位分布，并与上面计算结果比较。

（1）用正方形网格，每边取4个网格计算。

图3每边4个网格分割

网格分割如图3所示：

图44个网格分割计算结果

取α=1，计算得结果如图4：

图54个网格分割计算结果数值表示

具体数值表示结果如图5所示：

表1不同超松驰因子α下的迭代次数

取不同的超松驰因子α，得不同的收敛速度如下表：

超松驰因子α

1

1.1

1.2

1.3

1.4

1.5

1.6

1.7

1.8

1.9

循环次数

14

11

8

10

12

16

19

27

44

82

可见随着α从1到2的增加，收敛速度先变快后变慢。

在α=1.2时收敛速度最快。

图68个网格分割计算结果

（2）划分网格加倍，计算结果如下图所示

具体数值表示为：

与划分为4个网格的结果做比较有下表：

X轴

Y轴

0.00

0.25a

0.5a

0.75a

a

8格

4格

8格

4格

8格

4格

8格

4格

8格

4格

a

10.00

10.00

10.00

10.00

10.00

10.00

10.00

10.00

10.00

10.00

0.75a

0.00

0.00

4.31

4.29

5.36

5.27

4.31

4.29

0.00

0.00

0.5a

0.00

0.00

1.84

1.87

2.50

2.50

1.84

1.87

0.00

0.00

0.25a

0.00

0.00

0.69

0.71

0.96

0.98

0.69

0.71

0.00

0.00

0.00

0.00

0.00

0.00

0.00

0.00

0.00

0.00

0.00

0.00

0.00

表2相同节点位置不同网格划分计算数值对比

可见在拟合过程中，划分网格较少时，距离y=a较近处，数据偏小，而距离y=0较近处，数据偏大。

可见划分网格较多时可以使得到的电势分布趋于平滑。

三、第二题

图8题2示意图

计算电机电枢槽气隙磁位分布：

如图所示，忽略曲度效应，采用直角坐标系，计算二维恒定场。

令定子和转子磁铁中μ=∞，EF是对称线，考虑对称性有：

求G1（ABCDEF）内部和边界上全部节点磁位，计算精度取ω≤1/1000。

设g=8，h=66，BC=6，DE=5。

已知磁位

也满足拉普拉斯方程，在所求区域内部，利用超松弛迭代法给出二维场拉普拉斯方程等距剖分差分格式公式为：

在AB处，有

，转化为差分形式可得

，所以在AB边界处的差分形式方程为：

同理在EF处满足的差分形式方程为：

图9题2计算结果

按正方形网格划分，网格大小为1×1，根据对称性，得整个区域的磁位如图9：

在G1区域中，在每个节点用数值表示为（放大可看具体数值）：

图10题2计算结果数值表示

至此求得全部磁位。

附：

matlab代码

第一题：

clc;

x=5;

y=5;

foralpha_test=10:

19;

Fai_Last=zeros（x,y）;

Fai_Now=zeros（x,y）;

Fai_Last（:

y）=10;

Fai_Now（:

y）=10;

a=1;

alpha=alpha_test/10;

fori=2:

（x-1）

forj=2:

（y-1）

Fai_Now（i,j）=Fai_Last（i,j）+alpha/4*（Fai_Last（i,j+1）+Fai_Last（i+1,j）...

+Fai_Now（i,j-1）+Fai_Now（i-1,j）-4*Fai_Last（i,j））;

end

end

k=1;

while（max（max（abs（（Fai_Last-Fai_Now）./Fai_Last）））>0.001&&k<10000）

Fai_Last=Fai_Now;

fori=2:

（x-1）

forj=2:

（y-1）

Fai_Now（i,j）=Fai_Last（i,j）+alpha/4*（Fai_Last（i,j+1）+Fai_Last（i+1,j）...

+Fai_Now（i,j-1）+Fai_Now（i-1,j）-4*Fai_Last（i,j））;

end

end

k=k+1;

end

k_test（alpha_test-9）=k;

end

lx=linspace（0,a,x）;

ly=linspace（0,a,y）;

figure

（1）;

surf（lx,ly,Fai_Now）;

xlabel（'x轴/a','FontName','微软雅黑','FontSize',20）;

ylabel（'y轴/a','FontName','微软雅黑','FontSize',20）;

zlabel（'z轴/V','FontName','微软雅黑','FontSize',20）;

title（'电势φ分布','FontName','微软雅黑','FontSize',20）;

set（gca,'XTick',0:

1/（x-1）:

1,'FontSize',12）;

set（gca,'YTick',0:

1/（y-1）:

1,'FontSize',12）;

holdon

figure

（2）;

[x_s,y_s]=meshgrid（0:

x-1）;

x_s=x_s./（x-1）;

y_s=y_s./（y-1）;

plot（x_s,y_s,'k',x_s',y_s','k'）;

set（gca,'XTick',0:

1/（x-1）:

1,'FontSize',12）;

set（gca,'YTick',0:

1/（y-1）:

1,'FontSize',12）;

xlabel（'x轴/a','FontName','微软雅黑','FontSize',20）;

ylabel（'y轴/a','FontName','微软雅黑','FontSize',20）;

title（'电势φ分布','FontName','微软雅黑','FontSize',20）;

holdon

fori=1:

x

forj=1:

y

text_show=strcat（'φ=',num2str（Fai_Now（i,j）））;

text（（i-1）/（x-1）,（j-1）/（y-1）,text_show,'color','r','FontSize',13）;

end

end

holdoff

xlswrite（'homework.xls',rot90（Fai_Now））;

第二题：

clc;

h=66;

g=8;

BC=6;

DE=5;

length_x=BC+DE;

length_y=h+g;

div_x=length_x+1;

div_y=length_y+1;

alpha_test=10;

Fai_Last=zeros（div_x,div_y）;

Fai_Now=zeros（div_x,div_y）;

Fai_Last（:

div_y）=100;

Fai_Now（:

div_y）=100;

a=1;

alpha=alpha_test/10;

fori=1:

（div_x）

forj=1:

（div_y）

if（i>1&&ih+1&&j

Fai_Now（i,j）=Fai_Last（i,j）+alpha/4*（Fai_Last（i,j+1）+Fai_Last（i+1,j）...

+Fai_Now（i,j-1）+Fai_Now（i-1,j）-4*Fai_Last（i,j））;

end

if（i>BC+1&&i1&&j<=h+1）

Fai_Now（i,j）=Fai_Last（i,j）+alpha/4*（Fai_Last（i,j+1）+Fai_Last（i+1,j）...

+Fai_Now（i,j-1）+Fai_Now（i-1,j）-4*Fai_Last（i,j））;

end

if（i==div_x&&j>1&&j

Fai_Now（i,j）=Fai_Last（i,j）+alpha/4*（Fai_Last（i,j+1）+Fai_Last（i-1,j）...

+Fai_Now（i,j-1）+Fai_Now（i-1,j）-4*Fai_Last（i,j））;

end

if（i==1&&j>h+1&&j

Fai_Now（i,j）=Fai_Last（i,j）+alpha/4*（Fai_Last（i,j+1）+Fai_Last（i+1,j）...

+Fai_Now（i,j-1）+Fai_Now（i+1,j）-4*Fai_Last（i,j））;

end

end

end

k=1;

while（max（max（abs（（Fai_Last-Fai_Now）./Fai_Last）））>0.001&&k<10000）

Fai_Last=Fai_Now;

fori=1:

（div_x）

forj=1:

（div_y）

if（i>1&&ih+1&&j

Fai_Now（i,j）=Fai_Last（i,j）+alpha/4*（Fai_Last（i,j+1）+Fai_Last（i+1,j）...

+Fai_Now（i,j-1）+Fai_Now（i-1,j）-4*Fai_Last（i,j））;

end

if（i>BC+1&&i1&&j<=h+1）

Fai_Now（i,j）=Fai_Last（i,j）+alpha/4*（Fai_Last（i,j+1）+Fai_Last（i+1,j）...

+Fai_Now（i,j-1）+Fai_Now（i-1,j）-4*Fai_Last（i,j））;

end

if（i==div_x&&j>1&&j

Fai_Now（i,j）=Fai_Last（i,j）+alpha/4*（Fai_Last（i,j+1）+Fai_Last（i-1,j）...

+Fai_Now（i,j-1）+Fai_Now（i-1,j）-4*Fai_Last（i,j））;

end

if（i==1&&j>h+1&&j

Fai_Now（i,j）=Fai_Last（i,j）+alpha/4*（Fai_Last（i,j+1）+Fai_Last（i+1,j）...

+Fai_Now（i,j-1）+Fai_Now（i+1,j）-4*Fai_Last（i,j））;

end

end

end

k=k+1;

end

k_test（alpha_test-9）=k;

Fai_Final=zeros（（BC+DE）*2+1,h+g+1）;

Fai_Final（1:

BC+DE+1,1:

h+g+1）=Fai_Now;

Fai_temp=fliplr（Fai_Now'）';

Fai_Final（（BC+DE+1）:

（BC+DE）*2+1,1:

h+g+1）=Fai_temp;

lx=linspace（0,（BC+DE）*2+1,（BC+DE）*2+1）;

ly=linspace（0,h+g+1,h+g+1）;

figure

（1）;

colormap（'jet'）;

surf（lx,ly,Fai_Final'）;

set（gca,'XTick',1:

（BC+DE）*2+1,'FontSize',3）;

set（gca,'YTick',1:

5:

h+g+1,'FontSize',3）;

xlabel（'x轴','FontName','微软雅黑','FontSize',20）;

ylabel（'y轴','FontName','微软雅黑','FontSize',20）;

zlabel（'z轴','FontName','微软雅黑','FontSize',20）;

title（'磁标势φ分布','FontName','微软雅黑','FontSize',20）;

holdon

figure

（2）;

[x_s,y_s]=meshgrid（0:

（BC+DE）*2+1,0:

h+g+1）;

plot（x_s,y_s,'k',x_s',y_s','k'）;

xlabel（'x轴','FontName','微软雅黑','FontSize',20）;

ylabel（'y轴','FontName','微软雅黑','FontSize',20）;

title（'磁标势φ分布','FontName','微软雅黑','FontSize',20）;

axis（[0（BC+DE）*2+10h+g+1]）

holdon

fori=1:

（BC+DE）*2+1

forj=1:

h+g+1

text_show=strcat（'φ=',num2str（Fai_Final（i,j）））;

text（（i-1）,（j-1）,text_show,'color','r','FontSize',5）;

end

end

holdoff