Lab 08 2011 2012

Informatyka I Lab 08, r.a. 2011/2012
prow. Sławomir Czarnecki
Zadania na laboratorium nr. 8
Po utworzeniu nowego projektu, dołącz bibliotekę bibs.h.
1. Zdefiniuj funkcję void abx double** a,double* b,double* x,int n implementującą
( )
Algorytm Gaussa znajdowania rozwiązania x " !n liniowego układu równań
a x = b (1)
gdzie a " Mn�n jest daną, nieosobliwą macierzą rzędu n ( rank a = n ), b "!n jest danym
wektorem prawej strony.
Algorytm Gaussa szkic metody
Zakładamy, \e a0,0 `" 0 , modyfikujemy wiersze macierzy a " Mn�n oraz składowe wektora
ak ,0 1 ak ,0
b "!n o numerach k =1,2,...,n -1 obliczając a1 = ak , j - a0, j , bk = bk - b0 dla
k , j
a0,0 a0,0
j = k -1,k,..., n -1, tzn. macierz a i wektor b
a0,0 a0,1 ... a0,n-1 �ł łł
�ł łł b0
�ł
a1,0 a1,1 ... a1,n-1 śł �ł śł
b1
�ł śł
�ł śł
, (2)
�ł śł
... ... ... ... �ł śł
...
�ł śł
�ł śł
n-1,0
�łbn-1�ł
�ła an-1,1 ... an-1,n-1�ł
przekształcamy do macierzy a1 i wektora b1 równych odpowiednio
a0,0 a0,1 ... a0,n-1 b0
�ł łł �ł łł
�ł
a1,0 a1,0 a1,0 śł �ł a1,0 śł
�ł
a1,0 - a0,0 a1,1 - a0,1 ... a1,n-1 - a0,n-1 śł �ł b1 - b0 śł
�ł śł �ł śł
a0,0 a0,0 a0,0 a0,0
, . (3)
�ł śł �ł śł
... ... ... ... ...
�ł śł �ł śł
�ł
an-1,0 an-1,0 an-1,0 śł �ł an-1,0 śł
�ł - a0,0 an-1,1 - a0,1 ... an-1,n-1 - a0,n-1 bn-1 - b0
śł �ł śł
an-1,0
a0,0 a0,0 a0,0 �ł �ł a0,0 �ł
�ł śł �ł śł
�ł
Otrzymany w ten sposób, zmodyfikowany układ równań
a1 x = b1 (4)
charakteryzuje macierz a1 " Mn�n , w której wszystkie składowe w pierwszej kolumnie
poni\ej elementu a0,0 są równe 0:
a0,0 a0,1 ... a0,n-1
�ł łł b0
�ł łł
�ł 1 1 �ł 1 śł
0 a1,1 ... a1,n-1 śł b1
�ł śł
�ł śł
a1 = , b1 = (5)
�ł śł �ł śł
... ... ... ... ...
�ł śł
�ł śł
1
0 a1 ... a1 �ł
n-1,1 n-1,n-1
�łbn-1�ł
�ł
Opisane powy\ej przekształcenia, powtarzamy dla podmacierzy i wektora przedstawionych
poni\ej:
... ... ... ... ...
�ł łł ...
�ł łł
�ł 1 1 1
... a1,1 a1,2 ... a1,n-1 śł �ł 1 śł
b1
�ł śł
�ł śł
1
�ł śł �ł śł
... a1 a1 ... a1 , b2 . (6)
2,1 2,2 2,n-1
�ł śł
�ł śł
... ... ... ... ...
...
�ł śł
�ł śł
1
�ł
... a1 a1 ... a1 śł �łbn-1śł
n-1,1 n-1,2 n-1,n-1
�ł �ł
�ł �ł
1
Zakładając zatem, \e składowa a1,1 `" 0 , modyfikujemy wiersze macierzy a1 " Mn�n oraz
składowe wektora b1 "!n o numerach k = 2,3,...,n -1 obliczając
a1 a1
2 k ,1 1 2 1 k ,1 1
ak , j = a1 - a1, j , bk = bk - b1 dla j = k -1,k,..., n -1
k , j
1 1
a1,1 a1,1
1 1 1 1
�ł a1,1 a1,2 ... a1,n-1 łł �ł b1 łł
�ł śł �ł śł
a1 a1 a1 a1
�ł 2,1 1 2,1 1 2,1 1 1 1
2,1
a1 - a1,1 a1 - a1,2 ... a1 - a1,n-1 śł �ł b2 - b1 śł
2,1 2,2 2,n-1
1 1 1 1
�ł
a1,1 a1,1 a1,1 śł �ł a1,1 śł
�ł śł �ł śł
, . (7)
�ł ... ... ... ... śł �ł ... śł
�ł śł �ł śł
a1 a1 a1 a1
�ła1 n-1,1 1 n-1,1 1 n-1,1 1 n-1,1 1śł
- a1,1 a1 - a1,2 ... a1 - a1,n-1śł �łb1 - b1 śł
n-1,1 n-1,2 n-1,n-1
1 1 1 1
�ł
a1,1 a1,1 a1,1 śł �ł n-1 a1,1 �ł
�ł �ł �ł
Otrzymany w ten sposób, zmodyfikowany układ równań
a2 x = b2 (8)
charakteryzuje macierz a2 " Mn�n , w której dodatkowo, oprócz wszystkich składowych w
pierwszej kolumnie poni\ej elementu a0,0 , wszystkie składowe w drugiej kolumnie poni\ej
1
elementu a1,1 są tak\e równe 0
a0,0 a0,1 a0,2 ... a0,n-1
�ł łł b0
�ł łł
�ł 1 1 1 �ł 1 śł
0 a1,1 a1,2 ... a1,n-1 śł b1
�ł śł
�ł śł
2 2 2
�ł śł �ł śł
a2 = 0 0 a2,2 ... a2,n-1 , b2 = b2 (9)
�ł śł
�ł śł
... ... ... ... ...
...
�ł śł
�ł śł
2 2 2
�ł
�łbn-1śł
0 0 an-1,2 ... an-1,n-1śł
�ł �ł
�ł �ł
Kontynuacja opisanych powy\ej obliczeń prowadzi do modyfikacji układu równań z
początkowej postaci a x = b do postaci
an-1 x = bn-1 (10)
(z tzw. trójkątną macierzą układu równań), gdzie
a0,0 a0,1 a0,2 ... a0,n-1
�ł łł b0
�ł łł
�ł 1 1 1 �ł 1 śł
0 a1,1 a1,2 ... a1,n-1 śł b1
�ł śł
�ł śł
2 2 2
�ł śł �ł śł
an-1 = 0 0 a2,2 ... a2,n-1 , bn-1 = b2 . (11)
�ł śł
�ł śł
... ... ... ... ...
...
�ł śł
�ł śł
n-1 n-1
�ł
�łbn-1 śł
0 0 0 ... an-1,n-1śł
�ł �ł
�ł �ł
1 2
W przypadku gdyby któraś ze składowych a0,0 , a1,1 , a2,2 ,& była równa 0, konieczna jest
dodatkowa zmiana lub zmiany polegające na odszukaniu (w najprostszym wariancie)
pierwszego wiersza k ( poni\ej wiersza i -1 ), w którym składowa ak ,0 k > 0 , a1 k > 1 ,
( ) ( )
k ,1
2
ak ,2 k > 2 ,& jest ró\na od 0 i zamianie wiersza, w którym jest zerowa składowa ze
( )
znalezionym wierszem o numerze k. W praktyce, znajdowanie odpowiednich składowych jest
nieco bardziej zło\one i wynika z konieczności zapewnienia minimalizacji błędów
numerycznych, które sukcesywnie pojawiają się w trakcie przeprowadzania opisanych
powy\ej obliczeń.
Oznaczając na nowo macierz an-1 " Mn�n jaki i wszystkie jej składowe ain-1
, j
i, j = 0,1,..., n -1 , a tak\e wektor bn-1 "!n jak i wszystkie jego składowe bin-1
( )
i = 0,1,..., n -1 ponownie symbolami odpowiednio a = ai, j " Mn�n oraz b = bi " !n ,
( ) ( )
( )
stwierdzamy, \e zadanie poszukiwania rozwiązania układu równań liniowych (1)
sprowadzone zostało do zadania
a x = b (12)
z trójkątną macierzą a " Mn�n , gdzie
a0,0 a0,1 ... a0,n-1
�ł łł b0
�ł łł
�ł
�ł śł
0 a1,1 ... a1,n-1 śł b1
�ł śł
�ł śł
a = " Mn�n , b = " !n . (13)
�ł śł �ł śł
... ... ... ... ...
�ł śł
�ł śł
0 0 ... an-1,n-1�ł
�łbn-1�ł
�ł
Rozwiązanie x = xi " !n układu równań (12) mo\na w bardzo prosty sposób znalezć
( )
wykorzystując trójkątną postać macierzy a. Najpierw znajdujemy ostatnią składową xn-1 z
bn-1
ostatniego równania (13), xn-1 = , a kolejne składowe obliczamy rekurencyjnie ze
an-1,n-1
n-1
bi - ai, j xj
"
j=i+1
wzoru xi = i = n
( - 2, n - 3,...,1,0 .
)
ai,i
Algorytm Gaussa (najprostsza wersja) pseudo-kod
start
dla ka\dego wiersza i =1,2,...,n -1
{
jeśli ai-1,i-1 = 0 , to znajdz pierwszy wiersz k i d" k d" n -1 , w którym ak ,i-1 `" 0 , a
( )
następnie zamień wiersz i -1 z wierszem k aktualnej macierzy a oraz
zamień składową i -1 ze składową k aktualnego wektora b, tzn.
ai-1, j "! ak , j j = i -1,i,..., n -1 , bi-1 "! bk
( )
(mo\na ograniczyć się do zamiany składowych macierzy a jedynie od i -1 kolumny)
dla ka\dego wiersza k i d" k d" n -1
( )
{
ak ,i-1
c =
ai-1,i-1
ak ,i-1 = 0
dla ka\dej kolumny j i d" j d" n -1
( )
{
ak , j = ak , j - c ai-1, j
}
bk = bk - c bi-1
}
}
bn-1
xn-1 =
an-1,n-1
dla ka\dej składowej i = n - 2,n - 3,...,1,0
{
n-1
bi - ai, j xj
"
j=i+1
xi =
ai,i
}
koniec
2. Przetestuj działanie funkcji void abx double** a,double* b,double* x,int n dla
( )
macierzy A " M i wektora B " !N o losowo wygenerowanych składowych z przedziału
N�N
[-10,10 oraz N = 10 . Wyświetl na ekranie rozwiązanie X " !N układu równań A X = B , a
]
tak\e wygenerowaną macierz A " M oraz wektor B " !N przed i po wywołaniu funkcji
N�N
void abx ... .
( )
3. Dla wygenerowanej w zadaniu 2 macierzy A " M i wektora B " !N , wywołaj funkcję
N�N
int gauss int n,double** a,double* b,double* x z biblioteki bibs.h. Porównaj wynik z
( )
rozwiązaniem numerycznym X " !N z zadania 2.
4. Przetestuj działanie funkcji void abx double** a,double* b,double* x,int n dla
( )
następująco zdefiniowanej macierzy A " M2�2 i wektora B " !2 (N = 2)
� 1 1
�ł łł �ł łł
A = " !2 , (14)
�ł1 1śł " M2�2, B = �ł2śł
�ł �ł �ł �ł
gdzie � jest bardzo małą liczbą dodatnią, np.: � = 1.0�"10-20 .
X0
�ł łł
Komentarz. Rozwiązanie ścisłe (dokładne) X = " !2 układu równań
�ł śł
X1
�ł �ł
A X = B (15)
tzn. układu
� 1 X0 1
�ł łł �ł łł �ł łł
= (16)
�ł1 1śł �ł śł �ł2śł
X1
�ł �ł �ł �ł �ł �ł
jest następujące
1
�ł łł
2 -
�ł śł
�
1- X1
�ł łł
�ł1- 1 śł
�ł śł
�ł 1
1- śł
�ł łł
�
X0 �ł 1 śł �ł � śł �ł 1-� śł
�ł łł
�ł śł
X = = = ... = �ł śł = ... = . (17)
� �ł śł
�ł śł -
X1 2 �
�ł śł
�ł śł
�ł �ł �ł1- 2� śł
1
�ł śł
�ł - śł
1 2 �ł śł
�ł 1-� �ł
1-
�ł śł �
�ł śł
�ł � �ł
1
�ł śł
1-
�ł śł
�ł � �ł
Algorytm Gaussa (w najprostszej postaci) doprowadza układ równań (16) do układu równań
następującego:
� 1
�ł łł
X0 �ł 1 łł
�ł łł
�ł śł �ł śł
= , (18)
�ł śł
�ł0 1- 1 śł �ł2 1 śł
X1 -
�ł �ł
�ł � �ł �ł � �ł
który ma rozwiązanie ścisłe (dokładne)
1
2 -
1- X1
�
X1 = , X0 = . (19)
1
�
1-
�
Wa\ne ! Kolejność X1, X0 zapisanych w (19) składowych odpowiada kolejności w jakiej
obliczane są one w algorytmie Gaussa patrz pseudo-kod algorytmu Gaussa wy\ej.
1 1
Dla bardzo małego � > 0 , np. � = 1.0�"10-20 , wartości obu wyra\eń 2 - oraz 1- w (19)1
� �
1 1 1
obliczanych na komputerze są w przybli\eniu równe - , tzn. zarówno 2 - H" - jak i
� � �
1 1
1- H" - . A zatem ze wzoru (19)1 otrzymamy następujący wynik numeryczny
� �
1 1
2 - -
� �
X1 = H" = 1 , (20)
1 1
1- -
� �
i w konsekwencji ze wzoru (19)2 obliczymy
1- X1 1-1
X0 = H" = 0 . (21)
� �
Jednak, po prostym przekształceniu wzorów (19) (porównaj tak\e (17)) otrzymujemy
następujące, poprawne numerycznie wartości:
1- 2� 1
X1 = H"1, X0 = H"1. (22)
1-� 1- �
Otrzymane rozwiązanie numeryczne otrzymane najprostszą metodą Gaussa jest zatem nie do
zaakceptowania. Przyczyną, tak du\ego, numerycznego błędu jest fakt, \e składowa a0,0 = �
jest bardzo mała w porównaniu z innymi składowymi pierwszego wiersza (u nas w
porównaniu ze składową a0,1 =1), co nie mo\e być w prosty sposób poprawnie
zaimplementowane na komputerach o stałej długości słowa, w którym musi być
przechowywana reprezentacja liczb zmienno przecinkowych w postaci zarówno cechy jak i
mantysy.
Zauwa\my jednak, \e proste przestawienie kolejności równań (16)
1 1 X0 2
�ł łł �ł łł �ł łł
= (23)
�ł� 1śł �ł śł �ł1śł
X1
�ł �ł �ł �ł �ł �ł
a następnie zastosowanie algorytmu eliminacji Gaussa (w najprostszej wersji) prowadzi do
równowa\nego układu równań
1 1 X0 2
�ł łł �ł łł �ł łł
= , (24)
�ł0 1- � śł �ł śł �ł1- 2� śł
X1
�ł �ł �ł �ł �ł �ł
który tym razem, ma poprawne numerycznie rozwiązanie
1- 2�
X1 = H"1 (25)
1-�
oraz
X0 = 2 - x1 H" 1 . (26)
Wynika stąd, \e dobry pod względem numerycznym algorytm rozwiązywania układu
równań liniowych powinien w ogólnym przypadku brać pod uwagę mo\liwość m.in.
odpowiedniego przestawienia wierszy w celu eliminacji przedstawionych wy\ej
problemów numerycznych.
5. Raz jeszcze przetestuj działanie funkcji void abx ... , ale tym razem dla następująco
( )
zdefiniowanej macierzy A " M2�2 i wektora B " !2 (N = 2)
1 1 2
�ł łł �ł łł
A = " !2 (27)
�ł� 1śł " M2�2, B = �ł1śł
�ł �ł �ł �ł
(prosta zamiana kolejności równań (16)). Wyświetl na ekranie rozwiązanie X " !N .
Dla zdefiniowanej w (27) macierzy A " M i wektora B " !N , wywołaj raz jeszcze
N�N
funkcję int gauss int n,double** a,double* b,double* x z biblioteki bibs.h i wyświetl
( )
na ekranie rozwiązanie X " !N . Porównaj oba otrzymane wyniki z (niepoprawnym)
rozwiązaniem numerycznym (20), (21) oraz z (poprawnym) rozwiązaniem (25), (26).

Wyszukiwarka

Podobne podstrony:
Lab 11 12
Lab 11 12
Lab 11 12
Lab 11 12
Lab ME MI instrukcja 11 12 E
dach (11 12)
zjazdy 11 12
Quas primas Pius XI (11 12 1925)
Konsultacje sem letnim 11 12 I16# 12
Hydrologia cwiczenia 11 i 12
11 (12)

więcej podobnych podstron