哪位做过王佳斌的无槽直线电机推力仿真？？？

suifengsan

附上论文，还有自己写的径向充磁永磁体内置时的matlab程序，仿真出来结果很乱，同样的方法径向充磁永磁体外置的B就仿真出来了（这里只是求的气隙处磁密），很奇怪，帮**看一下吧！！！

1. clc;clear all;
   
2. %永磁体内置式
   
3. num = 11; %11次谐波
   
4. g = 0.0025; %气隙宽度
   
5. taup = 0.06;%永磁体宽度
   
6. hm = 0.005; %永磁体径向厚度
   
7. taum = 0.044; %永磁体高度
   
8. hw = 0.0325; %绕组厚度
   
9. Rr = 0.078; %永磁体外径
   
10.                                     
    
11. miur = 1.05;%回复磁导率
    
12. Brem = 1.2;%剩磁         
    
13. Rm = 0.083; %永磁体内径                     
    
14. Rs = 0.118; %电枢内径
    
15. r = 0.08425;%气隙半径
    
16. alphap = taum/taup;
    
17. %利用卡特系数修正Rs
    
18. % taut=0.02;%齿距
    
19. % b0=0.0012;%槽宽
    
20. % g1=g+hm/miur;
    
21. % b=b0/(2*g1);
    
22. % v=(b0/(2*g1*tan(b))-0.5*log(1+b*2))*(4/pi);
    
23. % kc=taut/(taut-v*g1);%卡特系数
    
24. % ge=g+(kc-1)*g1;
    
25. % Rs=Rm+ge;
    
26. 
    
27. for n = 1:num
    
28. 
    
29.         
    
30.         Mn = (2*n-1)*pi/taup;
    
31.         Pn = (4/taup)*Brem*sin(2*n-1)*(pi/2)*alphap;
    
32.         a = 1e-3;
    
33.         fa = 0;
    
34.         fb = 0;
    
35.      for x = Mn*Rr:a:Mn*Rm
    
36.         f1 = besselk(1,x)/(besseli(1,x)*besselk(0,x)+besselk(1,x)*besseli(0,x));
    
37.         f2 = besseli(1,x)/(besseli(1,x)*besselk(0,x)+besselk(1,x)*besseli(0,x));
    
38.         fa = fa+f1;
    
39.         fb = fb+f2;
    
40.       end
    
41.         FAn = Pn/Mn*a*fa;
    
42.         FBn = Pn/Mn*a*fb;
    
43.         C1n = besseli(0,Mn*Rs);
    
44.         C2n=besselk(0,Mn*Rs);
    
45.         C3n = besseli(0,Mn*Rr);
    
46.         C4n=besselk(0,Mn*Rr);
    
47.         C5n = besseli(0,Mn*Rm);
    
48.         C6n=besselk(0,Mn*Rm);
    
49.         C7n = besseli(1,Mn*Rm);
    
50.         C8n=besselk(1,Mn*Rm);
    
51.         a1n = (FAn*C5n/C6n+FBn)/(miur*(C5n/C6n-C1n/C2n)-(C5n/C6n-C3n/C4n));
    
52.         a2n = (FAn*C7n/C8n-FBn)/((C7n/C8n+C1n/C2n)-(C7n/C8n+C3n/C4n));
    
53.         b1n = C1n/C2n*a1n;
    
54.         b2n = C3n/C4n*a2n;
    
55.         zn = 0;
    
56.         for z = 0:1e-3:2*taup
    
57.             zn = zn+1;
    
58.             B1z(n,zn) = (a1n*besseli(0,Mn*r)-b1n*besselk(0,Mn*r))*sin(Mn*z);
    
59.             B1r(n,zn) = -((a1n*besseli(1,Mn*r)+b1n*besselk(1,Mn*r))*cos(Mn*z));
    
60.         end
    
61.         B1zn=sum(B1z);B1rn=sum(B1r);
    
62. end
    
63. 
    
64. 
    
65. z = 0:1e-3:2*taup;
    
66. figure(1)
    
67. plot(z,B1zn,'g',z,B1rn,'b')
    
68. title('气隙处沿z方向与r方向磁密大小');
    
69. xlabel('Z');
    
70. ylabel('B');
    
71. legend('z','r');

复制代码

suifengsan

论文忘记发了 题目是：A General Framework for the Analysis and Design of Tubular Linear Permanent Magnet Machines

suifengsan

比较经典的一篇论文

沙隆巴斯

内置式的没有考虑硅钢片磁化曲线？

suifengsan

论文里应该有考虑这个，我就想按照论文里的公式复现一遍，发现现在高词谐波的幅值越来越大 ，正常高次谐波应该幅值很小才对，看了很多天没有找到原因，拜托大佬帮我看看再。
可有偿帮忙！！