Wstęp do modelowania komputerowego w fizyce - ebook

WSTĘP

Modelowanie numeryczne to stosunkowo nowe i potężne narzędzie badawcze. Jego imponujący rozwój, stymulowany rozwojem technologii komputerowych, można obserwować od połowy XX wieku. Obecnie nie istnieje dziedzina nauki czy inżynierii, w której modelowanie numeryczne nie odgrywałoby kluczowej roli. Lista korzyści jest długa. Znajdują się na niej nowe możliwości, takie jak: przewidywanie właściwości układów, badanie właściwości niedostępnych bezpośrednio eksperymentalnie, uzyskiwanie danych ilościowych wymagających ogromnej liczby operacji matematycznych lub przetwarzania wielkich zbiorów danych (Big Data), uczenie maszynowe (Machine Learning) oraz wiele innych. Tym niezwykłym możliwościom towarzyszy łatwość stosowania, szybkie uzyskiwanie wyników oraz relatywnie niskie koszty. Dzięki wymienionym zaletom, badania obliczeniowe stały się niezależnym narzędziem badawczym (obok zawsze najważniejszego eksperymentu i idącej za nim teorii), które z jednej strony jest rozszerzeniem narzędzi teoretycznych, a z drugiej – tworzy nową metodykę, często przypominającą bardziej eksperyment niż teorię (eksperymentowanie na układach wirtualnych). Narzędzie to staje się nieocenionym wsparciem dla badań doświadczalnych, technologii oraz fizyki teoretycznej, gdzie modelowanie matematyczne często na ostatnim etapie wymaga złożonych obliczeń prowadzących do wyników ilościowych.

W świetle tych uwag staje się rzeczą oczywistą, że różne aspekty modelowania numerycznego powinny być ważnym elementem wykształcenia, szczególnie w jednostkach ukierunkowanych na nauki ścisłe i inżynieryjne. Tego rodzaju kursy pojawiają się w ofertach dydaktycznych na całym świecie, nawet w formie niezależnych specjalności na uczelniach wyższych, jak Informatyka Stosowana czy Big Data. Koncentrują się one na różnych aspektach zastosowań komputerów w nauce i technice, w szczególności na modelowaniu numerycznym. Niniejszy podręcznik jest odpowiedzią na potrzeby w tym ostatnim aspekcie, czyli modelowaniu komputerowym. Poziom został dobrany tak, aby kurs mogli realizować studenci pierwszych lat fizyki lub studiów inżynieryjnych. Jest to istotne, ponieważ ze względu na powszechność narzędzi kształcenie tego rodzaju powinno się rozpocząć już na wczesnych etapach studiów w szkole wyższej. Konstrukcja kursu oparta jest na idei nauczania problemowego (problem-based learning) tak, jak sugeruje to tytuł podręcznika: Wstęp do modelowania komputerowego w fizyce. Oznacza to, że nie jest to systematyczny kurs metod numerycznych, lecz pewien zbiór projektów-problemów, a metody numeryczne pojawiają się jako narzędzie ich rozwiązania. W ten sposób poznawanie tych metod staje się funkcją celu, jakim jest rozwiązanie konkretnego problemu, co stanowi czynnik motywujący. Oprócz rozdziału Pierwsze kroki, gdzie studenci poznają metody numeryczne dla podstawowych operacji matematycznych (znajdowanie miejsca zerowego oraz ekstremum funkcji), kurs jest podzielony na projekty: osiem podstawowych i sześć zaawansowanych. Każdy podstawowy projekt rozpoczyna się przedstawieniem podstawowych informacji z zakresu fizyki, koniecznych do realizacji zadania (Podstawy fizyczne), kolejny etap to prezentacja i dyskusja konkretnego zagadnienia wraz z opisem matematycznym (Problem), następnie omawiane są metody numeryczne i algorytmy służące do rozwiązania problemu (Metody numeryczne) oraz proponowane podstawowe wersje programów (repozytorium internetowe), projekt kończy się zestawem ćwiczeń do wykonania (Ćwiczenia). Projekty zaawansowane oparte są na wcześniej wprowadzonych podstawach fizycznych i otwierają możliwość pierwszych badań obliczeniowych, gdzie student może tworzyć swoje własne scenariusze badawcze, uruchamiając w ten sposób indywidualną kreatywność. Podręcznik zawiera materiał wystarczający do realizacji kursu i w zasadzie nie ma potrzeby sięgania do dodatkowych źródeł. Oczywiście rekomendowane jest poszerzanie wiedzy z wykorzystaniem literatury.

Metody numeryczne omawiane w tym podręczniku, oprócz wspomnianych wyżej prostych operacji matematycznych oraz pochodnych i kwadratury funkcji, koncentrują się głównie na równaniach różniczkowych zwyczajnych (ODE – Ordinary Differential Equations) i cząstkowych (PDE – Partial Differential Equations): zagadnienie wartości początkowej, zagadnienie wartości brzegowych, zagadnienie własne. Większość wprowadzanych metod to metody rekurencyjne, oparte na dyskretyzacji zmiennej niezależnej (dla ODE). Dla PDE przedstawiono idee metod różnic skończonych (FD – Finite Difference) i elementów skończonych (FD – Finite Elements), by przy okazji zaprezentować wybraną metodę macierzową (FD) rozwiązania układu równań oraz przykładową metodę iteracyjną optymalizacji (FE), przy czym w obu przypadkach wykorzystano symetrię, by zredukować zagadnienie D do jednego wymiaru, co umożliwia łatwą implementację i realizację projektu w niezbyt obszernych ramach czasowych.

Wartością dodaną kursu jest nauka fizyki, gdyż projekty obejmują wiele jej obszarów, takich jak: mechanika, elementy mechaniki kwantowej, elektrodynamika czy termodynamika. Stanowi on dobry pretekst do efektywnej powtórki oraz rozszerzania i pogłębiania wiedzy z zakresu fizyki. Każdy projekt jest poprzedzony krótkim omówieniem podstaw fizycznych. Omówienia te nie wyczerpują tematów, jednak są skonstruowane tak, żeby były spójne, zawierały informacje niezbędne do realizacji projektów oraz inspirowały do głębszych studiów.

Autorzy dziękują Annie Latosińskiej za twórczą weryfikację manuskryptu. Do tworzenia kodów oraz ich weryfikacji w językach C++ i Python przyczynili się w dużym stopniu studenci Fizyki Technicznej oraz Inżynierii Kwantowej Wydziału Podstawowych Problemów Techniki Politechniki Wrocławskiej.JAK KORZYSTAĆ Z PODRĘCZNIKA

Oczekuje się, że studenci przystępujący do kursu posiadają wiedzę z zakresu fizyki ogólnej na poziomie szkoły wyższej, podstawową znajomość obsługi komputera, wybranego programu graficznego oraz podstawowe umiejętności programowania w wybranym języku.

Przed rozpoczęciem realizacji projektów zaleca się przeczytanie rozdziału Pierwsze kroki oraz wykonanie zamieszczonych ćwiczeń. Opisane tam operacje matematyczne oraz metody numeryczne (wyznaczanie miejsca zerowego oraz mini- mum/maksimum funkcji) są stosunkowo proste i bardzo przydatne. Wykonanie ćwiczeń pozwoli na zapoznanie się ze środowiskiem pracy (komputer, system operacyjny, narzędzia programistyczne i graficzne), co umożliwi skupienie się tylko na aspektach merytorycznych podczas realizacji projektów.

Pierwsze projektów to projekty podstawowe, zawierające sekcje: Podstawy fizyczne, Problem, Metody numeryczne, Ćwiczenia. W projekcie dodatkowo znajduje się sekcja Redukcja ruchu pojedynczej planety do jednego wymiaru. Projekty zaawansowane ( – ) opierają się na wcześniej opisanych podstawach fizycznych, dlatego zawierają już tylko opis problemu, niekiedy opis metod numerycznych (jeśli wcześniej nie zostały wprowadzone) oraz zadania. Realizację projektów podstawowych rozpoczynamy poprzez zapoznanie się z podstawami fizycznymi i charakterystyką układu fizycznego wraz z opisem matematycznym (sekcje i ). W razie potrzeby (np. trudności w zrozumieniu) należy sięgnąć do dodatkowych źródeł. Sekcja trzecia poświęcona jest metodom numerycznym i jest ona kluczowa dla celu kursu. Prezentowane są tu pełne lub skrócone wyprowadzenia algorytmów, ponieważ ich świadome wykorzystanie jest zasadnicze do późniejszej analizy wyników, często obarczonych niefizycznymi cechami związanymi z cyfrowym charakterem analizy. Rekomenduje się również, aby czytelnik samodzielnie wyprowadził wszystkie formuły związane z metodami numerycznymi (jeśli nie są wyprowadzone) na podstawie opisanej metody (jest to zwykle jedno z ćwiczeń), przy czym ważniejsze niż samo wyprowadzenie jest zrozumienie zasady, na jakiej zbudowana jest metoda. Ostatnia część każdego projektu to ćwiczenia, jakie należy wykonać w programie.

Ćwiczenia podzielone są na kategorie: podstawowe (które powinny być wykonane), uzupełniające (których celem jest pogłębienie i utrwalenie wiedzy) oraz zaawansowane (zawierające trudniejsze zadania dla osób szczególnie zainteresowanych przedmiotem). W ramach ćwiczeń podstawowych (obowiązkowych) duży nacisk położony jest na podstawowe etapy pracy obliczeniowej, takie jak:

Testowanie programu. Jest to sprawdzenie dla szczególnego przypadku, dla którego wyniki znane są z innych źródeł (np. analityczne), czy program odtwarza dany wynik. Ten etap zwykle występuje w pierwszych uruchomieniach programu, po eliminacji błędów w kodzie, które nieuchronnie się pojawiają.

Testowanie wpływu parametrów kontrolnych oraz ustalenie ich prawidłowych wartości. Jest to bardzo ważny punkt, nieodłączny w pracy fizyka czy inżyniera zajmującego się obliczeniami. Tworzenie i implementacja algorytmów numerycznych zawsze wiąże się z cyfryzacją formuł analitycznych, gdzie pojawiają się pewne parametry techniczne (kontrolne) regulujące ten proces, takie jak krok czasowy czy parametr siatki zmiennej przestrzennej. Nieprawidłowe wartości tych parametrów mogą zaburzać wyniki i sprawiać, że są one nieprawdziwe (niefizyczne). Podstawową metodą badania wpływu parametrów kontrolnych jest test zbieżności, polegający na analizie wyniku w funkcji danego parametru kontrolnego. Pozwala to zwykle ustalić zakresy parametrów, przy których spodziewamy się prawidłowego wyniku. Problem bywa skomplikowany, ponieważ przy wielu parametrach kontrolnych mogą występować korelacje między nimi, jednak celem tego podręcznika jest jedynie zwrócenie uwagi na ten problem. Istnieją również inne metody badania wpływu parametrów kontrolnych, o których będzie można dowiedzieć się podczas kursu.

Eksperymenty wirtualne. Jest to najważniejszy etap – cel, do którego prowadziła cała wcześniejsza praca. Jeśli mamy pewność, że program działa poprawnie i znamy właściwe wartości parametrów kontrolnych, możemy przystąpić do analizy problemu. Odbywa się to zwykle poprzez eksperymentowanie na wirtualnym układzie fizycznym, czyli badanie zachowania układu dla różnych parametrów fizycznych (np. ruch planet dla różnych warunków początkowych).

Wersje podstawowe programów znajdują się w repozytorium internetowym (rys. 0.1).

Rysunek 0.1: Repozytorium internetowe; https://wppt.pwr.edu.pl/PhysModelCodeshttps://wppt.pwr.edu.pl/PhysModelCodes

Zaleca się, aby pierwsze – kody były napisane samodzielnie przez studentów w wybranym języku, aby sprawdzić umiejętność programowania i ewentualnie uzupełnić znajomość języka programowania. Wersje wzorcowe to kody w Pythonie, wraz z dokładnymi opisami. Dla zainteresowanych dostępne są również programy w Fortranie oraz C++. Programy nie są optymalizowane pod względem wykorzystania zasobów obliczeniowych, aby zachować przejrzystość implementacji algorytmów. Niektóre ćwiczenia wymagają modyfikacji programów, co sprawia, że stają się okazją do nauki programowania. Oprócz podstawowych projektów (pierwszych ośmiu), książka zawiera propozycje projektów zaawansowanych. Mają one na celu zademonstrowanie możliwych zastosowań poznanych algorytmów do analizy bardziej złożonych zagadnień fizycznych. Zaleca się, aby zaawansowane projekty były realizowane indywidualnie przez studentów zainteresowanych pracą badawczą o charakterze obliczeniowym. Wszystkie zadania proponowane w tym podręczniku można wykonać, korzystając z darmowych narzędzi dostępnych w Internecie (np. Spider dla Pythona, Force dla Fortrana czy Gnuplot do opracowań graficznych).PIERWSZE KROKI

Jest to pierwsza sesja w pracowni komputerowej, podczas której uczestnik poznaje podstawowe narzędzia, które będą wykorzystywane podczas nauki. Są to: elementy wybranego języka programowania, środowisko programistyczne oraz wybrany program graficzny. Opisane zostaną algorytmy dla podstawowych operacji matematycznych: wyznaczanie miejsca zerowego oraz ekstremum funkcji D, a także przećwiczone podstawowe operacje programistyczne: instrukcja pętli, tablicowanie funkcji D i D.

0.1 Znajdowanie miejsca zerowego funkcji 1d

0.1.1 Metoda bisekcji

Rysunek 0.2: Metoda bisekcji

Metoda bisekcji pozwala na znalezienie miejsca zerowego funkcji przy wykorzystaniu faktu, że funkcja, przechodząc przez miejsce zerowe, zmienia znak. Rozpoczynamy od sytuacji, w której jesteśmy pewni, że jedno i tylko jedno miejsce zerowe funkcji znajduje się w przedziale (rys. 0.2). Oznacza to, że znamy położenie miejsca zerowego z niepewnością . Strategia dalszego postępowania polega na redukcji przedziału niepewności w kolejnych krokach. Najefektywniejszym sposobem jest podzielenie przedziału z poprzedniego kroku (także wyjściowego) na pół poprzez znalezienie punktu środkowego . Otrzymujemy dwa przedziały o dwukrotnie mniejszej długości, jednak miejsce zerowe może znajdować się tylko w jednym z nich. W celu jego identyfikacji sprawdzamy warunek , czyli czy funkcja zmienia znak w przedziale . Jeśli tak, to miejsce zerowe znajduje się właśnie w tym przedziale; jeśli nie, to znajduje się ono w drugim przedziale. Operację powtarzamy, podstawiając za (lub ); po krokach niepewność redukuje się -krotnie. Całą procedurę kończymy, gdy niepewność będzie mniejsza od założonej wartości , czyli gdy spełniony będzie warunek . Można także uprzednio wyznaczyć konieczną liczbę kroków, korzystając z początkowej długości przedziału oraz faktu, że długość przedziału po krokach zmniejsza się -krotnie i wykonać tylko tyle kroków.

0.1.2 Metoda Newtona–Rhapsona oraz siecznych

Rysunek 0.3: Metoda Newtona–Raphsona

Metoda bisekcji jest bardzo skuteczna, jednak nie jest najbardziej efektywna. Popularnymi metodami alternatywnymi (zwykle bardziej efektywnymi) są metoda siecznych i metoda Newtona–Raphsona. Ta druga jest szczególnie użyteczna, jeśli wraz z funkcją znana jest jawna postać jej pochodnej oraz gdy wiadomo (lub z dużym prawdopodobieństwem przypuszcza się), że w otoczeniu miejsca zerowego funkcja ma człon liniowy w rozwinięciu potęgowym (rys. 0.3). Wówczas w punkcie będącym n-tym punktem iteracji tworzymy liniową aproksymację funkcji, co jest prostą operacją, jeśli znamy pochodną. Miejsce zerowe otrzymanej funkcji liniowej jest kolejnym przybliżeniem miejsca zerowego . Kryterium zakończenia oparte jest na zmianach położeń kolejnych punktów . Jeśli te zmiany w kilku kolejnych krokach nie przekraczają założonej wartości , wówczas procedurę kończymy.

Gdy nie znamy jawnej postaci pochodnej funkcji, możemy ją obliczyć numerycznie (np. metodą -punktową, patrz podrozdz. 2.3). Prowadzi to do algorytmu nazywanego metodą siecznych.

0.2 Znajdowanie minimum (maksimum) funkcji 1d

0.2.1 Metoda złotego podziału

Zagadnienie znajdowania minimum (maksimum) funkcji jednej zmiennej jest kluczowe, ponieważ wiele metod wyznaczania minimów w przestrzeniach wielowymiarowych (metody optymalizacji) opiera się na minimalizacjach kierunkowych. Należy zwrócić uwagę, że odwrócenie znaku funkcji powoduje, że maksimum staje się minimum, dlatego omawiane algorytmy dotyczyć będą tylko znajdowania minimum.

Aby zidentyfikować przedział dziedziny , w którym znajduje się minimum, zakładając, że jest tylko jedno, potrzebny jest dodatkowy punkt kontrolny położony wewnątrz przedziału (rys. 0.4). Możemy być pewni, że w przedziale znajduje się minimum, jeśli spełnione są warunki i .

Rysunek 0.4: Metoda złotego podziału

Konstruowanie algorytmu polega na wprowadzeniu czwartego punktu wewnątrz , aby powstały dwa -punktowe pokrywające się obszary, żeby można było dla każdego z nich sprawdzić powyższe warunki. Kluczowe jednak staje się położenie czwartego punktu, by redukcja przedziału niepewności w każdym kroku była taka sama. Cel ten osiąga się, jeśli czwarty punkt dzieli dłuższy z przedziałów i w proporcji złotego podziału, tzn. ( – dłuższy segment, – krótszy segment). Ten sam warunek musi spełniać podział startowy. W ten sposób dwa nowe pokrywające się przedziały mają taką samą długość. Przez podstawienie w powyższej proporcji otrzymujemy równanie , którego rozwiązanie jest proporcją złotego podziału, znanego w architekturze i sztuce już w starożytności. Po krokach przedział niepewności redukuje się -krotnie, czyli nieco wolniej niż w algorytmie bisekcji wyznaczania miejsca zerowego. Metoda złotego podziału jest bezpieczną i skuteczną metodą wyznaczania ekstremum, choć nie jest ona najbardziej efektywna.

0.2.2 Inne metody

Jeśli funkcja ma znaną jawną postać pochodnej, wówczas wyznaczanie minimum (lub maksimum) sprowadza się do znalezienia miejsca zerowego pochodnej z wykorzystaniem jednej z metod omówionych w poprzedniej sekcji.

Z kolei pewną analogią metody siecznych znajdowania miejsca zerowego funkcji jest metoda parabol do wyznaczania minimum. Mianowicie, korzystając z punktów , można jednoznacznie lokalnie aproksymować funkcję parabolą. Położenie minimum tej paraboli jest przybliżeniem minimum funkcji, jednak można go potraktować jako czwarty punkt kontrolny w opisanej wyżej procedurze iteracyjnej złotego podziału. Sprawdzenie warunków i pozwala na identyfikację przedziału, w którym znajduje się minimum. Ponownie aproksymujemy funkcję parabolą w nowym przedziale i procedurę powtarzamy aż do osiągnięcia założonej niepewności (w odniesieniu do długości przedziału w ostatnim kroku). Trzeba zauważyć, że metoda ta będzie skuteczna, jeśli funkcja ma człon kwadratowy w potęgowym rozwinięciu w otoczeniu minimum.

Rysunek 0.5: Metoda simplex

Metoda simplex w D (metoda kładki) jest najprostszą (jednak najmniej efektywną) metodą. Wykorzystuje ona okienko testowe „kładkę” (bez trzeciego punktu kontrolnego) o długości . Począwszy od pewnego punktu okienko testowe przesuwane jest w kierunku, w którym funkcja maleje (rys. 0.5). Na przykład przy przesuwaniu okienka w prawo sygnałem pojawienia się minimum jest spełnienie warunku . W takiej sytuacji przedział testowy jest zmniejszany dwukrotnie i przeszukiwanie zaczyna się od punktu w kierunku przeciwnym (w tym przypadku w lewo), aż do spełnienia warunku , co jest sygnałem do kolejnego zwrotu. Procedurę kończymy, gdy osiągnięta zostaje założona niepewność . Dużą zaletą metody jest to, że pozwala ona na znalezienie przedziału, gdzie występuje minimum. Po jego znalezieniu można stosować bardziej efektywne, opisane wcześniej metody.

0.3 Ćwiczenia

Obowiązkowe

1. Korzystając z programu FTABLE, wypisz do pliku tablicę wybranej przez siebie funkcji i zwizualizuj ją przy użyciu programu graficznego.
2. Zmodyfikuj program FTABLE tak, aby przedstawiał w postaci tabeli funkcję 2D. Utwórz tablicę zaproponowanej przez siebie funkcji i narysuj ją programem graficznym.
3. Przetestuj program BISEC, znajdując pierwiastki wybranej funkcji kwadratowej i porównując je z wartościami analitycznymi.
4. Korzystając z programu BISEC, wyznacz liczbę jako miejsce zerowe funkcji sinus. Jaką precyzję można osiągnąć? Wyjaśnij dlaczego.

Zaawansowane

1. Napisz program DMINIMUM, który znajdowałby minimum funkcji D z wykorzystaniem jednej z metod opisanych w tym rozdziale (złotego podziału, parabol, simplex). Przetestuj program, znajdując minimum zadanej funkcji kwadratowej. Wyznacz wartość liczby jako położenie minimum funkcji cosinus. Jaką precyzję można osiągnąć? Wyjaśnij dlaczego.1 PROJEKT: PROSTOKĄTNA SKOŃCZONA STUDNIA KWANTOWA – STACJONARNE RÓWNANIE SCHRÖDINGERA W JEDNYM WYMIARZE

W tym projekcie wykorzystuje się procedurę wyznaczania miejsc zerowych do znajdowania wartości własnych stanów związanych w prostokątnej studni potencjału (program QWELL). Traktując studnię jako najprostszy możliwy model atomu, można dopasować otrzymywane energie własne do pierwszych – energii własnych atomu wodoru, co jest przedmiotem głównego zadania w tym projekcie. Zagadnienie własne pojawia się w kilku obszarach fizyki (mechanika drgań, optyka falowa, mechanika kwantowa) i jest przedmiotem oddzielnego rozdziału (Projekt 6), a także projektu zaawansowanego (Projekt 12), gdzie rozważana jest studnia prostokątna częściowo wypełniona elektronami. Paradoksalnie, mimo że w omawianym projekcie stosowane są metody numeryczne dla najprostszych operacji matematycznych (znajdowanie miejsca zerowego), koncepcje fizyczne, na których jest on oparty, są zaawansowane. Jest to mechanika kwantowa, zwykle jeszcze nieznana na pierwszych latach studiów. Nauka podstaw mechaniki kwantowej wymaga około - -godzinnego kursu oraz znajomości zaawansowanej matematyki, dlatego tutaj ograniczymy się tylko do wprowadzenia wybranych i najprostszych jej idei, w zakresie niezbędnym do świadomego wykonania projektu.