Równania różniczkowe cząstkowe są jednym z podstawowych narzędzi matematycznych do opisu rzeczywistego układu i zachodzących w nim zmian w czasie. Stosuje się je począwszy od fizyki, gdzie wzięły swój początek, przez chemię, biologię, aż po socjologię i ekonomię. Niestety dla większości opisywanych zjawisk tworzą one skomplikowane układy równań. Ich rozwiązanie analityczne jest możliwe tylko w bardzo niewielu przypadkach. Ogromną większość można rozwiązać tylko numerycznie na bardzo silnych komputerach, a i tak zajmuje to bardzo dużo czasu: dni, tygodnie, a choćby miesiące. Próby wykorzystania sieci neuronowych i głębokiego uczenia bezpośrednio do modelowania natury na podstawie danych natrafiły na bardzo poważne problemy, o czym możecie przeczytać w osobnym artykule – Fizyka zmienia uczenie maszynowe. Na szczęście te narzędzia uczenia maszynowego fizycy byli wstanie wykorzystać w zupełnie nowy sposób. DeepONet Dzięki uniwersalnej teorii aproksymacji z teorii sieci neuronowych udało się przybliżać nie tylko ciągłe funkcje matematyczne, ale i nieliniowe operatory różniczkowe. Zespół naukowców Lu Lu, Pengzhan Jin i George Em Karniadakis z Chińskiej Akademii Nauk oraz Uniwersytetu Brown w Providence w USA stworzył architekturę sieci neuronowych do aproksymacji operatorów różniczkowych o nazwie DeepONet – deep operator networks. Celem jest nauczenie sieci neuronowej przybliżania operatora $G$ (oznaczenie za autorami) działającego na funkcje $u(x)$. Sieć DeepONet jest zbudowana z dwóch podsieci neuronowych i przyjmuje na wejściu dwa zestawy danych. Jedna podsieć jest odpowiedzialna za przybliżanie wartości funkcji $u(x)$ ze stałej listy sensorów w punktach $\{x_1, x_2, x_3, … , x_n\}.$ Druga podsieć ma zadanie lokalizować wynik działania operatora $G(u)$. Rezultatem jest złożenie tych dwóch wyników razem. Sieć DeepONet w dwóch różnych konfiguracjach podsieci odpowiedzialnej za przybliżanie funkcji $u(x)$; źródło. Po szczegóły naukowe i techniczne odsyłam do pracy naukowej „DeepONet: Learning nonlinear operators for identifying differentia equations based on the universal approximation theorem of operators”. Równania różniczkowe w Python Jako przykład zastosowania wybrałem równanie Burgersa, które jest jednym z podstawowych równań w mechanice płynów, zajmującą się zachowaniem płynów, gazów oraz wywieranych przez nie siłach. Dzięki mechanice płynów konstruuje się napęd odrzutowy, turbiny wiatrowe czy hamulce hydrauliczne, a choćby projektuje zastawki mechaniczne serca. Równanie Burgersa zastosowanie znalazło również w modelach dynamiki ruchu ulicznego. Wykorzystuje się je także w meteorologii do przewidywania ruchu atmosfery wirującej wokół Ziemi przy prognozowaniu pogody. Równanie Burgersa służy także do modelowania ruchu ulicznego Równanie Burgersa wiąże prędkość płynu $u$ z położeniem $x$ i czasem $t$ oraz lepkością płynu $\nu$. Człon po prawej stronie równania opisuje dyssypację, czyli rozpraszanie energii. W ruchu ulicznym lepkość odwzorowuje zachowania kierowców, mających różny czas reakcji przy ruszaniu i hamowaniu, czy zmianie pasa ruchu. $$ \frac{\partial u}{\partial t} + u \frac{\partial u}{\partial x} = \nu \frac{\partial^2 u}{\partial x^2} $$ gdzie $x$ – położenie, $t$ – czas, $u$ – prędkość płynu, $\nu$ – stała lepkości. jeżeli pominiemy lepkość, to równanie upraszcza się do nielepkiego równania Burgersa: $$ \frac{\partial u}{\partial t} + u \frac{\partial u}{\partial x} = 0 $$ Dla porównania zbadamy również niejednorodne lepkie równanie Burgersa, które posiada dodatkowy człon z funkcją nieróżniczkowalną $F(x,t)$. Niejednorodność pozwala odwzorować wpływ adekwatności środowiska, w którym modelujemy zjawisko. $$\frac{\partial u}{\partial t} + u \frac{\partial u}{\partial x} = \nu \frac{\partial^2 u}{\partial x^2} + F(x,t)$$ Biblioteka DeepXDE Autorzy architektury DeepONet stworzyli bibliotekę DeepXDE dla języka Python, która w przystępny sposób pozwala na wykorzystanie tej specjalnej sieci neuronowej do rozwiązywania równań różniczkowych. Instalacja DeepXDE przez PyPi: Instalacja DeepXDE przez Conda: Do budowy sieci neuronowych wykorzystuje ona biblioteki TensorFlow (instalacja) lub PyTorch (instalacja). Ma dość dobrze opisany interfejs oraz przykłady dla różnych zagadnień wraz z danymi testowymi (obliczenia analityczne lub numeryczne równań). W poniższym przykładzie wybrałem bibliotekę PyTorch jako podstawę dla budowy sieci neuronowych przez DeepXDE. Import koniecznych bibliotek: Nielepkie równanie Burgersa Przykład zaczniemy od najprostszej nielepkiej postaci równania Burgersa, które jak już wiemy ma następująca postać: $$\frac{\partial u}{\partial t} + u \frac{\partial u}{\partial x} = 0$$ Definicja przestrzeni i czasu Potrzebujemy zdefiniować geometrię czasoprzestrzenną, w której to będziemy modelować zachowanie się układu. Dla uproszczenia będzie to model jednowymiarowy, ale dzięki temu na wykresie 3D będzie można zobaczyć ewolucję w czasie. Położenie będzie w zakresie $x \in [-1,1]$ , natomiast czas $t \in [0,1]$ : Definicja równania różniczkowego W równaniu występuje pochodna cząstkowa prędkości płynu pierwszego rzędu po czasie oraz pochodna cząstkowa pierwszego rzędu po położeniu. Reprezentacją pochodnych cząstkowych dla wspomnianej listy sensorów jest macierz Jacobiego zwana jakobianem. Choć biblioteki PyTorch i TensorFlow mają swoje definicje jakobianów, to autorzy polecają stosować wersję „lazy” z DeepXDE, która wylicza elementy macierzy tylko wtedy, gdy jest to konieczne i dodatkowo zapamiętuje wyniki. Funkcja pde definiująca równanie różniczkowe ma dwa argumenty. Pierwszym jest dwuwymiarowy wektor czasoprzestrzenny $x$, który zawiera położenie $x[:0]$ i czas $x[:1]$ (odpowiednio zmienna $j=0, j=1$), który jest daną wejściową dla sieci. Drugim jest funkcja prędkości $u(x,t)$, która jest wynikiem działania sieci. Warunki brzegowe i początkowe Warunki brzegowe określają jak funkcja ma się zachowywać na brzegu określonego wcześniej obszaru. W mechanice płynów, czyli w niniejszym przykładzie często przyjmuje się, iż cząsteczki płynu nie poruszają się. Dlatego zastosowano tu warunek brzegowy Dirichleta i dla zmiennej położenia $x = 0$. Model musi mieć też określony swój stan początkowy, czyli w chwili $t = 0$. W tym przykładzie dane dla położenia zostały wygenerowane dzięki rozkładu normalnego. Budowa i trenowanie sieci neuronowej z głębokim uczeniem Gdy mamy już określoną dziedzinę (czasoprzestrzeń), równanie różniczkowe i warunki brzegowe oraz początkowe, to można zbudować rozwiązanie zależne od czasu, w którym określamy odpowiednio liczbę punktów do trenowania sieci próbkowanych w dziedzinie, liczbę punktów do próbkowania na brzegu oraz dodatkowo dla warunków początkowych. Następnie wybieramy sieć FNN (ang. fully-connected neural network) i w definicji sieci podajemy ile ukrytych warstw neuronów (4) i jak mają być szerokie (20), a także wybieramy funkcję aktywatora $tanh$ i inicjator sieci „Glorot normal”. W ten sposób definiujemy model podając mu zdefiniowane wcześniej rozwiązanie oraz sieć neuronową, optymalizator (Adam) i współczynnik szybkości uczenia się (ang. learning rate, $lr=1e-3$). Po 15 tys. iteracji, w których to model się uczył możemy zmienić optymalizator by osiągnąć lepszą dokładność w kolejnych krokach (L-BFGS). Równanie Burgersa (lepkie) W równaniu uwzględniającym lepkość środowiska mamy dodatkowy człon z drugą pochodną prędkości po położeniu i stałym współczynnikiem lepkości odpowiadający za rozpraszanie energii. $$\frac{\partial u}{\partial t} + u \frac{\partial u}{\partial x} = \nu \frac{\partial^2 u}{\partial x^2}$$ Aby uwzględnić wpływ lepkości wystarczy dodać do definicji funkcji pde macierz Hessego zwaną hesjanem, która reprezentuje właśnie różniczkę drugiego rzędu oraz wybrać współczynnik lepkości – np. $\nu = 0.5$. W wersji niejednorodnej dodajemy człon z funkcją sinus. Poniżej pełny kod skryptu wraz z dodatkową funkcją plot_3d do wizualizacji ewolucji w czasie (oś $t$). Ponieważ jest sporo parametrów, którymi można sterować uczenie modeli, to dla celów dydaktycznych pozostawiłem pełne funkcje bez próby wydzielania części wspólnych i innego refactoringu kodu. Porównanie wyników Poniżej zestawiłem wyniki jednowymiarowych ewolucji w czasie trzech wspomnianych wcześniej wersji równania Burgera. Widać na nich falę uderzeniową przemieszczającą się i z biegiem czasu (oś $t$) spłaszczającą się. Tabela 1. Nielepkie równanie Burgersa $\frac{\partial u}{\partial t} + u \frac{\partial u}{\partial x} = 0$ Po wprowadzeniu dużego czynnika lepkości i członu dyssypacji energii (pochodna drugiego rzędu) fala uderzeniowa gwałtownie się „rozpłynęła”. Tabela 2. Lepkie równanie Burgersa $\frac{\partial u}{\partial t} + u \frac{\partial u}{\partial x} = 0.5 \frac{\partial^2 u}{\partial x^2}$ Wprowadzona niejednorodność do równania w postaci sinusoidy wpłynęła na bazowy kształt całego rozwiązania. Tabela 3. Niejednorodne lepkie równanie Burgersa $\frac{\partial u}{\partial t} + u \frac{\partial u}{\partial x} = 0.5 \frac{\partial^2 u}{\partial x^2} + 5 sin(x)$ Testowanie i poprawność modeli Zespół naukowców przetestował sieć DeepONet na czterech różnych problemach związanych z różniczkami zupełnymi i cząstkowymi, z których wynika, iż sieć ta może dawać małe błędy aproksymacji. Przy okazji zauważono interesujące efekty w trakcie uczenia, których do tej pory nie obserwowano w głębokich sieciach neuronowych. Nie wszystko jest też jasne z teoretycznego punktu widzenia i potrzebne są dalsze badania nad sieciami neuronowymi, by lepiej zrozumieć ich zachowanie. Gdyby chcieć modele fizyczne opisane równaniami różniczkowymi i aproksymowane głębokimi sieciami neuronowymi w tradycyjny sposób testować, to można by zrobić to tylko jeżeli istnieje rozwiązanie analityczne lub udało się je rozwiązać numerycznie. Wówczas można by bezpośrednio obliczyć błąd aproksymacji dla konkretnego modelu. Z drugiej strony przewagą modeli naukowych jest możliwość ich uogólniania na przypadki, dla których nie mamy wystarczających lub w ogóle dostępnych danych. Bibliografia DeepONet: Learning nonlinear operators for identifying differential equations based on the universal approximation theorem of operators DeepXDE – instrukcja i przykłady DeepXDE – GitHub Regularization by noise and stochastic Burgers equations A physics-informed variational DeepONet for predicting the crack path in brittle materials Solving Partial Differential Equations Using Deep Learning and Physical Constraints Full Discretization of Stochastic Burgers Equation with Correlated Noise Diffusion approximation for self-similarity of stochastic advection in Burgers’ equation On the Numerical Solution of One-Dimensional Nonlinear Nonhomogeneous Burgers’ Equation Adam: A Method for Stochastic Optimization