طريقة الحجم المحدود للتدفق غير المستقر

في ميكانيكا الموائع ، يُعرف التدفق غير المستقر بأنه انتقال مائع تتغير خصائصه مع الزمن . ينتمي هذا النوع من الظواهر إلى نموذج ميكانيكا الأوساط المتصلة ، مما يعني أن المعادلات الحاكمة تتضمن فروقًا متناهية الصغر ، أي تفاضلات دوال في متغيرات مكانية وزمنية متصلة.

تُعد طريقة الحجم المحدود نهجًا في ديناميكا الموائع الحسابية حيث يمكن حل النظام عدديًا من خلال تجزئة الهندسة المحاكاة ، مما يسمح للحساب التدريجي للحل بمراعاة المشتقات الزمنية لخصائص المائع أيضًا . ^{[ 1 ]}

تأخذ معادلة الحفظ لنقل كمية قياسية في التدفق غير المستقر الشكل العام ^{[ 2 ].}

${\displaystyle {\frac {\partial \rho \phi }{\partial t}}+\operatorname {div} \left(\rho \phi \upsilon \right)=\operatorname {div} \left(\Gamma \operatorname {grad} \phi \right)+S_{\phi },}$

أين${\displaystyle \rho }$الكثافة و​${\displaystyle \phi }$هو شكل محافظ لجميع تدفق السوائل؛${\displaystyle \Gamma }$هو معامل الانتشار و${\displaystyle S}$هو مصطلح المصدر. معدل التدفق الصافي لـ${\displaystyle \phi }$عنصر خارج المائع ( الحمل الحراري )${\displaystyle \operatorname {div} \left(\rho \phi \upsilon \right),}$و${\displaystyle \operatorname {div} \left(\Gamma \operatorname {grad} \phi \right)}$معدل الزيادة${\displaystyle \phi }$بسبب الانتشار ؛${\displaystyle S_{\phi }}$معدل الزيادة${\displaystyle \phi }$بسبب المصادر و${\displaystyle {\frac {\partial \rho \phi }{\partial t}}}$معدل الزيادة${\displaystyle \phi }$عنصر سائل (عابر).

يعكس الحد الأول من المعادلة عدم استقرار التدفق، وهو غائب في حالة التدفقات المستقرة. يتم إجراء التكامل باستخدام طريقة الحجم المحدود للمعادلة الحاكمة على حجم تحكم، وكذلك على خطوة زمنية محدودة ∆t.

$\int \limits _{cv} \int _{t}^{t+\Delta t} \left( \frac {\partial \rho \phi }{\partial t}} \,\mathrm {d} t \right) \,\mathrm {d} V + \int _{t}^{t+\Delta t} \int \limits _{A} \left( n \cdot {\rho \phi u} \,\mathrm {d} A \right) \,\mathrm {d} t = \int _{t}^{t+\Delta t} \int \limits _{A} \left( n \cdot \left(\Gamma \operatorname {grad} \phi \right) \,\mathrm {d} A \right) \,\mathrm {d} t+\int _{t}^{t+\Delta t}\!\!\!\int \limits _{cv}S_{\phi }\,\mathrm {d} V\,\mathrm {d} t}$

يُشابه تكامل حجم التحكم للجزء المستقر من المعادلة تكامل معادلة الحالة المستقرة . لذا، ينبغي التركيز على تكامل المكون غير المستقر من المعادلة. لفهم أسلوب التكامل، نرجع إلى معادلة التوصيل الحراري غير المستقرة أحادية البعد . ^{[ 3 ]}

${\displaystyle \rho c{\frac {\partial T}{\partial t}}={\frac {\partial }{\partial x}}{\frac {k\partial T}{\partial x}}+S}$

${\displaystyle \int _{t}^{t+\Delta t}\!\!\!\int \limits _{cv}\rho c{\frac {\partial T}{\partial t}}\,\mathrm {d} V\,\mathrm {d} t=\int _{t}^{t+\Delta t}\!\!\!\int \limits _{cv}{\frac {\partial }{\partial x}}{\frac {k\partial T}{\partial x}}\,\mathrm {d} V\,\mathrm {d} t+\int _{t}^{t+\Delta t}\!\!\!\int \limits _{cv}S\,\mathrm {d} V\,\mathrm {d} t}$

${\displaystyle \int _{e}^{w}\!\!\!\int _{t}^{t+\Delta t}\left(\rho c{\frac {\partial T}{\partial t}}\,\mathrm {d} t\right)\,\mathrm {d} V=\int _{t}^{t+\Delta t}\left[\left(kA{\frac {\partial T}{\partial x}}\right)_{e}-\left(kA{\frac {\partial T}{\partial x}}\right)_{w}\right]\,\mathrm {d} t+\int _{t}^{t+\Delta t}{\bar {S}}\Delta V\,\mathrm {d} t}$

الآن، بافتراض أن درجة الحرارة عند العقدة سائدة في حجم التحكم بأكمله، يمكن كتابة الجانب الأيسر من المعادلة على النحو التالي ^{[ 4 ].}

${\displaystyle \int \limits _{cv}\!\!\!\int _{t}^{t+\Delta t}\left(\rho c{\frac {\partial T}{\partial t}}\,\mathrm {d} t\right)\,\mathrm {d} V=\rho c\left(T_{P}-{T_{P}}^{O}\right)\Delta V}$

باستخدام طريقة الفروق الخلفية من الدرجة الأولى ، يمكننا كتابة الطرف الأيمن من المعادلة على النحو التالي:

${\displaystyle \rho c\left(T_{P}-{T_{P}}^{0}\right)\Delta V=\int _{t}^{t+\Delta t}\left[\left(K_{e}A{\frac {T_{E}-T_{P}}{\delta x_{PE}}}\right)-\left(K_{w}A{\frac {T_{P}-T_{W}}{\delta x_{WP}}}\right)\right]\,\mathrm {d} t+\int _{t}^{t+\Delta t}{\bar {S}}\Delta V\,\mathrm {d} t}$

لتقييم الجانب الأيمن من المعادلة، نستخدم الآن معامل ترجيح.${\displaystyle \theta }$بين 0 و1، ونكتب تكامل${\displaystyle T_{P}}$

${\displaystyle I_{T}=\int _{t}^{t+\Delta t}T_{P}\,\mathrm {d} t=\left[\theta T_{P}+\left(1-\theta \right){T_{P}}^{0}\right]\Delta t}$

الآن، يعتمد الشكل الدقيق للمعادلة المتقطعة النهائية على قيمة${\displaystyle \Theta }$. باعتباره تباين${\displaystyle \Theta }$يكون${\displaystyle 0<\Theta <1}$، المخطط الذي سيتم استخدامه للحساب${\displaystyle T_{P}}$يعتمد ذلك على قيمة${\displaystyle \Theta .}$هكذا

${\displaystyle \rho c{\frac {\left(T_{P}-{T_{P}}^{0}\right)}{\Delta t}}\Delta x=\theta \left[\left(K_{e}{\frac {T_{E}-T_{P}}{\delta x_{PE}}}\right)-\left(K_{w}{\frac {T_{P}-T_{W}}{\delta x_{WP}}}\right)\right]+(1-\theta )\left[\left(K_{e}{\frac {T_{E}-T_{P}}{\delta x_{PE}}}\right)-\left(K_{w}{\frac {T_{P}-T_{W}}{\delta x_{WP}}}\right)\right]+{\bar {S}}\Delta x}$

• المخطط الصريح : في المخطط الصريح، يتم تحويل حد المصدر إلى خطي كما يلي:${\displaystyle b=S_{u}+{S_{P}}{T_{P}}^{0}}$نستبدل${\displaystyle \theta =0}$للحصول على التقطيع الصريح، أي: ^{[ 5 ]}

${\displaystyle a_{P}T_{P}=a_{w}{T_{w}}^{0}+a_{e}{T_{e}}^{0}+\left[{a_{P}}^{0}-\left(a_{w}+a_{e}-S_{P}\right)\right]{T_{P}}^{0}+S_{u}}$

أين${\displaystyle a_{P}={a_{P}}^{0}}$تجدر الإشارة إلى أن الطرف الأيمن يحتوي على قيم من الخطوة الزمنية السابقة، وبالتالي يمكن حساب الطرف الأيسر عن طريق المطابقة الأمامية في الزمن. تعتمد هذه الطريقة على التفاضل الخلفي، وخطأ اقتطاع متسلسلة تايلور فيها من الدرجة الأولى بالنسبة للزمن. يجب أن تكون جميع المعاملات موجبة. بالنسبة لقيمة k ثابتة وتباعد شبكي منتظم،${\displaystyle \delta x_{PE}=\delta x_{WP}=\Delta x}$يمكن كتابة هذا الشرط على النحو التالي:

${\displaystyle \rho c{\frac {\Delta x}{\Delta t}}>{\frac {2K}{\Delta x}}}$

تفرض هذه المتباينة شرطًا صارمًا على الحد الأقصى للخطوة الزمنية التي يمكن استخدامها، وتمثل قيدًا خطيرًا على المخطط. يصبح تحسين الدقة المكانية مكلفًا للغاية لأن الحد الأقصى للخطوة الزمنية الممكنة يجب تقليله بمربع${\displaystyle \Delta x}$^{[ 6 ]}

• مخطط كرانك-نيكلسون  : تنتج طريقة كرانك-نيكلسون من تحديد${\displaystyle \theta ={\frac {1}{2}}}$تصبح معادلة التوصيل الحراري غير المستقرة المتقطعة

${\displaystyle a_{P}T_{P}=a_{E}\left[{\frac {T_{E}+{T_{E}}^{0}}{2}}\right]+a_{W}\left[{\frac {T_{W}+{T_{W}}^{0}}{2}}\right]+\left[{a_{P}}^{0}-{\frac {a_{E}}{2}}-{\frac {a_{W}}{2}}\right]{T_{P}}^{0}+b,}$

أين${\displaystyle a_{P}={\frac {a_{W}+a_{E}}{2}}+{a_{P}}^{0}-{\frac {S_{P}}{2}}.}$

بما أن المعادلة تحتوي على أكثر من قيمة مجهولة لـ T عند المستوى الزمني الجديد، فإن الطريقة ضمنية، ويجب حل المعادلات الآنية لجميع نقاط العقد في كل خطوة زمنية. على الرغم من أن المخططات ذات${\displaystyle {\frac {1}{2}}<\theta <1}$بما في ذلك مخطط كرانك-نيكلسون، فإنّ الاستقرار المطلق لجميع قيم خطوة الزمن أمرٌ بالغ الأهمية، ومن الأهمية بمكان ضمان أن تكون جميع المعاملات موجبة للحصول على نتائج واقعية ومحدودة من الناحية الفيزيائية. ويتحقق ذلك إذا كان معامل${\displaystyle {T_{P}}^{0}}$يستوفي الشرط التالي

${\displaystyle {a_{P}}^{0}=\left[{\frac {a_{E}+a_{W}}{2}}\right],}$

مما يؤدي إلى

${\displaystyle \Delta t<\rho c{\frac {\Delta x^{2}}{K}},}$

وبالتالي، تعتمد طريقة كرانك-نيكلسون على التفاضل المركزي، ومن ثم فهي دقيقة من الدرجة الثانية زمنيًا. كما تعتمد الدقة الإجمالية للحساب على ممارسة التفاضل المكاني، لذا تُستخدم طريقة كرانك-نيكلسون عادةً بالتزامن مع التفاضل المركزي المكاني.

3. المخطط الضمني الكامل : عندما تُضبط قيمة Ѳ على 1، نحصل على المخطط الضمني الكامل. المعادلة المتقطعة هي: ^{[ 7 ]}

${\displaystyle a_{P}T_{P}=a_{W}T_{W}+a_{E}T_{E}+{a_{P}}^{0}{T_{P}}^{0}+S_{u}}$

${\displaystyle a_{P}={a_{P}}^{0}+a_{W}+a_{E}-S_{P}}$

يحتوي كلا طرفي المعادلة على درجات الحرارة عند الخطوة الزمنية الجديدة، ويجب حل نظام من المعادلات الجبرية عند كل مستوى زمني. تبدأ عملية التقدم الزمني بحقل أولي معين من درجات الحرارة. ${\displaystyle T^{0}}$يتم حل نظام المعادلات بعد اختيار خطوة زمنية${\displaystyle \Delta t}$ثم الحل${\displaystyle T}$يُسند إلى${\displaystyle T^{0}}$وتُكرر هذه العملية لتقدم الحل بخطوة زمنية إضافية. يتضح أن جميع المعاملات موجبة، مما يجعل الطريقة الضمنية مستقرة بشكل مطلق لأي حجم خطوة زمنية. ولأن دقة هذه الطريقة من الدرجة الأولى فقط بالنسبة للزمن، فإن الخطوات الزمنية الصغيرة ضرورية لضمان دقة النتائج. يُوصى باستخدام الطريقة الضمنية في الحسابات العابرة للأغراض العامة نظرًا لمتانتها واستقرارها المطلق.