|
楼主 |
发表于 2007-6-20 17:07:19
|
显示全部楼层
更详细的求助了
请帮忙用有限元法做个同步电机电机磁场的模型相关程序(基于MATLAB)如下:% 1:air-gap, 2:culasse stator, 3: fer rotor, 8 H$ @* X4 @' w- ?2 H2 n
% 4-7: stator conductor, 28 - 29: rotor conductor
8 ~, x% A- a, L9 {3 _ \# g8 `) _6 p# }% H, e% L
tic; clear;clc; close all;
! U6 y. P) D& y, n1 L/ a) ^3 R1 w. {0 R! l+ @
. c# [, y T9 Y& L8 m. } N% D5 H& \
. h1 z. w, E( M# C2 w* d) v
% ------------- Draw the geometry of stator, rotor et air-gap --------------) Z) w5 \- q$ i% T
- o, W4 n. M C2 s2 u! Aentref=1e-3; % air-gap depth
; u/ h3 Y/ z0 z X4 w9 [' W; PRs=39.385e-3; % radius of (rotor+airgap)
/ b5 ^3 Y/ z' |* N! J8 {Rr=Rs-entref; % radius of rotor
& J/ e/ j. w* R9 fRc=71.75e-3; % radius of stator
8 I6 `/ V; {9 D" H% g
8 Q' ~" [! g$ q+ H/ ^Nb=140; % number of turns per phase
/ R$ F4 z, a5 w* G" ]; lLong=0.125; % longitude of machine in rotation axis direction
$ U$ h5 o9 Z& q% s8 d, c1 s! [7 J/ Op=1; % ??
7 z0 d4 D6 {+ ef=50; % ??
% v8 F, V! U/ X
( C) y+ a& g: a) y$ u2 D2 gbeta=35*pi/180; beta=beta/p; % culasse rotor
6 N" S e+ D7 a8 ]( Qbeta1=50*pi/180; beta1=beta1/p; % epanouissement polaire5 |# G6 b; `* N
0 ^+ P' ]& y: w/ X0 I) e, \' XdessinP1; % a sub-programme to draw the geometry of the machine section with the parameters defined above
$ p, t& H/ O. v, L% --------------* K# p, N' S# d T3 D$ w3 q
) O2 U7 O4 }2 f2 s5 o0 \% k' W; v# w- z' l* o: u" u8 T
% ------------- Rotor positioning and initial mesh --------------
, Y" m# {/ _+ X2 m* a; g
$ _0 e8 b3 [. i1 [9 k. Q3 Yaa=-30*pi/180; % initial angle of rotor
: t. I5 b, ~* lgr=groto(gr,(aa)); % rotating the rotor to the angle specified above using the sub-programme 'groto'
% s8 b9 Y4 Q u8 j, `1 Qg=[gsf gr]; % stator + rotor given the whole machine geometry
2 m& ], K( w5 _# ^limites; % a sub-programme defining the bondary condition
( B5 y6 m! B- `& e: K) p0 q[p,e,t]=initmesh(g); % initial mesh for the partial differential equation resolution+ x* v5 {" W; o9 V/ e! c
$ x$ q/ Y$ Q T3 ^8 E; m % ------------- mesh refine if necessary ------------
5 a y; x% Y( m. M% [p,e,t]=refinemesh(g,p,e,t,[1 30]);
2 H4 e" L* A6 n% pdegplot(g); axis equal, hold on ; pdemesh(p,e,t), hold off) a: g+ G7 D- t6 a
% pause
; V4 f/ Z4 x) I. Z% -------------------------------------8 l% S' d0 S& f7 c0 x% {9 \) H
$ p g, z+ D5 i! |9 z7 ]0 w' z- k3 S* M$ O6 D, D
, k: v: Z# i2 f( O( S/ J+ {
% ------------- specification of parameters required for calculation -------" ]( Z# K: t/ n6 ~
kexc=0; 5 ?$ |+ j5 h7 @% l; N; z
sigmaexc=kexc*5e6; % conductivity of rotor coil conductor0 C3 r$ x. y1 P1 Y4 J2 U
nuo=1/(4e-7*pi); % inverse of permeability of vacuum ; L7 k8 @" i9 E' O
murs=500; murr=500; % relative permeability of stator and rotor
. Z; w2 f$ }+ M. F
5 t; I8 P; o. V7 C" A ks=1; kr=0; % ks=1 or 0: stator excited or not; kr=1 or 0: rotor excited or not.
# D$ ]$ v! U$ y6 e" h1 I
; l$ s3 D) h. D0 ^9 z Jex1=0; Jex11=0;( z2 L/ b, T, e! i" r
Jex2 = 10e6*ks; Jex21 = -Jex2;
( v; @/ o$ A" F/ M6 f6 z/ w/ l Jex3 = 0; Jex31 = 0;
) {9 r2 ~6 `2 B% Jex3=10e6*ks; Jex31=-Jex3; % current density of phase a (stator)+ ]" u; W6 a% X
% Jex1=-Jex3/2; Jex11=-Jex1; % current density of phase b (stator)
% ]" `4 d# R: S/ L) }4 e% Jex2=-Jex3/2; Jex21=-Jex2; % current density of phase c (stator)2 g) @6 \" C5 M/ F1 G
Jexr1=10e6*kr; Jexr2=-Jexr1; % current density of rotor excitation
; q% T. J" A% W
: f3 q D6 s+ x/ X) T garnissageP1; % a sub-programme for the giving the current density and permeability at each mesh element
& r' _% y A& k8 w% -------------
+ p5 A: Z4 ?# r6 f3 W' }) X, b; K3 U7 B( p( A2 O3 }2 X
( D+ z+ D* o1 J/ l6 W$ z- C' ^9 I9 d6 |+ J, L
% ------------- resolution of the partial differential equation and figure ploting ---------
4 h+ I: _2 M, @! @5 n4 w8 _2 S4 P6 M
u=assempde(bond,p,e,t,nu,0,J); % resolution of the partial differential equation
0 m5 W' N9 N2 v: o% U=ASSEMPDE(B,P,E,T,C,A,F) assembles and solves the PDE problem -div(c*grad(u))+a*u=f2 j% n, c0 l1 p8 m0 C
5 A& D6 H+ o% @" cfigure (1), pdemesh(p,e,t), axis equal, title ('Mesh of the machine section')
5 |; i% c1 K) zfigure (2), pdegplot(g); axis equal; hold on; pdecont(p,t,full(real(u)),30), hold off, title('Magnetic field distribution'), % plots using 30 levels.1 s8 o7 o D- o8 q& ?
% ------------------------------* O- H9 s* ^; @+ f6 h0 i
4 T! k! Z4 B/ `% l1 D' q$ @9 |% ], U
6 { q0 F6 v: f1 r, V6 y1 W6 F# ^% ]. B% u
% -------------- calculation of the energy stocked in the machine (energy sum in the magnetic field) --------------------
9 m9 j' Q" F, G% o4 r. I, w
: s! t1 c9 x" k1 x$ o% [ux,uy]=pdegrad(p,t,u);
, O6 K0 K7 R( f$ y! l% bx=uy;
# h$ s1 `# ?( U8 X* j% by=-ux;4 n/ ^% f4 I; Y* t& y5 U; d
% for i=1:length(bx);
! u* G2 ^1 p6 ] S4 O% b(i)=(bx(i)^2+by(i)^2)^.5;; F7 h% C6 c5 c6 J
% end
6 W/ X! c& C. T3 l9 W% nu_e=1/(4e-7*pi);% U/ Y5 i0 @% J9 g% r1 S
% nu_c=nu_e;
7 F" f$ ^' x& C% nu_s=1/500*nu_e;
! w. c' H f$ `. Q5 k% nu_r=1/500*nu_e;
1 M' X& y$ q1 s% M& m% h=sparse(1,ntrg);
( i8 k) K( i5 r$ @% ind_e=find(t(4,:)==2);
6 n+ h9 Q2 R2 F* X# ?% ind_s=find(t(4,:)==1);
2 n- h2 t' B2 x5 e% ind_r=find(t(4,:)==3);* I) a; }6 i/ k" Y
% ind_c=find((t(4,:)>=4) & (t(4,:)<=27));8 c5 K+ L7 N% K3 e4 Y% G+ ]
% h(ind_e)=nu_e*b(ind_e)';: z. ] B+ e; l5 g: a6 r4 t* a
% h(ind_s)=nu_s*b(ind_s);
& `4 O# C& A G. }% h(ind_r)=nu_r*b(ind_r);
5 p. L( u9 _/ |7 o, b% h(ind_c)=nu_c*b(ind_c)'; A$ ~5 \+ R% W% x
% aire=pdetrg(p,t);5 f' a& M2 d- i! D9 N/ q+ x, o
% vol_trg=Long*aire;& w0 L x* Z8 i8 z( P- w- T1 H( v
% energ_elem_e=1/2*(h(ind_e).*b(ind_e)').*vol_trg(ind_e);
) j0 d( H2 H K4 t% energ_tot_e=sum(energ_elem_e);
u: U# J# b+ I: f' b% energ_elem_s=1/2*(h(ind_s).*b(ind_s)').*vol_trg(ind_s);
( {1 a" x$ V2 n+ E: g% energ_tot_s=sum(energ_elem_s);6 a& f! ?6 E' S$ x6 ]0 W- p
% energ_elem_r=1/2*(h(ind_r).*b(ind_r)').*vol_trg(ind_r);2 F$ P/ y: N* i" V" D. ?
% energ_tot_r=sum(energ_elem_r);
- P$ z- R, h# h3 g$ w% energ_elem_c=1/2*(h(ind_c).*b(ind_c)').*vol_trg(ind_c);
1 Q4 f; a6 [: @0 S% energ_tot_c=sum(energ_elem_c);
, M# @2 Y! r u% N% energ_total=energ_tot_e+energ_tot_s+energ_tot_r+energ_tot_c;/ A' m6 H6 H- r$ u, C" A: d5 h D
% enr=energ_total;
- ]0 s" @ ^- S! n) [) B1 c; b, e
toc; |
|