Jakiś czas temu była dyskusja na temat sił więzów na blogu profesora Jadczyka. W książkach siły te są opisane jako siły centralne, ja twierdziłem że siły te są doosiowe. Czytelnik Bjab podał jednak nieco inny wzór
dv/dt = (dω/dt) x r + ω x (dr/dt) = ε x r + ω x v (1)
Argumentacja Bjab wydaje się poprawna ale zawsze trzeba to jeszcze sprawdzić w praktyce.
Tutaj można zobaczyć animacje na której bazuje większość testów. Bryła sztywna z trzema różnymi momentami bezwładności
Ix=1, Iy=2, Iz=3
Startowa prędkość kątowa ω(0, 1, 0.1) blisko pośredniej osi głównej momentu bezwładności Iy.
Wektory v1, a1 są dla punktu na osi X` czyli punkty których ramiona są prawie prostopadłe do startowej prędkości kątowej, wektory v2, a2 są dla punktu na osi Y` czyli punkty których ramiona są prawie równoległe do startowej prędkości kątowej.
Metoda krokowa. Na podstawie wartości startowych t liczone są przyspieszenia w odcinku czasu dt. Bryła jest obracana o wyliczone przyspieszenie i zaczyna się następny krok t+1. Dla dt = 0,001 jedno fiknięcie to około 10 tysięcy kroków.
Długości wektorów dla tej konfiguracji
Liczenie błędu
Liczyłem wektor prędkości liniowej dwoma metodami, licząc różnice wektorów położenia
v=dr/dt (2)
oraz ze wzoru
v = ω x r (3)
wektory przyspieszenia punktów z różnicy wektorów prędkości v liczone metodą (3)
a=dv/dt (4)
oraz wzoru Bjab
dv/dt = (dω/dt) x r + ω x (dr/dt) = ε x r + ω x v (1)
Następnie zbierałem dane o długości wektora, będącym różnicą tych samych wektorów liczonymi różnymi metodami. Błędy dla dt=0,001
Wektory prędkości v są punktem odniesienia. Wzory (2) i (3) nie budzą żadnych kontrowersji a błąd w wynikach pokazuje skale błędu obliczeniowego, wynikający z zastosowanej metody i algorytmu. Jak widać na wykresach błąd jest rzędu jednej tysięcznej. Co to oznacza? Że jeżeli wektor miałby metr długości to różnica między wektorami liczonymi dwoma różnymi metodami wynosi jeden minimetr. Czy to dużo? Jest to dużo więcej niż oczekiwałem, chce tylko nadmienić że sporo testów jakie przeprowadziłem podczas studiowania tego zagadnienia, udało mi się doprowadzić do sytuacji gdzie błędy były rzędu 10 -9 lub mniejsze.
Sprawdziłem więc jak wygląda błąd przy zmniejszaniu odcinka czasu dt, co powinno dawać większą dokładność wyliczeń. Test już jedynie na v2 i a2 aby uzyskać większość przejrzystość wyników.
Jak widać zmniejszenie odcinka czasu dt dziesięciokrotnie w przedziale 10-3 do 10-5 skutkuje zmniejszeniem się błędu dziesięciokrotnie. Dla odcinka czasu dt=0,00001 i dla wektora długości metra, błąd jest zaledwie setną minimetru.
Nie mam żadnych przesłanek by błąd ten wynikał z błędnych wzorów i na podstawie tych wyników mogę stwierdzić z dużą dozą prawdopodobieństwa, że wzór Bjaba (1) jest skuteczny i prawidłowo zinterpretowany. Oczywiście nie daję to jeszcze 100% pewności, to nadal wymaga sporo pracy i czasu (którego nie mam) ale nie jest to wzór nowy, tak jak i moje rozwiązania nie opierają się na nowych wzorach i udowadnianie że są prawdziwe to tak jak udowadnianie że ziemia nie jest płaska. Oczywiście można to udowadniać ale nie ma się co spodziewać wyniku innego niż pozytywny.
Jednak może trafiliśmy w jakiś przypadek szczególny? Dlatego przeprowadziłem testy dla ruchu obrotowego ciała sztywnego w różnych konfiguracjach
Jak widać największe błędy daje obrót wokoło pośredniego momentu bezwładności, często opisywany jako ruch "niestabilny" przy którym dochodzi do efektu Dzanibekowa. Pozostałe konfiguracje dają podobne błędy lub mniejsze.
Ja wiem że te wyliczenia nie przekonają nikogo, tak samo jak w XVw czy XVIw nie dało się nikogo przekonać że ziemia nie jest płaska ale dla mnie wyniki są na tyle zadowalające by uznać wzór (1) za prawdziwy, do momentu kiedy ktoś nie pokaże mi wyników wyliczeń które temu przeczą. Moje wyniki oparte są na wieloletniej pracy i ogromnej ilości wyliczeń a często są krytykowane przez osoby mające problemy ze zrozumieniem prostych wektorów.
