D2 Identyfikacja modelu dynamicznego


Politechnika Warszawska
Wydział Samochodów i Maszyn Roboczych
Instytut Maszyn Roboczych Ciężkich
Laboratorium Dzwignic
Ćwiczenie D2
Identyfikacja modelu dynamicznego
żurawia naściennego.
Wersja robocza
Tylko do użytku wewnętrznego SiMR PW
Opracowanie:
Dr inż. Andrzej Buczyński
Dr inż. Artur Jankowiak
Warszawa 2011
Wszelkie prawa zastrzeżone
1
1. CEL ĆWICZENIA
Celem ćwiczenia jest wprowadzenie studentów w dziedzinę zagadnień
modelowania zachowania się urządzeń dzwignicowych pod wpływem działania
obciążeń eksploatacyjnych. W części praktycznej ćwiczenie obejmuje identyfikację
modelu żurawia naściennego poprzez pomiar sił obciążających urządzenie podczas
podnoszenia i opuszczania ładunku. W części obliczeniowej, na podstawie
uzyskanych wyników, studenci wyznaczają wartości parametrów modelu
i przeprowadzają teoretyczne oszacowania jego działania. Sformułowanie wniosków,
jakie wynikają z porównania wyników pomiarów i przewidywań teoretycznych
stanowi ostatni etap ćwiczenia.
2. WPROWADZENIE
Urządzenia dzwignicowe, w swoim szerokim spektrum, obejmują też ważną grupę żurawi
stałych, w tym naściennych. Ich zaletą jest to, że montowane są na ścianach bądz podporach
hal, przez co nie zajmują powierzchni użytkowej podłogi. Są one produktami, które
doskonale nadajÄ… siÄ™ do optymalizacji stanowisk pracy, stÄ…d znajdujÄ… szeroki zakres
zastosowań jako urządzenia udzwigowienia stanowisk pracy. Na rys. 1. przedstawiono szkic
oraz zdjęcie dostępnego w laboratorium MRC żurawia naściennego DEMAG WSK-KBK.
Rys. 1. Schemat konstrukcyjny [4] i zdjęcie żurawia naściennego WSK-KBK [5]
Żuraw ten o udzwigu, Q = 400 [kg] będzie obiektem badań wykonywanych podczas
ćwiczenia. W urządzeniu należy wyszczególnić następujące elementy: konstrukcję nośną o
masie zredukowanej, mK, sztywności, KN, tłumieniu, CN; wciągnik o masie, mW; cięgno nośne
(łańcuch) o sztywności, KS, tłumieniu, CS; zblocze o masie, mZ; podnoszony ładunek o masie,
mH. Masa konstrukcji nośnej i masa wciągnika zazwyczaj traktowane są łącznie, jako masa
układu nośnego, mN = mK+ mW, co wynika z charakteru pracy urządzenia oraz jednakowego
położenia środków masy ładunku i wciągnika. Podobnie można określić masę obciążenia:
mG = mZ+ mH, jako sumÄ™ mas zblocza i Å‚adunku.
Proponowane rozwiązanie konstrukcji nośnej (rys.2a) i łańcuchowe zawieszenie zblocza
wskazuje na możliwość wprowadzenia uproszczeń w formułowaniu odpowiedniego modelu
urządzenia. W praktyce inżynierskiej [1], celem dokonania stosownych analiz,
np. koniecznych w fazie projektowania, zazwyczaj przyjmuje się zastępczy dwumasowy
model z tłumieniem, jaki zilustrowano (przy założeniu, że wciągnik znajduje się na końcu
wysięgnika) na rys. 2b.
2
a)
b)
kN, CN
mN
kN CN
mN
kS, CS
kS
CS
mG mG
Rys. 2. Dwumasowy model żurawia naściennego: a) schemat układu konstrukcyjnego,
b) zastępczy model dwumasowy
Postępowanie to znajduje uzasadnienie, gdy sztywności cięgna i konstrukcji nośnej
przyjmują zbliżone wartości. Doświadczenie wskazuje, że wartości współczynników
sztywności i tłumienia (KN, CN) konstrukcji nośnej oraz (KS, CS) cięgien nośnych mogą
znacznie różnić się od siebie. Wtedy dominującą rolę odgrywa element mniej sztywny. Stąd,
w przewidywaniach chwilowych obciążeń bądz zachowania się żurawia podczas podnoszenia
(opuszczania) ładunku, słuszne jest wykorzystanie uproszczonego, jednomasowego modelu.
3. UPROSZCZONY MODEL DYNAMICZNY ŻURAWIA NAŚCIENNEGO
Wyniki badań doświadczalnych i analitycznych [1] ujawniają, trzy fazy obciążenia
żurawia, jakie powstają podczas podnoszenia lub opuszczania ładunku.
Analizując proces podnoszenia, należy zauważyć, że stosowane w żurawiach naściennych
mechanizmy podnoszenia umożliwiają poderwanie ładunku, praktycznie z nominalną
prędkością podnoszenia, VP. Wówczas bezwładność ładunku generuje dodatkowe,
oscylacyjne, tłumione obciążenia żurawia. Ich wygaśnięcie do nieistotnego poziomu kończy
pierwszą i jednocześnie rozpoczyna drugą fazę. Tutaj, ładunek porusza się ze stałą
prędkością (VP), zaś stan obciążenia żurawia odpowiada warunkom równowagi statycznej.
Ostatnia faza obciążenia, jakim podlegają elementy konstrukcyjne urządzenia, następuje w
momencie hamowania ruchu ładunku. Tak jak działo się to podczas pierwszej fazy, również
obecnie, zródło dodatkowych obciążeń konstrukcji tworzy powiązanie dynamicznych
właściwości urządzenia z bezwładnością ładunku. Opuszczanie ładunku wywołuje podobne
efekty obciążenia żurawia. Stosowne analogie wiążą wszystkie fazy opuszczania
i podnoszenia.
Charakter działania żurawia, zarówno podczas podnoszenia, czy też opuszczania ładunku
prowadzi do konieczności rozważenia układu, który przedstawiono na rys.3. Założono tu, że
dominujące dla walorów dynamicznych układu będzie tu działanie cięgien nośnych1.
Przyjmiemy, że sprężyna o sztywności, KS i tłumik z oporem, CS, proporcjonalnym do
prędkości, wiążą masę, mH, ze sztywną płytą, która odwzorowuje zblocze żurawia.
Początkowe położenie masy i płyty, na osiach X, Z, w swobodnym stanie sprężyny,
oznaczono punktami OX i OX.
Przesunięcie płyty do punktu o współrzędnej, Z, traktowane jako kinematyczne
wymuszenie ruchu masy, wywołuje reakcje sprężyny i tłumika:
1
W rzeczywistych urządzeniach nie musi tak być.
3
& &
S = KS (Z - X ) ; T = CS (Z - X ) . (1a, b)
Obok wymienionych reakcji na masę działa siła bezwładności i ciężkości:
&&
B = -mH X ; G = -mH g . (2a, b)
y
y = Vp · t
t
Oy
z0
z
OZ
S T
kS
CS
x
mG
G
B
OX
Rys.3. Schemat obliczeniowy uproszczonego modelu żurawia
Warunek równowagi wymienionych sił: S+T+B+G=0, określa równanie ruchu masy:
& & &
KS (Z - X ) + CS (Z - X ) - mH X& - mH g = 0 . (3)
Przyjmiemy nową oś odniesienia ruchu płyty:
Y = Z - ZO . (4a)
Stąd, odpowiednie pochodne można zapisać w postaci:
& & & &
Y = Z ; Y& = Z& . (4b, c)
Wartość ZO, odpowiada położeniu płyty, jakie narzuca warunek równowagi statycznej masy,
gdy jej ciężar równoważy napięcie sprężyny:
mH g
ZO = . (5)
KS
Równanie ruchu (3), po wprowadzeniu zależności (4), (5) i wykonaniu odpowiednich
przekształceń, przyjmie następującą postać:
& & &
KS (Y - X ) + CS (Y - X ) - mH X& = 0 . (6)
Warunek kinematycznego wymuszenia, nałożony jako przemieszczenia płyty ze stałą
prędkością, VP, prowadzi równanie (6) do postaci:
& &
KS (VPt - X ) + CS (VP - X ) - mH X& = 0 , (6a)
&& &
mH X + CS X + KS X = KSVP t + CSVP . (6b)
Odnosząc się do jednostkowej masy, równanie ruchu zapisujemy w następującej formie:
4
2 2
& &
X& + 2hX + ÉO X = ÉOVP t + 2hVP , (7)
gdzie h = CS / 2m , oznacza współczynnik tÅ‚umienia, zaÅ› ÉO = KS / mH , czÄ™stotliwość
drgań własnych układu bez tłumienia.
Znane sformułowania teorii równań różniczkowych [2] kwalifikują równanie (7) do
zbioru niejednorodnych równań liniowych drugiego rzędu o stałych współczynnikach.
Dyskusję metod ich rozwiązań podano w pracy [2]. Skrócone uwagi, w wymiarze
koniecznym do wykonania ćwiczenia, zamieszczono w dodatku do instrukcji. Równanie (7),
jak wiadomo z teorii drgań [3], opisuje zachowanie się masy pod wpływem kinematycznego
wymuszenia jej ruchu. W przypadku maÅ‚ego tÅ‚umienia, tzn. gdy wartoÅ›ci h i ÉO, speÅ‚niajÄ…
warunek: h < ÉO, ogólne rozwiÄ…zanie równania (7) zapisuje siÄ™ w formie nastÄ™pujÄ…cej sumy:
X = e-ht (C1 cosÉ t + C2 sinÉ t) + VPt , (8)
2
gdzie É = ÉO - h , oznacza czÄ™stotliwość tÅ‚umionych drgaÅ„ masy.
Dwa nieznane współczynniki, C1, C2, traktowane są jako stałe całkowania równania (7),
2
gdy jego prawa strona speÅ‚nia warunek: ÉOVP t + 2hVP = 0 . Poszukiwanie ich wartoÅ›ci
wymaga sformułowania dwóch równań, które opisują stan początkowy ruchu masy.
W tym celu, wykorzystując dyskusję wstępnej fazy, przyjmiemy, że w chwili, t=0,
przemieszczenie płyty o stałej prędkości, VP, inicjuje przemieszczenie nieruchomej
&
masy, x(0) = 0 , odnotowane od jej początkowego położenia, x(0) = 0 ,. Założone warunki
wiążą ze sobą równanie(8) i jego pochodną względem czasu:
Å„Å‚
ôÅ‚X = e-ht (C1 cosÉ t + C2 sinÉ t) + VP t,
(9)
òÅ‚
&
ôÅ‚
ółX = e-ht[(C2É - C1h) cosÉ t - (C1É + C2h)sinÉ t] + VP.
&
Stosowne przekształcenia układu (9) i wprowadzeniu, x(0) = 0 i x(0) = 0 prowadzą do
poszukiwanych wartości:
VP
C1 = 0 ; C2 = - . (10a,b)
É
& &&
Przemieszczenie, X, masy, mH, jej prędkość, X i przyspieszenie, X , przy wykorzystaniu
(10) wynikają z różniczkowania równania (8):
Å„Å‚
îÅ‚ Å‚Å‚
e-ht
ôÅ‚X = VP - sinÉ tśł ,
ïÅ‚t
É
ðÅ‚ ûÅ‚
ôÅ‚
ôÅ‚
îÅ‚ h Å‚Å‚
ôÅ‚X = VP e-ht ëÅ‚ öÅ‚
&
òÅ‚
ïÅ‚1+ ìÅ‚ sinÉ t - cosÉ t ÷łśł , (11)
É
íÅ‚ Å‚Å‚
ðÅ‚ ûÅ‚
ôÅ‚
ôÅ‚ 2
Å‚Å‚
&
ôÅ‚X& = VPe-ht îÅ‚É - h2
sinÉ t + 2hcosÉ tśł .
ïÅ‚
ôÅ‚ É
ðÅ‚ ûÅ‚
ół
W dziedzinie badań drgań tłumionych, gdy poszukiwania dotyczą wartości
współczynnika tÅ‚umienia, h i czÄ™stotliwoÅ›ci drgaÅ„ wÅ‚asnych, É, wykorzystywana jest
prędkość ich wygaszania. Posługujemy się tutaj logarytmicznym dekrementem tłumienia,
definiowanym jako logarytm naturalny stosunku wartości kolejnych amplitud:
5
X
n
´ = ln = hT , (12)
X
n+1
gdzie X = X e-h(t+nT ) , X = X e-h[t+(n+1)T ] , jeżeli XO oznacza wartość początkowej
n O n+1 O
amplitudy, zaÅ› T=2Ä„/É opisuje okres drgaÅ„ tÅ‚umionych.
Wynika stąd, że drgania masy, mH, podlegają stopniowemu wygaszaniu, z prędkością
zależnÄ… od wartoÅ›ci parametrów, h i É. Wymieniony efekt eliminuje wzglÄ™dnÄ… prÄ™dkość
& &
masy, jaką opisuje następująca różnica, Y - X = 0 . Zamyka to pierwszą fazę podnoszenia.
Następnie, masa i płyta poruszają się z jednakową prędkością podnoszenia, VP, co odpowiada
warunkowi równowagi statycznej układu.
Ostatnia faza podnoszenia wiąże się z hamowaniem ruchu masy. Tutaj, przewidywanie
jej zachowania wymaga wprowadzenia następujących założeń. Przyjmiemy, że zatrzymanie
masy jest wynikiem zatrzymania płyty, co następuje w chwili t=0, gdy nieruchoma płyta
znajduje się w położeniu, Y=0. Przyjęty stan eliminuje, więc działanie wymuszenia
zewnętrznego. W tej samej chwili, t=0, położenie masy zwiążemy ze współrzędną, X=0
i narzucimy prędkość, VP, jej przemieszczenia w kierunku płyty. Wymienione warunki
prowadzą równanie (7) do uproszczonej postaci:
2
& &
X& + 2hX + ÉO X = 0 . (13)
Jak wiadomo [2], całkę ogólną równania (13) zapisujemy jako sumę:
X = e-ht (C1 cosÉ t + C2 sinÉ t) , (14)
gdzie, tak jak poprzednio, stałe całkowania wynikają z warunków początkowych.
& &&
Dalszą dyskusję równania (14), jego pochodnych, X , X oraz wyznaczenie wartości
stałych C1 i C2 pozostawiono jako jedno z zadań części obliczeniowej ćwiczenia. Podobne
zadanie stanowi analiza odpowiednich faz opuszczania Å‚adunku.
4. WYKONANIE ĆWICZENIA
W ćwiczeniu wykonywane są zadania doświadczalne, odnoszące się do identyfikacji
i wyznaczenia wartości parametrów dynamicznego modelu żurawia naściennego, jaki
zainstalowano w laboratorium MRC. Zadania obliczeniowe obejmujÄ… teoretyczne
przewidywania zachowania siÄ™ Å‚adunku podczas podnoszenia i opuszczania.
" część praktyczna
- wyznaczyć masę ładunku i prędkość jego podnoszenia.
- przeprowadzić pomiar przebiegu siły działającej na żuraw podczas
podnoszenia i opuszczania ładunku. Odpowiednie czynności należy
przeprowadzić dla położenia ładunku w połowie i przy pełnym zasięgu
żurawia.
" część obliczeniowa
- wyznaczyć okres drgań, T, wartość współczynnika tłumienia, h, częstotliwości
drgaÅ„ tÅ‚umionych, É.
- uzyskane wyniki wprowadzić do odpowiednich równań i pokazać
& &&
przewidywane przebiegi przemieszczenia, X, prędkości, X , przyspieszenia, X
ruchu Å‚adunku w pierwszej fazie podnoszenia.
6
- wyznaczyć stałe, C1 i C2 równania, jakie opisuje ruch masy w fazie jej
hamowania. Sformułować odpowiednie równania i tak jak poprzednio, pokazać
& &&
przebiegi przemieszczenia, X, prędkości, X , przyspieszenia, X ruchu
Å‚adunku.
&&
- przewidywane przebiegi przyspieszenia, X , porównać z wynikami pomiaru.
" sprawozdanie
- w sprawozdaniu należy przedstawić obliczenia i uzyskane wartoÅ›ci, , T, h, É.
- przedstawić dyskusję przewidywanych przebiegów przemieszczenia, X,
& &&
prędkości, X , przyspieszenia, X ,
&&
- we wnioskach ocenić porównanie przewidywanych przebiegów, X , wynikami
pomiarów.
5. WYMAGANY ZAKRES WIADOMOŚCI OGÓLNYCH
- podstawowe pojęcia z teorii drgań,
- rozwiązywanie niejednorodnych równań różniczkowych liniowych drugiego rzędu
ze stałymi współczynnikami,
- podstawy technik pomiarowych.
6. LITERATURA
[1] PiÄ…tkiewicz, A., Sobolski, R.,  Dzwignice , WNT, Warszawa, 1977.
[2] Aubowicz, H., Wieprzkowicz, B.,  Matematyka. Podstawowe wiadomości teoretyczne
i ćwiczenia dla studiów inżynierskich , Oficyna Wydawnicza Politechniki
Warszawskiej, Warszawa, 1996.
[3] Osiński, Z.,  Teoria drgań , PWN, Warszawa, 1980.
[4] Materiały firmy Demag Cranes & Components.
[5] Zbiory własne.
7
DODATEK
Celem dodatkowego rozdziału instrukcji jest prezentacja podstawowych
wiadomości dotyczących równań różniczkowych liniowych drugiego rzędu. Zakres
dyskusji obejmuje zagadnienia konieczne do wykonania ćwiczenia.
Równanie różniczkowe liniowe drugiego rzędu zapisywane w postaci:
& &
X& + q1(t)X + q2 (t)X = f (t) , (d1)
określane jest jako równanie niejednorodne. Jeżeli f(t)=0, wówczas równanie
to nazywa siÄ™ jednorodnym bÄ…dz uproszczonym:
& &
X& + q1(t)X + q2 (t)X = 0 . (d2)
Rozwiązanie równania niejednorodnego (d1), w pierwszej kolejności, wymaga
wyznaczenia całek szczególnych równania jednorodnego i jego całki ogólnej.
Następnie, gdy znana jest jedna całka szczególna równania niejednorodnego,
wówczas poszukiwana całka ogólna równania niejednorodnego jest sumą całki
ogólnej równania jednorodnego i całki szczególnej równania niejednorodnego
Przyjmiemy funkcje X1(t), X2(t), jako całki szczególne równania (d2). Jeśli
wybrane funkcje są ciągłe i liniowo niezależne, wówczas całkę ogólną tego równania
wyraża się w postaci sumy:
X (t) = C1X1(t) + C2 X (t) , (d3)
j 2
gdzie C1, C2 reprezentują dowolne stałe. Warunek liniowej niezależności założonych
całek, X1(t), X2(t), jak wiadomo, determinuje wyznacznik:
X1 X
2
W = `" 0 . (d4)
& &
X1 X
2
Założymy obecnie, że funkcja, Xn(t), reprezentuje całkę szczególną równania
niejednorodnego. Można więc poszukiwaną całkę ogólną równania niejednorodnego
zapisać jako sumę całki (d3) i założonej całki szczególnej:
X (t) = X + X = C1X1(t) + C2 X (t) + X (t) . (d5)
j n 2 n
Praktyka inżynierska mieści wiele zagadnień, w których funkcje q1(t), q2(t),
przyjmują stałe wartości. Badania takich zagadnień wiążą się z zastosowaniami
liniowych równań różniczkowych drugiego rzędu ze stałymi współczynnikami.
Równanie różniczkowe liniowe drugiego rzędu ze stałymi współczynnikami
wynika z równania (d1) i jako niejednorodne przyjmuje postać:
& &
X& + a1X + a2 X = f (t) , (d6)
jeżeli funkcje q1,2 reprezentują stałe wartości:q1(t)=a1 oraz q2(t)=a2.
Równanie jednorodne otrzymujemy z równania (d6), w wyniku założenia, f(t)=0,
po jego prawej stronie:
& &
X& + a1X + a2 X = 0 . (d7)
8
Rozwiązanie równania (d6) uzyskujemy na drodze omówionego powyżej
postępowania. Rozpoczynamy od rozwiązania równania jednorodnego (d7). W tym
celu przewidujemy rozwiÄ…zanie w postaci:
X = er t , (d8)
2
& &
gdzie r, dobieramy tak, aby funkcje: X = er t , X = r er t , X& = r er t , spełniały
równanie (d7), skąd otrzymujemy równanie charakterystyczne:
2
r + a1r + a2 = 0 , (d9)
równania (d7). Funkcja (d8) jest więc całką szczególną równania (d7), jeśli r stanowi
pierwiastek równania charakterystycznego (d9).
2
Wyróżnik równania kwadratowego (d9), " = a1 - 4a2 , określa trzy przypadki na
drodze poszukiwania wartości, r, jego rozwiązań.
1. " > 0 , r przyjmuje dwie rzeczywiste wartości, r1,2 = (-a1 ą ") / 2 . Istnieją
więc, zgodnie z relacją (d4),W<>0, dwie liniowo niezależne całki szczególne,
1 2
X1 = er t , X = er t , które, na podstawie (d3), prowadzą do całki ogólnej:
2
1 2
X (t) = C1er t + C2er t . (d10)
j
2. " = 0 , r stanowi pierwiastek podwójny, r = r1 = r2 = a1/2, (W=0). Dlatego całki
1 2
szczególne mają postać: X1 = er t , X = t er t . Całkę ogólną tworzy więc
2
następująca suma:
X (t) = C1er t + C2t er t . (d11)
j
3. " < 0 , równanie charakterystyczne ma pierwiastki zespolone, r1,2 = Ä… Ä… ² i ,
gdzie Ä… = -a1 / 2 , ² = - " . W tym przypadku funkcje X1 = ea t cos ² t ,
X1 = ea t sin ² t , opisujÄ… rzeczywiste caÅ‚ki szczególne, jakie w postaci poniższej
sumy tworzą całkę ogólną:
X (t) = eÄ… t (C1 cos ² t + C2 sin ² t) . (d12)
j
Równanie różniczkowe liniowe drugiego rzędu niejednorodne ze stałymi
współczynnikami rozwiązujemy wykorzystując metodę uzmienniania stałych, bądz
metodą przewidywania całki szczególnej. Druga z wymienionych metod jest bardzo
prosta i efektywna. Wykorzystanie jej wymaga jednak, aby funkcja f(t), pośród kilku
k k -1
postaci, była opisana wielomianem: Wk (t) = bkt + bk -1t + .. + b1t + b0 . Jak można to
łatwo zauważyć, odwołując się do równania (7), narzucony warunek jest spełniony,
bowiem wielomian o stopniu, k, jest tu ograniczony do formy liniowej:
W (t) = b1t + b0 = f (t) . (d13)
Metoda przewidywania całki szczególnej, niezależnie od właściwej postaci
funkcji, f(t), dopuszcza tÄ™ samÄ… strukturÄ™ przewidywanej funkcji, Xn(t), jakÄ…
reprezentuje prawa strona niejednorodnego równania (d6).
Niech wielomian (d13) reprezentuje poszukiwaną całkę szczególną:
9
X (t) = At + B . (d14)
n
& &
Podstawienie X = At + B , X = A , X& = 0 , do równania (d6) prowadzi do zależności:
n n n
0 + a1A + a2 (At + B) = b1t + b0 . (d15)
Porównując wartości współczynników przy odpowiednich potęgach po obu stronach
równania (d15), otrzymujemy:
a2A = b1,
Å„Å‚
(d16)
òÅ‚
A + a2B = b0.
óła1
Stąd, przewidywaną całkę szczególną reprezentuje funkcja:
b1 a2b0 - a1b1
X (t) = t + . (d17)
n
2
a2 a2
Odnosząc się obecnie do trzeciego przypadku rozwiązania równania jednorodnego,
jaki odpowiada warunkom obowiązującym w równaniu (7), rozwiązanie ogólne,
zapisane zależnością (d5), ostatecznie, przyjmie postać:
b1 a2b0 + a1b1
X (t) + X (t) = eÄ… t (C1 cos ² t + C2 sin ² t) + t + . (d18)
j n
2
a2 a2
10


Wyszukiwarka

Podobne podstrony:
Identyfikacja modelu matematycznego elementu
Identyfikacja Obiektów Dynamicznych
identyfikacja OBIEKTOW DYNAMICZNYCH
D2
2 Dynamika cz1
Sporządzenie modelu rozprzestrzeniania się zanieczyszczeń
,Modelowanie i symulacja systemów, Model dynamiczny
Kinematyka i Dynamika Układów Mechatronicznych
C w6 zmienne dynamiczne wskazniki funkcji
02 JTPWyjatki Klonowanie Identycznosc Rodzaje
3 dobór zmiennych do liniowego modelu ekonometrycznego
Identyfikacja leśnych siedlisk przyrodniczych NATURA 2000 na przykładzie Nadleśnictwa Oleśnica Śląsk
Identyfikacja zwiazkow organicznych
7 Dynamika ruchu obrotowego bryly sztywnej
antropologia identyfikacja płci
PHP6 i MySQL 5 Dynamiczne strony WWW Szybki start ph6ms5

więcej podobnych podstron