33, s.181-186, Gliwice 2007
GLOBALNE OBLICZANIE CAŁEK PO OBSZARZE W PURC DLA DWUWYMIAROWYCH ZAGADNIEŃ BRZEGOWYCH MODELOWANYCH RÓWNANIEM NAVIERA-LAMEGO I POISSONA
EUGENIUSZ ZIENIUK, KRZYSZTOF SZERSZEŃ, AGNIESZKA BOŁTUĆ
Instytut Informatyki, Uniwersytet w Białymstoku Sosnowa 64, 15-887 Białystok
e-mail: {ezieniuk, aboltuc, kszerszen}@ii.uwb.edu.pl
Streszczenie. W pracy przedstawiono globalny sposób numerycznego obliczania całek powierzchniowych w dwuwymiarowych zagadnieniach brzegowych.
Prezentowana technika opiera się na matematycznym zdefiniowaniu obszarów za pomocą parametrycznych płatów powierzchniowych oraz wykorzystaniu kwadratur całkowania numerycznego wyższych rzędów. Praktyczną realizację proponowanej procedury przedstawiono dla zagadnień brzegowych definiowanych równaniem Poissona.
1. WSTĘP
Do modelowania szerokiej klasy zagadnień z mechaniki ciała stałego jest stosowane między innymi tradycyjne równanie przemieszczeniowe Naviera-Lamego. Równanie to służy do modelowania liniowych zagadnień mechaniki i charakteryzuje się obecnością w nim sił masowych. Jest ono najczęściej rozwiązywane metodami numerycznymi, takimi jak metoda elementów skończonych (MES) oraz metoda elementów brzegowych (MEB). Siły masowe dosyć często w praktycznych obliczeniach są jednak pomijane, szczególnie jest to zauważalne przy zastosowaniu MEB [1]. Jest to podyktowane tym, że w MEB nie stosuje się dyskretyzacji obszaru, a jedynie dyskretyzację jego brzegu.
Uwzględnienie jednak sił masowych w MEB pociąga za sobą konieczność obliczania całek po obszarze. Praktycznie ich obliczenie tradycyjnymi metodami numerycznymi sprowadza się do podzielenia tego obszaru na tzw. komórki. Podzielenie obszaru na komórki z technicznego punktu widzenia jest bardzo zbliżone do podzielenia go na tzw. elementy skończone stosowane w MES. Różnica sprowadza się tylko do innego ich przeznaczenia. W MES elementy służą do ułatwienia aproksymacji rozwiązań w węzłach należących do elementów skończonych, modelujących w sposób dyskretny rozpatrywany obszar. Komórki natomiast w MEB służą tylko do ułatwienia numerycznego obliczania całek po obszarze. Całki te są obliczane na podstawie zsumowania całek z poszczególnych komórek, na które został podzielony rozpatrywany obszar.
We własnych pracach do rozwiązywania zagadnień brzegowych otrzymano parametryczny układ równań całkowych (PURC) [5]. Główną cechą PURC jest fakt, iż geometria brzegu jest uwzględniona w jego formalizmie matematycznym. Dlatego też w celu praktycznego
zdefiniowania obszarów zadawane są tylko: punkty narożne w przypadku segmentów liniowych oraz punkty brzegowe do wykreowania segmentów krzywoliniowych. Liczba tych punktów jest znacząco mniejsza niż liczba węzłów w przypadku MEB, ponadto zamodelowany w ten sposób brzeg jest brzegiem ciągłym. Zastosowana metoda do rozwiązania PURC pozwala na otrzymanie rozwiązań w dowolnym punkcie na brzegu z wysoką dokładnością. Dotychczas otrzymane PURC dotyczyły równania Laplace’a i Helmholtza, czyli równań, w których nie zachodziła potrzeba całkowania po obszarze.
Ostatnio rozwijany PURC dotyczył równania przemieszczeniowego Naviera-Lamego [5], ale w celu pozbycia się całkowania po obszarze pominięte zostały siły masowe.
Celem niniejszej pracy jest zaproponowanie i przetestowanie techniki obliczania całek powierzchniowych (występujących w PURC), polegającej na obliczaniu tych całek w sposób globalny, czyli bez dzielenia obszaru na komórki. Technika ta jest testowana głównie na zagadnieniach modelowanych równaniem Poissona, mających rozwiązania dokładne.
Zaproponowana technika, w odróżnieniu od techniki stosowanej w tradycyjnej MEB, charakteryzuje się tym, że nie wymaga dzielenia obszaru na komórki. W zaproponowanym sposobie obszar jest traktowany globalnie jako makroelement (czyli jedna komórka). Tylko w przypadkach bardzo skomplikowanych obszarów dopuszczalne jest jego podzielenie na niewielką liczbę podobszarów. Sposób ten dla rozpatrywanych przykładów okazał się bardzo efektywny, ponadto zamieszczone przykłady numeryczne potwierdzają wysoką dokładność zaproponowanej metody w porównaniu z rozwiązaniami analitycznymi.
2. GLOBALNE OBLICZANIE CAŁEK PO OBSZARZE W PURC
Dwuwymiarowe zagadnienie brzegowe modelowane równaniem Naviera-Lamego, w obszarach wielokątnych może być rozwiązywane za pomocą PURC, przedstawionym w [5], natomiast w obszarach krzywoliniowych za pomocą równania przedstawionego w [6]. W obu tych przypadkach w celu uproszczeń w trakcie otrzymywania PURC zostały pominięte siły masowe. Stosując sposób przedstawiony w tych pracach do modyfikacji BRC z uwzględnieniem sił masowych, otrzymano PURC w następującej postaci
}
{
( , ) ( ) ( , ) ( ) ( )( ) ( ) ( )
, ,) ( 5 .
0 * 1
1
1 1
1
1
∑ ∫ ∫
= Ω
∗
∗ − + Ω
=
−
y y y U
J s ds s b d
s s s s
s s
s p
n
r s
s
r r pr
r pr
p
r
r
u P
p U
u (1)
Funkcje podcałkowe Upr oraz Ppr∗ (w pierwszej całce) w jawnej postaci zostały przedstawione w [5,6], natomiast funkcja U dla całki po obszarze *p Ω ma następującą postać
,
) ln(
) 4 3 ( )
ln(
) 4 3 ( ) 1 ( 8 ) 1 , (
2 1
2 1
1
*
−
−
−
−
−
−
− −
=
2 2 2 2
2 2
2 1
η η η
η η η
U ηη ν η
η η ν η
µ ν y π
p s p=1,2,...n (2)
gdzie
5 . 0 2 2 2
1 ]
[η +η
=
η , η1 = y1−Γp(1)(s1) oraz η2 = y2−Γp(2)(s1).
W formule (1) obok całek przedziałowych zdefiniowanych na linii prostej w parametrycznym układzie odniesienia dodatkowo pojawia się druga całka po obszarze Ω. Całka ta też pojawiała się i w tradycyjnej MEB. Była ona obliczana na dwa różne sposoby, w niektórych szczególnych przypadkach w zależności od postaci funkcji podcałkowej można było ją sprowadzić do całki po brzegu [1].
W większości przypadków jednak bardziej ogólnym sposobem jest bezpośrednie całkowanie po obszarze Ω. W celu przeprowadzenia tego całkowania obszar ten był dzielony na komórki. Taki podział wizualnie jest bardzo podobny do podziału obszaru na elementy skończone stosowane w tradycyjnej MES, niezależnie jednak od tego podobieństwa merytorycznie podział ten ma zupełnie inne przeznaczenie. Obliczenie całki po obszarze praktycznie sprowadza się do obliczenia całek na poszczególnych komórkach, a następnie na zsumowaniu otrzymanych wartości. Całki na poszczególnych komórkach obliczane są na podstawie numerycznych kwadratur niższego rzędu, z niewielką liczbą współczynników.
W tym celu wykorzystywany jest następujący wzór [1]
∫ ∑ ∑
= =Ω
≈ Ω
= E
e e
N
n
n nF y w y
d y F B
1 1
) ( )
( )
( , (3)
przy czym w jest wagą, natomiast n F(yn) całkowaną funkcją obliczaną dla n -tego współczynnika zastosowanej kwadratury, E - liczbą komórek, na które podzielono obszar.
Celem niniejszej pracy jest przedstawienie nowej techniki obliczania niektórych całek po obszarze Ω, opartej na globalnym traktowaniu obszaru, czyli bez konieczności dzielenia jego na komórki. Realizacja tego celu stała się możliwa dzięki:
- zastosowaniu kwadratur wysokiego rzędu do całkowania numerycznego [4],
- wykorzystaniu płatów powierzchniowych (stosowanych w grafice komputerowej) do modelowania obszaru Ω [2].
W takim przypadku formuła (3) całkowania numerycznego może być zredukowania do następującego wyrażenia
∫ ∑
=Ω
≈ Ω
= M
m
m
mF y
w y
d y F B
1
) ( )
( )
( , (4)
przy czym w jest wartością wagową, zaś m F(ym) wartością całkowanej funkcji określoną dla m -tego współczynnika kwadratury całkowania numerycznego wyższego rzędu M .
Kolejnym problemem przy obliczaniu całki po obszarze pozostaje praktyczne uwzględnianie tego obszaru w całce. Na podstawie drugiej całki we wzorze (1) możemy stwierdzić, że obszar Ω, po którym obliczamy całkę w sposób bezpośredni, nie jest uwzględniony w formalizmie matematycznym tej całki. Na podstawie tej całki możemy tylko stwierdzić, że y∈Ω, ale rozpatrywany obszar musi być dodatkowo zdefiniowany, aby można było po nim obliczyć całkę. Dlatego też modelowanie obszarów za pomocą komórek należy traktować jako pośredni sposób uwzględnia obszarów w całkach po obszarze.
Do bezpośredniego uwzględniania obszaru Ω w całkach zaproponowano trójwymiarowe płaty powierzchniowe (Coonsa, Beziera) stosowane w grafice komputerowej [2]. Płaty te są bardzo efektywne, ponieważ za pomocą niewielkiej ilości punktów kontrolnych możemy kreować ich kształt w przestrzeni trójwymiarowej. W związku z tym, że w całkach obszary są płaskie, płaty powierzchniowe (Coonsa, Beziera) zostały sprowadzone do płaszczyzny przez wyzerowanie jednej współrzędnej. Dlatego też w płatach powierzchniowych Coonsa lub Beziera stosowanych do modelowania obszarów w całkach po obszarze należy przyjąć, że
3 =0
y , natomiast za y oraz 1 y należy podstawić wzory wynikające z płatów 2 powierzchniowych Coonsa lub Beziera przedstawianych w postaci parametrycznej.
Przykładowo dla płatów Coonsa w (2) należy podstawiać następujące wyrażenia [2]
( )( ) ( ) ( )
[ ]
0,1, 1,2w v,
), ( 1
) ( )
( 1
) ( 1
1 1 2 3 4
=
∈
− + +
− +
−
−
=
j
y P w v y vwP y
wP v y
P w v
yj j j j j
(5) gdzie Pi(yj),(i =1,2,3,4) są odpowiednio kolejnymi współrzędnymi punktów kontrolnych.
2.1. Modelowanie i modyfikowanie obszarów za pomocą płatów
W proponowanym sposobie całkowania obszar jest globalnie zdefiniowany w całce w postaci pojedynczej makrokomórki opisanej sparametryzowanym płatem powierzchni. Na rys. 1 przedstawiono trzy obszary modelowane płatami powierzchni Beziera i Coonsa.
Charakterystyczną cechą rozpatrywanych płatów jest możliwość efektywnego modelowania i modyfikowania kształtu obszarów w sposób globalny, za pomocą niewielkiego zbioru punktów kontrolnych. Kolejną zaletą wprowadzonych płatów jest też łatwość analitycznego obliczenia jakobianu. Płaty powierzchniowe dają możliwość zamodelowania obszaru w sposób globalny, czyli bez dzielenia na komórki. Całkowanie jednak po tak dużych obszarach wymaga zastosowania kwadratur wyższego rzędu. W dotychczasowej literaturze, między innymi w [1,3], można znaleźć głównie współczynniki w dla kwadratur niższego rzędu, jakie były i stosowane do całkowania po komórkach, na które został podzielony obszar.
W pracy z 2003 roku [4] przedstawiono zestawione tabelarycznie współczynniki w dla i kwadratur wyższego rzędu. Podano wartości 85, 126 oraz 175 współczynników dla kwadratur zdefiniowanych na powierzchni trójkątnej. Dodatkowo na podstawie zawartych w tej pracy wzorów możliwe jest bezpośrednie wygenerowanie współczynników dla kwadratury dowolnego rzędu M. Umożliwia to bezpośrednią ich generację w programie komputerowym oraz elastyczny dobór liczby niezbędnych współczynników.
3. PRZYKŁADY TESTOWE
Zaproponowana metoda była testowana na równaniu Poissona zdefiniowanym na różnych obszarach przedstawionych na rys.1, modelowanych za pomocą różnych płatów powierzchniowych. Zbadano też wpływ globalnego modelowania obszarów i obliczania całek po obszarze na wyniki rozwiązań uzyskiwane za pomocą PURC w porównaniu z rozwiązaniami dokładnymi.
a) b) c)
Rys. 1. Modelowanie obszaru: a) trójkątnym płatem Beziera, b) prostokątnym i trójkątnym płatem Coonsa, c) zmodyfikowanym prostokątnym płatem Beziera
Wyniki obliczeń PURC uzyskane w wybranych punktach na brzegu dla każdej z trzech geometrii z rys.1 zestawiono w przedstawionych poniżej tabelach. W każdym z przykładów zadawano warunki brzegowe typu Dirichleta w postaci jawnej. Dla obszaru trójkątnego (rys.1a) na brzegu zadano warunek w postaci elementarnej funkcji
( )
222 1 2
1,x x x
x
u = + ,
natomiast w postaci bardziej złożonej
( )
22 1 2 2 1 3 2 3 1 2
1,x x x 3x x 3x x
x
u = + + + dla obszarów
z rys.1b,c. Rozwiązania otrzymane na brzegu za pomocą PURC porównywano z warunkami Neumanna (kolumna 3), analitycznie otrzymanymi na podstawie zadanych warunków
Dirichleta. W kolumnach czwartej i piątej przedstawiono błąd rozwiązań PURC w stosunku do wartości dokładnych z kolumny 3 w zależności od podanej w nawiasie liczby rozwiązywanych równań algebraicznych.
Tabela 1. Wyniki na brzegu dla obszaru trójkątnego z rys.1a PURC
Nr. Punkt
obliczeniowy(x1, x2)
Rozwiązanie
dokładne Błąd rozwiązań [%]
(12 równań)
Błąd rozwiązań [%]
(24 równania)
1 (2.55, 0.45) 4.24264 0.00736 0.05729
2 (2.1, 0.9) 4.24264 0.00258 0.00876
3 (1.65, 1.35) 4.24264 0.00710 0.02227
4 (1.2, 1.8) 4.24264 0.00620 0.00884
5 (0.75, 2.25) 4.24264 0.00013 0.00429
6 (0.3, 2.7) 4.24264 0.01189 0.03144
Tabela 2. Wyniki na brzegu dla obszaru modelowanych płatami Coonsa z rys.1b PURC
Nr. Punkt
obliczeniowy(x1, x2)
Rozwiązanie
dokładne Błąd rozwiązań [%]
(16 równań)
Błąd rozwiązań [%]
(24 równania)
1 (2.55, 3) 92.40750 0.00002 0.00031
2 (2.1, 3) 78.03000 0.00005 0.00029
3 (1.65, 3) 64.86750 0.00009 0.00022
4 (1.2, 3) 52.92000 0.00008 0.00007
5 (0.75, 3) 42.18750 0.00002 0.00062
6 (0.3, 3) 32.67000 0.00006 0.00015
Tabela 3. Wyniki na brzegu dla obszaru modelowanych płatem Beziera z rys.1c PURC
Nr. Punkt
obliczeniowy(x1, x2)
Rozwiązanie
dokładne Błąd rozwiązań [%]
(16 równań)
Błąd rozwiązań [%]
(32 równania) 1 (0.62885, 4.88658) -119.29000 0.00265 0.00106
2 (2.34363, 2.34363) -93.19290 0.00020 0.00157
3 (4.88658, 0.62885) -119.29000 0.00265 0.00106
4 (8.0, 1.24672) 256.96100 0.00177 0.00020
5 (8.0, 2.49625) 330.56700 0.00015 0.00045
6 (8.0, 3.74766) 413.60300 0.00101 0.00021
Analiza uzyskanych rozwiązań na brzegu wskazuje na ich dużą zbieżność w stosunku do wartości dokładnych. Zostało to osiągnięte przy niewielkiej liczbie rozwiązywanych równań algebraicznych (od 12 do 32 równań) oraz zaproponowanym globalnym sposobie całkowania całek po obszarze.
4. WNIOSKI
Przedstawione w tabelach wyniki potwierdzają słuszność globalnego obliczania całek po obszarze w PURC modelowanych za pomocą płatów powierzchniowych. Takie ich obliczanie nie wymaga jego dzielenia na komórki. Wymaga jednak stosowania do numerycznego całkowania kwadratur wyższego rzędu. Ważną zaletą proponowanej techniki jest też łatwość
modyfikowania obszarów za pomocą punktów kontrolnych. Okazało się na podstawie przeprowadzonych testów, że nieznaczna modyfikacja brzegu globalnego obszaru trójkątnego (rys.1a) nie miała wpływu na pogorszenie wyników niezależnie od tego, że współczynniki do kwadratury były stosowane dla obszaru trójkątnego. W przypadkach bardziej skomplikowanych obszarów dopuszczalne jest podzielenie obszaru na niewielką liczbę podobszarów (rys. 1b). Bardzo interesujący jest przykład z obszarem pokazany na rys.1c. Do jego globalnego zamodelowania wstępnie był wykorzystany płat kwadratowy Beziera oraz współczynniki kwadratury również dla obszaru kwadratowego. Następnie ten płat został za pomocą punktów kontrolnych dość znacznie zmodyfikowany do kształtu pokazanego na rys.1c. Modyfikacja obszaru, jak się okazało na podstawie uzyskanych obliczeń, nie ma negatywnego wpływu w tym przykładzie na dokładność uzyskiwanych rozwiązań.
Praca naukowa finansowana ze środków budżetowych na naukę w latach 2005-2007 jako projekt badawczy 3T11F01528.
LITERATURA
1. Brebbia C. A., Telles J. C. F., Wrobel L. C.: Boundary element techniques, theory and applications in engineering. New York: Springer, 1984.
2. Foley J. D.: Wprowadzenie do grafiki komputerowej. Warszawa: WNT 2001.
3. Lyness J., Jespersen D.L.: Moderate degree symmetric quadrature rules for the trian- gle. “Journal of the Institute of Mathematics and its Applications”, Volume 15, Num- ber 1, February 1975, s. 19-32.
4. Wandzura S., Xiao H.: Symmetric quadrature rules on a triangle. “Computers and Mathematics with Applications”, Volume 45, s. 1829-1840, 2003.
5. Zieniuk E., Bołtuć A.: Non-element method of solving 2D boundary problems defined on polygonal domains modeled by Navier equation.” International Journal of Solids and Structures”, 2006, vol 43, s. 7939-7958.
6. Zieniuk E., Bołtuć A.: Krzywe Beziera w modelowaniu ciągłej geometrii brzegu w zagadnieniach brzegowych opisywanych równaniem Naviera. Prace Naukowe
„Transport" Politechniki Radomskiej nr 3(23), Radom 2005, s.561-566.
GLOBAL COMPUTATION OF DOMAIN INTEGRALS IN PIES FOR TWO-DIMENSIONAL BOUNDARY
PROBLEMS MODELLED
BY NAVIER-LAMÉ AND POISSON EQUATIONS
Summary. The paper presents a novel technique for global considerations and numerical integration of domains in 2D boundary problems. It base on computation of these integrals in global way, i.e. without division of the domain into cells. In proposed approach the domain is treated globally as single parametric surface and using numerical quadratures of high orders. Included numerical examples for boundary problems described by Poisson confirm high accuracy of proposed method compared with analytical results.
