[WikiDyd] [TitleIndex] [WordIndex

http://prpw.iem.pw.edu.pl/images/small_logo_prpw.jpg

Materiały zostały opracowane w ramach realizacji Programu Rozwojowego Politechniki Warszawskiej.

http://prpw.iem.pw.edu.pl/images/KAPITAL_LUDZKI.gif http://prpw.iem.pw.edu.pl/images/EU+EFS_P-kolor.gif

Laboratorium podstaw elektromagnetyzmu

FAT: Wyznaczenie impedancji charakterystycznej kabla jednożyłowego, ekranowanego

Prześledzimy teraz proces wyznaczenia impedancji charakterystycznej kabla, którego przekrój przedstawiony jest na rysunku:

Przykładowy kabel - przekrój poprzeczny

Nietypowy - prostokątny przekrój żyły wynika z chęci pełniejszego zaprezentowania sposobu opisu geometrii.

Impedancja charakterystyczna kabla określona jest wzorem:

wz1.png (1)

przy czym L i C są parametrami jednostkowymi (określonymi dla 1 metra długości) kabla.

Załóżmy, że chcemy wyznaczyć te parametry jednostkowe i następnie impedancję przy wykorzystaniu pierwszej równości wzoru (1). Musimy w tym celu wykonać analizę zarówno pola elektrycznego, jak i magnetycznego. Do wyznaczenia parametrów jednostkowych (odniesionych do jednego metra długości, przy założeniu, że kabel jest nieskończenie długi) wystarczające jest modelowanie pola w dwóch wymiarach (w przekroju poprzecznym kabla).

Pole elektrostatyczne

Przy modelowaniu pola elektrycznego w przekroju kabla założymy, że zarówno żyła, jak i pancerz są powierzchniami ekwipotencjalnymi. Interesuje nas jedynie pole w dielektryku, a elementy metalowe modelujemy przez zdefiniowanie na ich powierzchni warunków Dirichleta. Ze względu na symetrię wystarczy opisać jedynie 1/4 przekroju. Bazą dla każdego programu w języku FAT powinien być czytelny szkic geometrii z zaznaczonymi nazwami obiektów geometrycznych. Taki szkic dla naszego programu przedstawiony jest na rysunku:

Szkic modelu pola elektrostatycznego

Sam program wygląda następująco:

   1 problemp1
   2 "Analiza pola elektrostatycznego
   3 w~1/4 przekroju poprzecznego przykladowego kabla"
   4 
   5 eqtn esxy
   6 
   7 +mt1 [5*eps0]
   8 
   9 +bp1 0 .5m
  10 +bp2 0 4.5m
  11 +bp3 1m 0
  12 +bp4 1m .5m
  13 +bp5 4.5m 0
  14 
  15 +cr1 0 0 4.5m
  16 
  17 +sg1 bp2 bp5 -1
  18 
  19 +chzyla bp1 bp4 bp3
  20 +chpancerz sg1
  21 
  22 +mediel( mt1 0 .25m 0 ) zyla bp1 bp2 pancerz bp5 bp3
  23 
  24 bnd zyla diri 10
  25 bnd pancerz diri 0
  26 
  27 geng * [2]
  28 
  29 run
  30 
  31 vis
  32 
  33 ? enr() ()
  34 
  35 redir si

Jak już wspominaliśmy w poprzednim przykładzie program w języku FAT składa się z ciągu instrukcji, z których każda jest (w zasadzie) umieszczona w osobnej linii. Instrukcja to ciąg słów kluczowych (np. problem), identyfikatorów (np. mt1 lub chzyla), liczb, zbiorów liczb, stałych tekstowych i znaków specjalnych (np. + lub ?). W dowolnym miejscu programu może wystąpić komentarz.

Pierwsza linia naszego programu zawiera instrukcję problem, która inicjalizuje struktury danych. Każdy program powinien zaczynać się od tej instrukcji (od tej zasady są pewne wyjątki – omówimy je później). Opcjonalnym parametrem, jest nazwa zagadnienia, w postaci stałej tekstowej, czyli dowolnego ciągu znaków ujętych w pojedyncze apostrofy. Interpreter FAT gnoruje wszystko, co znajduje się przed instrukcją problem. Jeżeli pisząc program, pomylimy się w tej linii, to najczęstszym efektem jest przejście interpretera przez cały plik zawierający program, bez wykonywania jakichkolwiek operacji.

zawierający program, bez wykonywania jakichkolwiek operacji. Druga i trzecia linia programu to komentarz. Komentarzem jest dowolny ciąg znaków ujęty w cudzysłowy – może on ciągnąć się przez więcej niż jedną linię. FAT traktuje komentarze podobnie jak znaki odstępu (spacje). Komentarze możemy też umieszczać na końcu linii, po znaku '>'.

Czwarta linia programu jest pusta. FAT ignoruje puste linie – możemy ich używać do poprawienia czytelności programu.

Następna linia jest określeniem rodzaju symulowanego pola. Służy do tego instrukcja eqtn (od ang. equation czyli równanie), której parametr (esxy) oznacza pole elektrostatyczne we współrzędnych kartezjańskich.

W siódmej linii, rozpoczynającej się od znaku plus (+) jest zdefiniowany materiał dielektryka. Jak wynika z założeń (por. rys. A.2), wartość przenikalności dielektrycznej materiału zastosowanego w kablu wynosi 5ε0 . FAT „zna” wartość stałej dielektrycznej próżni (ε0 ). Definicja określa, że przez mt1 będziemy dalej rozumieć materiał, którego przenikalność względna jest pięciokrotnie większa niż ε0 . Oczywiście można zapisać to samo używając wartości liczbowych (5ε0 ≈ 44.2510−12 ), ale warto wiedzieć, że FAT zwalnia użytkownika z obowiązku obliczania wartości, które można wyrazić przez ε0.

W dziewiątej linii programu rozpoczyna się opis geometrii zagadnienia. Jest to zwykle najdłuższa część programu. Każda instrukcja rozpoczynająca się od znaku plus, to wprowadzenie do programu nowego obiektu. Będziemy to dalej nazywać definicją obiektu. Obiekty definiowalne przez użytkownika to: materiały, obiekty geometryczne (punkty bazowe, segmenty, łańcuchy i makroelementy) oraz obiekty pomocnicze (wielobieguny, okręgi, funkcje sklejane).

Z punktu widzenia języka modelowany obszar składa się z co najmniej jednego makroelementu. Makroelement jest dowolną figurą płaską, dla której zdefiniowane są:

Segment jest odcinkiem prostym, łukiem lub krzywą opisaną funkcją sklejaną. Definiuje się go przez określenie nazwy, podanie krańcowych punktów bazowych i rodzaju łączącej je krzywej. Jeżeli odcinek nie jest prostoliniowy, to wymagane jest wcześniejsze zdefiniowanie odpowiedniego okręgu lub funkcji sklejanej.

Segmenty prostoliniowe mogą być również definiowane niejawnie podczas deklarowania makroelementów. Jeżeli w opisie brzegu makroelementu podana zostanie sekwencja punktów bazowych, to segmenty łączące te punkty zostaną zdefiniowane jako odcinki proste.

Punkt bazowy jest charakterystycznym punktem geometrii obszaru. Definiuje się go przez podanie nazwy oraz współrzędnych.

Dodatkowym obiektem, który może znaleźć zastosowanie przy zadawaniu warunków brzegowych lub określenia rozkładu pola jest łańcuch, czyli dowolna ciągła sekwencja segmentów.

Składnia języka wymaga, aby używać tylko obiektów uprzednio zdefiniowanych. Ponieważ do definicji makroelementów wymagane jest zdefiniowanie segmentów, a do definicji segmentów wymagane jest określenie punktów bazowych, a więc opis geometrii rozpoczynamy zwykle od zdefiniowania tych ostatnich. Przed definicją segmentów, które są fragmentami okręgów należy też zdefiniować te okręgi. Definicje segmentów (o ile są one niezbędne) powinny poprzedzać definicje makroelementów i łańcuchów.

Linie 9–13 zawierają definicje punktów bazowych. Parametrami instrukcji definiującej punkt są, obok nazwy, jego współrzędne. Przy podawaniu dowolnych liczb (a w szczególności współrzędnych) możliwe jest używanie tzw. mnożników skalujących. W tym wypadku użyto przyrostka „m” oznaczającego „mili” (tę samą liczbę można zapisać na kilka różnych sposobów, np. 5m = 5e-3 = 0.005). FAT „wie”, że współrzędne węzłów są podawane w metrach, a więc nie ma potrzeby określania jednostek. Na marginesie warto dodać, że w języku FAT używany jest konsekwentnie układ SI, a więc współrzędne podawane są w metrach, ładunki w kulombach itd.

W każdym miejscu programu, w którym może wystąpić liczba, można zamiast niej wstawić proste wyrażenie składające się z liczb i czterech podstawowych operatorów. Należy pamiętać, że wszystkie operatory mają ten sam priorytet, np. wyrażenie 2+1*3 ma wartość 9.

Jedynym segmentem, który trzeba w naszym programie zdefiniować jawnie jest sg1 – fragment okręgu cr1. Zgodnie z tym, co napisaliśmy powyżej, najpierw trzeba zdefiniować okrąg – robimy to podając nazwę, współrzędne środka i promień (linia 15).

Przy definiowaniu segmentu sg1 (linia nr 17) podajemy jego nazwę, punkty krańcowe i rodzaj krzywej. Liczba ujemna, użyta w naszym programie oznacza, że segment jest fragmentem okręgu – wartość bezwzględna tej liczby określa o który okrąg chodzi. W następnym programie pokażemy inne możliwe parametry opisu odcinka. Należy zwrócić uwagę na fakt, że przy opisie odcinka będącego fragmentem okręgu istotna jest kolejność podawania punktów krańcowych. Dwa punkty „wycinają” z okręgu dwa różne łuki. FAT „wybiera” ten właściwy, zakładając, że przy opisie segmentu obiegamy okrąg zgodnie z ruchem wskazówek zegara. W naszym przypadku „od bp2 do bp5” oznacza łuk o kącie 90 stopni (od osi y do osi x). Jeżeli podalibyśmy punkty krańcowe w odwrotnej kolejności, to instrukcja „wybrałaby” z okręgu łuk 270 stopni.

Następny fragment programu jest nieco nietypowy. Zwykle definiujemy najpierw makroelementy, a potem łańcuchy. Chcieliśmy jednak zademonstrować, że do określania nowych obiektów możemy używać obiektów już zdefiniowanych.

Definicja łańcucha polega na określeniu jego nazwy i następnie ciągłej sekwencji odcinków. W przypadku pierwszego łańcucha (chzyla, linia nr 19) sekwencja odcinków zdefiniowana jest przez podanie trzech punktów bazowych. FAT sprawdzi, czy prostoliniowe segmenty między punktami bp1 i bp4 oraz bp4 i bp3 zostały wcześniej zdefiniowane. Jeżeli nie, to zostaną automatyczne utworzone te dwa segmenty o nazwach, w tym wypadku, sg_autg_001 oraz sg_autg_002.

Drugi łańcuch określony jest za pomocą segmentu sg1. Nie mamy tu możliwości zdefiniowania go jako sekwencji punktów bp2 i bp1, gdyż określiłoby to prosty segment (a nie łuk) między tymi punktami.

Dodajmy jeszcze, że każdy łańcuch ma kierunek określony przy jego definicji (np. dla chzyla jest to kierunek od bp1 przez bp4 do bp3). Kierunek ten jest istotny przy postprocesingu – pokażemy to w dalszych rozważaniach.

Definicja makroelementu jest nieco bardziej skomplikowana. Po nazwie następują parametry (cztery lub pięć) ujęte w nawiasy, a następnie określenie brzegu. Parametrami są kolejno:

Raster i waga określają sposób podziału makroelementu na elementy skończone. Waga jest używana do wyboru kolejności podziału makroelementów na trójkąty. Jeżeli mamy do wyboru kilka makroelementów, to jako pierwszy zostanie zdyskretyzowany ten, którego waga jest największa. Jeżeli dwa elementy mają taką samą wagę, to w pierwszej kolejności zostanie wygenerowana sieć w tym, którego raster jest mniejszy. Raster określa jednocześnie średnią długość boku trójkątów, które będą generowane w makroelemencie.

Przy opisie brzegu makroelementu wykorzystano zdefiniowane poprzednio łańcuchy. W naszym przypadku brzeg mediel opisany jest jako:

. . . zyla bp1 bp2 . . .

Uważny czytelnik zwrócił zapewne uwagę, że bp1 jest krańcowym punktem chzyla – nasuwa się pytanie „po co podawać go jeszcze raz?”. Jest to konieczne, gdyż sekwencja:

. . . zyla bp2 . . .

nie określa dokładnie, z którym końcem chzyla należy połączyć bp2. Zawsze, gdy przy opisie łańcucha lub brzegu makroelementu (jest to też łańcuch) występuje przejście od segmentu do sekwencji punktów (lub odwrotnie), należy podać ten koniec segmentu, z którym ma być połączona sekwencja punktów.

Następne linie naszego programu zawierają definicje warunków brzegowych. Każda definicja rozpoczyna się od słowa kluczowego bnd (od ang. boundary). Po nim następuje nazwa obiektu (łańcucha, punktu bazowego lub węzła sieci), w którym ten warunek zadajemy. W naszym programie są to: zyla dla pierwszego, a pancerz dla drugiego warunku. Zwróćmy uwagę, że przy podawaniu nazw uprzednio zdefiniowanych łańcuchów (lub makroelementów) możemy pominąć przyrostek ch (lub me). Po znalezieniu identyfikatora, który nie jest słowem kluczowym i nie rozpoczyna się od żadnego przedrostka, interpreter FAT sprawdza, czy nie jest to nazwa makroelementu, a następnie czy nie jest to nazwa łańcucha.

Po określeniu obiektu, na którym warunek ma być zadany, następuje określenie rodzaju (i wartości) warunku brzegowego. W naszym programie jest to w obu wypadkach warunek pierwszego rodzaju (Dirichleta), oznaczany słowem kluczowym diri. Po nim następuje wartość potencjału, który zostanie zadany na brzegu (10 V na chzyla i 0 V na chpancerz).

Kolejną operacją, która powinna być wykonana przez interpreter, jest generacja sieci elementów skończonych. Nastąpi to po wykonaniu instrukcji geng. Jej parametrami są gwiazdka (*) oznaczająca żądanie „wygładzenia” wygenerowanej siatki i zbiór ([1]) parametrów określających sposób wygładzenia.

Generacja sieci ES jest przeprowadzana przez automatyczny generator, którego algorytm „stara się” podzielić makroelement na trójkąty możliwie zbliżone do równobocznych i o podobnej wielkości (zależnej od zadanego rastra). Nie zawsze jest to możliwe, bądź ze względu na różnice w podziale brzegu makroelementu (część brzegu może przylegać do innego makroelementu o dużym rastrze, a inna część do jeszcze innego makroelementu o małym rastrze), bądź ze względu na złożony kształt makroelementu. Wygładzanie sieci polega na poprawieniu położenia już wygenerowanych węzłów, tak aby zbliżyć kształt trójkątów do równobocznych. FAT stosuje przede wszystkim bardzo prosty (ale skuteczny) algorytm polegający na przesunięciu każdego węzła do środka wielokąta złożonego z trójkątów otaczających węzeł.

Drugi algorytm wygładzania, który jest znacznie bardziej pracochłonny i nie przynosi zwykle tak dobrych rezultatów jak pierwszy, polega na regeneracji wszystkich możliwych par trójkątów występujących w siatce, tak aby uzyskać „lepsze” (bardziej równoboczne) trójkąty.

Jeżeli generator otrzymuje polecenie wygładzenia sieci, to standardowo używa jedynie pierwszego algorytmu. Użycie tylko drugiego lub obu możemy wymusić przez podanie parametru o wartości 2 lub 3. Podsumowując:

Efektem wygenerowania sieci jest pojawienie się nowych obiektów, których nazw możemy odtąd używać – węzłów i elementów. Nazwa węzła (ang. node) składa się z przyrostka nd i numeru – np. nd22. Nazwa elementu składa się z przyrostka el i numeru. Węzłów można używać do zadawania warunków brzegowych. Elementy mogą być przydatne tylko w postprocesingu. Posługując się węzłami i elementami, należy pamietać o tym, że ich parametry mogą się zmienić – kolejne instrukcje generacji mogą zmienić numer trójkąta lub współrzędne węzła.

Generacja sieci kończy przygotowanie danych niezbędnych do utworzenia modelu zagadnienia. Kolejna instrukcja programu – run – powoduje wykonanie analizy zagadnienia brzegowego. Po obliczeniu rozkładu potencjału przełączamy interpreter w tryb graficzny (vis) i polecamy mu narysowanie rozkładu gęstości energii i konturów dielektryka (? enr() ()). Wynik działania tej instrukcji przedstawiony jest na rysunku:

Rozkład energii pola elektrycznego w dielektryku

Ostatnia instrukcja programu (redir si) przełącza interpreter w interakcyjny tryb pracy.

Jak wynika z rysunku efektem wykonania polecenia ? enr() jest nie tylko mapa rozkładu gęstości energii w analizowanym obszarze, ale też i całkowita energia zgromadzona w obszarze (TOTAL ENERGY).

Gęstość energii obliczana jest zgodnie ze wzorem

wz2.png (2)

natomiast całkowita energia:

wz3.png (3)

Na rysunku występuje jeszcze pozycja SOURCE ENERGY – jest to wartość obliczana dla pola elektrostatycznego zgodnie ze wzorem:

wz4.png (4)

gdzie ϕ to potencjał, a – objętościowa gęstość ładunku. W tym przypadku wartość ta wynosi zero, gdyż w całym obszarze = 0.

Jak pamiętamy, naszym celem jest wyznaczenie pojemności jednostkowej kabla. Możemy wykorzystać w tym celu następujący wzór:

wz5.png (5)

w którym W to energia zgromadzona w kondensatorze, C – jego pojemność, a U – napięcie między okładkami. Nasz „kondensator” to odcinek kabla o długości jednego metra. Wartość 1 W (analizowaliśmy 1/4 przekroju) została obliczona przez nasz program jako TOTAL ENERGY, a napięcie między żyłą a pancerzem wynosi 10 V (różnica warunków brzegowych). Ostatecznie:

wz6.png (6)

Opisany sposób wyznaczenia pojemności C nie jest jedynym możliwym. Inna droga to obliczenie C na podstawie zależności definicyjnej:

wz7.png (7)

przy czym Q to ładunek zgromadzony na jednej z okładek kondensatora, a U – tak jak poprzednio – napięcie między okładkami. Zgodnie z prawem Gaussa, ładunek zgromadzony na okładce może być wyznaczony jako strumień wektora indukcji elektrycznej przez pewną zamkniętą powierzchnię otaczającą okładkę:

wz8.png (8)

Język FAT pozwala wyznaczyć strumień wektora D przez dowolną powierzchnię (w przypadku pola dwuwymiarowego powierzchnia redukuje się do linii). W naszym przypadku moglibyśmy wyznaczyć strumień przez powierzchnię żyły lub powierzchnię pancerza. Wystarczy w tym celu podać interpreterowi następujące instrukcje:

? flx(zyla)
? flx(pancerz)

– możemy dokonać tego w trybie pracy interakcyjnej. W pierwszym przypadku otrzymamy flx(zyla) = 0.34875e-9 C, zaś w drugim flx(pancerz) = 0.44268e-9 C.

Strumienie te odpowiadają 1/4 (symetria) ładunku zgromadzonego na żyle i pancerzu. Teoretycznie obie wartości powinny być równe – nie są ze względu na błędy wynikające z dyskretyzacji MES. Podstawienie ich do wzoru (7) daje następujące wyniki:

wz9.png (9)

Otrzymaliśmy trzy wartości C, C' , C" , z których jedna (C' ) różni się znacznie od pozostałych dwóch. Powstaje pytanie – która z nich jest najdokładniejsza? Aby rozwiązać ten problem, musimy powiedzieć nieco o błędach wynikających z zastosowania MES do aproksymacji zagadnienia brzegowego. MES poszukuje rozwiązania w postaci pewnych funkcji wielomianowych. W naszym przypadku, dla elementów trójkątnych, są to wielomiany pierwszego stopnia. Nie spełniają one dokładnie ani równania opisującego pole (równania Laplace’a) ani warunków brzegowych. Co więcej, pochodne tych wielomianów nie są ciągłe na granicy elementów, a to w praktyce oznacza, że nie są spełnione warunki ciągłości dla wektorów pola na granicach elementów. W efekcie nasze rozwiązanie różni się od modelowanego potencjału. Możemy (i musimy) zaakceptować ten fakt pod dwoma warunkami:

  1. jesteśmy w stanie wyznaczyć (przynajmniej w przybliżeniu) dokładność aproksymacji,
  2. w miarę zagęszczania sieci nasza aproksymacja powinna zbliżać się do rozwiązania dokładnego (zbieżność metody).

Drugi warunek jest wypełniany przez MES – jest to zresztą niezbędne dla jej praktycznej przydatności. Jeżeli chodzi o ocenę błędu aproksymacji, to można go wyznaczać różnymi metodami. Na poziomie składni FAT zakłada, że jakiś algorytm oceny rozwiązania jest wbudowany w interpreter. Implementacja FATpcv stosuje metodę polegającą na oszacowaniu nieciągłości wektorów pola na granicy elementów. Metoda ta pozwala oszacować względny i bezwzględny błąd dla globalnej energii zgromadzonej w obszarze. W naszym programie błąd ten wynosi około 1.5% i z taką właśnie dokładnością wyznaczona jest wartość C, która jest obliczona na podstawie globalnej energii. Lokalne błędy mogą być nawet o rząd wielkości większe od błędu globalnego. W naszym modelu tak właśnie jest w otoczeniu żyły, gdzie ze względu na geometrię i rodzaj zmienności poszukiwanej funkcji należałoby wygenerować gęstszą sieć.

Pole magnetostatyczne

Zajmiemy się teraz wyznaczeniem indukcyjności jednostkowej kabla. Obliczymy ją analogicznie jak pojemność, wykorzystując wzór:

wz11.png (10)

W przypadku pola magnetycznego musimy przeanalizować nie tylko pole w izolacji kabla, ale też uwzględnić przekrój żyły i pancerza, w których płynie prąd będący źródłem pola magnetycznego. Podobnie jak dla pola elektrostatycznego wystarczy zamodelować 1/4 przekroju kabla. Szkic będący bazą dla naszego drugiego programu jest przedstawiony na rysunku:

Szkic geometrii do modelu pola magnetycznego

Program niewiele różni się od modelu pola elektrycznego

   1 problemp2
   2 "Analiza pola magnetostatycznego
   3 w~1/4 przekroju poprzecznego przykladowego kabla"
   4 
   5 eqtn msxy
   6 +mt1 [mi0inv]
   7 
   8 +bp0 0 0
   9 +bp1 0 .5m
  10 +bp2 0 4.5m
  11 +bpa 0 5m
  12 +bp3 1m 0
  13 +bp4 1m .5m
  14 +bp5 4.5m 0
  15 +bpb 5m 0
  16 
  17 +cr1 0 0 4.5m
  18 +cr2 0 0 5m
  19 
  20 +sg1 bp2 bp5 -1
  21 +sg2 bpa bpb -2
  22 +mez( mt1 7.4612826e6 .25m 0 ) bp0 bp1 bp4 bp3 bp0
  23 +mediel( mt1 0 .25m 0 ) bp1 bp2 sg1 bp5 bp3 bp4 bp1
  24 +mep( mt1 -1e6 .25m 0 ) bp2 bpa sg2 bpb bp5 sg1
  25 
  26 +chzew sg2
  27 
  28 bnd zew diri 0
  29 
  30 geng * [3]
  31 run
  32 
  33 
  34 vis
  35 ? enr() ()
  36 
  37 redir si

Porównując programy „p1” i „p2” zauważymy następujące zmiany:

Parametr instrukcji eqtn przyjmuje tym razem wartość msxy, co oznacza pole magnetostatyczne w kartezjańskim układzie współrzędnych.

Zmienił się opis materiału mt1, gdyż musi on teraz odzwierciedlać magnetyczne właściwości środowiska. Upraszczając nasz model, pomijamy nieistotne dla dokładności odwzorowania różnice między przenikalnością magnetyczną miedzi i dielektryka, zakładając, że jedna i druga wynosi µ0 . Współczynnik materiałowy mt1 określony jest jako mi0inv. Jest to predefiniowana stała, równa 1/µ0. Podobnie jak w polu elektrostatycznym, możemy używać jej wielokrotności do określania współczynników innych materiałów. Należy pamiętać, że zawsze podajemy odwrotność przenikalności magnetycznej.

Jeżeli chodzi o opis geometrii, to różnice między programami odpowiadają różnicom między rysunkami przedstawiającymi szkice modeli. Musieliśmy dodać punkty bazowe bp0, bpa i bpb, segment sg2 (i okrąg cr2 oraz dwa makroelementy (mez i mep).

Nowe makroelementy mają różną od zera wartość funkcji źródła – odpowiada ona gęstości prądu wyrażonej w amperach na metr kwadratowy. Dla naszych rozważań nie jest istotna wartość prądu płynącego przez przewód. Ponieważ układ jest liniowy (nie występują w nim nieliniowe materiały), więc wartość indukcyjności nie zależy od wartości prądu. Istotne jest tylko, aby całkowity prąd płynący przez żyłę był równy całkowitemu prądowi płynącemu przez pancerz. Zakładamy, że prąd płynie żyłą w kierunku „dodatnim” (zgodnie z orientacją układu współrzędnych), a „wraca” pancerzem. Ponieważ stosunek pola powierzchni przekroju żyły do pola powierzchni przekroju pancerza wynosi 1/7.4612826, a więc przyjęto gęstość prądu 7.4612826 · 106 A/m2 w żyle i 106 A/m2 w pancerzu. Całkowity prąd płynący w przewodzie możemy obliczyć mnożąc pole powierzchni żyły przez zadaną gęstość prądu (pamiętamy ciągle o tym, że modelujemy 1/4 przekroju kabla): I = 4·1·10−3 ·0.5·10−3 ·7.4612826·106 = 14.922565 A.

Analogiczną wartość (ale o przeciwnym znaku) otrzymamy mnożąc pole przekroju pancerza przez gęstość prądu w pancerzu.

Warunki brzegowe dla pola magnetycznego sprowadzają się do zadania zerowej wartości potencjału na zewnętrznym brzegu kabla. FAT używa do modelowania pola magnetycznego wektorowego potencjału A, zdefiniowanego następującą zależnością:

wz12.png (11)

Przyjęcie stałej wartości potencjału A na jakiejś części brzegu oznacza, że składowa normalna indukcji magnetycznej jest na tym brzegu równa zero. Zgodnie z prawem Ampere’a, pole magnetyczne na zewnątrz przewodu powinno być równe zero. Warunek ten określa jednocześnie potencjał odniesienia dla całego obszaru.

Efekt uruchomienia naszego drugiego programu jest przedstawiony na rysunku:

rozkład gęstości energii pola magnetycznego

wz13.png (12)

Nie możemy skorzystać z zależności definiującej indukcyjność jako:

wz14.png (13)

gdyż trudno jest obliczyć wartość strumienia Ψ sprzężonego z prądem I. FAT pozwala łatwo obliczyć strumień indukcji magnetycznej przez dowolny łańcuch lub odcinek, ale da to jedynie górne i dolne oszacowanie Ψ. Dolne oszacowanie to strumień przez odcinek między punktami bp3 i bp5 (patrz rysunek). Możemy wyznaczyć je za pomocą następującej sekwencji instrukcji:

+chzz bp3 bp5
? flx(zz)

której wykonanie da nam wartość 0.48577 · 10−5 Wb.

Oszacowanie górne uzyskamy w następujący sposób:

+chzzz bp0 bp3 bp5 bpb
? flx(zz)

uzyskując 0.65581 · 10−5 Wb.

uzyskując 0.65581 · 10Wb.

Po podstawieniu powyższych wartości do wzoru (13) otrzymamy:

wz15.png

Jak widać, poprzednio obliczona wartość mieści się w powyższym zakresie.

Po podstawieniu parametrów jednostkowych kabla do wzoru (1) otrzymamy wartość jego impedancji charakterystycznej:

wz16.png

Dokładność obliczeń

Jak już wspomnieliśmy, rozwiązania uzyskane w poprzednich podrozdziałach obarczone są pewnym błędem wynikającym z samej istoty MES. Aby uzyskać dokładniejszy model, powinniśmy wygenerować w badanym obszarze „gęstszą” sieć elementów. Wraz ze zwiększeniem liczby elementów (i węzłów) rośnie jednak złożoność modelu – zwiększa się czas obliczeń. Pogodzenie wymagania dużej dokładności i jak największej prostoty modelu jest w ogólnym przypadku niemożliwe. Docelowy model i docelowe rozwiązanie uzyskuje się zwykle w wyniku powtarzanych iteracyjnie obliczeń, w których stopniowo poprawia się model na podstawie poprzednich wyników.

FAT umożliwia pewną automatyzację tego procesu. Aby zilustrować adaptacyjną MES, dodajmy do programu „p1”, po instrukcji run, następującą sekwencję poleceń:

anl ni 10 err .0001
geng +
run

Oznaczają one kolejno:

Po wykonaniu analizy adaptacyjnej uzyskamy rozkład gęstości energii przedstawiony na rysunku:

Rozkład gęstości energii pola elektrycznego po adaptacyjnym zagęszczeniu siatki

Całkowita energia wynosi teraz 0.21262 · 10−8 J,

W analizowanym przypadku FATpcv wykonuje nie 10, a 6 iteracji. Po szóstym zagęszczeniu sieci liczba elementów wynosi 1910, a liczba węzłów 1000 (w pierwotnej sieci odpowiednio 630 i 349), co stanowi górne ograniczenie na wielkość sieci w tej implementacji. Procentowy błąd względny dla energii zgromadzonej w obszarze maleje z 1.6% do 0.03%.


Materiały do Laboratorium zostały opracowane w ramach realizacji Programu Rozwojowego Politechniki Warszawskiej.

http://prpw.iem.pw.edu.pl/images/KAPITAL_LUDZKI.gifhttp://prpw.iem.pw.edu.pl/images/EU+EFS_P-kolor.gif

http://www.pr.pw.edu.pl/ jest projektem współfinansowanym przez Unię Europejską w ramach Europejskiego Funduszu społecznego (działanie 4.1.1 Programu Operacyjnego Kapitał Ludzki) i ma na celu poprawę jakości kształcenia oraz dostosowanie oferty dydaktycznej Politechniki Warszawskiej do potrzeb rynku pracy. Będzie on realizowany przez Uczelnię w latach 2008-2015.


2015-09-23 06:44