ثلاث طرق، ثلاث مناطق
بالنسبة لمعادلة الاختبار dy/dx = λy، هناك ثلاث طرق صريحة لمحاكاة المعادلات التفاضلية بمناطق استقرار محددة في المستوى المعقد hλ:
طريقة أويلر (من الدرجة الأولى): منطقة الاستقرار هي القرص |1 + hλ| ≤ 1، دائرة بنصف قطر 1 مركزها عند (-1، 0). يجب أن يقع hλ الحقيقي السالب في [-2، 0].
Runge-Kutta 2 (طريقة نقطة المنتصف) (من الدرجة الثانية): منطقة الاستقرار هي |1 + hλ + (hλ)²/2| ≤ 1. أكبر من قرص أويلر، لكنها محدودة دائماً.
Runge-Kutta 4 (من الدرجة الرابعة): منطقة الاستقرار تحقق |1 + hλ + (hλ)²/2 + (hλ)³/6 + (hλ)⁴/24| ≤ 1. يمتد hλ الحقيقي السالب إلى تقريباً -2.785. المنطقة أكبر بشكل كبير من قرص أويلر.
Backward Euler (ضمني): منطقة الاستقرار هي المستوى المعقد كاملاً باستثناء القرص |1 - hλ|⁻¹ > 1، بشكل مكافئ |1/(1-hλ)| ≤ 1. بالنسبة لـ λ في النصف الأيسر من المستوى (Re(λ) < 0)، هذا مستقر بشكل غير مشروط — لا قيد على h من ناحية الاستقرار.
دالة التضخيم
بالنسبة لأي طريقة Runge-Kutta، معامل التضخيم لكل خطوة R(hλ) هو تقريب متعدد الحدود لـ e^(hλ):
- أويلر: R(z) = 1 + z (مختزل عند الدرجة 1)
- RK2: R(z) = 1 + z + z²/2 (مختزل عند الدرجة 2)
- RK4: R(z) = 1 + z + z²/2 + z³/6 + z⁴/24 (مختزل عند الدرجة 4)
منطقة الاستقرار هي {z ∈ ℂ : |R(z)| ≤ 1}. تضخيم الحل الحقيقي: |e^z| = e^(Re(z)). بالنسبة لـ Re(z) < 0 (معادلة تفاضلية مستقرة)، الحل الحقيقي يتضاءل. الطريقة الرقمية مستقرة إذا |R(z)| ≤ 1 — مطابقة السلوك المتضاءل.
القيم الذاتية النقية الخيالية: الأنظمة التذبذبية
العديد من الأنظمة الفيزيائية لها قيم ذاتية نقية خيالية: λ = iω (تذبذبات بدون تخميد). نظام الكتلة و النابض، ميكانيكا الحركة المدارية، ديناميكا البندول.
بالنسبة لـ λ = iω: hλ = ihω يقع على المحور الخيالي.
استقرار أويلر على المحور الخيالي: |1 + ihω|² = 1 + (hω)² > 1 لأي h > 0. أويلر غير مستقر لأي حجم خطوة على القيم الذاتية النقية الخيالية. التذبذب 'المحسوب' ينمو بدون حد.
استقرار RK4 على المحور الخيالي: منطقة الاستقرار تمتد إلى تقريباً |hω| ≤ 2.83 على المحور الخيالي. بالنسبة لـ h صغيرة بما فيه الكفاية، يتعامل RK4 مع التذبذبات غير المخمدة. لا يمكن لأويلر أن يفعل ذلك.
هذه الهندسة هي السبب في أن أويلر يفشل على الأنظمة المحافظة (كتلة-نابض، مدارات، معادلات الموجة) حتى مع h صغيرة، بينما يتعامل RK4 معها بشكل جيد.
الهندسة الجامدة للمشاكل
نظام معادلات تفاضلية جامد له قيم ذاتية بحجم مختلف جداً. نسبة الصلابة: κ = max|Re(λᵢ)| / min|Re(λᵢ)| >> 1.
لماذا الصلابة مكلفة للمحاكيات الصريحة:
الاستقرار يتطلب h·max|λᵢ| ≤ C (حيث C يعتمد على الطريقة). القيمة الذاتية الأكثر سلبية تحدد الحد.
الدقة لديناميات بطيئة تتطلب h·min|λᵢ| ≥ ε (حل الوضع الأبطأ بشكل كافٍ).
إذا كانت κ كبيرة، فإن هذين المتطلبين يفرضان h صغيراً: صغيراً بما يكفي لاستقرار الوضع السريع، كبيراً بما يكفي لأخذ عينات الوضع البطيء. عدد الخطوات يتناسب مع κ.
الصورة الهندسية في طيف القيم الذاتية: القيم الذاتية لمصفوفة جاكوبيان ∂f/∂y تشكل مجموعة من النقاط في المستوى المعقد. منطقة الاستقرار في المحاكاة الصريحة يجب أن تحتوي على جميع النقاط h·λᵢ. إذا امتدت القيم الذاتية من -1 إلى -1000، يجب على منطقة الاستقرار أن تغطي نطاقاً من 1000 على طول المحور الحقيقي — مما يتطلب h صغيراً جداً.
المحاكيات الضمنية: منطقة استقرار Backward Euler تغطي النصف الأيسر من المستوى كاملاً. جميع القيم الذاتية بـ Re(λ) < 0 داخل منطقة الاستقرار تلقائياً بغض النظر عن h. القيد على h يأتي فقط من الدقة، و ليس من الاستقرار.
نسبة الصلابة و التكلفة
فكر في شبكة تفاعل كيميائي بتفاعلات سريعة (مقياس زمني 10⁻⁶ s) و تفاعلات بطيئة (مقياس زمني 1 s).
نسبة الصلابة: κ = 10⁶ / 1 = 10⁶.
مع RK4 (حد الاستقرار h·|λ_max| ≤ 2.785): h_max = 2.785 / 10⁶ ≈ 2.8 × 10⁻⁶ s.
للمكاملة على 10 s من وقت التفاعل: الخطوات = 10 / (2.8 × 10⁻⁶) ≈ 3.6 × 10⁶.
مع Backward Euler (مستقر بشكل غير مشروط): يمكن اختيار h لدقة التفاعلات البطيئة. h = 10⁻² s (100 عينة على مقياس 1 s). الخطوات = 10 / 10⁻² = 1000.
نسبة التكلفة: 3.6 مليون خطوة صريحة مقابل 1000 خطوة ضمنية — عامل 3600. كل خطوة ضمنية تتطلب حل نظام خطي (التكلفة لكل خطوة أعلى)، لكن التكلفة الإجمالية أقل بكثير للمشاكل الجامدة جداً.
لماذا أنابيب الأبعاد n ليست ما تعتقد
في 2D، 'الأنبوب' بنصف قطر ε حول منحنى C هو مجموعة النقاط على بعد ε من C. القطع العرضي هو دائرة بنصف قطر ε. يزداد حجم الأنبوب بالتناسب مع طوله.
في n أبعاد، تتغير هندسة الأنبوب بشكل أساسي، بسبب الظاهرة من الفصل 9:
مفارقة الزاوية n-الأبعادية: في فضاء n-الأبعادي، يقع معظم حجم مكعب n-الأبعادي في الزوايا — و ليس في المنطقة المركزية. مع زيادة n، يذهب جزء الحجم ضمن مسافة ε من المركز إلى الصفر لأي ε ثابت.
تطبيقاً على أنابيب حل المعادلات التفاضلية:
في 2D: إذا كان الحل الحقيقي يمر عبر مركز الأنبوب، معظم النقاط القريبة قريبة من المنحنى. الاضطرابات الصغيرة تبقيك قرب الحل الحقيقي.
في أبعاد عالية: معظم النقاط ضمن صندوق الأنبوب الحدودي بعيدة فعلاً من منحنى الحل الحقيقي. 'حجم' الأنبوب يهيمن عليه الزوايا — مناطق بعيدة عن المركز في عدة أبعاد في نفس الوقت.
النتيجة للمحاكاة: مع 28 معادلة تفاضلية مقترنة (مسألة اعتراض البحرية في هامينج)، اضطراب بحجم ε في كل بُعد يمكن أن ينتج إزاحة إجمالية بقدر ε√28 ≈ 5.3ε من الحل الحقيقي. يجب فهم الأنبوب من حيث معيار L2 عبر جميع الأبعاد، و ليس فقط الحد الأقصى للإزاحة في أي بُعد واحد.
الاستقرار في أبعاد عالية: نظام حيث كل مكون يتضاءل بشكل مستقل (كل قيمة ذاتية لها جزء حقيقي سالب) قد لا تزال تظهر إزاحات مدمجة كبيرة لأن أخطاء المكونات تضيف في معيار L2. الأنبوب 28-الأبعادي ليس فقط 28 أنبوب 1-الأبعادي مستقل — الهندسة تقترنهم.
من الهندسة إلى التصميم
الرؤى الهندسية من الفصول 18-20 تتجمع معاً كمجموعة من مبادئ التصميم لمحاكاة رقمية:
اختيار حجم الخطوة: h يجب أن يضع h·λ داخل منطقة الاستقرار لكل قيمة ذاتية. بالنسبة للأنظمة الجامدة، تزيل الطرق الضمنية قيد الاستقرار، مما يترك فقط متطلبات الدقة.
تراكم الخطأ في أبعاد عالية: الخطأ العام هو متجه في فضاء n-الأبعادي. ينمو معياره كـ √n مرات الخطأ لكل مكون. محاكاة عالية الأبعاد تحتاج متطلبات دقة أكثر صرامة لكل خطوة.
التغذية الراجعة كمثبت: إذا كانت المحاكاة تتضمن تغذية راجعة (الناتج المحسوب يؤثر على المدخلات اللاحقة، كما في نظام التوجيه)، التغذية الراجعة المتقاربة تخفف الأخطاء. يمكن للمحاكاة أن تتسامح مع مدخلات غير دقيقة للكميات داخل حلقة التغذية الراجعة.
عدم الاستقرار كإشارة: للمشاكل بحقول اتجاهات متباعدة، عدم الاستقرار يمكن استغلاله: اتجاه التباعد يحمل معلومات حول خطأ الشرط الأولي، مما يمكّن التعديل التصحيحي.