Pożar lasu w Monte Carlo ;-)

Oczywiście w Monte Carlo czyli stolicy Monako lasów już nie ma. Na 1 km kwadratowym ledwo znalazło się miejsce dla "paru drzewek", a parkingi dla samochodów wydrążone są w skale. Ale ja dziś nie o tym chciałem... ;-)

Dziś będzie o modelu pożaru lasu ("forest fire"), którego klasyczna wersja będąca typowym automatem komórkowym pochodzi z końca lat 80tych XX wieku i jest jednym z przykładów dla

self-organized criticality oraz dla perkolacji.

Nasz model tym będzie się różnić od klasycznego, że zamiast synchronicznie stosować reguły, użyjemy algorytmu uaktualnienia Monte Carlo.

Powoduje to że w porównaniu z modelem klasycznym nasze drzewa muszą palić się nieco dłużej niż przez 1 krok czasu. Uznajemy że czas ten jest PROPORCJONALNY do wielkości czyli wieku drzewa.

Obok więc klasycznego parametru N (linia 6) czyli długości boku "macierzy świata" (World) mamy parametr FireTimeDiv , który określa czas "płonięcia" drzewa przez podzielenie jego wieku. Domyślnie ma on wartość 10, co możemy rozumieć tak, że drzewo stuletnie płonie 10 godzin (czyli 10 kroków M C modelu).
Ponadto mamy parametr InitT (linia 9) wskazujący jak gęsty jest las w porównaniu z lasem wypełnionym maksymalnie. Domyślnie 0.75 oznaczający że w lesie jest 75% możliwej maksymalnie liczby drzew.
Ostatni parametr to IgnitionP (linia 8) czyli prawdopodobieństwo zapalenia w danym kroku jednego z sąsiednich drzew przez drzewo, które już płonie. Wartości tego parametru nie przekładają się bezpośrednio na model klasyczny, w którym drzewo płonęło zawsze tylko jeden krok czasu.
W jaki sposób budujemy las? Ponieważ nie chcemy (na razie) modelować jego wzrostu zasiewamy go wg. jakiegoś rozkładu (linie 26-33). Możemy stworzyć monokulturę w której wszystkie drzewa mają po 100 lat (linia 29), ale możemy też użyć rozkładu płaskiego w którym tylko najstarsze drzewa mają 100 lat, średnia wynosi 50, i każdy wiek pomiędzy 0 a 100 lat jest równie prawdopodobny.
No prawie ;-) To co napisałem nie do końca zgadza się z tym co jest w programie. Co należy zrobić żeby linia 30 robiła dokładnie to?
Linia 31 tworzy nam las w którym wiek drzew ma rozkład zbliżony do normalnego, który najbardziej lubią statystycy i fizycy. Czy to lepiej odzwierciedla leśną rzeczywistość? Trochę wątpię, ale musimy zostawić to zagadnienie na później.

Procedura draw() służy nam tylko do wywołania dwóch typowych dla symulacji zadań - wizualizacji i uaktualnienia świata modelu.

Wizualizacja klasycznie już rysuje kwadraciki (choć łatwo by je było zastąpić kółkami), a może trójkątami czy nawet "choinkami") różniące się kolorem:

• Jeśli w tablicy World jest 0 to oznacza że jest ona pusta. Albo nigdy nie było tam drzewa, albo się ono już spaliło. 
• Liczba dodatnia oznacza żywe drzewo, którego intensywność zielonego koloru jest proporcjonalna do wieku. Choć nasze drzewa mają najwyżej 100 lat, zabezpieczamy się na przyszłość i nigdy drzewo nie jest bardziej zielone niż 255 (linia 55).
• Jeśli w komórce macierzy jest liczba ujemna to oznacza że drzewo właśnie płonie i jego kolor jest losową mieszanką składowej czerwonej i zielonej co daje nam głównie różne odcienie koloru żółtego.

A skąd się bierze pożar i jak przebiega? 

Cóż, pożar wywoła sam operator programu. Naciskając dowolny przycisk na klawiaturze aktywuje obsługę zdarzenia keyPressed(), a w tej obsłudze losowana jest jedna z komórek świata (linie 69 i 70) i jej wartość jest zmieniana na ujemną (linia 71), wynikającą z podzielenia wieku drzewa przez FireTimeDiv. Ponieważ jest to operacja przebiegająca na liczbach całkowitych, to drzewa młodsze niż FireTimeDiv spalałyby się w czasie 0 kroków. Odejmujemy więc jeszcze 1, czyli czas pożaru pojedynczego drzewa nie może być krótszy niż 1 krok M C.
Sam krok Monte Carlo przebiega w sposób następujący:

1. Ustalamy liczbę losowań M jako N do kwadratu.
2. Wykonujemy M nawrotów pętli (linia 78-100). Licznikiem pętli jest małe m, co jak najbardziej działa, jako że język Processing, tak jak praktycznie wszystkie języki C-podobne jest wrażliwy na wielkość liter (case sensitive) czyli zmienne M i m są różne. Ale nie polecam wam na przyszłość tej sztuczki. Pochodzi z arsenału konkursów na "obfuscated code" :-)
3. W pętli oczywiście losujemy kolejne komórki.
4. Jeśli komórka jest ujemna (test w linii 83) to 

• losujemy jednego z jej sąsiadów stosując znaną już metodę "zapinania" świata symulacji w torus za pomocą operacji modulo N (linie 86-87)
• jeśli wylosowany sąsiad jest żywym drzewem to z prawdopodobieństwem IgnitionP próbujemy go podpalić w taki sam sposób jak w przypadku pierwszego zapłonu (linie 89-94)
• Natomiast wartość komórki aktualnej powiększamy o jednostkę (linia 97) co gwarantuje nam że w końcu osiągnie 0 i przestanie się palić. Jeśli wszystko dobrze zrobimy to nigdy nie ujrzymy na konsoli znaku ? pochodzącego z testu kontrolnego (linia 98) .

Jak już zadziała zaprezentowany powyżej kod to proponuje się pobawić parametrami gęstości lasu (InitT) oraz prawdopodobieństwa zapłonu, które możemy rozumieć jako odwrotność wilgotności.
Zauważyliście już zapewne, że nasz las nie odrasta. Cóż... Las rośnie lata, a pożar trwa godziny. Pod tym względem klasyczny model "forest fire" jest BARDZO ODLEGŁY OD RZECZYWISTOŚCI.