W ojciech SZELĄG 1*, Piotr SUJKA21
ANALIZA NIEUSTALONEGO POLA MAGNETYCZNEGO W PRZETWORNIKACH ELEKTROMAGNETYCZNYCH Z UWZGLĘDNIENIEM ZJAWISKA HISTEREZY
Streszczenie. W artykule przedstawiono algorytm analizy dwuwymiarowego nieustalonego pola magnetycznego w przetwornikach elektromagnetycznych z uwzględnieniem zjawiska histerezy. Do wy
znaczania rozkładu i przebiegu pola wykorzystano metodę elementów skończonych oraz metodę kolej
nych kroków czasowych. Zjawisko histerezy odwzorowano za pomocą modelu Prelsacha. Opracowany program komputerowy wyznaczania rozkładu pola z uwzględnieniem zjawiska histerezy wykorzystano do symulacji wymuszanego napięciowo pola magnetycznego w osiowosymetrycznym dławiku kubkowym.
ANALYSIS OF ELECTROMAGNETIC FIELDS IN ELECTROMAGNETIC CONVERTERS TAKING INTO ACC O U NT HYSTERESIS PHENOMENON
Sum m ary. In the paper the algorithm of the electromagnetic field analysis in electromagnetic con
verters has been presented. The phenomenon of hysteresis has been taken into account by means of the B oriented Preisach model. For analysis of the electromagnetic field the finite element method and step by step algorithm have been used. The elaborated program has been used for simulation voltage- excited electromagnetic field in the axisymmetrical reactor with the ferrite core
Key w o rd s : electromagnetic field, Boriented Preisach model, finite element method
1. W S TĘP
Do dokładnej analizy zarówno ustalonych, ja k i nieustalonych stanów pracy przetworników e lektrom agnetycznych i elektrom echanicznych oraz do wyznaczania ich param etrów funkcjonalnych potrzebny je s t precyzyjny m odel zjaw isk elektrom agnetycznych. Rozkład i przebieg pola m agne
tycznego w przetw ornikach zależy m .in. od w łaściw ości m agnetycznych użytych do ich budowy m ateriałów ferrom agnetycznych. W łaściw ości m agnetyczne ferrom agnetyków w ynikają z rodziny pętli histerezy. Zjawisko histerezy je st w przypadku m ateriałów magnetycznie miękkich zazwyczaj niekorzystne. Związane są z nim histerezowe straty m ocy w rdzeniach przemagnesowywanych polem przemiennym. Histereza m agnetyczna materiału wpływa ponadto na odkształcenie przebiegów indu
kowanych w uzwojeniach napięć oraz płynących w nich prądów. Pozostałość magnetyczna, w prze
twornikach działających na zasadzie zm iany lepkości cieczy magnetoreologicznej pod wpływem pola magnetycznego, powoduje natom iast powstanie pasożytniczego momentu oporowego [6],
Zjaw isko histerezy m agnetycznej wykorzystuje się w m agnesach trwałych. Magnesy wykonuje się z m ateriałów charakteryzujących się pętlą histerezy o dużym natężeniu powściągającym i dużej pozostałości m agnetycznej. Dla współczesnych m agnesów natężenie powściągające przekracza 1000kA/m, a pozostałość m agnetyczna 1,4T. M ateriały charakteryzujące się histerezą magnetyczną s ą w ykorzystyw ane także do zapisu obrazu i dźwięku oraz stosuje się je jako elem enty pamięciowe w technice kom puterowej.
Przy projektowaniu i przy analizie stanów pracy przetworników elektromechanicznych i elektroma
gnetycznych właściwości m agnetyczne ferromagnetyków uwzględniane są zazwyczaj w sposób przy
bliżony. Pomija się zjawisko histerezy i korzysta z jednowartościowej krzywej magnesowania. Podej
ście takie z uwagi na m ałą dokładność je st niewystarczające np. przy analizie stanów nieustalonych w obwodach magnetycznych z podmagnesowaniem czy przy wyznaczaniu strat m ocy w rdzeniu.
11 Dr hab. inż., Politechnika Poznańska, Instytut Elektrotechniki Przemysłowej, ul. Piotrowo 3a, 60-965 Poznań, tel.
061 6652116, fax 061 6652381, szelag@solput.poznan.pl
21 Mgr inż., Politechnika Poznańska, Instytut Elektrotechniki Przemysłowej, ul. Piotrowo 3a, 60-965 Poznań, Piotr.Sujka@domdata.depfa-it.com
168 Szeląg W., Sujka P.
W artykule do odwzorow ania zjaw iska histerezy w ykorzystano m odel Preisacha. M odel ten je st przydatny, np. w procesie w yznaczania rozkładu pola m agnetycznego za p om ocą m etody elem en
tó w skończonych, do określania m agnetyzacji m ateriału ferrom agnetycznego w elem encie dyskre- tyzującym na podstaw ie indukcji m agnetycznej B [1, 7].
2. M O D EL PR E ISA C H A
Z ale żn o ść w ektora indukcji m agnetycznej B od wektora natężenia pola m agnetycznego H oraz w ektora m agnetyzacji H, w ferrom agnetyku opisuje równanie:
B = p 0(H + H |). (1)
P rzy rozpatryw aniu zjaw isk elektrom agnetycznych z uwzględnieniem zjaw iska histerezy przyj
m uje się najczęściej, że w ektory te s ą rów noległe [2, 3] i korzysta się z rodziny pętli histerezy B = fBH(H) m ateriału. Przy powyższych założeniach do m atem atycznego opisu histerezy m ożna w ykorzystać skalarny m odel Preisacha [3, 4].
O dw zorow anie histerezy i „pam ięć globalną” w m odelu Preisacha uzyskuje się przez superpo
zycję nieskończonej liczby elem entarnych operatorów histerezy yap = Ya,p{u(t)}- O perator ta ki może być przedstaw iony za p o m o cą prostokątnej pętli histerezy - rys. 1. Przez a i p oznaczono w artość sygnału w ejściow ego u, przy którym następuje skokow a zm iana sygnału w yjściow ego elem entarne
go operatora. S ygnał w yjściow y m oże przyjm ow ać tylko dwie w artości: +1 i - 1 . Jeśli syg n a ł w ej
ściow y je s t m onotonicznie rosnący, to sygnał w yjściow y zm ienia się zgodnie z g a łę zią abcde, nato
m iast je śli je s t m onotonicznie m alejący, to zm ienia się on według gałęzi edfba (rys. 1). O perator cechuje się zatem p a m ię cią lokalną. Dla każdego operatora ya p określa się w m odelu funkcję wagi
p (a ,p ), z ja k ą dany operator w pływ a na wypadkow y sygnał w yjściow y f(t) z m odelu. S ygnał w yj
ściow y f(t) wyznacza się z równania
f ( t ) = J J n ( a ,p ) Yo, p { u ( t ) } d a d p . ( 2)
a 2(3
Parę liczb (a ,p ) m ożna traktow ać ja ko punkt na płaszczyźnie a, 3, leżący w obszarze trójkąta zaw artego m iędzy p ro stą a = p oraz prostym i p = p0 i a = a 0 , przy czym p0 = - a 0 (rys. 2). Dla u ( t) ^ p 0 układ znajduje się w stanie nasycenia ujemnego, a dla u ( t ) ź a 0 w stanie nasycenia do
datniego. Każdem u punktow i P(a,p) z obszaru T trójkąta A B C je s t przyporządkow any operator ya p . W artości sygnału w ejściow ego, przy którym następuje „przełączenie” tego operatora, s ą równe w spółrzędnym a i p tego punktu.
A (a 0,p0)
Rys. 2. Trójkąt Preisacha Fig. 2. The Preisach triangle
W zależności od przebiegu sygnału wejściowego u(t) trójką t zostaje podzielony na dwa obszary: a) S ~ ( t) , w którym sygnał w yjściow y każdego z opera
torów ya| Pi je s t równy -1 oraz b) S+( t ) , w którym sygnał w yjściow y każdego z operatorów y0| p je st rów ny 1. O bszary S~ oraz S+ oddzielone s ą linią schodkow ą. Położenie w ierzchołków tej linii zm ienia się w czasie. Ich odcięte i rzędne pokryw ają się odpow iednio z lokalnym i w artościam i m inim alnym i i m aksym alnym i sygnału wejściow ego w kolejnych chw ilach czasowych. U w zględniając podział obszaru trójkąta Preisacha, z rów nania (2) uzyskuje się:
f ( t) = [ f p ( a , p ) d a d p - JJp(<x,p)d<xdp. (3)
s-(t) S-(t)
W praktyce w celu w yznaczenia sygnału f(t) dogod
nie je s t rów nanie (3) przedstaw ić w postaci:
Rys. 3. Podział trójkąta Preisacha Fig. 3. The division of Preisach triangle
f (0 = 2 || p ( a , ß ) d a d ß - F ( a 0,ß0 ) s*(t)
gdzie: F(a0 , ß0) = | | p(a, ß) d a d ß - całka po obszarze trójkąta Preisacha T.
(4)
Bio rą c pod uwagę, że całka po obszarze S +(t) je s t równa sum ie całek po obszarach Q trapezów (rys. 3)
||p ( a . ß ) d a d ß = ^ | | p(a, ß )da dp s-(t) k=1 Q ,(t)
oraz że całka po obszarze k-tego trapezu je s t równa różnicy całek po obszarach trójkątów EGE' oraz F G F
| | |i(a ,ß )d a d ß = J J p (a ,ß )d a d ß - J | p (a,ß )dadß = F(Mk,m k_1) - F ( M k,m k)
Q ,(t) EGF' EGF'
z zależności (4), przy w zrastającym sygnale wejściowym , uzyskuje się:
f(t) = - F(a0,p0)+ 2 2 >’[F(Mn.mk. 1)-F (M k,mk)] + 2[F(Mn.mn_1) -F (M n,u(t)l
(5)
(6)
przy czym : n(t) - liczba trapezów.
N atom iast, je śli sygnał u(t) je s t sygnałem m alejącym to [3]:
f ( t ) = - F(a0,p0) + 2 x ' | F ( M k,m k. 1)-F (M k ,mk) ]+ 2F(u(t),m n_1) .
(7)
(8)
k=1
W pow yższych zależnościach przez Mk oraz mk oznaczono odpowiednio k-te maksimum i k-te m inim um przebiegu u(t).
W artości całek F(Mk ,m k_1),F (M k ,m k ) w zależności (6) wyznacza się doświadczalnie w sposób p rzedstaw iony szczegółow o w pracach [1, 3, 7],
W opracow anym algorytm ie m odelow ania histerezy, w celu wyznaczenia sygnału wyjściowego f(t) zgodnie z za le żn o ścią (7) I (8) korzysta się z w artości funkcji F (a,,Pj) wyznaczonych dla w y
branych punktów w obszarze trójkąta T. Przyjęto, że punkty te leżą w węzłach siatki dyskretyzują-
170 Szeląg W., Sujka P.
cych trójkąt Preisacha. Siatkę tw o rzą rów nom iernie rozłożone linie równoległe do osi a i do osi p.
W artości funkcji F w węzłach siatki zapam iętuje się w tablicy. Na ich podstawie, w ykorzystując liniow ą aproksym ację funkcji, oblicza się w artość funkcji F dla dowolnego punku w obszarze trójkąta Preisacha. Ze w zrostem gęstości siatki zw iększa się dokładność odwzorowania histerezy.
W opracow anym algorytm ie w yznacza się i zapam iętuje w artości ekstrem alne (M i,m i), ... , (Mn,mn) sygnału w ejściow ego u(t) dla kolejnych chwil czasowych. Z zapam iętanego ciągu m inimów i m aksim ów elim inuje się przy tym „nieaktualne” w artości ekstrem ów [3], S ygnał w yjściow y f(t) oblicza się w program ie z zależności (7) i (8).
Przy m odelow aniu zjaw iska histerezy w ferrom agnetyku za p om ocą klasycznego m odelu Prei
sacha sygnałem w ejściow ym je s t natężenie pola m agnetycznego u(t) = H(t), a sygnałem w yjścio
w ym indukcja m agnetyczna f(t) = B(t). T aki m odel zjaw iska histerezy je s t niewygodny np. przy wy
znaczaniu rozkładu pola m eto d ą elem entów skończonych. W yn ika to z tego, że przy znanej indukcji B odpow iadające je j natężenie pola H należy do b ie ra ć iteracyjnie. Dogodniej je s t posługiw ać się w takich przypadkach m odelem odw rotnym Preisacha, tj. m odelem , w którym sygnałem w ejściow ym je s t indukcja m agnetyczna B, a sygnałem wyjściow ym natężenie pola m agnetycznego H.
Zasada działania m odelu odw rotnego je s t identyczna ja k m odelu klasycznego. Różnica polega na tym, że do ta b licy z w artościam i funkcji F wprow adza się w ielkości proporcjo
nalne do natężenia pola lub m agnetyzacji, a nie do w artości indukcji m agnetycznej. Przy ich obliczaniu korzysta się z rodziny pętli histerezy H = f(B ) (rys. 4) [1, 7].
Rys. 4. Rodzina pętli histerezy H = f(B) Fig. 4. The hysteresis loops H = f(B)
3. ALG O R Y TM W Y ZN A C ZA N IA PO LA M AG N ETYC ZN EG O
W rozw ażaniach przyjęto, że rozkład pola elektrom agnetycznego w przestrzeni zaw ierającej fer- rom agnetyk m iękki m agnesy trw ałe oraz podobszary o stałej przenikalności m agnetycznej opisują rów nania [5]:
rot(v rotA ) = J + J n (9)
J = y g ra d V B - (10)
w których: v - reluktyw ność m agnetyczna środowiska, A - w ektorow y potencjał m agnetyczny, J - w e kto r gęstości prądu przew odnictw a w podobszarach o konduktywności elektrycznej y, V 9 - skalarny potencjał elektryczny, J m - w ektor gęstości prądów A m pera w obszarze ferrom agnetyka, zależny od w ektora m agnetyzacji Hi:
J m = r o tH j. (11)
W ystępujący w rów naniu (9) w ektorow y potencjał m agnetyczny A je st powiązany z w ektorem indukcji m agnetycznej relacją B = rotA. Z ałożono ponadto, że w łaściw ości m agnetyczne ferrom a
gnetyka m iękkiego w y n ik a ją z zależności:
(12) W ogólnym przypadku w celu wyznaczenia rozkładu i przebiegu pola m agnetycznego w przetw orniku elektrom agnetycznym przy w ym uszeniu napięciowym należy równania (9) I (10)
rozw iązać łącznie z rów naniam i opisującym i rozpływ prądu w obwodach elektrycznych przetw ornika. Przykładowo, dla uzwojeń przetw ornika wykonanych z cienkich przew odów przyjmują one postać:
u = Ri + — V (13)
dt gdzie:
u - w e kto r napięć zasilających, i - w e kto r prądów w uzwojeniach,
R - diagonalna m acierz rezystancji uzwojeń, 41 - w e kto r strum ieni skojarzonych z uzwojeniami.
T a k zdefiniow any m odel nazywam y m odelem polowo-obwodowym. Jeśli m odel połowy obw odow y je s t w ykorzystyw any do rozpatrywania stanów dynam icznych przetworników elektrom echanicznych, to w rozw ażaniach należy dodatkow o uwzględnić równanie opisujące dynam ikę przetwornika [5], Szczegółow y algorytm rozwiązywania równań takiego polowo- obw odow ego m odelu zjaw isk dla m aszyn elektrycznych o strukturze walcowej zam ieszczono w [5], W przedstaw ionych tam rozw ażaniach przyjęto, że pole w części elektrom agnetycznie czynnej m aszyny je s t dw uwym iarowe, a trójw ym iarow ość zjawisk w obszarze połączeń czołowych ujęto w sposób przybliżony przez w prow adzenie do równań m odelu ich rezystancji Ra i indukcyjności U Do w yznaczania rozkładu pola w ykorzystano m etodę elem entów skończonych.
W przedstaw ionych w pracy [5] rozw ażaniach przyjęto, że w łaściw ości magnetyczne ferrom agnetyka m iękkiego w yn ika ją z równania H = vB i opisane są jednow artościow ą krzywą m agnesow ania. Podobnie, korzystając z krzywej odmagnesow ania m odelowano m agnesy trwałe.
Z akładano przy tym, że ich w łaściw ości w yn ika ją z zależności (12).
Dla indukcji m agnetycznej B = 0 wynikająca z pętli histerezy reluktywność magnetyczna v ferrom agnetyka je st nieokreślona. Ujęcie wykorzystujące tak zdefiniowaną reluktywność środowiska je st zatem nieprzydatne w algorytmie symulacji pola elektromagnetycznego z uwzględnieniem zjawiska histerezy. W celu przystosowania przedstawionego algorytmu do symulacji pola elektrom agnetycznego z uwzględnieniem histerezy magnetycznej przyjęto, że natężenie pola H zarówno w ferromagnetyku m iękkim, jak i twardym opisane je s t zależnością (12). Do określania magnetyzacji Hi z rodziny pętli histerezy wykorzystano klasyczny model Preisacha [3].
Zgodnie z algorytm em m etody elem entów skończonych w wyniku dyskretyzacji przestrzeni z rów nań (9), (10) i (13) uzyskuje się [5]
S + K ( l- C ) p - z ę - z Tp - ( R + p L e) I - u
przy czym:
S - m acierz współczynników ; <p = <p(t) - w ektor zastępczych m agnetycznych potencjałów wekto
row ych <p w węzłach siatki dyskretyzującej; 0m - w ektor przepływów w otoczeniu węzłów siatki odw zorow ujący m agnetyzację Hi ferrom agnetyka; K - w ektor konduktywności podobszarów w otoczeniu w ęzłów siatki; z - m acierz transform ująca w ektor prądów Oczkowych w w ektor przepły
w ó w zw iązanych z podobszaram i w otoczeniu węzłów; La - m acierz indukcyjności połączeń czoło
wych; R - m acierz złożona z rezystancji elem entów przewodzących w części elektrom agnetycznie czynnej i rezystancji połączeń czołowych; p = d/dt - operator różniczkowania. W obszarze uzwojeń w ykonanych z cienkich przewodów C = 1, natom iast dla uzwojenia w ykonanego z elementów m a
syw nych przew odzących m acierz C transform uje potencjały cp w węzłach siatki w strumienie skoja
rzone z każdym z elem entów m asywnych.
Dla przetw orników charakteryzujących się sym etrią płaszczyznow ą zastępczy wektorowy po
te n cja ł m agnetyczny <p rów ny iloczynowi składowej A z wektorowego potencjału magnetycznego w kierunku osi z i długości I rdzenia: tp = IAZ . Powyższe rów nania m ożna także w ykorzystać do opisu pola w układach charakteryzujących się sym etrią o s io w ą tj. pola, dla którego potencjał wektorowy
172 Szeląg W., Sujka P.
ma tylko składow ą o b w odow ą A s . Należy w tym celu przyjąć, że rezystancje Re i Indukcyjności U są równe zeru oraz 1 = 2%t, gdzie r - prom ień.
W celu rozw iązania rów nania (14) dyskretyzuje się czas. Przyjęcie schem atu różnicowego C ranka-N icholsona prow adzi do układu równań różnicowych:
9 S k + — K ( l - C ) At
—U T
A t
- 9 z - S R — —L
At
-f
' (i-o)ek m - ' ' (i - a)sk_I K(i
- S u k
r (l-o )u k-'
— Z TA t
- ( l- S ) z 1
(15) .
k-1- ( l - S ) R - f - — L At
W powyższym równaniu indeksem k oznaczono w ielkości dla chwili czasow ej t = tk, a indeksem k-1 w ielkości zw iązane z c h w ilą t = tk -1, w spółczynnik w agow y S = 0,5, At - długość kroku czasowego.
Zgodnie z przyjętym sposobem odw zorowania ferrom agnetyka elem enty m acierzy współczyn
ników po lewej i po praw ej stronie równania (15) nie za le żą od indukcji m agnetycznej. Funkcją indukcji są tylko elem enty w ektora przepływu 6m . Szczegółow y algorytm wyznaczania elem entów tego wektora przedstaw iono w [5], W artykule potrzebną do obliczenia elem entów w ektora em m agnetyzację Hi (B) określa się za pom ocą odwrotnego m odelu Preisacha.
T rudność w odw zorow yw aniu pętli histerezy polega na tym, że w artość wektora m agnetyzacji w przypadku pętli histerezy zależy nie tylko od aktualnej w artości indukcji B w elem encie dyskretyzu- jącym , lecz także od „historii” je g o m agnesowania. Ta zakodow ana je s t w ciągu w ystępujących po sobie m aksym alnych i m inim alnych w artości indukcji. Proces m agnesowania m oże przebiegać w inny sposób w każdym z elem entów dyskretyzujących obejm ujących ferrom agnetyk. N ależy zatem
•dla każdego z nich p a m iętać historię całego cyklu przem agnesowania.
N ieliniow y układ równań (15) rozw iązuje się dla każdej chwili czasow ej m eto d ą Newtona- Raphsona [5], Przy gęstej dyskretyzacji przestrzeni m oże zaw ierać on do kilkudziesięciu tysięcy nieliniow ych rów nań opisujących rozkład potencjału <p oraz do kilkudziesięciu równań dla obw odów elektrycznych przetwornika.
4. B A D AN IA S Y M U LA C YJN E
Na podstaw ie przedstaw ionego algorytm u wyznaczania pola z uw zględnieniem histerezy opra
cowano program kom puterow y do sym ulacji i w izualizacji stanów pracy przetw orników elektrom a
gnetycznych. Do odw zorow ania histerezy w ykorzystano m odel odw rotny Preisacha. Program napi
sano w ję zyku program ow ania Delphi 5 i wdrożono do obliczeń na sprzęcie kom puterow ym klasy PC Pentium.
W artykule przedstaw iono w ybrane w yniki obliczeń przeprowadzonych dla dławika kubkowego.
Założono, że m ateriał ferrytow y, z którego w ykonano rdzeń, m a pętle histerezy ja k na rys. 4 [3].
Sym ulow ano pracę dław ika przy w ym uszeniu napięciowym.
C e ch ą ch a rakterystyczną zaim plem entow anego m odelu Preisacha je s t to, że przy starcie obli
czeń, tj. przy braku „historii” m agnesow ania m ateriału, sygnał w yjściow y z m odelu zgodnie z zależ
nościam i (7) i (8) odpow iada ujem nem u nasyceniu. Dlatego każdorazowo przed przystąpieniem do obliczeń sym ulow ano proces częściow ego rozm agnesow ania ferrom agnetyka. Polegał on na m o- notonicznym zm niejszaniu do ok. 0,02 T am plitudy sinusoidalnie zm ieniającej się indukcji m agne
tycznej B.
R ozpatrywano stan załączenia na zaciski dławika napięcia sinusoidalnie zm iennego o często
tliw o ści 50 Hz i am plitudach 100 V, 250 V i 400 V. Fazę napięcia dobrano w taki sposób, by skła
dow a aperiodyczna prądu w uzw ojeniu była stosunkow o m ała. Uzyskane pętle histerezy, w tym sam ym punkcie rdzenia, dla napięć 250 V i 400 V, przedstawiono na rys. 5. M ała dokładność od
w zorow ania kształtu histerezy w ynika z przyjętej rzadkiej dyskteryzacji trójkąta Preisacha. Na rys. 6
pokazano natom iast, w jednostkach względnych, obliczone przebiegi prądu w uzwojeniach dla napięć zasilających 100 V i 250 V.
b j ; ; ; W lUfiZ---
: ; : r
• i i
M Í
)2\
; H|A/m]
-2A6E4 ¡-1.545E4i-T.03E4 i-5143
1 1 1 514
0.21
1.03E411.545E4I 2.06E4
0.41
■0.62 -0.82
a ) : : : B >T i0 .4 2 .
: ¡ i
Ó.32
; J
i / ; A !
¡ c i ...
HlA/m -355.1 1-2 6 6 .3 i-1 7 7 .6 ¡-0 8 .7 8 8 8 .7 8 i A r?7 .6 ! 2 66 .3 3 5 5 .
I : / L A l l j
I / i /
y r 1 1
L ______ ________
-0 ,3 2 I :
•0 ,4 2 I í
R y s . 5 . P ę t l a h i s t e r e z y d l a : a ) U m = 2 5 0 V , b ) U „ = 4 0 0 V F i g . 5 . T h e h y s t e r e s i s l o o p f o r : a ) U m = 2 5 0 V ; b ) U m = 4 0 0 V
5. W N IO SK I
R ezultaty obliczeń w ska zu ją na przydatność opracowanego algorytm u i programu obliczenio
w ego do w yznaczania rozkładu pola m agnetycznego z uwzględnieniem zjawiska histerezy. D okład
ność przedstaw ionej m etody zależy m.in. od gęstości dyskretyzacji trójkąta Preisacha. Z agęszcze
nie dyskretyzacji powoduje je d n a k w zrost czasu obliczeń. Porów nując czasy obliczeń sym ulacji pola m agnetycznego w ym uszanego napięciowo uzyskanych przy przyjęciu jednowartościowej krzywej m agnesow ania przy posługiw aniu się: a) m odelem , w którym ferrom agnetyk odwzorowany je s t przez w prow adzenie prądów m agnetyzacji m agnetycznej (prądów Am pera) oraz b) m odelem kla
sycznym o zm ieniającej się przenikalności m agnetycznej środowiska, stwierdzono, że dla m etody klasycznej obliczenia przebiegają ok. 9-krotnie krócej. M etody klasycznej nie można niestety w yko
rzystać do uwzględnienia zjaw iska histerezy w ferrom agnetykach. W ynika to z występowania nie
oznaczonej w artości przenikalności m agnetycznej dla natężenia pola H = 0.
LITERATURA
1. G yim dthy Sz., Iv in y i A.: Im plem entation and application o f the 2D vectorial Preisach m odel fo r field calculation with the finite elem ent m ethod, Proc. 6-th International I GTE Sym posium , Graz, Austria, Septem ber 1994, pp. 89-94.
174 Szeląg W., Sujka P.
2. H enrotte F., N icolet A., Deliance F., G enon A., Legros W .: M odeling of ferro m a g n e tic m aterials in 2D finite e lem ent problem s using Preisach m odel, IEEE Trans. On M agn., Vol. 28, No.5, S eptem ber 1991, pp 2614-2616.
3. M ayergoyz I.D.: M athem atical m odels o f hysteresis, Springer Verlag, New Y ork 1991.
4. N inet O., P eccolo M.A., Fraisse H., M asson J.P.: 2D field calculations w ith hysteresis fo r the characterisation o f m agnetic circuits, E lectrim acs’96, Saint-Nazaire, Septem ber 1996.
5. Szeląg W .: Analiza stanów pracy i synteza silników synchronicznych m agnetoelęktrycznych, ujęcie polowe, W ydaw nictw a P olitechniki Poznańskiej, Poznań 1998.
6. Szeląg W ., N owak L, M yszkow ski A.: H am ulec elektrom agnetyczny z cie czą m agnetorolo- giczną, Prace Naukow e Instytutu M aszyn, Napędów i Pom iarów Elektrycznych Politechniki W rocław skiej Nr 48, SM E 2000, s.206-213.
7. Takahashi N., M iyabara S., Fujiwara K.: Problem s in practical finite elem ent analysis using Preisach hysteresis m odel, IEEE Trans. On Magn., Vol. 35, No. 3, M ay 1999, pp1243-1246.
Recenzent: Dr hab. inż. A dam Jagiełło Profesor Politechniki Krakowskiej
W płynęło do R edakcji dnia 2 m arca 2001 r.
Abstract
In the m agnetic field analysis increasing attention is paid to m agnetic properties representation o f fe rrom agnetic m aterial. M agnetic m aterials are usually described in term s o f a single valued B = f(H ) curve in finite e lem ent analysis. However, in som e applications, such as in loss evaluations, the hysteresis b ehaviour o f the m aterial is im portant since hysteresis loss can be a significant co m p onent o f the total loss. In these cases, usually the Preisach m odel is used to describe the hystere
sis phenom enon.
The classical Preisach m odel describes the non-linear relation between the field strength H and m agnetic flux density B or the m agnetisation Hi. In considerations it has been assum ed th a t B and H are parallel. In scalar Preisach m odel, it is assum ed that the m agnetic m aterial consist o f m any e lem entary interacting particles and each o f them can be represented by a rectangular elem entary hysteresis loop, as show n in Fig. 1. The m ain properties and geom etric interpretation (Fig. 2, Fig. 3) o f the Preisach m odel are presented. The Preisach m odel can be num erically im plem ented by using the form ula (8). T his m odel is not suitable fo r B oriented finite elem ent method, because the addi
tional iteration is required to find field H and m agnetisation Hi = B /p 0 - H from the calculated B.
T herefore in the paper, in order to calculate H, the Preisach m odel with inverse distribution function p(a, p) is introduced (B-oriented m odel).
T he algorithm o f the electrom agnetic field analysis in electrom agnetic converters has been pre
sented. T he elaborated m ethod bases on com bined solutions o f the m agnetic field equation (9), the e lectric circuits equation (13). For analysis o f the electrom agnetic field the finite e lem ent m ethod and step by step algorithm have been used. The hysteresis has been taken into account by m eans of the B-oriented Preisach m odel. The field problem has been considered as tw o-dim ensional. In order to describe the field distribution the m odified m agnetic potential <p has been applied. A trian
g ular grid has been constructed and the field and circuit equations have been approxim ated by the system o f ordinary differential equations (14). In order to solve these equations the time discretization has to be carried out. The derivatives have been considered by using the Crank-N icholson formula.
Thus, the differential equation (14) have been substituted with the system of algebraic equation (15). In order to solve nonlinearity the Newton-Raphson iterative method has been used.
T he elaborated algorithm and program has been used fo r sim ulation voltage-exciting electro
m agnetic field in the axisym m etrical reactor.