Strona domowa ćwiczeń z Metody Numerycznych II A
w roku akademickim 2001/02.

Metody Numeryczne II A

Zasady zaliczania

• 57% (57 pkt.) - zadania na ćwiczeniach (programy)
• 30% (30 pkt.) - test z wiadomości z wykładu
• 13% (13 pkt.) - obecność na ćwiczeniach

Ćwiczenia odbywają się w czwartki w godzinach 9:15 - 12:00 w sali OKWF B.

Wiadomości

Tutaj będą pojawiały się różne informacje, np. upłynięcie terminu oddawania zadania domowego, termin testu, informacje o odwołaniu zajęć czy informacje o wynikach (czy też informacje o zmianach na stronie).

• Po uwzględnieniu zadań przysłanych za pomocą e-mail i oddanych osobiście do dnia 14-06-2002 punktacja uległa zmianie. Poprawione wyniki zostały uwzględnione w tabeli wyników oraz Liscie obecności i oddanych zadań. Zmianie uległy także oceny niektórych osób. Wyniki wysłałem do dr Wernera.
  
• dr hab. Tomasz Werner (wykładowca) wyjeżdża w okolicach 16. czerwca,
  więc jak ktoś odda zadania na ostatnią chwile może czekać na wpisanie zaliczenia.
• Zadania domowe można oddawać do 16. czerwca (tydzień przed końcem sesji letniej):
  muszę mieć czas na umówienie się/sprawdzenie zadań.
  Po upływie tego terminu nie przyjmuję zadań domowych!.
  Przypominam, że zadania obowiązkowe będą oceniane na co najwyżej 0.8 punkta "zadaniowego".
• Tutaj znaleźć można skróconą listę wyników, (grupa 3.) z poprawkami.
• Tutaj znaleźć można listę obecności i oddanych zadań, (grupa 3.) stan na dzień 14. czerwca.
• Ostateczne progi ocen oraz zaliczenia
  

  progi	44	54	64	74	84	94
  oceny	3	3+	4	4+	5	5+

  
  
• Za test można dostać 30 pkt. (10 pytań * 3 punkty)
• Zostały uaktualnione zasady zaliczania (punktacja) (patrz wyżej)
• Pojawiła się lista zadań domowych z linkiem do pełnego opisu każdego zadania
• Zostały poprawione skrypty sortowanie.ps.gz (tex), gauss_LU.ps.gz (tex), eigen.ps.gz (tex), eigenvalues.ps.gz (tex), ode-RK.ps.gz (tex) (dawniej ode.ps.gz)
• W informacjach o oddawaniu zadań za pomocą e-mail będą pojawiać się informacje o przysłanych w ten sposób zadaniach
• Pojawiła się nowa propozycja zadania zaliczeniowego (LIC, line integral convolution)
• Proszę o zgłaszanie propozycji zajęcia się zadaniami zaliczeniowymi (jak i zgłaszanie własnych zadań)
• W poradach pojawił się wstęp do programowania modularnego i programu make
• Pojawiły się pierwsze propozycje zadań zaliczeniowych

Spis treści

• Zasady zaliczania
• Wiadomości
• Spis treści
• Zadania domowe
• Oddawanie zadań za pomocą e-mail
• Zadanie zaliczeniowe
• Spis zadań domowych
• Zagadnienia przedstawiane na ćwiczeniach
  • 21-02: Sortowanie.
  • 28-02: Rozwiązywanie układu równań liniowych. Wstęp.
  • 07-03: Rozwiązywanie układu równań liniowych. Metoda Gaussa.
  • 14-03: Rozwiązywanie układu równań liniowych. Metody iteracyjne.
  • 21-03: Znajdowanie wartości i wektorów własnych macierzy.
  • 11-04: Równania różniczkowe zwyczajne. Metody jednokrokowe.
  • 25-04: Równania różniczkowe zwyczajne. Metody wielokrokowe.
  • 09-05: Równania różniczkowe zwyczajne. Metody predyktor-korektor.
  • 16-05: Transformata Fouriera.
• Planowane przyszłe zadania domowe
• Porady
• Linki
• Literatura

Zadania (programy) domowe

Programy będą oceniane w skali ciągłej [0,1]. Za rozszerzenie zadania domowego można dostać od 0 do 1/3 punkta. Na koniec zajęć ilość punktów uzyskanych za oddane programy będzie dzielona przez ilość obowiązkowych zadań domowych i przeliczana do 57% (57 pkt.). Dodatkowo za oddanie przynajmniej jednego zadania domowego osobiście będą przyznawane dodatkowe punkty.

Jest 6. zadań domowych i jedno rozszerzenie za 1/3 punkta. Oby otrzymać ilość punktów za zadania domowe należy przemnożyć sumę otrzymanych ocen przez 8, oraz dodać 6 1/3 punkta jeśli oddało się przynajmniej jeden program osobiście (uwaga: mogę wymagać osobistego oddania programu!). (6 1/3)*8 + 6 1/3 = 57 pkt. Za zadania dodatkowe (nieobowiązkowe) przyznaję odpowiednio: za zadanie 4. (wszyskie wartości własne) 10 pkt, za zadanie 7. (metoda RK ze zmiennym krokiem) 10 pkt, zaś za zadanie 8. (predyktor-korektor) 8 pkt.

Za szczególnie "ładne" rozwiązania można dostać dodatkowe punkty. Możliwe jest zatem uzyskanie za programy więcej niż 57% (57pkt). Każde zadanie domowe ma termin oddania (celem tego rozwiązania jest uniknięcie oddawania wszystkich zadań na ostatnich i przedostatnich zajęciach). Po przekroczeniu tego terminu punktacja obniżana będzie do maksymalnie 0.8 pkt (punkt o szczególnie "ładnych" rozwiązaniach ma zastosowanie także tutaj). Dokładniej ocena jaką by się dostało za zadanie oddane w terminie jest mnożona przez 0.8. Uwaga: rozszerzenia zadań domowych oraz zadania nieobowiązkowe nie mają terminów.

W poradach znajduje się opis preferowanego sposobu (stylu) kodowania (zapisywania programu). Nie jest on obowiązkowy. Służy on głównie jako wskazówka, jak należy pisać programy.

Celem każdego zadania domowego powinno być napisanie zadanego algorytmu (metody) oraz jego przetestowanie.

Oddawanie zadań domowych za pomocą poczty elektronicznej

Programy można wysyłać pocztą (mój adres e-mail to jnareb@fuw.edu.pl), ale preferowane jest oddawanie programów osobiście. Uwaga: proszę zaznaczać które zadanie domowe rozwiązuje przysłany program. Przynajmniej połowa programów powinna (ale nie musi) być oddana osobiście (czy to na ćwiczeniach, czy też umówiwszy się ze mną na konkretny termin). Oceniane jest działanie programu oraz sposób jego zapisania (t.j. ocena bazuje także na źródłach programu). W większości przypadków styl programu nie ma wpływu na ocenę.

Wadą oddawania programów za pomocą poczty elektronicznej jest to, że mogę nie mieć czasu na sprawdzenie przysłanych programów. Zawsze gdy program zostanie sprawdzony i oceniony, wysyłam propozycję oceny i uwagi co do programu na adres zwrotny, tak więc proszę sprawdzać pocztę. Chwilę przysłania programu który wykonuje zlecone działanie uznaję za czas oddania. W razie potrzeby mogę wysyłać potwierdzenie odbioru (proszę podać w treści listu, że się chce otrzymać potwierdzenie odbioru).

Wolałbym, by programy napisane przez więcej niż jedną osobę były oddawane na ćwiczeniach. Przysłanie programu za pomocą e-mail może służyć jako oddanie programu w sensie terminu jego oddania (napisania).

Jeśli czyjś list nie dotarł, to proszę go wysłać ponownie z dopiskiem kiedy został wysłany po raz pierwszy. Można prosić o potwierdzenie odbioru. Na każdych ćwiczeniach dostępna jest lista oddanych i sprawdzonych oraz przysłanych za pomocą e-mail zadań domowych.

Docierają do mnie listy wysyłane z adresów: astrouw.edu.pl, tempac.okwf.fuw.edu.pl, ds2.uw.edu.pl, go2.pl, konto.pl, poczta.fm. Najprawdopodobniej nie docierają listy wysyłane bezpośrednio z komputerów w pracowni (bpcX, apcX) ani z primusa.

Zadanie zaliczeniowe

Zaliczyć ćwiczenia można będzie także oddając program zaliczeniowy. Program na zaliczenie powinien obejmować zakresem przynajmniej dwa zagadnienia (lub jedno mocno rozszerzone) oraz rozwiązywać jakiś fizyczny problem. Przykłady poniżej (więcej informacji udzielam osobiście lub przez e-mail). Można także zgłaszać własne propozycje. Program (kod źródłowy) powinien zostać udostępniony do celów dydaktycznych (polecam wybór którejś z otwartych licencji np. GPL, lub uczynienie programu jawnie public domain).

Programy zaliczeniowe można oddawać zarówno u mnie (w mojej grupie), jak i u mgr Tomasza Pawłowskiego (w grupie 2.). Proszę zgłaszać zamiar oddania programu zaliczeniowego!

• Symulacja ruchu planet (najlepiej w 3D).
  Metody: rozwiązywanie układu równań różniczkowych zwyczajnych metodą ze zmiennym krokiem (kontrola dokładności metody). Przynajmniej jedna metoda ze zmiennym krokiem powinna być napisana osobiście. W wersji uproszczonej (trzeba to jakoś skompensować) można skorzystać z kilku metod ze stałym krokiem, w tym metody LeapFrog.
  Grafika: Najlepiej 3D (z rzutowaniem prostym/perspektywicznym do wyboru) z możliwością obrotów, przesuwania, centrowania na obiekcie. Powinna istnieć możliwość włączenia [zanikających] śladów.
  Interfejs: Powinna istnieć możliwość zadania (np. z pliku) układu początkowego. Podczas symulacji powinna być możliwe zmienianie kroku dla metod stałokrokowych i dokładności dla metod ze zmiennym krokiem. Jeśli zaimplementowanych zostało kilka metod powinno być możliwa ich zmiana w trakcie symulacji.
  Fizyka: Program powinien wyświetlać (umożliwiać wyświetlanie) całkowitą energię i moment pędu (które powinny być zachowane).
  Przykłady:
  1. Układ podwójny (dwie gwiazdy i planeta): badanie własności. Orbita ciasna wokół jednej gwiazdy, orbita szeroka wokół obu gwiazd.
  2. Układ planetarny w którym nie wszystkie planety mają tę samą płaszczyznę ekliptyki; planeta posiada księżyc o płaszczyźnie orbity prostopadłej do ekliptyki. Tutaj przydałaby się możliwość oglądania symulacji w układzie nieinercjalnym związanym ze słońcem, planetą, księżycem.
  3. Zagadnienie trzech ciał (lub zredukowane zagadnienie trzech ciał, z ciałem próbnym nie wpływającym na dwa pozostałe) oglądane w układzie w którym dwa z nich się nie poruszają. Uwaga na niestabilności numeryczne!
  
  
• Rozwiązanie stacjonarnego jednowymiarowego równania Schroedingera (stany związane). Dla chętnych dwuwymiarowe/trójwymiarowe.
  Metody: Adaptacyjne (z kontrolą dokładności) całkowanie numeryczne dla znalezienia elementów macierzowych. Znajdowanie wszystkich (lub kilku najniższych) wartości własnych macierzy hermitowskiej (symetrycznej w większości przypadków). Znajdowanie wektorów własnych. Przynajmniej jedna z metod powinna być napisana osobiście (być może jako jedna do wyboru).
  Grafika: Wystarczy wykres potencjału (ew. z poziomami energetycznymi) oraz kilku funkcji falowych stanów własnych (uwaga: funkcja zespolona). Może to być zaimplementowane za pomocą wypisywania wyników do pliku i użycia zewnętrznego programu (np. gnuplot).
  Interfejs: Powinna być możliwość wyboru dokładności przez użytkownika, jakiś element wyboru potencjału (zadanie parametru, wybranie funkcji z listy, wybranie poprawki z listy, ew. podanie wzoru lub zadanie potencjału przez podanie wartości w punktach w pliku), ew. wybór bazy w której liczymy elementy macierzowe hamiltonianu.
  Fizyka: Funkcje własne powinno się normować.
  Przykłady:
  
  1. Potencjał anharmoniczny (a x^2 + b x^3 + c x^4, a >> b,c lub inny podobny), z możliwością wyboru anharmoniczności. Jako bazę można wziąć funkcje własne potencjału harmonicznego (wielomiany Hermite'a). Porównanie z potencjałem harmonicznym.
  2. Potencjał z dwiema studniami potencjału (a x^2 + b x^4), ew. lekko przechylony (a x^2 + b x^3 + c x^4). Baza jak wyżej.
  3. Gładka studnia potencjału. Jako bazę można wziąć fale płaskie (baza Fouriera), cosinusy lub bazę dla prostokątnej funkcji potencjału. Porównanie z prostokątną studnią.
  4. Trójkątna (piłokształtna) studnia potencjału, skończona bądź nieskończona. Reszta jak wyżej.
  
  
• Rozwiązywanie równania Schroedingera zależnego od czasu (jednowymiarowego lub lepiej (ale trudniej) dwuwymiarowego).
  Metody: Metoda SOD, z przybliżeniem laplasjanu za pomocą różnicy centralnej (metoda różnic skończonych), zapewniająca zachowanie normy (szczegóły na życzenie). Metoda Split FFT, korzystająca z operatora ewolucji (szczegóły na życzenie). Nie potrzeba samemu implementować szybkiej transformaty Fouriera.
  Grafika: Zobrazować ewolucję funkcji falowej w czasie (uwaga: funkcja zespolona); wykres funkcji w czasie (można to ewentualnie zrobić za pomocą innego programu, czy to zapisując wyniki do pliku jako osobne serie, czy to periodycznie zapisując wyniki do pliku/kolejki FIFO (named pipe), czy to wyrzucając wyniki na standardowe wyjście). Zagadnienie dwuwymiarowe jest trudniejsze do zwizualizowania, ale ciekawsze.
  Interfejs: Powinna istnieć możliwość wyboru metody rozwiązywania jak i kroku/dokładności metody.
  Fizyka: funkcje falowe powinny być unormowane. Program powinien kontrolować normę funkcji falowej (np. pokazując jaka jest bądź jaka byłaby bez normalizacji norma funkcji falowej; w postaci bieżącej wartości jak i w postaci wykresu (na żądanie)). Przydało by się, aby jako jeden z testów użyć ewolucji swobodnej paczki gaussowskiej, jako drugiego ewolucji stanu własnego [podstawowego] dla jakiegoś potencjału.
  Przykłady:
  1. Rozpraszanie paczki falowej (gaussowskiej) lub innego stanu (np. fali płaskiej) na zadanym potencjale (np. prostokątnej studni potencjału czy garbie potencjału; uwaga: mogą być problemy z niegładkością). Można porównać z wynikiem analitycznym.
  2. Ewolucja funkcji falowej (np. paczki gaussowskiej) nie będącej stanem własnym w stacjonarnym potencjale, np. dla potencjału dwujamowego jako warunek początkowy funkcja gaussowska (stan podstawowy potencjału harmonicznego) w jednym z dołków.
  3. Ewolucja funkcji falowej w potencjale zmiennym w czasie (być może stanu podstawowego potencjału stacjonarnego). Przykładem może być okresowo przechylany potencjał dwujamowy bądź harmoniczny (a x^2 + (b + b' sin(w t)) x^3 + c x^4) bądź studnia potencjału o zmiennych parametrach (szerokości, głębokości) bądź przechylanym dnie. Zależność potencjału od czasu może dla przykładu być okresowa bądź wykładnicza (zanikająca wykładniczo).
  4. Zderzenie ze sobą dwu oddziałujących cząstek (trudne).
  
  
• Wizualizacja dwuwymiarowego pola wektorowego metodą splotu wzdłuż linii pola (LIC, line integral convolution)
  Metody: Rozwiązywanie równań różniczkowych zwyczajnych metodą Runge-Kutty 4. rzędu (w celu znalezienia linii pola). Całkowanie metodą trapezów (do obliczenia splotu). Interpolacja liniowa, ewentualnie interpolacja za pomocą spline'ów parametrycznych (krzywe Beziera) do znalezienia wartości w punktach całkowania (równoodległe).
  Grafika: Dla każdego punktu pewnego obszaru należy wyliczyć całkę:
  
  LIC(P,K)(x) = 1/N \int_{-e}^{e} K(t) P(g_x(t)) dt
  gdzie K(t) oznacza funkcję odpowiedzi (splatany sygnał; funkcję rozmazywania), g_x(t) oznacza punkt znajdujący się na linii pola przechodzącej przez x, odległy o czas t od niego, zaś P(x) jest pewną wartością (kolorem losowego obrazu). N oznacza normę funkcji odpowiedzi K(t). Potrzebne będzie zatem zrobienie (kolorowej lub w odcieniach szarości) bitmapy; może nie być wyrzucana na ekran a tylko zapisywana do pliku, np. za pomocą biblioteki gd jako plik gif. Interfejs: Użytkownik powinien móc wybrać jedną z predefiniowanych funkcji odpowiedzi (np. K(t) = 1 lub K(t) = x < 0 ? x+1 : 0, obie zadane dla e = 1). Przy zapisywaniu do pliku użytkownik powinien móc podać jego nazwę. Przy pokazywaniu wyników na ekranie mile widziana byłaby możliwość interaktywnej pracy: zwiększania/zmniejszania długości rozmycia, zmiany oglądanego obszaru, zwiększenia/zmniejszenia ilości punktów w obrazie P(x) itp.
  Komputerowe: Należy pamiętać o obcinaniu kolorów (wyników) do odpowiedniego zakresu.
  Przykłady:
  1. Wizualizacja pola elektrostatycznego między dwoma równoimiennymi (różnoimiennymi) ładunkami. Jako obraz P(x) można wziąć losowo rozmieszczone punkty, kolorowane za pomocą wartości potencjału (za pomocą jakiejś mapy kolorów).
  2. Wizualizacja różnych rodzajów stabilnych punktów granicznych: zlewu, źródła, siodła i ich wersji zdegenerowanych.
  
  

Własne propozycje też powinny mieć podobną strukturę, tzn. jako podstawę mieć jakąś (jakieś) metodę numeryczną, rozwiązywać (prezentować) jakieś zagadnienie fizyczne i do tego mieć dobry interfejs użytkownika. Program na zaliczenie powinien obejmować zakresem przynajmniej dwa zagadnienia (lub jedno mocno rozszerzone). Przynajmniej jedna metoda numeryczna powinna być napisana samodzielnie (propozycja powinna przedstawiać minimum pracy własnej). Przynajmniej jedna metoda powinna być z zagadnień poruszanych w obecnym semestrze, tzn. sortowania, równań liniowych, wartości i wektorów własnych, transformaty Fouriera, filtrowania, splotu, równań różniczkowych zwyczajnych (lub cząstkowych). Zagadnienie fizyczne może być zastąpione czymś innym, ale zawsze powinno to być zastosowanie danej metody numerycznej, a nie tylko pokazanie jej własności. Interfejs użytkownika powinien być w miarę interakcyjny; jeśli jednak umiejętności nie pozwalają na implementację grafiki można się (jednak nie we wszystkich zadaniach) obejść.

Uwaga: zrobienie zadania zaliczeniowego nie daje punktów za obecności i test.

Spis zadań domowych

• Zad. 1: Porównanie sortowania o koszcie O}(N^2) (np. sortowania przez wstawianie) i o koszcie O(N log N) (np. sortowania stogowego).
  Pełny opis zadania: 21-02 , 1.0*8 = 8 pkt.
• Zad. 2: Rozwiązywanie układu równań liniowych za pomocą rozkładu L U z częściowym (dla chętnych: pełnym) wyborem elementu podstawowego.
  Pełny opis zadania: 07-03 , 1.0*8 = 8 pkt.
• Zad. 2': Iteracyjne poprawianie rozwiązania układu równań liniowych do uzyskania zadanej dokładności (lub przekroczenia maksymalnej liczby iteracji) (uzupełnienie).
  Pełny opis zadania: 14-03 , 1/3*8 = 2 2/3 pkt.
• Zad. 3: Znajdowanie wektorów i wartości własnych macierzy za pomocą metody potęgowej i odwrotnej metody potęgowej.
  Pełny opis zadania: 21-03 , 1.0*8 = 8 pkt.
• Zad. 4^*: Znajdowanie wszystkich wartości własnych macierzy symetrycznej za pomocą transformacji podobieństwa (nieobowiązkowe).
  Pełny opis zadania: 21-03 , 1.0*10 = 10 pkt.
• Zad. 5: Rozwiązanie zagadnienia początkowego dla równania ruchu tłumionego oraz oscylatora harmonicznego metodą Eulera, ulepszoną metodą Eulera (midpoint) oraz metodą Heuna.
  Pełny opis zadania: 11-04 , 1.0*8 = 8 pkt.
• Zad. 6: Rozwiązanie zagadnienia początkowego dla oscylatora harmonicznego i wahadła matematycznego metodą Runge-Kutty i metodą wielokrokową.
  Pełny opis zadania: 25-04 , 1.0*8 = 8 pkt.
• Zad. 7^*: Metoda Runge-Kutty (lub dowolna inna) ze zmiennym (adaptacyjnym) krokiem (ew. dodatkowo ze zmiennym stopniem) (nieobowiązkowe).
  Pełny opis zadania: 09-05 , 1.0*10 = 8 pkt.
• Zad. 8^*: Metoda typu predyktor-korektor (nieobowiązkowe).
  Pełny opis zadania: 09-05 , 1.0*8 = 8 pkt.
• Zad. 9: Badanie własności dyskretnej zespolonej transformaty Fouriera.
  Pełny opis zadania: 16-05 , 1.0*8 = 8 pkt.

Zagadnienia przedstawiane na ćwiczeniach

1. Sortowanie tablicy

21-02-2002

1.1. Opis teoretyczny

• W pliku sortowanie.ps.gz (tex) opisane są metody sortowania przez wstawianie, selekcję oraz sortowanie stogowe. Podany jest też sposób wykorzystania bibliotecznej funkcji qsort do sortowania szybkiego (quicksort).

1.2. Zadania

1. Napisać procedurę i przykładowy program sortujące tablicę liczb zmiennoprzecinkowych (float) w porządku malejącym dowolnie wybraną metodą.
2. Napisać program pokazujący sortowanie krok po kroku (metoda sortowania dowolna); rozmiar tablicy może być mały. Można oznaczyć w jakiś sposób bieżący element (bieżący kawałek tablicy).

1.3. Zadanie domowe

Napisać dwie procedury sortujące, jedną sortującą w czasie kwadratowym O(N), drugą sortująca w czasie O(N log N). Porównać szybkości dwu zaimplementowanych algorytmów dla różnych rozmiarów tablicy. W przypadku korzystania z grafiki (może być zapis wyników do pliku i wykorzystanie programu gnuplot) zrobić wykres zależności czasu realizacji od rozmiaru tablicy (tablica losowa) w odpowiedniej skali. W przypadku nie korzystania z grafiki porównać czasy dla kilku znaczących wielkości tablicy oraz porównać zachowanie dla tablicy [prawie] uporządkowanej, uporządkowanej [prawie] odwrotnie i ułożonej w porządku losowym.

Jako oszacowanie szybkości można przyjąć ilość wykonanych porównań, ilość wykonanych zamian (jeśli można), albo zmierzyć czas spędzony na wykonywaniu sortowania np. za pomocą funkcji clock z time.h lub times z sys/times.h (na ogół nie daje dobrych wyników ze względu na małą rozdzielczość).

Czas oddania: do połowy marca. Termin już minął, za zadanie można otrzymać maksymalnie 0.8 pkt.

2. Programowanie modularne. Wstęp do rozwiązywania układu równań liniowych metodą Gaussa.

28-02-2002

2.1. Opis teoretyczny

2.2. Opis teoretyczny

• Opis metody Gaussa rozwiązywania układu równań liniowych znaleźć można w gauss_LU.ps.gz (tex).
• Opis programowania modularnego w C/C++ (w tym używanie programu make znaleźć można w poradach.

2.3. Zadania

1. Napisać procedurę (funkcję) mnożącą macierz kwadratową przez wektor, która będzie służyć do sprawdzania poprawności wyników.
2. Napisać procedurę (funkcję), która rozwiązuje układ równań liniowych z macierzą trójkątną górną. Sprawdzić otrzymany wynik. Kod tej procedury przyda się w następnym zadaniu; krok ten można pominąć.
3. Napisać procedurę (funkcję), która mając daną macierz w postaci rozkładu na macierz trójkątną dolną z jednynkami na diagonali i macierz trójkątną górną (rozkład LU) rozwiązuje układ równań z tą macierzą. Napisać procedurę sprawdzającą wynik (czy to przez wymnożenie macierzy składowych, czy przez mnożenie przez wektor za pomocą rozkładu).
4. Opakować procedury operujące na rozkładzie LU w niezależnie kompilowany moduł.

3. Rozwiązywanie układów równań liniowych. Rozkład LU metodą Gaussa-Crouta i Dollitle'a.

07-03-2002

3.1. Opis teoretyczny

• W pliku gauss_LU.ps.gz (tex) opisane są metody rozwiązywania układów równań lioniowych. Opisana została metoda rozwiązywania układów z macierzą trójkątną, metoda Gaussa oraz metoda Dollitle'a. Przedstawiony został także problem wyboru elementu podstawowego.

3.2. Zadania

1. Rozwiązać układ równań (klasyczną) metodą Gaussa dla zadanego wektora prawych stron. Można założyć, że nie trzeba dokonywać przestawień wierszy.

3.3. Zadanie domowe

Opis zadania domowego znaleźć można w gauss_LU.ps.gz (tex). Poniżej skrócony opis.

Należy napisać funkcję dokonującą rozkładu L U macierzy z częściowym wyborem elementu podstawowego oraz funkcję rozwiązująca układ równań A x = b na podstawie rozkładu obliczanego przez poprzednią funkcję.

Funkcja licząca rozkład L U macierzy powinna go zapisywać w jednej kwadratowej tablicy (w postaci "spakowanej"). Funkcja może nadpisywać macierz A. Jako jeden z parametrów powinien byc przekazywany wektor na którym będzie zapisywana permutacja wierszy. Funkcja powinna zwracać 0 jeśli rozkład się udał, w przeciwnym wypadku powinna zwrócić numer wiersza (liczony od 1) macierzy w którym nastąpił błąd.

Funkcja rozwiązująca układ równań powinna jako argumenty przyjmować wartości obliczone przez poprzednią funkcję. Funkcja może nadpisywać wektor b prawej strony. Funkcja powinna zwracać 0 jeśli udało się rozwiązać układ równań, w przeciwnym wypadku powinna zwrócić numer wiersza (liczony od 1) macierzy w którym nastąpił błąd.

Funkcje mogą używac jako zmiennych (wektorów) roboczych dodatkowych parametrów. Sposób przekazywania macierzy do funkcji pozostawia się piszącemu program (patrz "Porady dotyczące programowania" w poradach). Nie należy do tego używać zmiennych globalnych.

Obie funkcje powinny zostać przetestowane w pisanym programie, tzn. powinno zostać sprawdzone, czy rozkład L U jest rozkładem dla macierzy A o odpowiednio przestawionych wierszach, tzn. czy A = P L U oraz czy rozwiązanie układu spełnia go tzn. czy A x = b. Dla chętnych: wypisać błędy obliczeń.

Czas oddania: do drugich zajęć po świętach (11 kwietnia). Termin już minął, za zadanie można otrzymać maksymalnie 0.8 pkt.

4. Rozwiązywanie układów równań liniowych. Metody iteracyjne.

14-03-2002

4.1. Opis teoretyczny

• W pliku iterate_lin.ps.gz (tex) opisana została metoda iteracyjnego poprawiania rozwiązania układu równań liniowych. Przedstawione zostały także metody iteracyjne rozwiązywania układów równań (stosowane m. in. do rozwiązywania układów z macierzami rzadkimi) oparte o podział macierzy.

4.2. Zadania

1. Rozwiązać układ równań z dużą macierzą rzadką, diagonalnie dominującą, symetryczną, dodatnie określoną za pomocą którejś z metod iteracyjnych, np. metody Jacobiego lub Gaussa-Seidla.

4.3. Zadanie domowe

Opis rozszerzenia poprzedniego zadania domowego znaleźć można w iterate_lin.ps.gz. Poniżej skrócony opis.

Napisać funkcję dokonującą jednego kroku iteracyjnego poprawienia rozwiązania równania A x = b oraz funkcję poprawiającą uzyskane rozwiązanie do odsiągnięcia zadanej dokładności (lub przekroczenia maksymalnej liczby iteracji). Jako warunek stopu można przyjąć ||r|| < eps ||A x||. Funkcje mogą używac jako zmiennych (wektorów) roboczych dodatkowych parametrów. Funkcje mogą nadpisywać swoje parametry ale powinno zostać to opisane.

Maksymalna ocena za to uzupełnienie wynosi 1/3 pkt.

5. Znajdowanie wartości i wektorów własnych macierzy.

21-03-2002 - 04-04-2002

5.1. Opis teoretyczny

• W pliku eigen.ps.gz (tex) przedstawione zostały metoda potęgowa znajdowania największej co do modułu wartości własnej i odpowiadającego jej wektora oraz odwrotna metoda potęgowa, za pomocą której możemy znaleźć wartość własną najbliższą zadanej wartości i odpowiadający jej wektor własny. Obie metody wymagają, by znajdowana wartość własna była jednokrotna aby znaleziony wektor był wektorem własnym.
• W pliku eigenvalues.ps.gz (tex) przedstawiony został sposób diagonalizacji macierzy za pomocą transformacji podobieństwa. Przedstawiona została transformacja Householdera. Opisano sprowadzanie macierzy symetrycznej za pomocą transformacji Householdera do macierzy trójdiagonalnej symetrycznej oraz metoda diagonalizacji macierzy trójdiagonalnej przy pomocy obrotów (metoda Jacobiego).

5.2. Zadania

1. Rozwiązać układ równań z dużą macierzą rzadką, diagonalnie dominującą, symetryczną, dodatnie określoną za pomocą którejś z metod iteracyjnych, np. metody Jacobiego lub Gaussa-Seidla.
2. Znaleźć największą co do modułu wartość własną (i odpowiadający jej wektor własny) zadanej macierzy za pomocą metody potęgowej z zadaną ilością iteracji.

5.3. Zadanie domowe

Napisać program implementujący metodę potęgową i odwrotną metodę potęgową. Program powinien umożliwiać co najmniej podanie ilości iteracji albo żądanej dokładności oraz wartości w pobliżu której będziemy szukać wartości własnej. Program powinien wypisywać co najmniej macierz której wartości i wektorów własnych poszukujemy, znaleziony wektor własny oraz wartość własną. Jako oszacowania błedu można użyć ||A v - lambda v||. Do odwrotnej metody potęgowej program powinien używać rozkładu LU lub podobnego do rozwiązywania układu równań.

Czas oddania: do 20 kwietnia. Termin już minął, za zadanie można otrzymać maksymalnie 0.8 pkt.

Dla chętnych: napisać program znajdujący wszystkie wartości własne macierzy symetrycznej za pomocą transformacji podobieństwa.

6. Równania różniczkowe zwyczajne. Metody jednokrokowe.

11-04-2002 - 18-04-2002

6.1. Opis teoretyczny

• W pliku ode-RK.ps.gz (tex) (jak Ordinary Differential Equations) znajduje się opis zagadnień początkowych dla równań różniczkowych zwyczajnych, przedstawione zostały metoda Eulera, midpoint, Heuna oraz metoda Runge-Kutty 4. rzędu.
• Skrypt do tworzenia "biegnących" wykresów przy pomocy programu gnuplot. Należy go wydobyć za pomocą komendy tar xvf wahad.tar; skrypt znajdzie się w pliku wahad1.scr. Wersja poprawiona i opatrzona kometarzami znajduje się tutaj. Należy pamiętać o daniu prawa do wykonywania skryptu za pomocą chmod a+x wykresik.sh. Skrypt uruchamiamy tak jak program, t.j. wpisując ./wykresik.sh nazwa_pliku_z_danymi.

6.2. Zadania

1. Rozwiązać za pomocą metody Eulera równanie ruchu tłumionego (ruchu lepkiego). Porównać rozwiązanie z rozwiązaniem analitycznym.

6.3. Zadanie domowe

Napisać program rozwiązujący zagadnienie początkowe dla równania ruchu tłumionego oraz oscylatora harmonicznego metodą Eulera, midpoint i Heuna. Porównać wynik z rozwiązaniem analitycznym. Można sprawdzić rząd zbieżności metody albo zrobić wykres.

Czas oddania: do połowy maja. Termin już minął, za zadanie można otrzymać maksymalnie 0.8 pkt.

7. Równania różniczkowe zwyczajne. Metody wielokrokowe.

25-04-2002

7.1. Opis teoretyczny

• W pliku rozn.ps znajduje się krótkie wprowadzenie do równań różniczkowych zwyczajnych i wstęp do metod różnicowych (numerycznych metod rozwiązywania równań różniczkowych w określonych punktach). Przedstawiona (i opisana) jest metoda Eulera i jej modyfikacje: ulepszona metoda Eulera (metoda Heuna) oraz zmodyfikowana metoda Eulera. Opisana jest ogólna postać metody Runge-Kutty i podana jest klasyczna metoda Runge-Kutty. Opisane są także liniowe metody wielokrokowe oraz metody typu predyktor-korektor.
• W pliku ode-wielokrokowe.ps.gz (tex) znajduje się poprawiona wersja powyższego skryptu. Dodatkowo przedstawiony został sposób wyznaczania wzorów na metody wielokrokowe za pomocą kwadratur Newtona-Cotesa.

7.2. Zadania

1. Rozwiązać za pomocą metody wielokrokowej równanie oscylatora harmonicznego. Jako kroki początkowe można wziąć rozwiązanie analityczne. Zobaczyć jaki jest rząd metody.

7.3. Zadanie domowe

Napisać program rozwiązujący zagadnienie początkowe dla oscylatora harmonicznego x'' = - k x oraz wahadła matematycznego x'' = - k sin(x) metodą Runge-Kutty 4. rzędu oraz jawną metodą wielokrokową 4. rzędu. Metodę wielokrokową można wystartowac za pomocą jawnego rozwiązania analitycznego (dla oscylatora harmonicznego) bądź metody Runge-Kutty. Przeanalizować zachowanie (np. za pomocą wykresu fazowego) dla różnych wartości prędkości początkowej.

Czas oddania: do 23. maja. Termin już minął, za zadanie można otrzymać maksymalnie 0.8 pkt.

8. Równania różniczkowe zwyczajne. Metody typu predyktor-korektor.

09-05-2002

8.1. Opis teoretyczny

• W pliku ode-przyklady.ps.gz (tex) znajdują się przykłady równań różniczkowych, które można rozwiązywać w zadaniach nieobowiązkowych 7 i 8. Przedstawione zostały oscylator (model) Duffinga, oscylator van der Pola oraz wahadło matematyczne, wszystkie z siłą wymuszającą. Został przedstawiony model Lorentza jako przykład autonomicznego (niezależnego jawnie od czasu) równania dającego zachowanie chaotyczne. Zredukowane zagadnienie trzech ciał (model satelity w polu grawitacyjnym Ziemi i Księżyca) jest przykładem zagadnienia wymagającego zastosowania metod ze zmiennym krokiem.
• W pliku ode-wielokrokowe.ps.gz (tex) znaleźć można opis metod wielokrokowych typu predyktor-korektor.
• Wkrótce powinien pojawić się opis metody Runge-Kutty ze zmiennym krokiem.

8.2. Zadania domowe (nieobowiązkowe)

Dla chętnych: napisać program rozwiązujący jakieś zagadnienie początkowe metodą typu Runge-Kutty (lub dowolną inną) ze zmiennym (adaptacyjnym) krokiem (ew. dodatkowo o zmiennym rzędzie), przy zadanej przez użytkownika dokładności, kroku minimalnym i maksymalnym, minimalnym i maksymalnym zwiększeniu kroku. Zadanie to może być jednocześnie częścią zaliczenia zadania z metodą Runge-Kutty i metodą wielokrokową, jeśli dodatkowo jako jedne z zadań do wyboru będzie oscylator harmoniczny i wahadło matematyczne. Uwaga: należy zbadać wybrane zagadnienie!

Dla chętnych: napisać program rozwiązujący jakieś zagadnienie początkową metodą wielokrokową typu predyktor-korektor. Do obliczenia kroku metody niejawnej (korektora) można zastosować metodę iteracji prostej bądź inną (np. metodę Newtona, siecznych czy wstecznego różniczkowania). Zadanie to może być jednocześnie częścią zaliczenia zadania z metodą Runge-Kutty i metodą wielokrokową, jeśli dodatkowo jako jedne z zadań do wyboru będzie oscylator harmoniczny i wahadło matematyczne. Uwaga: należy zbadać wybrane zagadnienie!

9. Dyskretna zespolona transformata Fouriera i jej własności.

16-05-2002

9.1. Opis teoretyczny

• W pliku four.ps znajduje się wstęp do transformaty Fouriera i opis dyskretnej transformaty Fouriera. Opisane są szeregi Fouriera, przejście do szeregów w zmiennych zespolonych, a następnie do transformaty Fouriera. Następnie opisana została dyskretna transformata Fouriera, którą potrafimy obliczać na komputerach (dla dowolnych funkcji/danych).
• W pliku fft.ps.gz (tex) znajduje się opis algorytmu szybkiej transformaty Fouriera. Przedstawione są tam również własności transformaty Fouriera i jej dyskretnej wersji. Podane są także najprostsze zastosowania transformaty Fouriera.
• FOURIER TRANSFORM - opis z Numerical Recipies (w języku angielskim)
  • 12.0 Introduction
  • 12.1 Fourier Transform of Discretely Sampled Data
  • 12.2 Fast Fourier Transform (FFT)

9.2. Zadania

• Sprawdzić, że zastosowanie wzoru na dyskretną zespoloną transformatę Fouriera, a następnie wzoru na odwrotną transformatę Fouriera daje nam wejściową funkcję.
• Sprawdzić jak wygląda transformata Fouriera funkcji czysto rzeczywistej parzystej i czysto rzeczywistej nieparzystej.

9.3. Zadanie domowe

Sprawdzić, że odwrotna transformata Fouriera powoduje odzyskanie (z jakim błędem?) oryginalnego sygnału. Sprawdzić jak wygląda transformata Fouriera funkcji czysto rzeczywistej parzystej, czysto rzeczywistej nieparzystej, czysto urojonej parzystej, czysto urojonej nieparzystej. Zobaczyć jak wyglada transformata Fouriera funkcji stałej, okresowej, okresowej z okresem niewspółmiernym do okresu próbkowania.

Można obliczać dyskretną zespoloną transformatę Fouriera bezpośrednio z definicji. Uwaga na normalizację transformaty odwrotnej (można stosowac dowolną normalizację).

Czas oddania: do 10. czerwca. Termin już minął, za zadanie można otrzymać maksymalnie 0.8 pkt.

Planowane przyszłe zadania domowe

6. Napisać program obliczający zespoloną dyskretną transformatę Fouriera. Nie trzeba implementowac algorytmu szybkiej transformaty Fouriera. Zbadać własności transformaty Fouriera: zachowanie dla funkcji parzystych i nieparzystych, czysto rzeczywistych i czysto urojonych. Obliczyć transformatę odwrotną.
   
7. Pokazać filtrowanie za pomocą dyskretnej transformaty Fouriera. Zaimplementować filtr górnoprzepustowy, dolnoprzepustowy i pasmowy. Odfiltrować składowe o małej amplitudzie. Pokazać zastosowania.
   
8. Zaimplementować splatanie za pomocą transformaty Fouriera. Porównać wynik z zastosowaniem splatania "z definicji". Należy pamiętać o uzupełnieniu sygnału zerami (jeśli sygnał traktujemy jako zwykłą funkcję) lub "zawijaniu" indeksów przy liczeniu z definicji (jeśli traktujemy sygnał jako funkcję okresową).
   
9. Rozwiązać numerycznie równanie dyfuzji, przynajmniej dwiema metodami jawnymi. Pokazać jak ewoluują wyniki w czasie.