Translate

piątek, 18 maja 2012

29. Ciągła dyskretność

    Zastanawia mnie problem lokalności zjawisk. Oddziaływania pomiędzy obiektami fizycznymi mają zasięg raczej nieskończony, ale same obiekty zdają się mieć ograniczone rozmiary. W mikroświecie atomy też zdają się mieć ograniczone rozmiary, ale ich opis na bazie prawdopodobieństwa mówi, że istnieje jednak niezerowe prawdopodobieństwo znalezienia elektronu w wielkich odległościach od atomu. Nic nie wiadomo, jak elektron miałby się przenieść tak daleko. Nic nie wiadomo, jak miałby wrócić i czy w ogóle mógłby wrócić. Na te pytania nie ma odpowiedzi, a powinny istnieć.
    Opis mikroświata, który jest kwantowy ze swojej natury, dokonywany jest za pomocą ciągłych funkcji matematycznych. Dziwi mnie taka sytuacja, gdy zjawiska dyskretne opisuje się funkcjami ciągłymi. Lubię się jednak dziwić. 
   Ciekawi mnie, co by się stało, gdyby można było zastąpić, póki co, oczywisty postulat, mówiący że prawdopodobieństwo znalezienia elektronu w całej przestrzeni jest równe jedności, do mniej ambitnego stwierdzenia, że to samo prawdopodobieństwo jest równe jedności w skończonym obszarze, a gdy chodzi o atomy, całkiem małym obszarze porównywalnym z rozmiarami tychże atomów.
    Jak by wyglądał opis matematyczny, gdyby wolno było założyć, że to prawdopodobieństwo nie jest opisywane ciągłymi funkcjami w tym skończonym obszarze, ale musi być opisane skończonym zbiorem punktów wystąpienia elektronu w przestrzeni wokółjądrowej. 
    Gdyby wolno było dokonać takich założeń, to nowy opis powinien składać się z następujących elementów:
  1. Ciągłej, nieskończonej przestrzeni, w której "działają" ciągłe funkcje matematyczne o nieskończonej dziedzinie.
  2. Skończonych obiektów fizycznych, które są opisane skończonym zasięgiem występowania ale w ramach określonych przez nieskończoną przestrzeń i nieskończone funkcje.
    Powstaje pytanie, co by miało ograniczać zasięg występowania elektronów, skoro nie zabraniają im tego funkcje falowe oraz prawo Coulomba. Właśnie, to coś musiałoby być pewną procedurą matematyczną, posiadającą trzy cechy:
  1. Niezależność od funkcji falowych
  2.  Dyskretność
  3. Ograniczony zasięg
   Być może żądam istnienia czegoś w rodzaju parametrów ukrytych Bohma, ale nie jestem pewny, czy właśnie tego. 
    Sądzę, że wyrażone wątpliwości i postawione pytania pozwalają zarysować program badawczy w obszarze teorii atomu. Program ten powinien się składać z następujących elementów:
  1. Dołączenie do istniejącej teorii atomu procedury ograniczającej zasięg i dyskretyzującej pojawianie się elektronu w przestrzeni wokójądowej, co pociąga za sobą,
  2. Zredukowanie znaczenia ciągłych funkcji falowych do swego rodzaju kanwy (matematycznej przestrzeni) występowania obiektów fizycznych.
  3. Zastąpienia całkowań sumami wynikającymi dyskretnego opisu pojawiania się elektronu.
    Program już jest, ale pozostają kwestie narzędzi jego realizacji oraz możliwości sukcesu lub porażki. Narzędzia zdają się już być używane, choć nie w tym miejscu, a to drugie jest sprawą otwartą. Będę się cieszył, gdy pojawi się sukces, ale nie będę cierpiał zbytnio, gdy program zakończy się komunikatem :Fatal error". Ostatecznie, praca naukowa jest niemal ciągłym pasmem porażek. Istotne jest tu słowo "niemal" - w nim mieszczą się sukcesy :)
    

niedziela, 13 maja 2012

28. Rysunki orbitali sigma oraz pi

Jesteśmy już przygotowani do narysowania prawdziwego orbitalu cząsteczkowego. Dla chemika sprawa to najważniejsza, żeby zrozumieć, co oznacza utworzenie wiązania chemicznego. Tu zajmiemy się stroną matematyczno-graficzną takiego procesu. Na początek warto przygotować rysunek orbitalu wiążącego stanu podstawowego cząsteczki wodoru. Przypominam, iż wszelkie stałe, które należy wprowadzić, aby uzyskać œ ścisły obraz orbitalu cząsteczkowego, zostaną tu pominięte. Nie czas na zajmowanie się tak dokładnymi postaciami funkcji orbitalowych. Na razie najważniejsza jest sama istota, samo sedno matematycznej procedury, która opiera się na dodawaniu (bądź odejmowaniu) orbitali atomowych. W najprostszym wypadku cząsteczki wodoru, wystarczy dodać orbitale 1s atomów wodoru, aby uzyskać zadowalający obraz wiązania sigma. Kluczowym pojęciem pozwalającym tworzyć orbitale cząsteczkowe jest operacja dodawania (lub odejmowania) funkcji orbitalowych. Musimy szczegółowo odpowiedzieć sobie na pytanie, co to znaczy dodać do siebie funkcje typu exp(-r), jeœli jedna z nich znajduje się w początku układu współrzędnych, a druga jest przesunięta na osi x (najlepiej w kierunku wartości dodatnich) o jakąś, chwilowo niezdefiniowaną wielkość. Matematyczna postać takiej sumy jest prosta: exp(-r1) + exp(-r2). Promienie wodzące r1oraz r2 są związane ze sobą; ich końce wskazują ten sam punkt w przestrzeni. Ponieważ "pracujemy" w przestrzeni dwuwymiarowej, promienie wskazują ten sam punkt na płaszczyźnie xy. Odpowiedni program przedstawiam poniżej.
'----------dane
skala = 50
ustawx = 200
ustawy = 230
kontur = .25
startx = -2
koniecx = 3.5
starty = -5
koniecy = 5
d = 1.5          'odlegloœć między "jšdrami"
boksiatki = .02
a1 = 1: b1 = 1: a2 = 1: b2 = 1: c1 =1: c2 = 2
screen 12
for x = startx to koniecx step boksiatki
  for y = starty to koniecy step boksiatki
    r1 = sqr(x ^ 2 + y ^ 2)
    x2 = d - x
    'wzór na r2 wynika z uproszczen geometrycznych.
    r2 = sqr(x2 ^ 2 + y ^ 2)
    orbital1 = a1*exp(-b1*r1)
    orbital2 = a2*exp(-b2*r2)
    orbital = c1*orbital1 + c2*orbital2
      xekran = x * skala + ustawx
      yekran = y * skala + ustawy
    if orbital >= kontur then pset (xekran, yekran), 3
    if orbital <= -kontur then pset (xekran, yekran), 5
  next y
next x
circle (ustawx, ustawy), 4, 4                     'jądro1
circle (ustawx + d * skala, ustawy), 4, 6         'jadro2
Poniżej przedstawiłem rysunki orbitalu wiążącego sigma ss i antywiążącego ss, otrzymane za pomocą tego programu.


Orbital sigma ss wiążący
Orbital sigma ss antywiążący



Analogicznie możemy otrzymać orbital wiążący sigmapp i jego antywiążący odpowiednik, co można uzyskać za pomocą poniższego programu.
'program rysuje przekrój orbitalu molekularnego typu sigma, utworzonego
'przez dwa orbitale typu px.
'---------------
skala = 50
ustawx = 200
ustawy = 230
kontur = .15
startx = -2.7
koniecx = 6
starty = -5
koniecy = 5
dlugosc = 2.5
boksiatki = .02
'----------------
screen 12
for x = startx to koniecx step boksiatki
  for y = starty to koniecy step boksiatki
    promien1 = sqr(x * x + y * y)
    x2 = dlugosc - x
    promien2 = sqr(x2 * x2 + y * y)
    orbital1 = x * exp(-promien1)
    orbital2 = x2 * exp(-promien2)
    orbital = .7 * orbital1 + .7 * orbital2
      xekran = x * skala + ustawx
      yekran = y * skala + ustawy
    if orbital >= kontur then pset (xekran, yekran), 3
    if orbital <= -kontur then pset (xekran, yekran), 5
  next y
next x
circle (ustawx, ustawy), 4, 4
circle (ustawx + dlugosc * skala, ustawy), 4, 6
Poniżej dwa rysunki przedstawiające wynik działania programu dla orbitalu wiążącego i antywiążącego (skala nie jest zachowana). 
Orbital sigma pp wiążący

Orbital sigma pp antywiążący

Na koniec wersja programu, która umożliwia zobrazowanie orbitali pi, utworzonych z dwóch orbitali py.
skala = 50
ustawx = 200
ustawy = 230
kontur = .1
startx = -2
koniecx = 4
starty = -5
koniecy = 5
dlugosc = 2
boksiatki = .02
screen 12
for x = startx to koniecx step boksiatki
  for y = starty to koniecy step boksiatki
    promien1 = sqr(x * x + y * y)
    x2 = dlugosc - x
    promien2 = sqr(x2 * x2 + y * y)
    orbital1 = y * exp(-promien1)
    orbital2 = y * exp(-promien2)
    orbital = .7 * orbital1 + .7 * orbital2
      xekran = x * skala + ustawx
      yekran = y * skala + ustawy
    if orbital >= kontur then pset (xekran, yekran), 3
    if orbital <= -kontur then pset (xekran, yekran), 5
  next y
next x
circle (ustawx, ustawy), 4, 4
circle (ustawx + dlugosc * skala, ustawy), 4, 6
W zależności od znaku we wzorze na zmienną orbital otrzymujemy przekrój orbitalu wiążącego lub antywiążącego.
Pi wiążący

Pi antywiążący

 Przypominam o stałych. Ich właściwe dobranie spowoduje, że funkcje nazywane tu dla uproszczenia orbitalami staną się orbitalami w sensie niemalże ścisłym.

26. Algorytmy rysowania orbitali

    Orbital to funkcja. Orbital, to ładna nazwa dla funkcji. Podobnie ładne nazwy, to np. sinus albo tangens. Te ostatnie funkcje opisują na przykład ruch punktu w obiekcie, który drga. Fizycy umówili się, że orbital opisuje zachowanie się elektronu w atomie lub cząsteczce związku chemicznego. Słowo orbital kojarzy się z orbitą, ale zachowanie się elektronu opisywane przez orbital jest o wiele bardziej niezwykłe, niż zachowanie się komety Halleya na orbicie okołosłonecznej. Jeżeli orbital jest funkcją, to powinniśmy umieć narysować jej wykres. Kłopot polega na tym, że wykres taki zajmuje CAŁĄ dostępną przestrzeń fizyczną. Nie możemy narysować wykresu orbitalu, bo o wartości tej funkcji należą do czwartego wymiaru geometrycznego.

    Na pytanie - co widzimy w książkach, gdy patrzymy na kształty orbitali - odpowiadam - że widzimy przekroje orbitali, a nie wykresy. Różnicę między pojęciem wykresu i przekroju najlepiej wyjaśnić na przykładzie atlasu geograficznego. Jeśli kształt góry albo rowu morskiego określimy, jako wykres, to warstwice na kartce papieru będą odpowiednikiem przekrojów przy ustalonej (my to ustalamy!) wartości wysokości nad poziomem morza. Czyli w celu narysowania przekroju orbitalu należy ustalić wartość liczbową, przy której będziemy stawiać punkty przekroju.

    Podstawową przeszkodą w upowszechnieniu umiejętności precyzyjnego rysowania przekrojów orbitali za pomocą komputerów, jest przekonanie, że odpowiednie rysunki można uzyskać jedynie za pomocą programów dysponujących zaawansowanymi procedurami graficznymi i obliczeniowymi. Jest tak zapewne, gdy żądamy narysowania trójwymiarowych przekrojów. W wypadku przekrojów na płaszczyźnie dwuwymiarowej sytuacja jest o wiele prostsza. Należy dysponować podstawowymi umiejętnościami stosowania języków programowania oraz właściwie zdefiniowanym zadaniem.

    Sięgniemy do języka programowania, jakim jest Small Basic oraz QBASIC. Mam nadzieję, że narzędzie to ma jeszcze zwolenników. Prostota struktury, całkiem sensowny interfejs użytkownika i niewyszukane instrukcje powodują, że QBASIC jest... wymagający. A tak! Wymaga takiego opracowania procedur obliczeniowych, aby dały się przedstawić w możliwie prosty sposób i liczyły się szybko. To jest zadanie dla dydaktyki.

    Procedura rysowania orbitalu jest stosunkowo prosta. Najpierw wybieramy funkcję, np.
psi = a*exp(-c*r).
Wielkości a oraz c to stałe, r to promień wodzący elektronu. Gdy stałe odpowiednio dobierzemy, to znaczy, że rozpatrujemy orbital 1s atomu wodoru. Gdy dobierzemy nieodpowiednio, to później możemy się oczywiście poprawić zmieniając pierwotne liczby na właściwe. Teraz ustalamy rozsądną wartość psi, obliczamy w każdym punkcie (może nie w każdym, ale w wielu) przestrzeni wartość orbitalu. Punkt, dla którego wartość funkcji jest większa od przyjętej liczby (zmienna kontur) malujemy na ekranie monitora. Tak powstaje kontur orbitalu. Dwuwymiarowy przekrój orbitalu 1s realizuje poniższy program, napisany w języku Small Basic. Program składa się ze sporej liczby linijek, ale wystarczy chwilę przyjrzeć się jego treści, aby dostrzec, że większą część wydruku zajmują dane i instrukcje rysowania punktów na ekranie. Fragment obliczeniowy znajduje się w środku zewnętrznej pętli For...EndFor. Wydaje mi się, że jego skomplikowanie nie jest zbyt duże.

'program orbital 1s, contour plotting
'Author, Wojciech Szczepankiewicz. Silesian University of Technology
'Gliwice, Poland
'Start of BlockData
    GraphicsWindow.BackgroundColor="Black"
    e = 2.7182818 'Natural log base
    scale = 50
    trans_x = 320
    trans_y = 200
    psi_value = .123
    start_x = -2.5
    end_x = 2.5
    start_y = -5
    end_y = 5
    net_step = .02
'End of BlockData
' Main part    
For x = start_x To end_x Step net_step
  For y = start_y To end_y Step net_step
    radius = math.SquareRoot(x * x + y * y)
    orbital_1s = math.Power(e, -radius)
     If orbital_1s >= psi_value Then
        x_screen = x * scale + trans_x
        y_screen = y * scale + trans_y
        GraphicsWindow.SetPixel(x_screen, y_screen, "Blue")
     EndIf
  EndFor
EndFor
GraphicsWindow.DrawText(10,400,"End of plot")
Działanie programu można zobaczyć na stronie Small Basica. Wynikiem działania programu jest po prostu koło. W jego wnętrzu wartość orbitalu jest większa od przyjętej w zmiennej orbital_1s.
    Proszę wykonać kilka eksperymentów z programem i odkryć, jak wpływa wartość zmiennej kontur na promień koła. Łatwo jest pokazać, że program po niewielkich modyfikacjach może namalować inne orbitale. Proszę przeanalizować poniższy program.

'-------------- poczatek danych
skala = 50
ustawx = 300
ustawy = 230
kontur = .125
startx = -4
koniecx = 4
starty = -4
koniecy = 4
boksiatki = .02
'-------------- koniec danych
screen 12
for x = startx to koniecx step boksiatki
   for y = starty to koniecy step boksiatki
      promien = sqr(x * x + y * y)
      orbital = x * exp(-promien)
      'orbital = y * exp(-promien)
      xekran = x * skala + ustawx
      yekran = y * skala + ustawy
      if orbital >= kontur then pset (xekran, yekran), 3
      if orbital <= -kontur then pset (xekran, yekran), 5
   next y
next x

    Program generuje rysunek, przedstawiony poniżej. Proszę poeksperymentować z wartościami zmiennych a, c oraz kontur i zbadać, jak wpływają na postać przekroju. Orbital 2py można uzyskać przez usunięcie apostrofu sprzed odpowiednio zdefiniowanej zmiennej 'orbital'.
    Proszę sprawdzić, co się stanie, gdy we wzorze na zmienną orbital wstawimy y zamiast x.
Teraz narysujemy wybrane orbitale d. Program jest prawie taki sam, jak poprzednie. Zmiany dotyczą jedynie parametrów wyświetlania na ekranie i oczywiście wzoru definiującego orbital.
'-------------- poczatek danych
skala = 50
ustawx = 300
ustawy = 230
kontur = .125
startx = -3.7
koniecx = 3.7
starty = -4
koniecy = 4
boksiatki = .02
a = 1: c = 1
'-------------- koniec danych
screen 12
for x = startx to koniecx step boksiatki
   for y = starty to koniecy step boksiatki
      promien = sqr(x * x + y * y)
      orbital = a * x * y * exp(-c * promien)
      xekran = x * skala + ustawx
      yekran = y * skala + ustawy
      if orbital >= kontur then pset (xekran, yekran), 3
      if orbital <= -kontur then pset (xekran, yekran), 5
   next y
next x

    Poniżej kopia ekranu z rysunkiem orbitalu.

Orbital 3dz2-r2 wymaga pewnego komentarza. Taki orbital powinien być prostopadły do ekranu, ale jeśli umówimy się, że oś y będzie tymczasowo "podstawiać" oś z, to otrzymamy odpowiedni rysunek. Program jest prawie taki, jak poprzednio.
'-------------- poczatek danych
skala = 27
ustawx = 300
ustawy = 230
kontur = .125
startx = -5.4
koniecx = 6
starty = -7
koniecy = 7
boksiatki = .06
a = 1: c = 1
'-------------- koniec danych
screen 12
for x = startx to koniecx step boksiatki
  for y = starty to koniecy step boksiatki
    promien = sqr(x * x + y * y)
    orbital = a * (3 * y ^ 2 - promien ^ 2) * exp(-c * promien)
      xekran = x * skala + ustawx
      yekran = y * skala + ustawy
    if orbital >= kontur then pset (xekran, yekran), 3
    if orbital <= -kontur then pset (xekran, yekran), 5
  next y
 next x
    A wynik jego działania pokazuje przekroje dwóch "hantli" otoczonych pierścieniem.

    Wydawałoby się, że o rysowaniu przekrojów obitali wiemy już sporo. Okazuje się jednak, orbitale skonstruowane są tak, iż możliwy jest jeszcze jeden sposób rysowania ich przekrojów. Sposób ten jest niezwykle szybki nawet dla bardzo wolnych komputerów, ale tę szybkość uzyskuje się kosztem straty części informacji o orbitalu. Chodzi o tę część postaci matematycznej funkcji falowej, która jest zależna od promienia wodzącego. Wystarczy nieco uważniej przyjrzeć się rozpatrywanym powyżej funkcjom p albo d, aby zauważyć, iż zawierają czynniki zależne od współrzędnych oraz od promienia. Generalnie rzecz biorąc, orbitale można przedstawić w postaci iloczynów obu części. Jeœli pominiemy część zależną od promienia, to pozostają tylko funkcje współrzędnej x (np. orbital 2px) lub y, albo obu jednocześnie (np. orbital 3dxy). Oczywiście, dla orbitali w przestrzeni trójwymiarowej należy dołączyć fragmenty pochodzące od współrzędnej z, ale przypominam, iż my zajmujemy się tylko przekrojami leżącymi w płaszczyŸnie xy.    
    Rysowanie części kątowych orbitali jest pomysłowe. Wyobraźmy sobie prostokątny układ współrzędnych na płaszczyźnie xy. Przez początek układu niech przechodzi prosta o pod kątem teta do osi x. Na tej prostej odłożymy odcinek o długości określonej przez wartość części kątowej orbitalu. Na płaszczyźnie xy część kątowa wybranego orbitalu zależy tylko od kąta teta (w przestrzeni trójwymiarowej od jeszcze jednego kąta oznaczanego zwyczajowo fi). To jest dość pogmatwany wywód, ale stanie się bardziej przyjazny, gdy rozważymy konkretny orbital (oczywiście pozbawiony stosowanych stałych). Niech to będzie orbital 2px. Przypomnijmy jego postać matematyczną: orbital = x*EXP(-promien).
Pominięcie części zależnej od promienia prowadzi do bardzo prostego wzoru:
orbital = x
Nie pomyliło się coś autorowi? Zapewniam, że nie! Tę zwykłą zmienną x można tak zastosować, żeby odtworzyła przybliżoną postać omawianego orbitalu. Potrzebujemy w tym celu nieco trygonometrii. Współrzędna x jakiegoœ punktu prostej przechodzącej pod kątem teta do osi OX wynosi cos(teta) (podobnie jak współrzędna = sin(teta)). Po tej drobnej modyfikacji część kątowa przyjmie postać:
orbital = cos(teta)
    Tę wartość orbitalu odłożymy na anonsowanej już prostej. Trzeba uważać, żeby się nie pomylić, bo na prostej wyznaczonej przez kąt teta odkładamy odcinek wyznaczony przez kosinus tego kąta. Koniec odłożonego odcinka wyznacza punkt części kątowej, który to punkt należy narysować na ekranie monitora. Ale skoro tak, to powstaje pytanie, jakie są współrzędne (x,y) końca odłożonego odcinka. Wypadałoby poinformować komputer, w którym miejscu ekranu ma zapalić punkt świetlny. Ponownie odwołujemy się do wprowadzonej przed chwilą odrobiny trygonometrii. Wzory na współrzędne x oraz y części kątowej orbitalu 2px mają ostatecznie postać:
x = abs(orbital)*cos(teta)
y = abs(orbital)*sin(teta)
Koniecznie należy wyjaśnić, dlaczego ze zmiennej "orbital" wyciąga się wartość bezwzględną. Ta zmienna reprezentuje tu długość odcinka, a długość odcinka nie może być wartością ujemną. W tym cała tajemnica. Czas na program, dzięki któremu uzyskamy odpowiedni rysunek. Zastosujemy ponownie "QBASIC". Poniższy, niewielki program jest tak prymitywny, że poznany uprzednio program musi się wydać niezwykle zaawansowanym narzędziem w porównaniu do wydruku pokazanego poniżej, który umożliwia narysowanie częœci kątowych orbitali 2p oraz 3d.

'program rysuje części kątowe orbitali
pix2 = 2 * 3.1415
skala = 100
ustawx = 300
ustawy = 220
a = 1: b = 2
screen 12
for teta = 0 to pix2 step .001
  'orbital = a * (sin(teta))            '2p
  'orbital = a * cos(teta)              '2px
  'orbital = a * sin(teta) * cos(teta)  '3dxy
   orbital = a * cos(teta) ^ 2 - b * sin(teta) ^ 2'3dx^2-y^2
  x = abs(orbital) * cos(teta)
  y = abs(orbital) * sin(teta)
  xe = skala * x + ustawx
  ye = skala * y + ustawy
  pset (xe, ye)
next teta

    W celu uzyskania rysunku innego orbitalu wystarczy usunąć apostrof sprzed definicji zmiennej orbital. Podobnie, jak w poprzednich przykładach proszę o sprawdzenie wpływu zmiennych oraz b na kształt obitali.

25. Dystrybucja zmiennej dla chaotycznych rozwiązań modułu funkcji logarytmicznej

W poprzednim wpisie określiłem obszar rozwiązań chaotycznych modułu funkcji logarytmicznej. W niniejszym wpisie uzupełniam informację o chaotycznym zachowaniu się tej funkcji przez podanie dystrybucji zmiennej x przy określonej wartości zmiennej c, leżącej w pobliżu środka obszaru chaotycznego. W tym celu napisałem program, którego treść przedstawiam poniżej:
'Dystrybucja x dla modulu f. logarytmicznej
CLS : SCREEN 12
vec = 1200    'rozmiar wektora zbierania danych
DIM d(vec)
x = .69
llos = 10000  'liczba iteracji
c = 1.7       'xmin=4E-6 xmax=11.9
OPEN "dystryb.txt" FOR APPEND AS #1
FOR i = 1 TO llos
  x = ABS(LOG(c * x))    'funkcja logarytmiczna
  yd = INT(100 * x)      'wskaznik w wektorze
  d(yd) = d(yd) + 1
NEXT i
FOR i = 1 TO vec
   PRINT #1, i, d(i)      'drukowanie wynikow
   'PSET (i / 4, 400 - 2 * d(i))
NEXT i
CLOSE #1
Program generuje plik danych, które posłużyły do wykonania wykresu (program EXCEL) przedstawionego poniżej:
Minimalna wartość zmiennej wynosi około 4*10-6 a maksymalna około 11.9. Uzyskana krzywa jest "poszarpana". Wydaje się, że jej wygładzenie jest możliwe poprzez zastosowanie zmiennych w podwójnej precyzji.
     (27.05.2012-22.17) Obliczenia wykonane w podwójnej precyzji nie zmieniły "poszarpania" wykresu.

24. Wykładnik Lapunowa dla modułu funkcji logarytmicznej

   Poprzednio pokazałem, że kwadrat funkcji logarytmicznej zachowuje się chaotycznie dla pewnego zakresu odpowiednio zdefinowanej zmiennej c. Obecnie pokażę, że również funkcja logarytmiczna (dokładnie jej moduł) zachowuje się podobnie. Tym razem jednak zmienna c jest czynnikiem modyfikującym wartość zmiennej niezależnej x:
   Należy zwrócić uwagę, że wartość bezwzględna jest integralną częścią algorytmu Lapunowa. Nie zachodziła więc potrzeba stosowania funkcji "modułu z modułu". Poniżej przedstawiłem treść programu wykonującego obliczenie współczynników Lapunowa:
'wykladnik Lapunowa
'dla rownania abs. logarytmicznego
CLS : SCREEN 12
N = 4000                          ' liczba iteracji
OPEN "lapunow4.txt" FOR APPEND AS #1
FOR c = .1 TO 4 STEP .01          'wspolcz. w rown. logarytmicznym
  suma = 0
  x = .7
  FOR i = 1 TO N
    suma = suma + LOG(ABS(1 / x)) 'suma log pochod.
    x = ABS(LOG(c * x))           'iterowane rown. logarytmiczne
  NEXT i
  Lapunow = suma / N              'wspoczynnik Lapunowa
  PRINT #1, c, Lapunow
  PSET (c * 150, 100 - Lapunow * 100)  'obraz na ekranie
NEXT c
CLOSE #1
   Program generuje pomocniczy obraz na ekranie oraz plik zawierający zbiór punktów. Na poniższym wykresie przedstawiłem zależność wartości wykładnika Lapunowa od wartości zmiennej c:
   Moduł funkcji logarytmicznej zachowuje się chaotycznie dla zakresu c(0.36, 2.71) przy założonej pojedynczej precyzji obliczeń w języku QBASIC 4.5 oraz kroku obliczeń 0.1.

niedziela, 6 maja 2012

23. Funkcje (zbyt) porządne

    Pamiętam ten moment, gdy mój Wykładowca chemii ogólnej zadał mi pytanie egzaminacyjne o warunki, jakie powinna spełniać funkcja falowa, aby można było ją użyć do obliczeń kwantowochemicznych. Nie pamiętam dokładnie swojej odpowiedzi, ale musiała być niezbyt zadowalająca, bo ogólna ocena egzaminacyjna sięgnęła oceny dobrej. Innych pytań nie pamiętam. Po latach mogę się zastanowić nad warunkami nadającymi funkcjom falowym klasę Q. Trzymając się ściśle reguł mechaniki kwantowej, nie ma się nad czym zastanawiać. Po prostu, należy stosować te reguły i basta! Tyle, że takie podejście, co prawda zapewnia sukces rachmistrzowski, ale zatrzymuje naukę w jej postępie. Tego chyba nie chcemy.
    Chciałbym się przyjrzeć warunkowi znikania funkcji falowej w nieskończoności. Warunek wydaje się być oczywisty i poddawanie go w wątpliwość może być objawem kompletnego niezrozumienia zagadnienia (ignorancji lub szaleństwa). Tylko ignorant lub szaleniec mógłby poddawać w wątpliwość sprawy oczywiste. W tym miejscu powołam się na poglądy twórców mechaniki kwantowej (Bohr i Pauli). Teorie muszą być dostatecznie szalone, aby mogły być prawdziwe. Myślę, że poddawanie w wątpliwość spraw kardynalnych jest dostatecznie szalone... Warunek ten jest konieczny do przeprowadzenia całkowania z kwadratem, bo tylko funkcje znikające (dążące do zera) w nieskończoności są całkowalne (czyli ich całki oznaczone są zbieżne do konkretnej wartości). Konsekwencją fizyczną warunku znikania funkcji falowej w nieskończoności jest fakt, że do policzenia energii elektronu w atomie musimy zaangażować całą dostępną przestrzeń fizyczną, a nawet więcej - całą trójwymiarową przestrzeń matematyczną:
    Opis każdego pojedynczego atomu wodoru angażuje przestrzeń całego wszechświata. Co za pazerność! Żartuję tu sobie swobodnie, ale fakt ten jest zadziwiający. Wydaje się, że atomy są jednak zlokalizowane. Może ich granice są nieostre, ale mieszczą się w skończonej przestrzeni. Czy nie budzi zdumienia fakt, że istnieje niezerowe prawdopodobieństwo, że elektron z któregoś z atomów wodoru, które tworzą na przykład tkanki mojego ciała, może się znaleźć w tym momencie w galaktyce Andromedy? Prawdopodobieństwo to jest maleńkie, ale istnieje tyle atomów wodoru, że ich elektrony co jakiś czas wybierają się na natychmiastową wycieczkę do sąsiedniego miasta lub na granice wszechświata. Mechanika kwantowa takiej natychmiastowej wędrówki, tam i z powrotem, nie wyklucza ba musi ją dopuścić i niczego nam nie powie o mechanizmie takiej wędrówki. Nie podoba mi się taka sytuacja i nie mam zamiaru bezradnie rozłożyć ręce i rzec - no cóż taka jest mechanika kwantowa - absurdalna. Myślę, że należy rozróżnić tu dwie rzeczy, piękny matematyczny opis zjawiska i jego fizyczną realizację. Zjawiska kwantowe chyba są lokalne, choć ich opis może być nielokalny. Wiem, że jest sporo danych świadczących o nielokalności zjawisk kwantowych, ale czy wszystkich? Problem ten trochę przypomina opis działania dźwigni, który pozwolił Archimedesowi na sformułowanie znanego zdania: "Dajcie mi punkt podparcia, a poruszę Ziemię!" Matematyczna konstrukcja pozwala podnieść Ziemię, tyle, że jest to fizycznie niewykonalne, nawet w zasadzie.
    Bazując na powyższych rozważaniach, przypuszczam,że funkcja falowa, która obejmuje całą przestrzeń jest tylko matematyczną kanwą dla zjawisk fizycznych. Czymś w rodzaju nieskończonej sieci, której węzły wypełnia fizyka. Tyle tylko, że powinien to być skończony zbiór takich węzłów. Owszem, wolno nam zaangażować całą tę sieć do opisu zjawisk fizycznych, dzięki czemu uzyskamy wyniki zgodne z doświadczeniem, ale przeczuwam, że ten sam efekt można uzyskać stosując jedynie skończony i to niewielki zbiór punktów sieci. Gdyby taka procedura istniała, mechanika kwantowa mogłaby zrezygnować z całek na korzyść prostszych sum, a tym samym stać się opisem, który jest naprawdę kwantowy. Zdaję sobie sprawę, że powyższe wywody są dość chaotyczne, ale może w chaosie właśnie tkwi sedno sprawy.

piątek, 15 kwietnia 2011

18. Energia stanu podstawowego kationu cząsteczki wodoru z wykorzystaniem metody Monte Carlo i algorytmu Metropolisa

Przygotowany uprzednio program do obliczania energii stanu podstawowego atomu wodoru z wykorzystaniem metody Monte Carlo i algorytmu Metropolisa zaadaptowałem do obliczenia energii stanu podstawowego cząsteczki H2+. Funkcja falowa, którą wybrałem do obliczeń była kombinacją liniową orbitali 1s atomów wodoru 1:
Parametrem obliczeniowym - odległość między protonami Rj w układzie. Rysunek poniżej przedstawia sytuację geometryczną użytą do obliczeń:
W wyniku przeprowadzonych obliczeń otrzymałem zbiór danych, które są przedstawione na poniższym wykresie:
Wartość minimalna energii wynosi około -0.567 j.at. przy odległości między protonami wynoszącej 2.5 angstrema.
Tekst programu, który wygenerował zbiór wyników jest przytoczony poniżej:
DECLARE FUNCTION fufal! (x!, y!, z!)
DECLARE FUNCTION GRNF! ()  'gaussowski generator pseudolosowy
                           'jednostki atomowe
llos = 500000            'liczba losowan
h = .01                  'krok rozniczkowania
h2 = h * h               'kwadrat kroku rozniczkowania
DIM SHARED Rj            'odl. miedzy protonami
OPEN "H2plus.txt" FOR APPEND AS #1
FOR Rj = .9 TO 6 STEP .1    'obl. energi w zal. od Rj

En = 0                   'zerowanie energii do sumowania
x = 0: y = 1.2: z = 0    'poczatkowy punkt lancucha Markowa
FOR i = 1 TO llos
  xzapas = x             'zabezpieczenie starych wspolrzednych
  yzapas = y
  zzapas = z

  psi = fufal(x, y, z)    'wartosc psi w punkcie (x,y,z)
  psikw = psi * psi       'kwadrat psi

           'Losowanie nowego punktu  (xn,yn,zn)
  xn = 2 * RND - 1        'zakres (-1, 1)
  yn = 2 * RND - 1
  zn = 2 * RND - 1

  pr = SQR(xn*xn + yn*yn + zn*zn)  'promien punktu (xn,yn,zn)
  skok = .6 * GRNF         'skok przypadkowy wg rozkladu Gaussa

  x = x + xn * skok / pr  'poprawione wspolrzedne nowego punktu
  y = y + yn * skok / pr
  z = z + zn * skok / pr
 
  psinowa = fufal(x, y, z)     'wartsc psi w nowym punkcie
  psinowakw = psinowa * psinowa

'obliczanie warunku przejscia lancucha do nowego punktu

 rel = psinowakw / psikw   'stosunek nowej i starej psi
  IF rel < 1 THEN          'przejscie bezwarunkowe gdy rel >= 1
    IF RND > rel THEN      'przejscie odrzucone
      x = xzapas
      y = yzapas
      z = zzapas
      psinowa = fufal(x, y, z)
    END IF
 END IF

 'obliczanie calki z energii lokalnych
  'obliczanie pochodnej w punkcie (x,y,z)
  fufal2 = 2 * fufal(x, y, z)
  d1 = (fufal(x - h, y, z) - fufal2 + fufal(x + h, y, z)) / h2
  d2 = (fufal(x, y - h, z) - fufal2 + fufal(x, y + h, z)) / h2
  d3 = (fufal(x, y, z - h) - fufal2 + fufal(x, y, z + h)) / h2
 
  Lappsi = -(d1 + d2 + d3) / 2'druga poch. psi w punkcie (x,y,z)

  V1 = -1/(SQR((x - Rj/2)^2 + y*y + z*z))'potencjal  psi/psi=1
  V2 = -1/(SQR((x + Rj/2)^2 + y*y + z*z))
  En = En + Lappsi/psinowa+V1+V2+1/Rj    'suma energii lok.
 NEXT i

   PRINT #1, Rj, En / llos   'druk zm. wariacyjnej i en. calk.
NEXT Rj
  CLOSE #1
END

FUNCTION fufal (x, y, z)
  x1 = x - Rj / 2
  x2 = x + Rj / 2
  r1 = SQR(x1 * x1 + y * y + z * z)
  r2 = SQR(x2 * x2 + y * y + z * z)
  fufal = EXP(-r1) + EXP(-r2)
END FUNCTION

FUNCTION GRNF
  r1 = SQR(-LOG(1 - RND))
  GRNF = r1 * SIN(2 * 3.14159265# * RND)
END FUNCTION