Przeczytaj

Modelowanie ruchów Browna

Napiszmy program, w którym zamodelujemy proces poruszania się cząsteczki w płyniepłynpłynie. Celem będzie utworzenie trasy cząsteczki. Wynikiem programu powinna być zatem lista punktów w układzie współrzędnych, w których była cząsteczka w danych odstępach czasowych. Listę taką będziemy mogli wykorzystać do utworzenia szkicu trasy (np. w arkuszu kalkulacyjnym).

Na początku cząsteczka umieszczana jest w określonym punkcie układu współrzędnych. Zasymulowane jest n następujących po sobie kolizji, w wyniku których cząsteczka porusza się o ustaloną wcześniej stałą odległość w losowym kierunku. W tym celu wyznaczany jest wektor o długości r, nachylony do osi OX pod losowym kątem. W wyniku przesunięcia punktu położenia początkowego o wyznaczony wektor modelowany jest ruch.

Ponieważ nie ma możliwości tworzenia wykresów toru ruchu w języku C++, nie będziemy tworzyć tablic, które przechowują kolejne pozycje cząsteczki. Zamiast tego zostaną one wypisane na ekranie.

Modelowanie ruchów Browna w języku C++

Napiszmy program, który symuluje poruszanie się cząsteczki zgodnie z ruchami Browna według specyfikacji:

Specyfikacja problemu:

Dane:

x — początkowe położenie cząsteczki w kartezjańskim układzie współrzędnych, współrzędna x
y — początkowe położenie cząsteczki w kartezjańskim układzie współrzędnych, współrzędna y

Wynik:

Program wypisuje pary współrzędnych — trasę, po której poruszała się cząsteczka.

Na początku w programie umieszczone zostaną dwie zmienne x oraz y, które będą stanowiły położenie początkowe cząsteczki. Ważne jest, aby dobrać odpowiedni typ danych do tego problemu. Tak złożony problem nie zostanie zrealizowany na liczbach całkowitych. Odpowiednim typem dla współrzędnych będzie typ zmiennoprzecinkowy double.

// początkowe położenie cząsteczki
double x = 0.0;
double y = 0.0;

Kolejnym krokiem będzie utworzenie pętli, która wykona się określoną liczbę razy i przeprowadzi szereg kolizji z interesującą nas cząsteczką. Wewnątrz pętli trzeba wyznaczyć losowy wektor o stałej długości, zgodnie z którym przemieszczona będzie cząsteczka w wyniku kolizji.

Przypomnijmy sobie, w jaki sposób wyznaczamy wyżej wymieniony wektor.

Na początku trzeba wygenerować losowy kąt z przedziału $\left\langle0,\ 2\pi\right)$ . W języku C++ można użyć do tego funkcji rand() z biblioteki <cstdlib>. Funkcja ta wybiera losową liczbę całkowitą od 0 (włącznie) do wartości RAND_MAX (włącznie). RAND_MAX to stała, która również jest częścią biblioteki <cstdlib>.

Czy jesteśmy w stanie, korzystając z funkcji rand(), wyznaczyć losową liczbę z przedziału $\left\langle0,\ 2\pi\right)$ ? SkalowanieskalowanieSkalowanie przedziałów może odbywać się poprzez wykonywanie działań na wylosowanej liczbie. Przykładowo jeżeli dodamy do wylosowanej liczby liczbę 5, przeskalujemy przedział na $⟨$ 5, RAND_MAX+5 $⟩$ . Tak samo sytuacja ma się z dzieleniem. Jeśli podzielimy otrzymaną liczbę przez stałą RAND_MAX, otrzymamy liczbę z przedziału $\left\langle0,\ 1\right\rangle$ .

Kolejnym krokiem będzie pomnożenie otrzymanej liczby przez $2\pi$ . Stała $\pi$ dostępna jest w bibliotece <cmath> jako M_PI. Powyższe kroki można zapisać za pomocą następującego kodu:

// losowy kąt z przedziału <0, 2pi>
double fi = double(rand()) / RAND_MAX * 2 * M_PI;

Jeżeli chcemy być bardzo dokładni, możemy wyeliminować dowolny kraniec przedziału, gdyż funkcje sinus oraz cosinus dla 0 oraz $2\pi$ dadzą taki sam wynik. Można to zrobić za pomocą pętli, w której będziemy wywoływali funkcję rand() do momentu, gdy wylosowana liczba będzie różniła się od 0 albo $2\pi$ .

Teoretycznie więc wylosowanie wektora nachylonego pod kątem 0 stopni jest dwa razy bardziej prawdopodobne niż wylosowanie wektora nachylonego pod każdym innym kątem, ponieważ wektor nachylony pod kątem 0 jest równy wektorowi nachylonemu pod kątem $2\pi$ . Musimy jednak pamiętać, że zarówno generatory liczb pseudolosowych, jak i reprezentacja liczb zmiennoprzecinkowych nie są dokładne. Taka zatem też będzie symulacja. Możemy więc pozwolić sobie na zignorowanie faktu, że przedział jest domknięty obustronnie.

Po wyznaczeniu odpowiedniego kąta możemy przejść do przekształcenia go na składowe wektora. W języku C++ dostępne są funkcje sin() oraz cos(), które możemy wykorzystać, jeśli tylko załączymy do programu bibliotekę <cmath>.

// długość wektora przesunięcia
double r = 0.1;
// składowe wektora przesunięcia
double dx = r * cos(fi);
double dy = r * sin(fi);

Wektor taki możemy wykorzystać przez dodawanie jego składowych do współrzędnych cząsteczki.

Utwórzmy pętlę, która n razy powtórzy proces przemieszczania cząsteczki. Możemy użyć zarówno pętli for, jak i while. W naszym przypadku będzie to pętla for.

// liczba kolizji
int n = 100;

// wykonuj kolejne kolizje
for(int i = 0; i < n; ++i) {
    // losowy kąt z przedziału <0, 2pi>
    double fi = double(rand()) / RAND_MAX * 2 * M_PI;
    // składowe wektora przesunięcia
    double dx = r * cos(fi);
    double dy = r * sin(fi);

    // przemieść cząsteczkę
    x += dx;
    y += dy;

    // wypisz aktualną pozycję
    cout << x << ' ' << y << endl;
}

Po połączeniu wszystkich części kodu w całość, otrzymamy następujący program:

Linia 1. kratka include otwórz nawias ostrokątny iostream zamknij nawias ostrokątny prawy ukośnik prawy ukośnik funkcja cout. Linia 2. kratka include otwórz nawias ostrokątny cstdlib zamknij nawias ostrokątny prawy ukośnik prawy ukośnik funkcja rand otwórz nawias okrągły zamknij nawias okrągły przecinek stała RAND podkreślnik MAX. Linia 3. kratka include otwórz nawias ostrokątny cmath zamknij nawias ostrokątny prawy ukośnik prawy ukośnik funkcje sin otwórz nawias okrągły zamknij nawias okrągły przecinek cos otwórz nawias okrągły zamknij nawias okrągły przecinek stała M podkreślnik PI. Linia 5. using namespace std średnik. Linia 7. int main otwórz nawias okrągły zamknij nawias okrągły otwórz nawias klamrowy. Linia 8. prawy ukośnik prawy ukośnik liczba kolizji. Linia 9. double n znak równości 100 średnik. Linia 10. prawy ukośnik prawy ukośnik początkowe położenie cząsteczki. Linia 11. double x znak równości 0 kropka 0 średnik. Linia 12. double y znak równości 0 kropka 0 średnik. Linia 13. prawy ukośnik prawy ukośnik długość wektora przesunięcia. Linia 14. double r znak równości 0 kropka 1 średnik. Linia 16. prawy ukośnik prawy ukośnik wypisz początkową pozycję. Linia 17. cout otwórz nawias ostrokątny otwórz nawias ostrokątny x otwórz nawias ostrokątny otwórz nawias ostrokątny apostrof apostrof otwórz nawias ostrokątny otwórz nawias ostrokątny y otwórz nawias ostrokątny otwórz nawias ostrokątny endl średnik. Linia 19. prawy ukośnik prawy ukośnik wykonuj kolejne kolizje. Linia 20. for otwórz nawias okrągły int i znak równości 0 średnik i otwórz nawias ostrokątny n średnik plus plus i zamknij nawias okrągły otwórz nawias klamrowy. Linia 21. prawy ukośnik prawy ukośnik losowy kąt z przedziału otwórz nawias ostrokątny 0 przecinek 2pi zamknij nawias ostrokątny. Linia 22. double fi znak równości double otwórz nawias okrągły rand otwórz nawias okrągły zamknij nawias okrągły zamknij nawias okrągły prawy ukośnik RAND podkreślnik MAX asterysk 2 asterysk M podkreślnik PI średnik. Linia 23. prawy ukośnik prawy ukośnik składowe wektora przesunięcia. Linia 24. double dx znak równości r asterysk cos otwórz nawias okrągły fi zamknij nawias okrągły średnik. Linia 25. double dy znak równości r asterysk sin otwórz nawias okrągły fi zamknij nawias okrągły średnik. Linia 27. prawy ukośnik prawy ukośnik przemieść cząsteczkę. Linia 28. x plus znak równości dx średnik. Linia 29. y plus znak równości dy średnik. Linia 31. prawy ukośnik prawy ukośnik wypisz aktualną pozycję. Linia 32. cout otwórz nawias ostrokątny otwórz nawias ostrokątny x otwórz nawias ostrokątny otwórz nawias ostrokątny apostrof lewy ukośnik t apostrof otwórz nawias ostrokątny otwórz nawias ostrokątny y otwórz nawias ostrokątny otwórz nawias ostrokątny endl średnik. Linia 33. zamknij nawias klamrowy. Linia 35. return 0 średnik. Linia 36. zamknij nawias klamrowy.

#include <iostream>		// funkcja cout
#include <cstdlib>		// funkcja rand(), stała RAND_MAX
#include <cmath>		// funkcje sin(), cos(), stała M_PI

using namespace std;

int main() {
	// liczba kolizji
    double n = 100;
	// początkowe położenie cząsteczki
    double x = 0.0;
    double y = 0.0;
	// długość wektora przesunięcia
    double r = 0.1;

	// wypisz początkową pozycję
    cout << x << ' ' << y << endl;

    // wykonuj kolejne kolizje
	for(int i = 0; i < n; ++i) {
		// losowy kąt z przedziału <0, 2pi>
        double fi = double(rand()) / RAND_MAX * 2 * M_PI;
		// składowe wektora przesunięcia
        double dx = r * cos(fi);
        double dy = r * sin(fi);

		// przemieść cząsteczkę
        x += dx;
        y += dy;

        // wypisz aktualną pozycję
        cout << x << '\t' << y << endl;
	}
    
    return 0;
}

Uruchomienie programu będzie skutkowało wypisaniem na ekranie podobnych danych:

Linia 1. 0 0. Linia 2. 0 kropka 0999969 0 kropka 000786181. Linia 3. 0 kropka 00787206 minus 0 kropka 0381113. Linia 4. 0 kropka 0427464 0 kropka 0556105. Linia 5. 0 kropka 0788219 minus 0 kropka 0376555. Linia 6. minus 0 kropka 00724931 minus 0 kropka 0885647. Linia 7. minus 0 kropka 106451 minus 0 kropka 0759522. Linia 8. minus 0 kropka 165377 0 kropka 0048417. Linia 9. minus 0 kropka 0859927 minus 0 kropka 0559701. Linia 10. minus 0 kropka 041807 minus 0 kropka 145679. Linia 11. minus 0 kropka 0439401 minus 0 kropka 245656. Linia 12. 0 kropka 0019576 minus 0 kropka 156811. Linia 13. 0 kropka 0651871 minus 0 kropka 234284. Linia 14. 0 kropka 0406234 minus 0 kropka 33122. Linia 15. minus 0 kropka 0590152 minus 0 kropka 339714. Linia 16. minus 0 kropka 0922941 minus 0 kropka 245414. Linia 17. 0 kropka 00726296 minus 0 kropka 236013. Linia 18. 0 kropka 0912202 minus 0 kropka 181688. Linia 19. 0 kropka 0253476 minus 0 kropka 10645. Linia 20. 0 kropka 0854836 minus 0 kropka 0265517. Linia 21. 0 kropka 135901 0 kropka 0598085. Linia 22. 0 kropka 235641 0 kropka 0526049. Linia 23. 0 kropka 141407 0 kropka 086069. Linia 24. 0 kropka 214697 0 kropka 154103. Linia 25. 0 kropka 314654 0 kropka 157036. Linia 26. 0 kropka 414497 0 kropka 162632. Linia 27. 0 kropka 342518 0 kropka 232052. Linia 28. 0 kropka 244491 0 kropka 212288. Linia 29. 0 kropka 154327 0 kropka 169038. Linia 30. 0 kropka 074082 0 kropka 109367. Linia 31. minus 0 kropka 00409218 0 kropka 0470066. Linia 32. 0 kropka 0461429 0 kropka 133473. Linia 33. minus 0 kropka 00581461 0 kropka 0480305. Linia 34. minus 0 kropka 101072 0 kropka 0784604. Linia 35. minus 0 kropka 160924 0 kropka 158571. Linia 36. minus 0 kropka 0672781 0 kropka 193647. Linia 37. minus 0 kropka 145249 0 kropka 131033. Linia 38. minus 0 kropka 124466 0 kropka 0332161. Linia 39. minus 0 kropka 0920116 minus 0 kropka 0613708. Linia 40. minus 0 kropka 191232 minus 0 kropka 0738311. Linia 41. minus 0 kropka 223297 0 kropka 0208888. Linia 42. minus 0 kropka 152155 minus 0 kropka 0493884. Linia 43. minus 0 kropka 166758 minus 0 kropka 148316. Linia 44. minus 0 kropka 070572 minus 0 kropka 175672. Linia 45. 0 kropka 0187325 minus 0 kropka 220668. Linia 46. minus 0 kropka 078226 minus 0 kropka 245144. Linia 47. minus 0 kropka 0156225 minus 0 kropka 167164. Linia 48. minus 0 kropka 112798 minus 0 kropka 143564. Linia 49. minus 0 kropka 103592 minus 0 kropka 0439882. Linia 50. minus 0 kropka 0387719 minus 0 kropka 120135. Linia 51. minus 0 kropka 0136602 minus 0 kropka 0233394. Linia 52. 0 kropka 00486575 minus 0 kropka 121608. Linia 53. 0 kropka 0603724 minus 0 kropka 204789. Linia 54. 0 kropka 160352 minus 0 kropka 206802. Linia 55. 0 kropka 260352 minus 0 kropka 206994. Linia 56. 0 kropka 183905 minus 0 kropka 271459. Linia 57. 0 kropka 105886 minus 0 kropka 208904. Linia 58. 0 kropka 0957168 minus 0 kropka 109423. Linia 59. 0 kropka 0664444 minus 0 kropka 0138032. Linia 60. 0 kropka 120104 minus 0 kropka 0981874. Linia 61. 0 kropka 218993 minus 0 kropka 0833243. Linia 62. 0 kropka 147898 minus 0 kropka 0129994. Linia 63. 0 kropka 231436 0 kropka 0419678. Linia 64. 0 kropka 187277 minus 0 kropka 0477535. Linia 65. 0 kropka 281103 minus 0 kropka 0131623. Linia 66. 0 kropka 380951 minus 0 kropka 00764263. Linia 67. 0 kropka 468213 minus 0 kropka 0564827. Linia 68. 0 kropka 452019 0 kropka 0421974. Linia 69. 0 kropka 437682 0 kropka 141164. Linia 70. 0 kropka 352553 0 kropka 0886958. Linia 71. 0 kropka 316432 minus 0 kropka 00455297. Linia 72. 0 kropka 368742 minus 0 kropka 0897805. Linia 73. 0 kropka 354026 minus 0 kropka 188692. Linia 74. 0 kropka 254473 minus 0 kropka 179243. Linia 75. 0 kropka 282155 minus 0 kropka 0831508. Linia 76. 0 kropka 278221 minus 0 kropka 183073. Linia 77. 0 kropka 180178 minus 0 kropka 163385. Linia 78. 0 kropka 0836462 minus 0 kropka 137277. Linia 79. 0 kropka 178587 minus 0 kropka 168683. Linia 80. 0 kropka 175093 minus 0 kropka 268622. Linia 81. 0 kropka 252828 minus 0 kropka 205716. Linia 82. 0 kropka 171577 minus 0 kropka 264009. Linia 83. 0 kropka 0964678 minus 0 kropka 197989. Linia 84. 0 kropka 0870619 minus 0 kropka 297546. Linia 85. 0 kropka 00959819 minus 0 kropka 360786. Linia 86. minus 0 kropka 0802308 minus 0 kropka 404727. Linia 87. minus 0 kropka 144619 minus 0 kropka 328214. Linia 88. minus 0 kropka 0866339 minus 0 kropka 246742. Linia 89. minus 0 kropka 0710553 minus 0 kropka 147963. Linia 90. minus 0 kropka 1602 minus 0 kropka 10265. Linia 91. minus 0 kropka 127582 minus 0 kropka 197181. Linia 92. minus 0 kropka 227005 minus 0 kropka 207908. Linia 93. minus 0 kropka 127202 minus 0 kropka 214193. Linia 94. minus 0 kropka 126229 minus 0 kropka 314188. Linia 95. minus 0 kropka 182729 minus 0 kropka 231679. Linia 96. minus 0 kropka 133993 minus 0 kropka 144359. Linia 97. minus 0 kropka 188997 minus 0 kropka 227873. Linia 98. minus 0 kropka 288867 minus 0 kropka 222784. Linia 99. minus 0 kropka 196731 minus 0 kropka 183913. Linia 100. minus 0 kropka 227777 minus 0 kropka 278972. Linia 101. minus 0 kropka 327731 minus 0 kropka 281992.

0 0
0.0999969 0.000786181
0.00787206 -0.0381113
0.0427464 0.0556105
0.0788219 -0.0376555
-0.00724931 -0.0885647
-0.106451 -0.0759522
-0.165377 0.0048417
-0.0859927 -0.0559701
-0.041807 -0.145679
-0.0439401 -0.245656
0.0019576 -0.156811
0.0651871 -0.234284
0.0406234 -0.33122
-0.0590152 -0.339714
-0.0922941 -0.245414
0.00726296 -0.236013
0.0912202 -0.181688
0.0253476 -0.10645
0.0854836 -0.0265517
0.135901 0.0598085
0.235641 0.0526049
0.141407 0.086069
0.214697 0.154103
0.314654 0.157036
0.414497 0.162632
0.342518 0.232052
0.244491 0.212288
0.154327 0.169038
0.074082 0.109367
-0.00409218 0.0470066
0.0461429 0.133473
-0.00581461 0.0480305
-0.101072 0.0784604
-0.160924 0.158571
-0.0672781 0.193647
-0.145249 0.131033
-0.124466 0.0332161
-0.0920116 -0.0613708
-0.191232 -0.0738311
-0.223297 0.0208888
-0.152155 -0.0493884
-0.166758 -0.148316
-0.070572 -0.175672
0.0187325 -0.220668
-0.078226 -0.245144
-0.0156225 -0.167164
-0.112798 -0.143564
-0.103592 -0.0439882
-0.0387719 -0.120135
-0.0136602 -0.0233394
0.00486575 -0.121608
0.0603724 -0.204789
0.160352 -0.206802
0.260352 -0.206994
0.183905 -0.271459
0.105886 -0.208904
0.0957168 -0.109423
0.0664444 -0.0138032
0.120104 -0.0981874
0.218993 -0.0833243
0.147898 -0.0129994
0.231436 0.0419678
0.187277 -0.0477535
0.281103 -0.0131623
0.380951 -0.00764263
0.468213 -0.0564827
0.452019 0.0421974
0.437682 0.141164
0.352553 0.0886958
0.316432 -0.00455297
0.368742 -0.0897805
0.354026 -0.188692
0.254473 -0.179243
0.282155 -0.0831508
0.278221 -0.183073
0.180178 -0.163385
0.0836462 -0.137277
0.178587 -0.168683
0.175093 -0.268622
0.252828 -0.205716
0.171577 -0.264009
0.0964678 -0.197989
0.0870619 -0.297546
0.00959819 -0.360786
-0.0802308 -0.404727
-0.144619 -0.328214
-0.0866339 -0.246742
-0.0710553 -0.147963
-0.1602 -0.10265
-0.127582 -0.197181
-0.227005 -0.207908
-0.127202 -0.214193
-0.126229 -0.314188
-0.182729 -0.231679
-0.133993 -0.144359
-0.188997 -0.227873
-0.288867 -0.222784
-0.196731 -0.183913
-0.227777 -0.278972
-0.327731 -0.281992

Dla zainteresowanych

Aby utworzyć wykres przemieszczania się naszej cząsteczki w arkuszu kalkulacyjnym, możemy skorzystać z opcji importu danych. Po utworzeniu odpowiedniego wykresu, ścieżka przebyta przez cząsteczkę będzie prezentować się następująco:

R1aOPOUnre1Q7

Ilustracja przedstawia przedstawia oś XY o tytule: Kolejne pozycje cząsteczki. Umieszczono na niej wiele punktów połączonych przerywanymi liniami. Tworzą one nieregularny kształt. — Źródło: Contentplus.pl sp. z o.o., licencja: CC BY-SA 3.0.

Słownik

dyfuzja

samorzutne mieszanie się substancji

płyn

substancja, która charakteryzuje się łatwością w zmianie kształtu; płynami są głównie ciecze oraz gazy, ale także plazma

radian

jednostka miary łukowej kąta płaskiego, $2\pi$ radianów to tyle samo, co $360^{\circ}$

skalowanie

powiększenie (zwiększenie skali działania) systemu

Wprowadzenie

Symulacja interaktywna

Modelowanie ruchów Browna

Modelowanie ruchów Browna w języku C++

Słownik