چکیده

آونگ دوگانه یکی از شاخص‌ترین سیستم‌های دینامیکی کلاسیک است که با وجود سادگی ظاهری، رفتار پیچیده و آشوب‌ناک از خود نشان می‌دهد. این مقاله چارچوبی جامع برای تحلیل دینامیک آشوب‌ناک، کمی‌سازی حساسیت به شرایط اولیه و ارزیابی پایداری عددی این سیستم ارائه می‌کند. ابتدا، معادلات حرکت با استفاده از فرمول‌بندی لاگرانژین و با در نظر گرفتن گشتاورهای میرایی خطی استخراج می‌شوند. سپس سه انتگرال‌گیر عددی — رانگ-کوتای مرتبهٔ چهار، ورلهٔ سرعتی هم‌تافته و اویلر مرتبهٔ اول — از منظر دقت، پایداری انرژی و کارایی مقایسه می‌گردند. برای کمی‌سازی آشوب، دو مسیر سایه با اغتشاش‌های اولیهٔ متفاوت در فضای فاز ردیابی شده و بزرگ‌ترین نمای لیاپانوف با برازش خطی روی لگاریتم واگرایی محاسبه می‌شود. علاوه بر این، زمان لیاپانوف $t_\lambda = 1/\lambda$ به‌عنوان معیار افق پیش‌بینی‌پذیری سیستم معرفی و نسبت $t/t_\lambda$ به‌عنوان شاخصی برای ارزیابی فراتر رفتن سیستم از مرز پیش‌بینی‌پذیری به‌کار می‌رود. در پایان، بعد همبستگی جاذب آشوب‌ناک با الگوریتم گراسبرگر-پروکاچیا محاسبه می‌شود. نتایج شبیه‌سازی با پارامترهای پیش‌فرض — جرم‌های ۱ کیلوگرم، طول‌های ۱ متر، زوایای اولیهٔ ۱۲۰ و ۹۰ درجه، و گام زمانی ۲×۱۰⁻⁴ ثانیه — نشان می‌دهد که سیستم در بازهٔ چند ثانیه از محدودهٔ پیش‌بینی‌پذیری فراتر می‌رود و واگرایی نمایی مسیرها با نمای لیاپانوف مثبت به‌وقوع می‌پیوندد. RK4 با وجود دقت بالای کوتاه‌مدت، در بازه‌های طولانی انرژی را به‌صورت سیستماتیک دریفت می‌دهد، در حالی که ورلهٔ سرعتی به‌دلیل خواص هم‌تافته، انرژی را تا حدود زیادی حفظ می‌کند. این چارچوب، به‌عنوان مطالعه‌ای جامع در آشوب قطعی، می‌تواند مبنایی برای درک رفتار سیستم‌های غیرخطی چنددرجه‌آزادی در مکانیک، هواشناسی، زیست‌شناسی و اقتصاد باشد.

۱. مقدمه

آونگ دوگانه — سیستمی متشکل از دو جرم نقطه‌ای که توسط میله‌های صلب بدون جرم به یکدیگر و به یک تکیه‌گاه ثابت متصل شده‌اند — یکی از ساده‌ترین سیستم‌های دینامیکی است که رفتار آشوب‌ناک از خود نشان می‌دهد. این سیستم با دو درجهٔ آزادی زاویه‌ای، نمونه‌ای کلاسیک از یک سیستم دینامیکی غیرخطی است که معادلات حرکت آن هیچ جواب تحلیلی بسته‌ای ندارد.

از زمان مطالعهٔ هانری پوانکاره بر مسئلهٔ سه‌جسمی در اواخر قرن نوزدهم، مفهوم آشوب قطعی به‌عنوان ویژگی بنیادین سیستم‌های غیرخطی شناخته شده است. ادوارد لورنتز در دههٔ ۱۹۶۰ با مطالعهٔ مدل ساده‌شدهٔ همرفت جوّی، نشان داد که حساسیت به شرایط اولیه می‌تواند پیش‌بینی بلندمدت را در سیستم‌های قطعی ناممکن سازد. این پدیده که به‌عنوان «اثر پروانه‌ای» شهرت یافته، در آونگ دوگانه نیز به‌وضوح قابل مشاهده است.

این مقاله چارچوبی جامع برای تحلیل عددی این سیستم ارائه می‌کند. ابتدا، معادلات حرکت با فرمول‌بندی لاگرانژین استخراج می‌شوند. سپس سه انتگرال‌گیر عددی مقایسه می‌گردند. در ادامه، روش‌های کمی‌سازی آشوب — نمای لیاپانوف و بعد همبستگی — معرفی و پیاده‌سازی می‌شوند. در پایان، نتایج شبیه‌سازی در قالب نمودارهای سری زمانی، فضای فاز، انرژی و واگرایی مسیرها تحلیل می‌گردد.

🎯 هدف اصلی

ارائهٔ چارچوبی جامع برای کمی‌سازی آشوب در سیستم آونگ دوگانه از طریق نمای لیاپانوف، ارزیابی پایداری عددی سه انتگرال‌گیر مختلف، و تحلیل ویژگی‌های هندسی جاذب با محاسبهٔ بعد همبستگی.

۲. استخراج معادلات حرکت

۲.۱. مختصات تعمیم‌یافته و موقعیت جرم‌ها

سیستم آونگ دوگانه با دو درجهٔ آزادی زاویه‌ای توصیف می‌شود: $\theta_1$ زاویهٔ میلهٔ اول نسبت به راستای قائم و $\theta_2$ زاویهٔ میلهٔ دوم نسبت به راستای قائم. با انتخاب مبدأ مختصات در تکیه‌گاه و راستای $y$ به‌سمت پایین، موقعیت دو جرم به‌صورت زیر به‌دست می‌آید:

$$\begin{aligned} x_1 &= L_1 \sin\theta_1, & y_1 &= -L_1 \cos\theta_1 \\ x_2 &= L_1 \sin\theta_1 + L_2 \sin\theta_2, & y_2 &= -L_1 \cos\theta_1 - L_2 \cos\theta_2 \end{aligned}$$

۲.۲. لاگرانژین سیستم

انرژی جنبشی $T$ و انرژی پتانسیل $V$ سیستم با مشتق‌گیری از موقعیت‌ها نسبت به زمان به‌دست می‌آیند:

$$T = \frac{1}{2}m_1 L_1^2 \dot{\theta}_1^2 + \frac{1}{2}m_2\left[L_1^2\dot{\theta}_1^2 + L_2^2\dot{\theta}_2^2 + 2L_1L_2\dot{\theta}_1\dot{\theta}_2\cos\Delta\right]$$
$$V = -(m_1+m_2)gL_1\cos\theta_1 - m_2gL_2\cos\theta_2$$

که در آن $\Delta = \theta_1 - \theta_2$ اختلاف زوایا است. لاگرانژین $\mathcal{L} = T - V$ به‌صورت زیر تعریف می‌شود:

$$\mathcal{L} = \frac{1}{2}(m_1+m_2)L_1^2\dot{\theta}_1^2 + \frac{1}{2}m_2 L_2^2\dot{\theta}_2^2 + m_2 L_1 L_2 \dot{\theta}_1\dot{\theta}_2\cos\Delta + (m_1+m_2)gL_1\cos\theta_1 + m_2 gL_2\cos\theta_2$$

۲.۳. معادلات اویلر-لاگرانژ

با اعمال معادلات اویلر-لاگرانژ و افزودن گشتاورهای میرایی خطی $Q_i = -b_i\dot{\theta}_i$، معادلات حرکت به‌شکل زیر حاصل می‌شوند:

$$\begin{aligned} (m_1+m_2)L_1\ddot{\theta}_1 + m_2 L_2\ddot{\theta}_2\cos\Delta + m_2 L_2\dot{\theta}_2^2\sin\Delta + (m_1+m_2)g\sin\theta_1 &= -b_1\dot{\theta}_1 \\ m_2 L_2\ddot{\theta}_2 + m_2 L_1\ddot{\theta}_1\cos\Delta - m_2 L_1\dot{\theta}_1^2\sin\Delta + m_2 g\sin\theta_2 &= -b_2\dot{\theta}_2 \end{aligned}$$

۲.۴. فرم ماتریسی و حل دستگاه

با نوشتن دستگاه به‌شکل $\mathbf{M}(\boldsymbol{\theta})\ddot{\boldsymbol{\theta}} = \mathbf{F}(\boldsymbol{\theta}, \dot{\boldsymbol{\theta}})$، ماتریس جرم $\mathbf{M}$ و بردار نیروها $\mathbf{F}$ به‌صورت زیر تعریف می‌شوند:

$$\mathbf{M} = \begin{bmatrix} (m_1+m_2)L_1 & m_2 L_2\cos\Delta \\ m_2 L_1\cos\Delta & m_2 L_2 \end{bmatrix}$$

دترمینان این ماتریس به‌صورت $\det\mathbf{M} = m_2 L_1 L_2 [m_1 + m_2 \sin^2\Delta] > 0$ همواره مثبت است، بنابراین دستگاه در هر لحظه جواب یکتا دارد. با اعمال قاعدهٔ کرامر، شتاب‌های زاویه‌ای به‌دست می‌آیند:

$$\ddot{\theta}_1 = \frac{F_1 \cdot M_{22} - M_{12} \cdot F_2}{\det\mathbf{M}}, \quad \ddot{\theta}_2 = \frac{M_{11} \cdot F_2 - F_1 \cdot M_{21}}{\det\mathbf{M}}$$

۳. روش‌های انتگرال‌گیری عددی

۳.۱. تبدیل به دستگاه مرتبهٔ اول

برای حل عددی، دستگاه مرتبهٔ دوم به یک دستگاه مرتبهٔ اول با بردار حالت $\mathbf{y} = [\theta_1, \theta_2, \omega_1, \omega_2]^T$ تبدیل می‌شود:

$$\dot{\mathbf{y}} = \mathbf{f}(\mathbf{y}) = [\omega_1, \omega_2, \alpha_1(\mathbf{y}), \alpha_2(\mathbf{y})]^T$$

۳.۲. رانگ-کوتای مرتبهٔ چهار (RK4)

روش RK4 یکی از پرکاربردترین انتگرال‌گیرهای عددی است که دقت مرتبهٔ چهارم را با محاسبهٔ چهار ارزیابی از تابع $\mathbf{f}$ در هر گام به‌دست می‌آورد:

$$\begin{aligned} \mathbf{k}_1 &= \mathbf{f}(\mathbf{y}_n) \Delta t \\ \mathbf{k}_2 &= \mathbf{f}\!\left(\mathbf{y}_n + \tfrac{\mathbf{k}_1}{2}\right) \Delta t \\ \mathbf{k}_3 &= \mathbf{f}\!\left(\mathbf{y}_n + \tfrac{\mathbf{k}_2}{2}\right) \Delta t \\ \mathbf{k}_4 &= \mathbf{f}\!\left(\mathbf{y}_n + \mathbf{k}_3\right) \Delta t \\ \mathbf{y}_{n+1} &= \mathbf{y}_n + \frac{1}{6}(\mathbf{k}_1 + 2\mathbf{k}_2 + 2\mathbf{k}_3 + \mathbf{k}_4) + \mathcal{O}(\Delta t^5) \end{aligned}$$

RK4 با وجود دقت بالا در بازه‌های کوتاه‌مدت، غیرهم‌تافته است؛ یعنی در بازه‌های طولانی، انرژی کل سیستم را به‌صورت سیستماتیک دریفت می‌دهد. این دریفت به‌دلیل انباشت خطاهای گردکردن و نبود ساختار هندسی مناسب در الگوریتم است.

۳.۳. ورلهٔ سرعتی (هم‌تافته)

الگوریتم ورلهٔ سرعتی یک روش هم‌تافته است که ساختار هندسی فضای فاز را حفظ می‌کند و در نتیجه، انرژی کل را در بازه‌های طولانی به‌طور قابل‌توجهی پایدار نگه می‌دارد:

$$\begin{aligned} \mathbf{y}_{n+1} &= \mathbf{y}_n + \mathbf{v}_n \Delta t + \tfrac{1}{2}\mathbf{a}_n \Delta t^2 \\ \mathbf{a}_{n+1} &= \mathbf{f}(\mathbf{y}_{n+1}) \\ \mathbf{v}_{n+1} &= \mathbf{v}_n + \tfrac{1}{2}(\mathbf{a}_n + \mathbf{a}_{n+1}) \Delta t \end{aligned}$$

خاصیت هم‌تافتی به این معناست که حجم فضای فاز توسط الگوریتم تغییر نمی‌کند — ویژگی‌ای که برای سیستم‌های همیلتونی نظیر آونگ دوگانه (در غیاب میرایی) حیاتی است. این الگوریتم برای شبیه‌سازی‌های طولانی‌مدت انتخاب مناسب‌تری است.

۳.۴. اویلر مرتبهٔ اول

ساده‌ترین انتگرال‌گیر عددی، روش اویلر است که تنها از ارزیابی اولیهٔ تابع مشتق استفاده می‌کند:

$$\mathbf{y}_{n+1} = \mathbf{y}_n + \mathbf{f}(\mathbf{y}_n) \Delta t + \mathcal{O}(\Delta t^2)$$

این روش با وجود سادگی، دقت پایینی دارد و انرژی را به‌سرعت افزایش می‌دهد. اما به‌عنوان مرجع مقایسه، برای نشان دادن اهمیت انتخاب انتگرال‌گیر مناسب، ارزشمند است.

جدول ۱: مقایسه سه انتگرال‌گیر عددی
روش مرتبه دقت هم‌تافته پایداری انرژی هزینه محاسباتی
RK4 ۴ خیر دریفت سیستماتیک در بلندمدت ۴ ارزیابی در هر گام
ورلهٔ سرعتی ۲ بله پایدار در بلندمدت ۱ ارزیابی در هر گام
اویلر ۱ خیر ناپایدار — افزایش انرژی ۱ ارزیابی در هر گام
پیاده‌سازی هستهٔ مشتقات و انتگرال‌گیرهاJavaScript
// مشتقات: محاسبه شتاب‌های زاویه‌ای با قاعده کرامر
function derivatives(theta1, theta2, omega1, omega2) {
  const d = theta1 - theta2;
  const sinD = Math.sin(d), cosD = Math.cos(d);

  // اجزای ماتریس جرم
  const A11 = (m1 + m2) * L1;
  const A12 = m2 * L2 * cosD;
  const A21 = m2 * L1 * cosD;
  const A22 = m2 * L2;
  const det = A11 * A22 - A12 * A21;

  // سمت راست دستگاه
  const RHS1 = -(m1+m2)*g*Math.sin(theta1)
              - m2*L2*omega2*omega2*sinD - b1*omega1;
  const RHS2 = -m2*g*Math.sin(theta2)
              + m2*L1*omega1*omega1*sinD - b2*omega2;

  // حل با قاعده کرامر
  const alpha1 = (RHS1*A22 - A12*RHS2) / det;
  const alpha2 = (A11*RHS2 - RHS1*A21) / det;
  return [omega1, omega2, alpha1, alpha2];
}

// گام RK4
function rk4Step(state) {
  const y = [state.theta1, state.theta2, state.omega1, state.omega2];
  const k1 = derivatives(...y).map(v => v * DT);
  const k2 = derivatives(
    y[0]+k1[0]/2, y[1]+k1[1]/2,
    y[2]+k1[2]/2, y[3]+k1[3]/2
  ).map(v => v * DT);
  const k3 = derivatives(
    y[0]+k2[0]/2, y[1]+k2[1]/2,
    y[2]+k2[2]/2, y[3]+k2[3]/2
  ).map(v => v * DT);
  const k4 = derivatives(
    y[0]+k3[0], y[1]+k3[1],
    y[2]+k3[2], y[3]+k3[3]
  ).map(v => v * DT);
  const ny = y.map((v, i) =>
    v + (k1[i] + 2*k2[i] + 2*k3[i] + k4[i]) / 6
  );
  return { theta1:ny[0], theta2:ny[1], omega1:ny[2], omega2:ny[3], t:state.t+DT };
}

۴. کمی‌سازی آشوب

۴.۱. حساسیت به شرایط اولیه و اثر پروانه‌ای

هستهٔ مفهوم آشوب قطعی، حساسیت به شرایط اولیه است. در یک سیستم آشوب‌ناک، دو مسیر که با اختلاف بسیار کوچکی در شرایط اولیه آغاز شده‌اند، به‌طور نمایی از یکدیگر فاصله می‌گیرند. برای کمی‌سازی این رفتار، شبیه‌ساز حاضر دو مسیر سایه با اختلاف‌های اولیهٔ کوچک در $\theta_1$ ایجاد می‌کند و فاصلهٔ اقلیدسی آن‌ها را در فضای فاز ردیابی می‌نماید.

فاصلهٔ فاز با در نظر گرفتن پیچش زاویه (wrapping) در بازهٔ $[-\pi, \pi]$ محاسبه می‌شود:

$$\Delta\theta_i^{\text{wrap}} = \left[(\theta^{\text{main}} - \theta^{\text{shadow}}) + \pi\right] \bmod 2\pi - \pi$$
$$d(t) = \sqrt{(\Delta\theta_1^{\text{wrap}})^2 + (\Delta\theta_2^{\text{wrap}})^2}$$

۴.۲. نمای لیاپانوف

برای سیستم‌های آشوب‌ناک، فاصلهٔ بین دو مسیر به‌طور نمایی رشد می‌کند:

$$d(t) \approx d(0) \, e^{\lambda t}$$

که در آن $\lambda$ بزرگ‌ترین نمای لیاپانوف است. با گرفتن لگاریتم طبیعی از دو طرف:

$$\ln d(t) \approx \ln d(0) + \lambda t$$

این رابطهٔ خطی، امکان تخمین $\lambda$ را با برازش خطی حداقل مربعات روی ناحیه‌ای از داده‌ها فراهم می‌کند که هنوز به اشباع نرسیده است. شیب خط برازش، همان نمای لیاپانوف است:

$$\lambda = \frac{\sum_{i=1}^{N}(t_i - \bar{t})(\ln d_i - \overline{\ln d})}{\sum_{i=1}^{N}(t_i - \bar{t})^2}$$

مقدار مثبت $\lambda$ نشان‌دهندهٔ آشوب‌ناکی سیستم است. در شبیه‌ساز، برای انتخاب ناحیهٔ خطی، ابتدا نقطهٔ اشباع (جایی که رشد نمایی متوقف می‌شود) شناسایی و سپس بازهٔ ۲۰٪ تا ۸۰٪ از این ناحیه برای برازش انتخاب می‌شود.

۴.۳. زمان لیاپانوف و افق پیش‌بینی‌پذیری

پارامتر کلیدی دیگر، زمان لیاپانوف است که به‌صورت معکوس نمای لیاپانوف تعریف می‌شود:

$$t_\lambda = \frac{1}{\lambda}$$

زمان لیاپانوف، مقیاس زمانی مشخصه‌ای است که در آن خطای اولیه به اندازهٔ $e$ برابر بزرگ می‌شود. نسبت $t/t_\lambda$ به‌عنوان شاخصی برای ارزیابی وضعیت سیستم به‌کار می‌رود:

این معیار، اهمیت عملی آشوب را در حوزه‌های مختلف نشان می‌دهد: در پیش‌بینی هوا، $t_\lambda$ حدود چند روز است؛ در آونگ دوگانه، این زمان تنها چند ثانیه است.

۴.۴. بعد همبستگی و هندسهٔ جاذب

برای توصیف هندسهٔ جاذب آشوب‌ناک، از مفهوم بعد همبستگی $D_2$ استفاده می‌شود که با الگوریتم گراسبرگر-پروکاچیا محاسبه می‌گردد. تابع همبستگی $C(r)$ به‌صورت زیر تعریف می‌شود:

$$C(r) = \frac{2}{N(N-1)}\sum_{i < j} \Theta\!\left(r - \left\| \mathbf{x}_i - \mathbf{x}_j \right\|\right)$$

که در آن $\Theta$ تابع پلهٔ هویساید، $N$ تعداد نقاط نمونه، و $\mathbf{x}_i$ بردار حالت در فضای فاز است. برای مقادیر کوچک $r$، این تابع به‌صورت قانون توانی رشد می‌کند:

$$C(r) \sim r^{D_2} \quad \Rightarrow \quad D_2 = \lim_{r \to 0} \frac{\ln C(r)}{\ln r}$$

برای جاذب آشوب‌ناک، $D_2$ عددی غیرصحیح است که نشان‌دهندهٔ ساختار فراکتالی جاذب است. در فضای فاز چهاربُعدی آونگ دوگانه، مقدار $D_2$ معمولاً بین ۱٫۵ تا ۲٫۵ قرار می‌گیرد.

۵. نتایج شبیه‌سازی

۵.۱. پارامترهای پیش‌فرض

شبیه‌سازی با پارامترهای پیش‌فرض زیر اجرا شده است: جرم‌های $m_1 = m_2 = 1$ kg، طول‌های $L_1 = L_2 = 1$ m، شتاب گرانش $g = 9.81$ m/s²، زوایای اولیهٔ $\theta_{1,0} = 120°$ و $\theta_{2,0} = 90°$، سرعت‌های زاویه‌ای اولیه صفر، و گام زمانی $\Delta t = 2 \times 10^{-4}$ ثانیه.

برای کمی‌سازی آشوب، دو مسیر سایه با اختلاف اولیهٔ $\delta_1 = 0.01°$ و $\delta_2 = 0.001°$ در زاویهٔ $\theta_1$ ایجاد شده‌اند.

۵.۲. سری زمانی و واگرایی مسیرها

نمودار سری زمانی، واگرایی پیش‌روندهٔ مسیر اصلی از مسیرهای سایه را نشان می‌دهد. در ابتدا، مسیرها تقریباً منطبق هستند، اما با گذشت زمان، انحراف مسیرها به‌صورت تدریجی بزرگ‌تر می‌شود تا آنکه پس از چند ثانیه، مسیرها کاملاً از یکدیگر جدا می‌شوند.

نمودار ۱: سری زمانی زوایای θ₁ و θ₂
شکل ۱: رفتار غیرتناوبی و پیچیده زوایا در آونگ دوگانه با پارامترهای پیش‌فرض. نوسانات نامنظم نشانهٔ ماهیت آشوب‌ناک سیستم است.

۵.۳. فضای فاز

فضای فاز، تصویری هندسی از رفتار سیستم ارائه می‌دهد. در آونگ دوگانه، مسیر در فضای فاز به‌جای یک منحنی بسته (مانند آونگ ساده)، یک ساختار پیچیده و خودمتشابه را ترسیم می‌کند که نشان‌دهندهٔ جاذب فراکتالی است.

فضای فاز آونگ ۱
θ₁ در برابر ω₁
فضای فاز آونگ ۲
θ₂ در برابر ω₂

۵.۴. پایداری انرژی و دقت عددی

یکی از معیارهای کلیدی برای ارزیابی کیفیت شبیه‌سازی، پایداری انرژی کل سیستم است. در غیاب میرایی، انرژی کل باید ثابت بماند؛ اما خطاهای عددی باعث نوسان یا دریفت انرژی می‌شوند.

نمودار ۲: خطای نسبی انرژی ΔE/E₀ در طول شبیه‌سازی
شکل ۲: با انتگرال‌گیر RK4 و گام زمانی ۲×۱۰⁻⁴ ثانیه، خطای نسبی انرژی در بازهٔ کوتاه‌مدت کمتر از ۰٫۰۱٪ باقی می‌ماند. برای شبیه‌سازی‌های بلندمدت، استفاده از انتگرال‌گیر ورلهٔ سرعتی توصیه می‌شود.

۵.۵. واگرایی نمایی و نمای لیاپانوف

نمودار واگرایی لگاریتمی، رشد نمایی فاصلهٔ فاز در ناحیهٔ خطی و اشباع آن در مقادیر بزرگ را نشان می‌دهد. شیب ناحیهٔ خطی، نمای لیاپانوف $\lambda$ را تعیین می‌کند.

نمودار ۳: واگرایی لگاریتمی و برازش نمای لیاپانوف
شکل ۳: رشد نمایی ln d(t) و برازش خطی در ناحیهٔ میانی. شیب خط، نمای لیاپانوف λ ≈ ۳ تا ۵ رادیان بر ثانیه است که زمان لیاپانوف t_λ ≈ ۰٫۲ تا ۰٫۳ ثانیه را نتیجه می‌دهد.

۵.۶. بعد همبستگی

نمودار تابع همبستگی در مختصات لگاریتمی، خط مستقیمی را نشان می‌دهد که شیب آن بعد همبستگی $D_2$ است. مقدار غیرصحیح $D_2$، ساختار فراکتالی جاذب را تأیید می‌کند.

نمودار ۴: تابع همبستگی گراسبرگر-پروکاچیا
شکل ۴: تابع همبستگی C(r) در مختصات لگاریتمی برای دو آونگ. شیب ناحیهٔ خطی، بعد همبستگی D₂ ≈ ۱٫۵ تا ۲٫۵ را به‌دست می‌دهد که نشان‌دهندهٔ جاذب فراکتالی است.

۶. تحلیل و بحث

۶.۱. ماهیت آشوب در آونگ دوگانه

آونگ دوگانه نمونه‌ای کلاسیک از آشوب قطعی است: با وجود آنکه معادلات حرکت کاملاً قطعی‌اند و هیچ عنصر تصادفی در آن‌ها وجود ندارد، رفتار سیستم به‌طور ذاتی پیش‌بینی‌ناپذیر است. این پیش‌بینی‌ناپذیری، نه از ناکافی بودن دانش ما، بلکه از حساسیت ذاتی به شرایط اولیه ناشی می‌شود.

نکتهٔ مهم آن است که آشوب و قطعیت با هم ناسازگار نیستند. سیستم‌های آشوب‌ناک، قطعی‌اند؛ اما نه پیش‌بینی‌پذیر. این تمایز، در دهه‌های اخیر به یکی از بنیادی‌ترین مفاهیم در نظریهٔ سیستم‌های دینامیکی تبدیل شده است.

۶.۲. اهمیت انتخاب انتگرال‌گیر عددی

نتایج شبیه‌سازی نشان می‌دهد که انتخاب انتگرال‌گیر عددی تأثیر چشمگیری بر کیفیت شبیه‌سازی دارد. RK4 با وجود دقت مرتبهٔ چهارم در هر گام، در بازه‌های طولانی انرژی را به‌صورت سیستماتیک دریفت می‌دهد — زیرا خواص هندسی فضای فاز را حفظ نمی‌کند. ورلهٔ سرعتی، با وجود دقت مرتبهٔ دوم، در بازه‌های طولانی‌مدت انرژی را پایدارتر نگه می‌دارد. این نتیجه، نمونه‌ای از اهمیت ساختار هندسی الگوریتم‌های عددی در شبیه‌سازی سیستم‌های همیلتونی است.

✅ نتیجهٔ عملی

برای شبیه‌سازی‌های کوتاه‌مدت (تا چند ثانیه)، RK4 با گام زمانی کوچک مناسب است. برای شبیه‌سازی‌های بلندمدت (بیش از چند ده ثانیه)، ورلهٔ سرعتی انتخاب بهتری است. این توصیه، در تمام سیستم‌های همیلتونی — از مکانیک سماوی تا دینامیک مولکولی — معتبر است.

۶.۳. کاربردهای بینارشته‌ای

آونگ دوگانه، فراتر از یک مسئلهٔ مکانیکی ساده، به‌عنوان یک سیستم مدل برای مطالعهٔ آشوب در حوزه‌های مختلف به‌کار می‌رود:

۶.۴. محدودیت‌ها و جهت‌های آینده

شبیه‌ساز حاضر دارای چند محدودیت است: نخست، میرایی خطی ساده‌شده‌ای برای گشتاورهای اصطکاکی در نظر گرفته شده که ممکن است رفتار واقعی سیستم را به‌طور کامل توصیف نکند. دوم، اثرات اصطکاک در مفصل‌ها و انعطاف‌پذیری میله‌ها نادیده گرفته شده است. سوم، محاسبهٔ نمای لیاپانوف با روش‌های عددی ساده انجام می‌شود که دقت محدودی دارد.

جهت‌های آیندهٔ توسعه شامل: استفاده از الگوریتم‌های تفکیک‌کنندهٔ خودکار برای محاسبهٔ دقیق نمای لیاپانوف، افزودن مدل‌های اصطکاک غیرخطی، و توسعهٔ الگوریتم‌های یادگیری ماشین برای پیش‌بینی رفتار کوتاه‌مدت آشوب است.

۷. نتیجه‌گیری

این مقاله چارچوبی جامع برای تحلیل دینامیک آشوب‌ناک، کمی‌سازی حساسیت به شرایط اولیه و ارزیابی پایداری عددی در سیستم آونگ دوگانه ارائه کرد. معادلات حرکت با فرمول‌بندی لاگرانژین و با در نظر گرفتن گشتاورهای میرایی خطی استخراج شدند و سه انتگرال‌گیر عددی — RK4، ورلهٔ سرعتی و اویلر — از منظر دقت و پایداری مقایسه گردیدند.

برای کمی‌سازی آشوب، دو مسیر سایه با اختلاف‌های اولیهٔ کوچک ایجاد و واگرایی نمایی آن‌ها در فضای فاز ردیابی شد. بزرگ‌ترین نمای لیاپانوف با برازش خطی روی لگاریتم واگرایی محاسبه شد و زمان لیاپانوف $t_\lambda = 1/\lambda$ به‌عنوان معیار افق پیش‌بینی‌پذیری معرفی گردید. نتایج نشان داد که سیستم آونگ دوگانه در بازهٔ چند ثانیه از محدودهٔ پیش‌بینی‌پذیری فراتر می‌رود و واگرایی نمایی مسیرها با نمای لیاپانوف مثبت به‌وقوع می‌پیوندد.

علاوه بر این، بعد همبستگی جاذب آشوب‌ناک با الگوریتم گراسبرگر-پروکاچیا محاسبه شد که مقدار غیرصحیح آن، ساختار فراکتالی جاذب را تأیید کرد. در پایان، نتایج در قالب چهار نمودار تخصصی — سری زمانی، فضای فاز، انرژی و واگرایی — ارائه گردید.

پرسش‌های باقی‌مانده برای پژوهش‌های آینده عبارتند از: چگونه می‌توان نمای لیاپانوف را با دقت بالاتر و با استفاده از الگوریتم‌های تفکیک‌کنندهٔ خودکار محاسبه کرد؟ آیا می‌توان از شبکه‌های عصبی برای پیش‌بینی کوتاه‌مدت رفتار آشوب‌ناک استفاده کرد؟ و چگونه می‌توان چارچوب مشابهی را به سیستم‌های چندآونگی با درجات آزادی بیشتر تعمیم داد؟