1. PENDAHULUAN

Kebutuhan peningkatan kecepatan mekanisme untuk peningkatan produktivitas menyebabkan pengaruh inersia menjadi tidak dapat diabaikan. Salah satu cara untuk mengurangi pengaruh ini adalah dengan mebuat komponen mekanisme seringan mungkin melalui reduksi dimensi. Namun, konstruksi yang ramping tersebut menyebabkan faktor elastisitas batang tidak lagi dapat diabaikan sehingga karakteristik dinamik batang berubah terhadap posisi.

Sebagai konsekuensi kekakuan batang yang berubah terhadap posisi, perilaku dinamik batang elastik diatur oleh persamaan diferensial parsial nonlinier dengan koefisien yang berubah-waktu (time-varying). Hal ini menyebabkan solusinya tidak dapat diperoleh secara analitik. Metode perturbasi dapat digunakan untuk melinierkan persamaan gerak. Solusi persamaan yang telah dilinierkan ini dapat diperoleh dengan menggunakan algoritma solusi tunak (steady-state solution algorithm) [1].

Salah satu ciri sistem yang mempunyai persamaan gerak berkoefisien berubah-waktu adalah ketidakstabilan parametrik yang dapat terjadi pada berbagai pita frekuensi [2]. Berbagai kasus menarik telah diteliti, misalnya Farhang [1], melakukan analisis terhadap mekanisme engkol peluncur dengan batang penghubung dan poros engkol fleksibel, Hsiao dan Yang [3], mengkaji pengaruh kurva awal terhadap kestabilan pada batang elastik. Beale dan Lee [4], melakukan penelitian tentang perilaku dinamik mekanisme engkol peluncur akibat batang penghubung fleksibel dengan metode elemen hingga. Mahyuddin dan Midha [5], menganalisis bahaya instabilitas parametrik pada mekanisme Cam-Follower.

Penelitian ini menyelidiki pengaruh beberapa parameter terhadap ketidakstabilan parametrik mekanisme empat batang dengan batang penghubung elastik. Parameter yang dikaji adalah perbandingan panjang batang engkol (crank) terhadap batang penghubung (coupler), \(\mu\), perbandingan massa batang penghubung terhadap pengikut (follower), S. Di samping itu, ditinjau juga pengaruh redaman, \(\eta\), terhadap kestabilan dinamik sistem.

Hasil analisis kestabilan parametrik terhadap parameter frekuensi dan amplitudo dipetakan pada bidang 2-dimensi. Amplitudo frekuensi merupakan perbandingan kecepatan engkol terhadap frekuensi pribadi ekivalen, sedangkan amplitudo eksitasi adalah perbandingan panjang engkol dan penghubung, \(\mu\). Amplitudo respon titik tengah batang penghubung diperoleh dengan algoritma solusi tunak dan kemudian dibandingkan dengan respon yang diperoleh melalui integrasi langsung dengan metode Runge-Kutta.

2. PERSAMAAN GERAK BATANG ELASTIK

Persamaan gerak batang diturunkan berdasarkan asumsi batang adalah balok Euler-Bernoulli.

2.1 Persamaan Keseimbangan

Keseimbangan diamati untuk posisi batang elastik sebelum dan sesudah terdeformasi seperti terlihat pada Gambar 1. Posisi satu elemen (AB) sepanjang dx untuk dua posisi tersebut disajikan secara lebih rinci pada Gambar 2, di mana u dan v menyatakan defleksi dalam arah aksial, x, dan transversal, y, sedangkan v₀ menyatakan kelengkuangan awal batang. Diagram benda bebas elemen AB disajikan pada Gambar 3.

Gambar 1 Batang elastik sebelum dan sesudah terdeformasi.

Gambar 2 Deformasi elemen di batang elastik.

Gambar 3 Diagram benda bebas batang elastik.

Diagram benda bebas elemen AB, menghasilkan persamaan keseimbangan gaya dan momen berikut:

\[\frac{\partial \mathbf{P}}{\partial \mathbf{x}} = \mathbf{y} \mathbf{A} \mathbf{a}_{\mathbf{x}} \tag{2.1}\]

\[\frac{\partial Q}{\partial x} + \frac{\partial}{\partial x} \left( P \frac{\partial (v + v_0)}{\partial x} \right) = \gamma A a_y\] (2.2)

\[Q + \frac{\partial M}{\partial x} = 0 \tag{2.3}\] di mana P, Q dan M berturut-turut adalah gaya aksial, gaya lintang dan momen lentur, sedangkan γ dan A adalah massa jenis dan luas penampang batang. Komponen percepatan dalam arah aksial dan transversal sembarang titik, T, pada batang elastik dalam koordinat lokal (oxy) seperti diperlihatkan pada Gambar 4 dinyatakan sebagai a_x dan a_y. Substitusi Persamaan(2.3) ke persamaan (2.2) memberikan:

\[\frac{\partial^{2} M}{\partial x^{2}} + \frac{\partial}{\partial x} \left( P \frac{\partial (v + v_{0})}{\partial x} \right) = \gamma A a_{y}\] (2.4)

Gambar 4 Mekanisme empat batang.

2.2 Persamaan Gaya dan Defleksi

Batang elastik dapat dimodelkan sebagai gabungan sistem yang memiliki kekakuan dan redaman. Redaman yang terjadi diasumsikan viskoelastik linear. Salah satu yang sering digunakan adalah model Kelvin-Voigt [1], yang memberikan hubungan tegangan dan regangan sebagai berikut:

\[\sigma_{n} = E\varepsilon_{n} + \eta \frac{\partial \varepsilon_{n}}{\partial t}\] (2.5)

di mana \(\sigma_n\) dan \(\epsilon_n\) menyatakan tegangan dan regangan aksial batang sehingga gaya aksial yang terjadi dapat dinyatakan sebagai:

\[P = A\sigma_n = EA\varepsilon_n + \eta A \frac{\partial \varepsilon_n}{\partial t}\] (2.6)

Persamaan momen model Kelvin-Voigt [1] adalah sebagai berikut:

\[\mathbf{M} = \mathbf{E}\mathbf{I}\frac{\partial^2 \mathbf{v}}{\partial \mathbf{x}^2} + \eta \mathbf{I}\frac{\partial}{\partial \mathbf{t}} \left(\frac{\partial^2 \mathbf{v}}{\partial \mathbf{x}^2}\right) \tag{2.7}\]

Regangan yang terjadi dapat dinyatakan dalam defleksi aksial dan transversal batang sebagai berikut:

\[\varepsilon_{n} = \frac{\partial u}{\partial x} + \frac{1}{2} \left( \frac{\partial v}{\partial x} \right)^{2} + \frac{\partial v}{\partial x} \frac{\partial v_{0}}{\partial x}\] (2.8)

2.3 Percepatan

Koordinat titik sembarang T di batang elastik dalam sistem koordinat global OXY (lihat Gambar 4) adalah:

\[X_T = r\cos\phi + (x + u)\cos\beta + (v + v_0)\sin\beta\] \[Y_T = r\sin\phi - (x + u)\sin\beta + (v + v_0)\cos\beta\] (2.9) (2.10)

Diferensiasi Persamaan (2.9) dan (2.10) dua kali terhadap waktu memberikan percepatan, a_x dan a_y, yang dalam koordinat lokal (0xy) dapat dinyatakan sebagai:

\[a_{x} = -r\cos(\phi + \beta)\left(\frac{d\phi}{dt}\right)^{2} - r\sin(\phi + \beta)\frac{\partial^{2}\phi}{\partial t^{2}} + \frac{\partial^{2}u}{\partial t^{2}} - \left(x + u\right)\left(\frac{\partial\beta}{\partial t}\right)^{2} + 2\frac{\partial v}{\partial t}\frac{\partial\beta}{\partial t} + \left(v + v_{o}\right)\frac{\partial^{2}\beta}{\partial t^{2}}\] \[a_{y} = -r\sin(\phi + \beta)\left(\frac{d\phi}{dt}\right)^{2} + r\cos(\phi + \beta)\frac{\partial^{2}\phi}{\partial t^{2}} + \frac{\partial^{2}v}{\partial t^{2}} - \left(x + u\right)\left(\frac{\partial\beta}{\partial t}\right)^{2} - 2\frac{\partial u}{\partial t}\frac{\partial\beta}{\partial t} - \left(v + v_{o}\right)\frac{\partial^{2}\beta}{\partial t^{2}}\] (2.11)

2.4 Persamaan Gerak Batang Elastik

Persamaan gerak batang elastik diperoleh dengan memasukkan Persamaan (2.6) dan (2.7) ke Persamaan (2.1) dan (2.4) sehingga menghasilkan:

\[EA \frac{\partial \varepsilon_{n}}{\partial x} + \eta A \frac{\partial}{\partial t} \left( \frac{\partial \varepsilon_{n}}{\partial x} \right) = \gamma A a_{x}\] \[EI \frac{\partial^{4} v}{\partial x^{4}} + \eta I \frac{\partial}{\partial t} \left( \frac{\partial^{4} v}{\partial x^{4}} \right) - \frac{\partial}{\partial x} \left[ \left( EA \varepsilon_{n} + \eta A \frac{\partial \varepsilon_{n}}{\partial t} \right) \frac{\partial (v_{0} + v)}{\partial x} \right] = \gamma A a_{y}\] \[(2.14)\]

Untuk memudahkan evaluasi yang dapat berlaku umum melalui variasi eksitasi parametrik, dilakukan proses nondimensionalisasi. Persamaan (2.13) dan (2.14) dapat dinyatakan dalam bentuk tuna-dimensi [6] berikut:

\[\frac{\partial \overline{\varepsilon}_{n}}{\partial \rho} + \frac{\Gamma \omega_{0}}{\omega_{T}} \frac{\partial}{\partial \tau} \left( \frac{\partial \overline{\varepsilon}_{n}}{\partial \rho} \right) = \frac{\omega_{0}^{2} \lambda^{2}}{\omega_{T}^{2}} \overline{a}_{x}\] (2.15)

\[\lambda^2 \frac{\partial^4 \overline{v}}{\partial \rho^4} + \frac{\Gamma \omega_0}{\omega_T} \lambda^2 \frac{\partial}{\partial \tau} \left( \frac{\partial^4 \overline{v}}{\partial \rho^4} \right) - \frac{\partial}{\partial \rho} \left( \overline{\epsilon}_n + \frac{\Gamma \omega_0}{\omega_T} \frac{\partial \overline{\epsilon}_n}{\partial \tau} \right) \frac{\partial \left( \overline{v}_0 + \overline{v} \right)}{\partial \rho} =\]

\[-\frac{\omega_0^2 \lambda^2}{\omega_T^2} \overline{a}_y\]

(2.16)

Parameter tuna-dimensi yang digunakan adalah:

\[\begin{split} &\overline{u} = \frac{u}{L} \ ; \ \overline{v} = \frac{v}{L} \quad ; \ \mu = \frac{L}{r} \\ &\tau = \omega_0 t \ ; \ \lambda^2 = \left(\frac{\omega_T}{\omega_L}\right)^2 = \frac{I}{AL^2} \ ; \ S \equiv \frac{mL}{M} \ ; \Gamma = \frac{\eta \omega_T}{E} \\ &\omega_T^2 = \frac{EI}{mL^4} \ ; \ \omega_L^2 = \frac{EA}{mL^2} \ ; \ m = \gamma A \quad ; \ \lambda = \frac{\omega_t}{\omega_T} \end{split}\] di mana \(\omega_T\), \(\omega_L\), m dan M berturut-turut menyatakan frekuensi pribadi transversal dan longitudinal, massa batang penghubung dan massa peluncur. Selanjutnya dengan menggunakan metode peturbasi di mana u, v dan \(v_0\) dinyatakan dalam deret \(\lambda\) dan menggunakan

parameter tuna dimensi \(\rho = \frac{x}{L}\), maka Persamaan (2.15) dan (2.16) dapat dituliskan sebagai:

\[\begin{split} \frac{\partial \hat{\epsilon}_n}{\partial \rho} + \frac{\Gamma \omega_o}{\omega_\tau} \frac{\partial}{\partial \tau} \left( \frac{\partial \hat{\epsilon}_n}{\partial \rho} \right) &= \frac{\omega_o^2 \lambda^2}{\omega_t^2} \Bigg[ -\frac{1}{\mu} \cos(\phi + \beta) \dot{\phi} - \frac{1}{\mu} \sin(\phi + \beta) \ddot{\phi} - \rho \dot{\beta} \Bigg] \\ \frac{\partial^2 \hat{v}}{\partial \tau^2} + \frac{\Gamma \omega_\tau}{\omega_o} \frac{\partial}{\partial \tau} \left( \frac{\partial^4 \hat{v}}{\partial \rho^4} \right) + \frac{\omega_\tau^2}{\omega_o^2} \frac{\partial^4 \hat{v}}{\partial \rho^4} + \left( g l(\tau) + \rho \dot{\beta}^2 \right) \frac{\partial \hat{v}}{\partial \rho} + \\ \left( g l(\tau) \rho + \frac{1}{2} \rho^2 \dot{\beta}^2 - \frac{\omega_\tau^2 h(\tau)}{\omega_o^2} \right) \frac{\partial^2 \hat{v}}{\partial \rho^2} - \dot{\beta}^2 \hat{v} + \left( g l(\tau) + \rho \dot{\beta}^2 \right) \frac{\partial \hat{v}_o}{\partial \rho} + \\ \left( g l(\tau) \rho + \frac{1}{2} \rho^2 \dot{\beta}^2 - \frac{\omega_\tau^2 h(\tau)}{\omega_o^2} \right) \frac{\partial^2 \hat{v}_o}{\partial \rho^2} - \dot{\beta}^2 \hat{v}_o - g 2(\tau) - \rho \ddot{\beta} = 0. \end{split}\]

Integrasi Persamaan (2.17) menghasilkan pesamaan gerak batang elastik berikut:

\[\begin{split} & \widehat{\epsilon}_n + \frac{\Gamma \omega_o}{\omega_\tau} \frac{\partial \widehat{\epsilon}_n}{\partial \tau} = \\ & \frac{\omega_o^2 \lambda^2}{\omega_t^2} \Bigg[ - \Bigg( \frac{1}{\mu} \cos(\phi + \beta) \dot{\phi}^2 + \frac{1}{\mu} \sin(\phi + \beta) \ddot{\phi} \Bigg) \rho - \frac{1}{2} \rho^2 \dot{\beta}^2 \Bigg] + \lambda^2 h(\tau) \end{split}\]

3 MEKANISME EMPAT BATANG

Analisis perilaku elastik pada mekanisme empat batang diperoleh dengan cara yang sama dengan mekanisme engkol peluncur seperti pada penelitian yang dilakukan Farhang [1]. Perbedaannya hanya pada kondisi batas batang elastik.

3.1 Kondisi Batas

Untuk sambungan engsel ideal tuna gesekan pada kedua ujung batang elastik dan engkol kaku sempurna diperoleh kondisi batas untuk batang penghubung mekanisme empat batang sebagai berikut:

Pada ujung x = 0

\[u(0,t) = 0 (3.1)\]

\[v(0,t) = 0 (3.2)\]

\[\frac{\partial^2 \mathbf{v}(0,\mathbf{t})}{\partial \mathbf{x}^2} = 0 \tag{3.3}\]

Pada ujung x = L, batang elastik bergerak mengikuti lintasan batang pengikut yang kaku sehingga diperoleh kondisi batas berikut:

\[u(L,t) = v(L,t)\tan(\theta_4 - \beta). \tag{3.4}\] di mana θ₄ adalah sudut batang pengikut, dan

\[\frac{\partial^2 v(L,t)}{\partial x^2} = 0 \tag{3.5}\]

Persamaan (3.1) – (3.5) memberikan lima kondisi batas sehingga dibutuhkan satu kondisi batas lagi untuk menyelesaikan persamaan gerak batang elastik. Untuk itu, dipilih kondisi keseimbangan momen di engsel batang pengikut seperti diperlihatkan pada Gambar 5.

Gambar 5 Diagram benda bebas batang pengikut.

\[\sum M_{04} = 0\] \[-F_{04} \cos \alpha \left(\frac{e}{\cos \alpha} + R_{g4}\right) - \hat{F}R_B + P \sin(\theta_L - (\beta + \theta_4))R_B + Q \sin\left(\frac{\pi}{2} - \theta_L + (\beta + \theta_4)\right)R_B = 0\] (3.6)

di mana \(F_{04} = ma_{g4}\) adalah gaya inersia batang pengikut, \(\alpha = (\theta_4 + \delta) - \frac{\pi}{2}\), \(R_{g4}\) adalah jarak engsel \(\theta_4\) ke pusat massa batang 4 (G₄), \(\theta_4\) adalah panjang batang pengikut, dan e menyatakan jarak \(\theta_4\) ke G₄ serta \(\theta_4\) adalah gaya batang pengikut. Persamaan (3.6) dapat disederhanakan

\[-F_{04}\cos\alpha\left(\frac{e}{\cos\alpha} + R_{g_4}\right) - \hat{F}R_B +\] \[P\{\sin\theta_L\cos(\beta + \theta_4) - \cos\theta_L\sin(\beta + \theta_4)\}R_B +\] \[Q\{\cos\theta_L\cos(\beta + \theta_4) - \sin\theta_L\sin(\beta + \theta_4)\}R_B = 0\] (3.7)

Persamaan (3.7) dapat disederhanakan lebih jauh dengan menyatakan

\[F = \frac{F_{04}}{m_3} \cos \alpha \left( \frac{e}{\cos \alpha} + R_{g_4} \right) - \hat{F}R_B\] (3.8)

Substitusi Persamaan (2.3), (2.6) dan (2.7) ke Persamaan (3.7) memberikan:

\[F + \left(EA\varepsilon_{n} + \eta A \frac{\partial \varepsilon_{n}}{\partial t}\right) \left(\frac{-\partial(v + v_{0})}{\partial x} \cos(\beta + \theta_{4}) - \sin(\beta + \theta_{4})\right) R_{B} + \left(EI \frac{\partial^{3} v}{\partial x^{3}} + \eta I \frac{\partial}{\partial t} \left(\frac{\partial^{3} v}{\partial x}\right)\right) \left(\frac{\partial(v + v_{0})}{\partial x} \sin(\beta + \theta_{4}) + \cos(\beta + \theta_{4})\right) R_{B} = 0\] (3.9)

Persamaan (3.9) dapat disederhanakan menjadi:

\[F + \left(EA\varepsilon_{n} + \eta A \frac{\partial \varepsilon_{n}}{\partial t}\right) \left(\frac{-\partial(v + v_{0})}{\partial x}\cos(\beta + \theta_{4}) - \sin(\beta + \theta_{4})\right) R_{B} + \left(EI \frac{\partial^{3} v}{\partial x^{3}} + \eta I \frac{\partial}{\partial t}\left(\frac{\partial^{3} v}{\partial x}\right)\right) \left(\cos(\beta + \theta_{4})\right) R_{B} = 0\] (3.10)

Dengan membandingkan Persamaan (2.18) dengan Persamaan (3.10), diperoleh besar konstanta integrasi \(h(\tau)\) berikut [6]:

\[h(\tau) = \left(\frac{\omega_0}{\omega_T}\right)^2 \left(gl(\tau) + \frac{1}{2}\dot{\beta}^2\right) - \frac{F}{AE} \frac{1}{\sin(\beta + \theta_4)} \frac{1}{R_B}\] (3.11)

dengan \[gl(\tau) = \frac{1}{\mu} \left[ \cos(\phi + \beta)\dot{\phi}^2 + \sin(\phi + \beta)\ddot{\phi} \right]\] (3.11a)

3.2 Persamaan Gerak

Persamaan deformasi transversal batang elastik (2.18) dapat diselesaikan dengan menggunakan deret modal. Komponen yang merupakan fungsi \(\tau\) dan \(\rho\) dapat dipisahkan dan selanjutnya penerapan kondisi batas pada Persamaan (3.1) – (3.6) memberikan:

\[\widehat{\mathbf{v}}(\rho,\tau) = \sum_{i=1}^{\infty} T_i \sin(i\pi\rho)\] (3.12)

di mana T_i menyatakan amplitudo modus ke-i batang elastis dengan fungsi kelengkungan awal

\[\hat{\mathbf{v}}_0(\rho) = \sum_{i=1}^{\infty} C_i \sin(i\pi\rho)\] (3.13)

Substitusi Persamaan (3.12) dan (3.13) ke Persamaan(2.18) dan penyelesaiannya menghasilkan:

\[\begin{split} &\sum_{i=1}^{\infty}\ddot{T}_{i}(\tau)\sin(i\pi\rho) + \sum_{i=1}^{\infty}(i\pi)^{4}\dot{T}_{i}(\tau)\sin(i\pi\rho) + \\ &\frac{\omega_{T}^{2}}{\omega_{0}^{2}}\sum_{i=1}^{\infty}(i\pi)^{4}T_{i}(\tau)\sin(i\pi\rho) + \left(gl(\tau) + \rho\dot{\beta}^{2}\right)\sum_{i=1}^{\infty}(i\pi)T_{i}(\tau)\cos(i\pi\rho) - \\ &\left(gl(\tau) + \frac{1}{2}\rho\dot{\beta}^{2} - \frac{\omega_{T}^{2}h(\tau)}{\omega_{0}^{2}}\right)\sum_{i=1}^{\infty}(i\pi)^{2}T_{i}(\tau)\sin(i\pi\rho) - \\ &\dot{\beta}^{2}\sum_{i=1}^{\infty}T_{i}(\tau)\sin(i\pi\rho) + \left(gl(\tau) + \rho\dot{\beta}^{2}\right)\sum_{i=1}^{\infty}(i\pi)C_{i}(\tau)\cos(i\pi\rho) - \\ &\left(gl(\tau) + \frac{1}{2}\rho\dot{\beta}^{2} - \frac{\omega_{T}^{2}h(\tau)}{\omega_{0}^{2}}\right)\sum_{i=1}^{\infty}(i\pi)^{2}C_{i}(\tau)\sin(i\pi\rho) - \\ &\otimes \dot{\beta}^{2}\sum_{i=1}^{\infty}C_{i}(\tau)\sin(i\pi\rho) - g2(\tau) - \rho\ddot{\beta} = 0 \end{split}\]

\[dengan g2(\tau) = \frac{1}{11} \left[ \sin(\phi + \beta)\dot{\phi}^2 + \cos(\phi + \beta)\ddot{\phi} \right]\] (3.14a)

Perkalian Persamaan (3.14) dengan \(\sin(j\pi\rho)\) untuk proses ortogonalisasi serta integrasi hasilnya dalam interval \(0 \le \rho \le 1\) memberikan:

\[\begin{split} \ddot{T}_{j}(\tau) + \frac{\Gamma \omega_{T}}{\omega_{0}} (j\pi)^{4} \dot{T}_{j}(\tau) + & \left[ \frac{\omega_{T}^{2}}{\omega_{0}^{2}} (j\pi)^{4} - \frac{1}{2} (j\pi)^{2} g l(\tau) - \left( \frac{5}{4} + \frac{(j\pi)^{2}}{6} \right) \dot{\beta} + \\ \frac{\omega_{T}^{2} h(\tau)}{\omega_{0}^{2}} (j\pi)^{2} \right] T_{j}(\tau) + & \sum_{i=1}^{\infty} D_{ij} T_{j}(\tau) = \frac{C_{j}}{2} (j\pi)^{2} g l(\tau) + C_{j} \left( \frac{5}{4} + \frac{(j\pi)^{2}}{6} \right) \dot{\beta} - \\ \frac{\omega_{T}^{2} h(\tau)}{\omega_{0}^{2}} (j\pi)^{2} C_{j} + \frac{4}{j\pi} O_{j} g 2(\tau) + \frac{2}{j\pi} (l - O_{j}) - \sum_{i=1}^{\infty} D_{ij} C_{i} \end{split}\] \[(3.15)\] dengan

\[\begin{split} D_{ij} &= \begin{cases} \left[ -4j \left( gl(\tau) - \frac{\omega_T^2}{\omega_0^2} h(\tau) \right) + 2i\dot{\beta}^2 (2O_{i+j} - 1) \right] \frac{i}{i^2 - j^2} + \\ \frac{4i^3j}{(i^2 - j^2)} \left( 2gl(\tau)O_{i+j} + \dot{\beta}^2 (2O_{i+j} - 1); i \neq j \\ 0; & i = j \end{cases} \end{split}\] \[(3.16)\]

\[O_{n} = \begin{cases} 1 \text{ jika ganjil} \\ 0 \text{ jika Genap} \end{cases}\] (3.17)

Persamaan (3.15) merupakan persamaan diferensial amplitudo modus getar batang elastik.

4. SOLUSI TUNAK DAN KESTABILAN PARAMETRIK

Persamaan (3.15) dapat diselesaikan dengan metode solusi tunak. Metode ini telah banyak digunakan untuk memperoleh respon dinamik sistem, di antaranya oleh Hsu dan Cheng [7], dan Midha et al. [8]. Keduanya menyatakan bahwa solusi periodik dapat dibagi menjadi beberapa subinterval di mana parameter dianggap konstan pada tiap subinterval. Selanjutnya, persamaan diferensial orde dua pada Persamaan (3.15) dapat dituliskan sebagai persamaan orde satu berikut:

\[\dot{\mathbf{x}}(\mathbf{t}) = \mathbf{A}(\mathbf{t}) \, \mathbf{x}(\mathbf{t}) + \mathbf{f}(\mathbf{t}) \tag{4.1}\] dengan A menyatakan koefisien matriks, x(t) simpangan dan f(t) gaya eksitasi. Untuk sistem periodik berlaku hubungan

\[A(t) = A(T_a + t) \tag{4.2}\]

\[f(t) = f(T_f + t) \tag{4.3}\] dan

\[x(t) = x(T+t) \tag{4.4}\] dengan T, T_a dan T_f masing-masing menyatakan periode solusi, koefisien matriks A dan gaya eksitasi.

Penyelesaian Persamaan (4.1) di atas dapat dilakukan dengan diskretisasi periode T atas N subinterval dimana elemen matriks A dianggap konstan pada tiap subinterval. Konstanta matriks A dihitung dengan merata-ratakan koefisien matriks tersebut pada kedua ujung interval. Diskretisasi matriks A diilustrasikan pada Gambar 6, di mana dalam subinterval ke-i, elemen matriks A_i adalah konstan.

Gambar 6 Diskretisasi persamaan state-space.

Persamaan (4.1) untuk subinterval ke-i menjadi:

\[\dot{x}(t) = A_i x(t) + f(t); t_{i-1} \le t \le t_i\] (4.5)

Solusi Persamaan (4.5) diberikan oleh:

\[x(t) = e^{A_i(t-t_{i-1})} x(t_{i-1}) + \int_{t_{i-1}}^{t} e^{A_i(t-\tau)} f(\tau) d\tau\] (4.6)

Solusi pada simpul subinterval diberikan oleh:

MESIN Vol. XVII No. 3

\[x_{i+1} = x(t_i)\]; \(i = 1, 2, ..., N-1\); \(t_{i-1} \le t \le t_i\) (4.7)

Persyaratan kompatibilitas pada simpul subinterval, memberikan hubungan berikut:

\[x_{i+1} = e^{A_i \Delta t} x_i - S_i\]; \(i = 1, 2, ..., N\) (4.8)

\(S_{i} = -\int_{0}^{t} e^{A_{i}(t_{i}-\tau)} f(\tau) d\tau\) (4.9)

Persyaratan periodisitas menyatakan solusi berulang setelah satu periode (N subinterval), sehingga hubungan solusi pada subinterval ke-1 dan ke-N+1 dapat dinyatakan sebagai:

\[x_1 = x_{N+1} = e^{A_N \Delta t} x_N - S_N\] (4.10)

Kombinasi Persamaan (4.8), (4.9) dan (4.10) untuk satu periode dapat dinyatakan dalam bentuk matriks berikut:

\[\text{[rumus tidak dapat ditampilkan dengan baik — lihat PDF asli]}\]

dimana I adalah matriks identitas dan \(B_i\) adalah pendekatan deret Taylor \(e^{A_i \Delta t}\) dengan jumlah suku n

\[e^{A_i \Delta t} = B_i = I + \sum_{i=1}^{n} \frac{(A_i \Delta t)^j}{j!}\] (4.12)

Persamaan (4.11) dapat disederhanakan dengan menyatakan \(x_i\) secara eksplisit dengan mengalikan baris pertama dengan \(B_2\) kemudian ditambahkan ke baris kedua, selanjutnya \(B_3\) dikalikan dengan baris kedua dan dijumlahkan dengan baris ketiga. Langkah ini dilakukan secara berurutan untuk baris lainnya, sehingga Persamaan (4.11) berubah menjadi:

\[\begin{bmatrix} M_1 & -I & 0 & 0 & \dots & 0 & 0 \\ M_2 & 0 & -I & 0 & \dots & \dots & \dots \\ \vdots & \vdots & \ddots & \ddots & \ddots & \ddots & \vdots \\ M_{N-2} & 0 & 0 & 0 & \dots & -I & 0 \\ M_{N-1} & 0 & 0 & 0 & \dots & 0 & -I \\ M_N & 0 & 0 & 0 & \dots & 0 & B_N \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \\ \vdots \\ x_{N-2} \\ X_{N-1} \\ x_N \end{bmatrix} = \begin{bmatrix} G_1 \\ G_2 \\ \vdots \\ G_{N-2} \\ G_{N-1} \\ G_N \end{bmatrix}\] \[(4.13)\] dengan

\[G_1 \equiv S_1 \tag{4.14a}\]

\[G_i \equiv S_i + B_i G_{i-1}\]; \(i = 1, 2, ..., N\) (4.14b)

\[M_i = \prod_{i=1}^{i} B_i = B_i B_{i-1} B_{i-2} \cdots B_2 B_1\] (4.14c)

Baris terakhir Persamaan (4.13) memberikan

\[x_1 = (M_N - I)^{-1} G_N (4.15)\]

sedangkan solusi pada selang ke-i diperoleh secara berulang berdasarkan hubungan berikut:

\[x_i = M_{i-1}x_i - G_{i-1} = M_{i-1}(M_N - I)^{-1}G_N - G_{i-1};\]
\(i = 2, 3, \dots, N\) (4.16)

\(M_N\) merupakan matriks monodromi yang digunakan untuk mengevaluasi kestabilan solusi dengan asumsi periodik. Persamaan (4.16) dapat diselesaikan jika \(det(M_N - I) \neq 0\) (4.17)

4.1 Kestabilan Parametrik

Banyak metode yang dapat digunakan untuk mengevaluasi kestabilan solusi tunak, di antaranya metode Hill's, metode Perturbation dan Floquet. Metode Hill's merupakan metode analitik yang sulit diterapkan untuk sistem dengan banyak derajat kebebasan, sebaliknya metode Floquet merupakan pendekatan numerik dan mudah diterapkan untuk sistem dengan banyak derajat kebebasan.

Metode Floquet berdasarkan pada evaluasi eigenvalue matriks monodromi [9]. Langkah awal metode tersebut adalah dengan menyatakan persamaan gerak dalam orde satu seperti berikut:

\[\dot{y}(t) = A(t) y(t)\] (4.18)

Solusi Persamaan (4.18) telah dikembangkan oleh Richards [2], yaitu:

\[\phi(t,0) = \begin{bmatrix} y_1(t) & y_2(t) & \cdots & y_n(t) \end{bmatrix} \begin{bmatrix} y_1(0) & y_2(0) & \cdots & y_n(0) \end{bmatrix}^{-1}\] (4.19)

A(t) dan \(\phi(t)\) adalah matriks dengan ukuran (n x n) dan \(y_i\) merupakan vektor yang berukuran (n x 1). Substitusi Persamaan (4.18) ke turunan pertama Persamaan (4.19) memberikan

\[\dot{\phi}(t,0) = A(t) \, \phi(t,0)\] (4.20)

di mana \(\phi(t)\) adalah periodik, sehingga \(\phi(t+T)\) adalah juga solusi dari Persamaan (4.18) dengan hubungan:

\[\phi(t+T,0) = C\phi(t,0) \tag{4.21}\]

C merupakan matriks konstan dan T periode \(\phi(t)\). Pada t=0, Persamaan (4.21) dapat dituliskan menjadi

\[\phi(T,0) = C \phi(0,0) \tag{4.22}\]

Dari Persamaan (4.19) untuk t = 0 diperoleh

\[\phi(0,0) = I \tag{4.23}\]

Dengan memasukkan persamaan (4.23) ke (4.22) didapatkan harga matriks C yaitu :

\[C = \phi(T,0) \tag{4.24}\]

Persamaan (4.19), (4.24) dan (4.24) memberikan:

\[y(t+T) = C y(t)\] (4.25)

di mana matriks C adalah matriks monodromi yang menghubungkan solusi setelah satu periode, T.

Suatu sistem dikatakan stabil jika respon dinamiknya berbatas untuk \(t \to \infty\). Untuk m periode, Persamaan (4.25) dapat dituliskan sebagai :

\[y(t+mT) = C^{m} y(t)\] (4.26)

Persamaan (4.26) menunjukkan bahwa karakteristik respon sangat dipengaruhi oleh matriks C.

4.2 Evaluasi Kestabilan Parametrik

Kondisi respon adalah ekivalen dengan kondisi matriks \(C^m\), sehingga dimungkinkan untuk menganalisis kestabilan melalui evaluasi matriks monodromi C. Matriks C dapat diubah menjadi bentuk Jordan Canonical berikut:

\[C = P^{-1}QP \tag{4.28}\] atau

\[C^{m} = P^{-1}Q^{m}P \tag{4.29}\] di mana karakteristik matriks \(C^m\) diwakili oleh \(Q^m\). Matriks Q adalah matriks diagonal yang mengandung nilai eigen matriks C, \(\lambda_i\), yang diperoleh dari

\[\det(\mathbf{C} - \lambda \mathbf{I}) = 0 \tag{4.30}\]

Matriks Q tersebut dapat dituliskan sebagai berikut;

\[Q = \begin{bmatrix} \lambda_{1} & 0 & \dots & 0 & 0 \\ \vdots & \lambda_{2} & \dots & \ddots & 0 \\ \vdots & \vdots & \ddots & \ddots & \vdots \\ \vdots & \vdots & \ddots & \ddots & \vdots \\ \vdots & \vdots & \ddots & \lambda_{n-1} & \vdots \\ 0 & 0 & \dots & \lambda_{n} \end{bmatrix}\](4.31)

Suatu sistem akan stabil jika \(|\lambda_i| < 1\) dan tidak stabil jika \(|\lambda_i| > 1\). Untuk harga \(|\lambda_i| = 1\) respon sistem bersifat periodik dengan periode T. Pada kondisi khusus yaitu \(\lambda_i = -1\) maka respon periodik dengan periode 2T.

Matriks monodromi C tersebut dapat dihitung dari persamaan (4.14)

\[C = M_N = \prod_{i=1}^{N} B_i \approx \prod_{i=1}^{N} \left( I + \sum_{j=1}^{n} \frac{(A_i \Delta t)^j}{J!} \right)\] (4.32)

di mana n adalah jumlah suku pendekatan Deret Taylor.

5 ANALISIS DINAMIK MEKANISME EMPAT BATANG

Analisis kestabilan dan respons batang penghubung elastik pada sebuah mekanisme empat batang seperti yang diperlihatkan pada Gambar 4 akan ditinjau dalam bagian ini. Material batang elastik adalah Aluminium (modulus elastisitas, E=68 GPa, dan massa jenis, \(\rho=2700\) kg/m³), sedangkan dimensi mekanisme dirangkum pada Tabel 1.

Tabel 1 Parameter geometri mekanisme empat batang.

5. 1 Peta Kestabilan Parametrik

Kaji kestabilan parametrik dilakukan untuk batang penghubung elastik dengan berbagai variasi parameter amplitudo dan frekuensi yang masing-masing diwakili panjang engkol nisbah terhadap oleh penghubung, μ, dan kecepatan putar tuna dimensi, Ω, dengan metode grid. Hasil analisis disajikan sebagai peta kestabilan parametrik pada bidang Ω-μ, dengan variasi \(\Omega\) antara 0 – 4 dan \(\mu\) antara 0.01 – 0.20. Harga \(\Omega\)dan u dimasukkan ke Persamaan (3.15) untuk kemudian eigenvalue matriks monodromi, λ_i, dievaluasi. Jika ada \(|\lambda_i| > 1\), maka pada koordinat \((\Omega, \mu)\) diplot tanda *,mengindikasikan ketidakstabilan parametrik untuk kombinasi parameter amplitudo dan frekuensi tersebut.

5. 1.1 Pengaruh Redaman

Peta kestabilan parametrik untuk kasus \(\hat{F}=0\) dan S=0.5 disajikan pada Gambar 7 dan 8, berturut-turut untuk faktor redaman, \(\eta=0.0001\) dan 0.001, yang dianggap mewakili redaman internal material. Terlihat bahwa ketidakstabilan dapat terjadi pada parameter frekuensi, \(\Omega<1\), terutama untuk parameter amplitudo yang relatif besar. Perbandingan Gambar 7 dan 8 juga menunjukkan daerah tak-stabil berkurang dengan meningkatnya redaman.

Gambar 7 Peta kestabilan untuk \(\hat{F} = 0\), S = 0.5, \(\eta = 0.0001\).

5. 1.2 Pengaruh Gaya Batang Pengikut

Pengaruh gaya batang pengikut terhadap kestabilan parametrik dipelajari dengan meninjau kasus dengan amplitudo gaya eksitasi sebesar 15 P_c, di mana P_c adalah beban kritis Euler. Peta kestabilan disajikan pada Gambar 9. Membandingkan Gambar 8 dan 9 terlihat bahwa adanya gaya tarik yang diberikan batang pengikut menjadikan daerah stabil respon batang elastik lebih luas.

Gambar & Peta kestabilan untuk \(\hat{\mathbf{F}} = 0\), S = 0.5, \(\eta = 0.001\).

5. 1.3 Pengaruh Rasio Massa Penghubung-Pengikut Pengaruh besarnya massa penghubung ditinjau dengan mengambil rasio, S=0.001 untuk kasus \(\hat{F}=0\) dan \(\eta=0.001\). Peta kestabilan yang diperoleh disajikan pada Gambar 10. Perbandingan Gambar 10 dengan Gambar 8 menunjukkan bahwa makin besar rasio S, makin luas daerah stabil batang elastik.

Gambar 9 Peta kestabilan untuk \(\hat{\mathbf{F}} = 15 \; P_c, \; S = 0.5, \; dan \; \eta = 0.001.\)

Gambar 10 Peta kestabilan untuk \(\hat{\mathbf{F}} = 0\), S = 0.001, dan \(\eta = 0.001\).

5. 2 Validasi Analisis Stabilitas dan Respon Dinamik. Untuk mengkaji keabsahan analisis stabilitas yang disajikan di atas, maka persamaan gerak batang elastik untuk kasus stabil dan tak-stabil diselesaikan dengan metode integrasi langsung. Runge-Kutta. Untuk kasus stabil ditinjau titik A pada Gambar 8 dengan \(\Omega = 0.6\) dan \(\mu = 0.03\), sedangkan untuk daerah tak stabil ditinjau titik B (\(\Omega = 2, \mu = 0.05\)).

Untuk titik A, amplitudo modus getar pertama yang diperoleh melalui algoritma solusi-tunak dan metode integrasi langsung disajikan pada Gambar 11. Terlihat bahwa respon metode integrasi langsung adalah tunak serta kedua respon bersesuaian.

Respon untuk titik B disajikan dalam Gambar 12. Terlihat bahwa amplitudo respon yang diperoleh melalui metode Runge-Kutta menunjukkan kecenderungan divergen. Hal ini mengindikasikan bahwa respon batang elastik dengan kondisi operasi \(\Omega=2, \mu=0.05\) adalah tidak stabil.

Solusi modus pertama yang disajikan pada Gambar 11 dan 12 ini mendukung keabsahan analisis stabilitas dinamik dan algorima solusi-tunak.

Gambar 11 Amplitudo modus pertama batang elastik mekanisme empat batang titik A.

Gambar 12 Amplitudo modus pertama batang elastik mekanisme empat batang kondisi takstabil, titik B.

6 KESIMPULAN

Elastisitas batang penghubung yang mengakibatkan variasi parameter dinamik terhadap posisi membuat sistem dapat mengalami ketidakstabilan parametrik (parametric instability). Dalam hal ini, respon sistem dapat menjadi tidak stabil pada berbagai pita frekuensi, bahkan pada kecepatan operasi yang relatif rendah.

Pengaruh variasi parameter \(\mu\), S dan \(\eta\) telah diinvestigasi dan dapat disimpulkan berberapa hal berikut:

• Perbandingan massa batang elastik terhadap peluncur atau batang pengikut sangat berpengaruh terhadap konstruksi peta kestabilan. Harga S yang lebih besar berkorelasi dengan bertambahnya luas daerah stabil pada peta kestabilan parametrik. Hal ini disebabkan karena batang penghubung dangan massa yang lebih besar akan lebih kaku sehingga amplituodo variasi kekakuan yang dihasilkan menjadi lebih kecil.
• Panjang engkol merepesentasikan parameter amplitudo eksitasi batang elastik. Parameter μ yang relatif besar akan dapat menimbulkan ketidakstabilan parametrik pada frekuensi rendah. Semakin besar harga μ, maka semakin mudah ketidakstabilan parametrik terjadi.
• 3. Konfigurasi gaya luar juga mempengaruhi kestabilan mekanisme. Gaya tarik cenderung memberikan luas daerah stabil yang lebih besar dibandingkan dengan gaya tekan.
• Redaman menyebabkan kombinasi parameter yang semula menghasilkan ketidakstabilan dapat menjadi stabil, terutama pada daerah μ kecil. Penambahan redaman mengakibatkan luas daerah tidak stabil berkurang.

Penelitian ini memberikan gambaran tentang perilaku dinamik batang penghubung elastik pada sebuah mekanisme empat-batang. Metodologi yang digunakan dapat diterapkan untuk meninjau perilaku dinamik mekanisme dengan batang elastik lainnya. Analisis kestabilan dapat memberikan gambaran tentang daerah operasi yang aman. Di samping itu, algorimta solusitunak dapat digunakan untuk memprediksi amplitudo getaran batang elastik, untuk kemudian dibandingkan dengan batasan yang diperbolehkan.