حل عددی معادلات دیفرانسیل جزئی حاکم بر میدان مغناطیسی در جریان روی هندسه نوک تیز در ماخ و ارتفاعات بالا و اثر آن بر ضرایب دراگ و لیفت
محورهای موضوعی : آمارسید میرالله حسینی متی کلایی 1 , محمد حاتمی 2 * , عزیز وظیفه شناس 3
1 - گروه مهندسی مکانیک، دانشکده مهندسی، دانشگاه آزاد اسلامی، واحد قائمشهر، قائمشهر، ایران
2 - گروه مهندسی مکانیک، دانشکده مهندسی، دانشگاه فردوسی مشهد، مشهد، ایران
3 - گروه مهندسی مکانیک، دانشکده مهندسی، دانشگاه آزاد اسلامی، واحد اسفراین، اسفراین، ایران
کلید واژه: Magnetohydrodynamic, Lift coefficient, Drag coefficient, Stall angle, Mach number,
چکیده مقاله :
در این مقاله به حل عددی اثر جریان مگنتو هیدرودینامیک بر روی دو هندسه آیرودینامیکی، بصورت دو بعدی و سه بعدی پرداخته شده است. نتایج برای دو هندسه دو و سه بعدی که تا حدودی سعی شده است شبیه به بال موشک باشد، در ارتفاعات بالا که فشار زیاد و دما کم می باشد و همچنین ماخ های 6 و 8 در 9 زاویه حمله مختلف، در دو حالت بدون مگنت و با مگنت ضرایب لیفت و دراگ بدست آمده است. در پایان مشاهده شد که با اضافه کردن مگنت به مسئله میزان ضریب لیفت افزایش پیدا می کند که بیشترین درصد افزایش ضریب لیفت برای هندسه سه بعدی و شرایط ماخ 8 و ارتفاع 50000 متر رخ داده است که 77.5% می باشد. همچنین با مقایسه دو هندسه سه بعدی و دو بعدی مشاهده شد که زاویه استال در هندسه دو بعدی در حالت ارتفاع 9000 و ماخ 6 رخ نداده اما در هندسه سه بعدی و همین شرایط حل زاویه استال برای حالت بدون MHD اتفاق افتاده، که علت این امر را می توان در خط جریانی بودن هندسه دو بعدی دانست. البته در ادامه مشاهده شد که با اضافه کردن مگنت در هندسه سه بعدی و شرایط حل مذکور این زاویه به تاخیر افتاده است.
In this paper, the effect of magneto-hydrodynamic flow on two aerodynamic geometries (2D & 3D) is investigated. The results (Lift and drag coefficients) for two and three-dimensional geometries, which have been tried to be similar to rocket wings, at high altitudes where the pressure is high and the temperature low, as well as Machs at 6 and 8 and at 9 different angles of attacks, are obtained in two modes with and without magnets. At the end, it was observed that adding a magnet to the problem increases the lift coefficient which maximum increment (77.5%) occurred for 3D geometry at Mach 8 and 50000 m height. Also, comparing the two-dimensional and three-dimensional geometries, it was observed that the stall angle did not occur in the two-dimensional geometry at 9000 altitude and Mach 6, but in the three-dimensional geometry and the same conditions, the stall angle was observed for the non-MHD mode, which is due to flow line of two-dimensional geometry. However, it was further observed that this angle was delayed by adding a magnet to the 3D geometry with the mentioned solution conditions.
[1] S. J. Liao, The proposed homotopy analysis technique for the solution of nonlinear problems. Ph.D. Thesis, Shanghai Jiao Tong University (1992)
[2] J. X. Li, Y. Yan, W. Q. Wang, Time-delay feedback control of a cantilever beam with concentrated mass based on the homotopy analysis method, Applied Mathematical Modelling 108:629-645(2022)
[3] P. Khaneh Masjedi,P. M. Weaver, Analytical solution for arbitrary large deflection of geometrically exact beams using the homotopy analysis method, Applied Mathematical Modelling 103: 516-542(2022)
[4] E. Botton, J.B. Greenberg, A. Arad, D. Katoshevski, V. Vaikuntanathan, M. Ibach, B. Weigand, An investigation of grouping of two falling dissimilar droplets using the homotopy analysis method, Applied Mathematical Modelling 104:486-498(2022)
[5] S. Abbasbandy, T. Hayat, On series solution for unsteady boundary layer equations in a special third grade fluid. Communications in Nonlinear Science and Numerical Simulation 16:3140-3146(2011)
[6] S. Abbasbandy, E. Shivanian, Predictor homotopy analysis method and its application to some nonlinear problems. Communications in Nonlinear Science and Numerical Simulation 16:2456-2468(2011)
[7] R. A. Van Gorder, K. Vajravelu, On the selection of auxiliary functions, operators, and convergence control parameters in the application of the Homotopy Analysis Method to nonlinear differential equations: A general approach. Communications in Nonlinear Science and Numerical Simulation 14:4078-4089(2009)
[8] Chuang Yu, Hui Wang, Dongfang Fang, Jianjun Ma, Xiaoqing Cai, Xiaoniu Yu, Semi-analytical solution to one-dimensional advective-dispersive-reactive transport equation using homotopy analysis method. Journal of Hydrology 565:422–428(2018)
[9] M.R. Shirkhani, H.A. Hoshyar, I. Rahimipetroudi, H.Akhavan, D.D.Ganji, Unsteady time-dependent incompressible Newtonian fluid flow between two parallel plates by homotopy analysis method (HAM), homotopy perturbation method (HPM) and collocation method (CM). Propulsion and Power Research 7:247-256(2018)
[10] S. Hussain, A. Shah, S. Ayub, A. Ullah, An approximate analytical solution of the Allen-Cahn equation using homotopy perturbation method and homotopy analysis method. Heliyon 5(12):1-9(2019)
[11] P.A. Naik, J. Zu and M. Ghoreishi, Estimating the approximate analytical solution of HIV viral dynamic model by using homotopy analysis method. Chaos, Solitons and Fractals, 131:1-21(2020)
[12] G. Zhang, Z. Wu, Homotopy analysis method for approximations of Duffing oscillator with dual frequency excitations. Chaos, Solitons and Fractals 127:342–353(2019)
[13] Y. Zhang, Y. Li, Nonlinear dynamic analysis of a double curvature honeycomb sandwich shell with simply supported boundaries by the homotopy analysis method. Composite Structures 221:1-12(2019)
[14] A. Jafarimoghaddam, On the Homotopy Analysis Method (HAM) and Homotopy Perturbation Method (HPM) for a nonlinearly stretching sheet flow of Eyring-Powell Fluids. Engineering Science and Technology, an International Journal 22:439-451(2019)
[15] Z. Odibat, On the optimal selection of the linear operator and the initial approximation in the application of the homotopy analysis method to nonlinear fractional differential equations. Applied Numerical Mathematics 137:203–212(2019)
[16] J. Rana, S. Liao, On time independent Schrödinger equations in quantum mechanics by the homotopy analysis method. Theoretical & Applied Mechanics Letters 9:376-381(2019)
[17] S.J. Liao, Beyond Perturbation: Introduction to the Homotopy Analysis Method, Chapman & Hall/CRC Press, Boca Raton (2003)
دسترسي در سايتِ http://jnrm.srbiau.ac.ir
سال نهم، شماره چهل و یکم، فروردین و اردیبهشت 1402
|
حل عددی معادلات دیفرانسیل جزئی حاکم بر میدان مغناطیسی در جریان روی هندسه نوک تیز در ماخ و ارتفاعات بالا و اثر آن بر ضرایب دراگ و لیفت
سید میرالله حسینی متی کلایی1، محمد حاتمی21، عزیز وظیفه شناس3
(1) گروه مهندسی مکانیک، دانشکده مهندسی، واحد قائمشهر، دانشگاه آزاد اسلامي، قائمشهر، ايران
(2) گروه مهندسی مکانیک، دانشکده مهندسی، دانشگاه فردوسی مشهد، مشهد، ایران
(3) گروه مهندسی مکانیک، دانشکده مهندسی، واحد اسفراین، دانشگاه آزاد اسلامي، اسفراین، ايران
تاريخ ارسال مقاله: 01/01/1400 تاريخ پذيرش مقاله: 16/06/1400
چکيده
در این مقاله به حل عددی اثر جریان مگنتو هیدرودینامیک بر روی دو هندسه آیرودینامیکی، بصورت دو بعدی و سه بعدی پرداخته شده است. نتایج برای دو هندسه دو و سه بعدی که تا حدودی سعی شده است شبیه به بال موشک باشد، در ارتفاعات بالا که فشار زیاد و دما کم می باشد و همچنین ماخ های 6 و 8 در 9 زاویه حمله مختلف، در دو حالت بدون مگنت و با مگنت ضرایب لیفت و دراگ بدست آمده است. در پایان مشاهده شد که با اضافه کردن مگنت به مسئله میزان ضریب لیفت افزایش پیدا می کند که بیشترین درصد افزایش ضریب لیفت برای هندسه سه بعدی و شرایط ماخ 8 و ارتفاع 50000 متر رخ داده است که 77.5% می باشد. همچنین با مقایسه دو هندسه سه بعدی و دو بعدی مشاهده شد که زاویه استال در هندسه دو بعدی در حالت ارتفاع 9000 و ماخ 6 رخ نداده اما در هندسه سه بعدی و همین شرایط حل زاویه استال برای حالت بدون MHD اتفاق افتاده، که علت این امر را می توان در خط جریانی بودن هندسه دو بعدی دانست. البته در ادامه مشاهده شد که با اضافه کردن مگنت در هندسه سه بعدی و شرایط حل مذکور این زاویه به تاخیر افتاده است.
واژههاي کليدي: مگنتوهیدرودینامیک، ضریب دراگ، ضریب لیفت، ماخ، زاویه استال.
[1] . عهدهدار مکاتبات: Email:m.hatami@xjtu.edu.cn, m-hatami@um.ac.ir
1- مقدمه
مطالعه برهم کنش میدان الکتریکی و مغناطیسی با جریان سیال رسانا منگتو هیدرودینامیک یا MHD نامیده می شود. تحت شرایط معینی می توان یک گاز را به سیال هادی الکتریسیته تبدیل نمود. حال اگر یک منبع قدرت خارجی برای تولید میدان الکتریکی به کار گرفته شود، بر هم کنش جریان الکتریکی عمود بر جریان گاز و یک میدان مغناطیسی منجر به نیروی شتاب دهندهای در جهت حرکت جریان گاز می شود. از این نیروی شتاب دهنده می توان برای کاهش اثرات منفی گرادیان فشار معکوس و در نتیجه کاهش ناحيه جدایش جریان و کاهش ضریب پسا بر روی سازه های هوافضایی مثل ایرفویل استفاده کرد. همچنین می توان از آن در سیستم پیشرانش موشک برای دستیابی به سرعت زیاد گاز خروجی استفاده نمود [۱] برای این منظور ابتدا بایستی گاز را هادی الکتریسیته کرد. یک روش ساده افزایش دمای گاز است. در دمای بالا اتم ها به یونها با بار مثبت و الکترون ها تقسیم میشوند. گاز حاصل پلاسما نام دارد. به عبارت دیگر میتوان گفت که واژه پلاسما به گاز يوتيزه شده ای اطلاق می شود که همه یا بخش قابل توجهی از اتمهای آن یکیاچند الکترون از دست داده و به یونهای مثبت تبدیل شده باشند[۲]. برای تولید گازهادی الکتریسیته میتوان از درصد کمی از مواد با انرژی یونش پایین به عنوان بذر استفاده کرد[۱]. البته با توجه به شرایط جریان در ماخهای بالا، دسترسی به دمای بالا و در نتیجه یونیزه شدن سيال عبوری از روی ایرفویل دور از تصور نیست. عملکرد یک سیستم ام اچ دی چه برای تولید قدرت و چه برای شتاب دادن به یک جریان گاز، بستگی به توانایی در تبدیل یک گاز خنثی به گاز با هدایت الکتریکی قابل ملاحظه دارد. بنابراین بافتن موادی که در ماخ های بالا دارای يون زایی بالاتر باشند، مورد توجه است. عامل دیگر نیز جنسی بالک با قابلیت رسانایی متفاوت بوده که بایستی مورد بررسی قرار گیرد. تحقیقات زیادی در زمینه کاریرد مگنتو هیدرو دینامیک در سازه های هوافضایی صورت گرفته است. اما تمرکز تحقیقات انجام شده بر روی تیروی پیشرانش موتور می باشد و لزوم بررسی بیشتر جریان هیدرو دینامیک و اثر آن بر روی جریان خارجی روی سطوح آیرو دینامیکی مانند ایرفویل احساس میشود.
در سال 1904، پرانت مفهوم لايه مرزي را بيان كرد و به اين وسيله ارتباط مهمي ميان جريان سيال ايدهآل و جريان سيال واقعي به وجود آمد. براي سيالاتي كه ويسكوزيته نسبتاً كمي دارند، اثر اصطكاك داخلي در سيال تنها در ناحيه باريكي از محيط كه مرز سيال را تشكيل ميدهد، قابل توجه است و لذا خارج از اين ناحيه باريك در نزديكي مرزهاي جامد،بايد جريانرا ايدهآلدرنظر گرفت ]3[.
ژائو و همکاران ]4[ بصورت عددی پدیده الکترو هیدرودینامیک را در یک کانال مسطح بررسی نمودند.الگورتیم حلعددی آنها ترکیبیاز روشهای المان محدود، روش مشخصهها، و حجم محدود بوده است. نتایج آنها شبیهسازی الگوی جریان تحت تاثیر پدیده الکتروهیدرودینامیک و پارامترهای موثر بر این پدیده را ارائه میدهد. آنها دریافتند با افزایش عدد رینولدز، تاثیر پدیده الکتروهیدرودینامیک بر جریان سیال کاهش میباشد. علاوه بر این بزرگی گردابه ها در عدد رینولدز ثابت به ولتاژ اعمال شده بستگی دارد و در سرعت های پایین، سرعت جریان تاثیری بر پارامترهای الکتریکی نمیگذارد، اما در سرعتهای بالا بطور موثری بر میدان الکتریکی تاثیر گذار میباشد.
تادا و همکاران ]5[ با استفاده از الکترود های سیمی به موازات جریان اصلی، به انجام یکسری آزمایشات به منظور تاثیر میدان الکتریکی بر جریان های داخل کانال پرداختند. در بررسی آنها الکترود های سیمی به فاصله مساوی از یکدیگر قرار گرفته اند و بوسیله منبع مولد ولتاژ بالای مثبت تغذیه میشدند. نتایج آنها افزایش آهنگ نرخ انتقال حرارت جابجایی بویژه در محدوده جریان های آرام، در اثر ترکیب جریانهای ثانویه ناشی از میدان های الکتریکی و جریان اولیه را نشان داده است. هیونو چون ]6[ بصورت تجربی تاثیر جریان سیال بر گردابههای پشت سیلندر در داخل یک کانال را مورد مطالعه قرار دادند. آنها در کار تجربی خود الکترود تزریق کننده را در زوایای مختلف 45، 90، 135 و 180 و فواصل 1، 1.5، 3 سانتی متر از سطح سیلندر قرار دادند. نتایج آنها با آشکار سازی خطوط جریان نشان داد جریان یونی تاثیر چشمگیری بر ساختار گردابه ای پشت سیلندر می گذارد و نیروي بازدارندگی فشاري تحت تاثیر بار یونی 4 به طور قابل ملاحظه اي کاهش می یابد.
لگر و همکاران ]7[ بر روي جریان هوا در امتداد یک DC نیز به بررسی تجربی تأثیر تخلیه کروناي صفحه تخت شیبدار پرداختند. آشکارسازي میدان جریان و نیز تعیین توزیع نشان داد که بدون اعمال میدان الکتریکی، PIV سرعت با استفاده از فنآوري جریان در لبه جلویی صفحه از آن جدا شده و دنباله هاي قابل توجهی در بالاي صفحه شکل میگیرد. اما با اعمال میدان الکتریکی، خطوط جریان به سمت صفحه منحرف میشود که با توجه به عدد رینولدز جریان و شیب صفحه، به برخورد یا عدم برخورد جریان با صفحه منجر شده و همزمان از میزان گردابه ها کاسته می شود. همچنین آنها دریافتند که تأثیر تخلیه یونها بر مشخصه هاي جریان در اعداد رینولدز پایین در قالب کاهش دنباله هاي جریان و کاهش نیروي بازدارندگی، بیشتر بوده و با افزایش سرعت جریان، تاثیر تخلیه کمتر میشود.
روبرتووگویلرمو ]8[ به صورت تجربی جریان روي یک بالواره متقارن در رینولدز پایین را تحت تاثیر محرك الکتروهیدرودینامیکی بررسی نمودند. نتایج آنها حاکی از آن بود که در رینولدز هاي پایین تاثیر محرك وابسته به ولتاژ اعمالی و فاصله نسبی بین محرك و خط جدایش است.
اسماعیل زاده و آقازینالی ]9[ جریان حول استوانه تحت تاثیر یک محرك الکتروهیدرودینامیکی میله صفحهاي چسبیده به سطح استوانه را با استفاده از روش اجزاء محدود به صورت عددي مورد مطالعه قرار دادند و دریافتند که با اعمال میدان الکتریکی، جریان سیال در مجاورت نیمه جلویی استوانه شتاب می گیرد و در نتیجه با افزایش اختلاف فشار سطوح جلویی و عقبی استوانه ضریب بازدارندگی افزایش می یابد. همچنین آنها نشان دادند که حجم ناحیه دنباله در پشت استوانه و ضریب بازدارندگی با افزایش سرعت ورودي به ازاي یک اختلاف ولتاژ اعمالی یکسان کاهش می یابند.
تطهیري و اسماعیل زاده ]10[ اثر پدیده EHD بر کنترل لایه مرزي ورقهاي روي استوانه به کمک دو الکترود سیمی و دو الکترود نواري چسبیده به سطح استوانه را به صورت عددي مورد بررسی قرار دادند. نتایج آنها بیان کننده این مطلب است که ضریب بازدارندگی افزایش یافته و در تمامی حالت ها افزایش سرعت در لایه مرزي مشاهده می شود که حاکی از افزایش مومنتم ذرات سیال در اثر اعمال میدان الکتریکی است. همچنین این تأثیرات در اعداد رینولدز پایین به نسبت بیشتر است. علاوه بر این آنها دریافتند که ضریب بازدارندگی تابع هندسه مسأله، ولتاژ اعمالی و رژیم جریان می باشد که در هندسه مورد نظر آنها ضریب بازدارندگی افزایش یافته است.
مطالعات زیادی در خصوص جریان هوا روی اجسام در کاربرد های مختلف از جمله انرژی باد و هوافضا انجام شده است. نان و همکاران ]11[ به مدلسازی جریان هوا روی ایرفویل توربین بادی پرداختند. ژو و همکاران ]12[ به مدلسازی جریان توربولانسی یک ایرفویل به کمک شبکه عصبی در رینولدزهای بالا پرداختند. ژی و همکاران ]13[ به مدلسازی جریان جت هوا روی ایرفویل در یک جسم آسان برا پرداختند. همچنین لی و همکاران ]14[ به مدلسازیLES جریان توربولانسی پیچیده روی ایرفویل یخ زده اشاره کردند و جدایش جریان را بررسی کردند. این روش همچنین توسط آلبرس و اشرودر ]15[ جهت مدلسازی، افزایش لیفت و کاهش دراگ در روی ایرفویلهای متحرک به کار گرفته شد. هولمن و فورست]16[ نیز از روش حجم محدود برای مدلسازی جریان لمینار و توربولانت روی ایرفویل بهره گرفتند. هدف از پروژه حاضر بررسی جریان سیال همراه با میدان مغناطیسی بر روی سطوح آیرودینامیکی، به منظور کاهش ضریب پسا و افزایش ضریب برا و همچنین کاهش گردابه های اتلاقی و به تاخیر انداختن زاویه استال می باشد. معادلات حاکم بر جریان شامل معادلات پیوستگی و مومنتوم همراه با عباراتی است که دربرگیرنده آثار میدانهای مغناطیسی می باشد که به کمک روش حجم محدود حل شده اند. معادلات ماکسول و قانون اهم هم بر مساله حاکم می باشد که از ترکیب آنها یک معادله انتقال برای محاسبه میدان مغناطیسی به دست می آید. در انتها اثر ارتفاع و ماخ و میدان مغناطیسی بر نیروهای دراگ و لیفت بررسی شده است.
2- مدلسازی و معادلات حاکم
جریان توربولانت تراکم پذیر به شکل برداری با در نظر گرفتن معادلات ناویر استوکس و انتقال به شکل زیر است ]16[:
(1) |
|
که در آن بردار متغیرهای میانگین، F شار لزجت مربوطه، R شار ویسکوز مربوطه و بردارQترم منبع را نشان میدهد.
معادلات حاکم در مدل استاندارد انرژي
سينتيک توربولانس (k) و ميزان پراکندگي آن از معادلات زير به دست مي آيند:
(2) |
|
(3) |
|
در اين معادلات، توليد انرژي سينتيک توربولانس، ناشي از گراديان سرعت است
توليد انرژي سينتيکي توربولانس، ناشي از نيروهاي شناوري،
تاثير نوسانات انبساطي در جريان هاي تراکم پذير بر روي ميزان پراکندگي هستند.
ثابت ها بوده،
اعداد پرانتل مغشوش براي
ميباشند.
ترم هاي تعريف شده توسط کاربر ميباشند.
مدل سازي لزجت مغشوش در مدل استاندارد : لزجت مغشوش يا لزجت ادي
از ترکيب
به صورت زير به دست ميآيد:
|
| (4) |
که عددي ثابت است.
2-1- ئوري مدل MHD
تاثيرات متقابل ميدان جريان سيال و ميدان مغناطيسي را ميتوان بر مبناي دو عامل بنيادي بررسي کرد: القاي جريان الکتريسيته ناشي از حرکت ماده رسانا در ميدان مغناطيسي و تاثير نيروي لورنتس در نتيجه برهم کنش جريان الکتريسيته و ميدان مغناطيسي. به طور کلي جريان الکتريکي القايي و نيروي لورنتس، تمايل به از بين بردن عامل ايجاد کننده ي خود دارند. بنابراين نيروي لورنتس ناشي از حرکات ماده رسانا در ميدان مغناطيسي در جهتي است که مانع آن حرکات گردد. القاي الکتريسيته همچنين ميتواند در حضور ميدان مغناطيسي متغير با زمان نيز، اتفاق بيافتد. به هر حال نتيجه اين امر بر هم زدن حرکت سيال توسط نيروي لورنتس است. ميدان هاي الکترومغناطيس توسط معادلات ماکسول1 توصيف مي گردند.
| (5) |
| (6) |
| (7) |
| (8) |
که در اين معادلات،شدت ميدان مغناطيسي و واحد آن تسلا،
ميدان الکتريکي و واحد آن ولت بر متر
است.
ميدان هاي القايي مغناطيسي و الكتريکي بوده، q و
به ترتيب چگالي بار و بردار چگالي جريان و واحدهاي آن ها
و
ميباشند.
توسط معادلات زير تعريف ميشوند.
| (9) |
| (10) |
که به ترتيب قابليت نفوذ مغناطيسي و ثابت دي الکتريک2ميباشند. براي موادي که به اندازه ي کافي هادي الکتريسيته باشند، نظير فلزات مايع، از چگالي بار الکتريکي q و جريان جابجايي
معمولاً صرفنظر ميگردد. در مطالعه تقابل ميان جريان سيال و ميدان الکترومغناطيس، دانستن چگالي جريان
ناشي از القا، حياتي است. به طور کلي دو طريق براي ارزيابي چگالي جريان وجود دارد. يکي از طريق حل معادله القاي مغناطيسي و ديگري حل معادله پتانسيل الکتريکي.
در اين روش براي اندازه گيري (چگالي جريان) معادله القاي مغناطيس، از قانون اهم3و معادلات ماکسول به دست ميآيد. اين معادله تداخل جريان سيال و ميدان مغناطيسي را ميسر ميسازد.
| (11) |
که هدايت الکتريکي ماده است. براي سيال با سرعت
در ميدان مغناطيسي
قانون اهم به شکل زير نوشته ميشود:
| (12) |
از قانون اهم و معادلات ماکسول معادلهي القا به شکل زير حاصل ميگردد.
| (13) |
با حل ميدان مغناطيسي چگالي جريان
را ميتوان از رابطه آمپر به دست آورد.
| (14) |
به طور کلي، ميدان مغناطيسي در يک مسئله MHD را ميتوان به دو ميدان خارجي
و يک ميدان القايي
ناشي از حرکت سيال تجزيه کرد. در اين صورت تنها ميدان القايي
نياز به حل دارد. در این مطالعه ضرایب دراگ و لیفت از روابط زیر محاسبه می شوند]15[.
(15) |
|
(16) |
|
3- تعریف مسئله، روش حل و استقلال از شبکه
در پژوهش حاضر از روش القای مغناطیسی استفاده می شود به این صورت که نيروي لورنتس با استفاده از ماژول فلوئنت به شکل يک ترم منبع به معادلات مومنتوم سيال افزوده شده است. معادلات ناشی از جریان حول دو هندسه (دو بعدی و سه بعدی) در ارتفاعات و ماخ بالا بررسی شده است. جريان سيال تراکم ناپذير و لزج فرض شده است و براي حل معادلات توربولانس از مدل استفاده شده است. روش حل عددی در فلوئنت بر مبنای روش حجم محدود است و مشخصات سیستم مورد استفاده مطابق شکل 1 میباشد. قابل ذکر است که میزان زمان دیتاگیری برای هندسه دو بعدی تقریبا 25 دقیقه و برای هندسه سه بعدی9ساعت میباشد.
براي رسیدن به نتایج قابل اطمینان باید شبکهي
رسم شده را مورد بررسی قرار داد تا این شبکه توانایی این را داشته باشد که داده ها را با سرعت بالا و همچنین با دقت کافی انتقال دهد تا نتایج قابل اطمینانی بدست آید. براي اینکه بتوان از حداکثر توان محاسباتی سیستم استفاده کرد و همچنین زمان اضافه اي صرف محاسبه نشود، باید با استفاده از مطالعه شبکه به یک شبکه معقول رسید. بدین منظور ابتدا از شبکه با تراکم پایین استفاده می شود. یک یا چند کمیت انتگرالی و یا توزیع یافته مورد بررسی قرار می گیرد. در هر مرحله تراکم افزایش می یابد تا جایی که پارامترهاي مورد بررسی ما تغییري نسبت به حالت قبل نکند. با ریز کردن تعداد تقسیمات برای دیواره های جسم و دیواره های مربع در نظر گرفته شده حول جسم، طبق شکل 2 و جدول1، شبکه بندی آغاز میشود.
قابل ذکر است که در شکل 1 برای قسمت cm 80 شبکه بندی شامل دو قسمت می باشد. قسمت اول مربوط به قسمت خمیده جلو و عقب جسم نوک تیز و قسمت دوم مربوط به دیواره افقی بالا و پایین جسممیباشد. شرایطحلنیز مطابق جدول2میباشد.
[1] -Maxwells equations
[2] - Magnetic permeability
[3] -Ohms law
شکل 1. طریقه قسمت بندی برای ایجاد شبکه مش.
شکل 2. هندسه دو بعدی شبکه بندی شده در گمبیت.
جدول 1. طریقه نامگذاری برای خطوط هندسه رسم شده.
L1 | 1500 cm |
L2 | قسمت خمیده جلو و عقب جسم نوک تیز |
L3 | به دیواره افقی بالا و پایین جسم (50 cm) |
جدول 2. شرایط حل برای مطالعه شبکه در هندسه 2 بعدی.
| ارتفاع | ماخ |
|
0 | 9000 m | 6 | 20 |
شکل 3. نمودار استقلال شبکه برای 5 شبکه بندی مختلف 2d.
جدول 3. مقادیر همگرا شده ضریب لیفت برای 5 شبکه بندی متفاوت 2d.
ضریب لیفت CL | نوع شبکه بندی | ||||||||||||
0.13842 | 90 -60 - 70 | ||||||||||||
0.13874 | 100 -70 - 80 | ||||||||||||
0.13911 | 110 -80 - 90 | ||||||||||||
0.13922 | 120 -90 - 100 | ||||||||||||
0.13926 | 125 -95 - 105 | ||||||||||||
| |||||||||||||
شکل 4. نمودار y+ برای دیواره بالا برای شبکه بندی 100-900-120.
جدول 4. شرایط حل برای مطالعه شبکه در هندسه 3 بعدی.
|
شکل 7. نمودار استقلال شبکه برای 4 شبکه بندی مختلف 3d.
جدول 5. مقادیر همگرا شده ضریب لیفت برای 4 شبکهبندی متفاوت 3d.
ضریب لیفت CL | نوع شبکه بندی |
0.85602 | 1 |
0.86011 | 2 |
0.86095 | 3 |
0.86114 | 4 |
با توجه به جدول بالا تمامی شبکه بندیهای صورت گرفته به ترتیب L1-L2-L3 می باشد. در شکل 3 ضریب لیفت برای زاویه 20 در جه برای حالت بدون MHD نشان داده شده است. همانطور که مشاهده می شود برای دو حالت آخر نتایج کاملا نزدیک به هم می باشند.
در جدول 3 مقادیر همگرا شده ضریب لیفت برای هر شبکه بندی نشان داده شده است. همانطور که از شکل و نمودار مشاهده می شود میزان تغییرات در دو شبکه آخر به میزان بسیار کمی با هم تفاوت دارند. پس برای کاهش زمان دیتاگیری از شبکه بندی 100-90-120 استفاده خواهد شد.
در شکل 4 نمودار y+ برای شبکه بندی 100-90-120 نشان داده شده است. که نشان می دهد شبکه انتخابی برای این هندسه مناسب است.
4-3- استقلال شبکه برای هندسه سه بعدی
برای هندسه سه بعدی نیز مطابق شکل 5 نتایج همگرا شده ضریب لیفت حاصل شد. شبکه بندی برای 4 نوع تعداد المان بررسی شد. همانطور که مشاهده می شود نتایج دو شبکه آخر در حد بسیار پایینی با یکدیگر اختلاف داشتند اما از نظر زمان محاسباتی 6 ساعت با یکدیگر اختلاف دارند. با صرف چنین زمانی و تغییر کم دقت و تعداد زیاد و اینکه این محاسبات باید براي موارد مختلف انجام پذیرد بهتر است براي صرفه جویی در زمان از شبکه با تعداد سلول 3824516 استفاده شود. شرایط حل مطابق جدول 4 می باشد. شکل 5 و 6 شماتیک هندسه سه بعدی و شبکه بندی شده آن در نرم افزار گمبیت نشان داده شده است. شکل 7 نیز صحت شبکه انتخابی را از نظر استقلال از شبکه نشان می دهد.
9000 متر برابر است با 29527.56 فوت بنابراین مقدار فشار و دمای هوا در این ارتفاع با استفاده از مراجع و با میان یابی بین 28000 و 30000 فوت به ترتیب برابر است با 30797.844 پاسکال و 229.747 درجه کلوین. مقدار B در ماژول فلوئنت در راستای Y و مقدار منفی در نظر گرفته می شود.
در این قسمت به بررسی صحت حل هندسه مورد مطالعه پرداخته می شود. برای این امر، در این پژوهش از مقاله ریاخووسکی]17[ استفاده می شود. هندسه مورد بررسی در این مقاله در شکل 8 نشان داده شده است. همانطور که مشاهده می شود هندسه این پژوهش نوک تیز بوده و از شرط تقارن برای کل هندسه استفاده شده است، به طوری که نصف هندسه در نرم افزار اوپن فوم تحت شرایط حل مطابق با جدول 6 شبیه سازی شده است.
شکل 8. شماتیک هندسه مورد مطالعه ریاخووسکی و همکاران]17[.
جدول 6. شرایط حل شبیه سازی ریاخووسکی و همکاران]17[ در نرم افزار اپون فوم.
هوا | گاز | |||||||||||||||||||||||||||||||||||||||
3000 |
| |||||||||||||||||||||||||||||||||||||||
216.5 |
| |||||||||||||||||||||||||||||||||||||||
3.37 | ماخ |
| ارتفاع | ماخ |
| |||||||||||||||||||||||||||||||||||||
0-0.7 | 9000 m | 6 | 40-0 | |||||||||||||||||||||||||||||||||||||
| ||||||||||||||||||||||||||||||||||||||||
شکل 10. اثر زاویه حمله و اضافه کردن مگنت بر ضریب لیفت در هندسه دو بعدی. | ||||||||||||||||||||||||||||||||||||||||
جدول 8. مقادیر عددی اثر زاویه حمله و اضافه کردن مگنت بر ضریب لیفت. | ||||||||||||||||||||||||||||||||||||||||
NON MHD | With MHD | |||||||||||||||||||||||||||||||||||||||
|
| |||||||||||||||||||||||||||||||||||||||
| ||||||||||||||||||||||||||||||||||||||||
شکل 11. اثر زاویه حمله و اضافه کردن مگنت بر ضریب دراگ در هندسه دو بعدی. | ||||||||||||||||||||||||||||||||||||||||
جدول 9. مقادیر عددی اثر زاویه حمله و اضافه کردن مگنت بر ضریب دراگ. | ||||||||||||||||||||||||||||||||||||||||
NON MHD | With MHD | |||||||||||||||||||||||||||||||||||||||
|
|
|
شکل 12. کانتور ماخ برای زاویه حمله 10 درجه بدون MHD برای هندسه دو بعدی. |
|
شکل 13. کانتور ماخ برای زاویه حمله 10 درجه با MHD برای هندسه دو بعدی. |
در شکل 9 فشار وارد شده از طرف سیال در طول هندسه مذکور، با نتایج روش حاضر مقایسه شده است. همانطور که مشاهده می شود نتایج شبیه سازی انجام شده با نرم افزار فلوئنت در کار حاضر با درصد قابلقبولی بانتایج کار ریاخووسکی و همکاران ]17[ تطابق دارد.
4- نتایج
در این قسمت به بررسی اثر جریان MHD و اثر زاویه حمله بر مقادیر ضریب دراگ و ضریب لیفت هندسه دو بعدی پرداخته می شود. در جدول 7 شرایط حل برای این هندسه نشان داده شده است. شکل 10 و 11 به ترتیب ضریب لیفت و ضریب دراگ را در 9 زاویه حمله متفاوت در دو حالت با مگنت و بدون مگنت نشان داده است. همانطور که مشاهده می شود در زوایای حمله پایین تفاوت چندانی در مقادیر ضرایب لیفت و دراگ مشاهده نمی شود. اما با افزایش زاویه حمله و با اضافه کردن مگنت افزایش چشمگیری در ضرایب مشاهده شده است. همچنین همانطور که مشاهده می شود در هیچکدام از حالت های بدون مگنت و با مگنت زاویه استال رخ نمیدهد که این امر ناشی از شرایط حل در فشار و ماخ بالاست که در پژوهشهای گذشته به آن اشارهای نشده است. جدول 8 و 9 مقادیر عددی نمودار را در زوایای حمله مختلف نشان میدهد.
در شکل های 12 و 13 کانتور ماخ حول هندسه دوبعدی نشان داده شده است. همانطور که مشاهده می شود در انتهای هندسه با جریان MHD جریان سکون کمتری به وجود آمده است.
4-5-2- نتایج هندسه سه بعدی
در این قسمت به بررسی اثر جریان MHD و اثر زاویه حمله بر مقادیر ضریب دراگ و ضریب لیفت هندسه سه بعدی پرداخته می شود. در این هندسه دو شرایط حل بررسی می شود که در جداول 10
و13 نشان داده شده است. برای حالت A یعنی فشار 9000 و ماخ 6، شکل 14 و 15 به ترتیب ضریب لیفت و ضریب دراگ را در 9 زاویه حمله متفاوت در دو حالت با مگنت و بدون مگنت نشان داده است. همانطور که در این هندسه نیز مشاهده می شود در زوایای حمله پایین تفاوت چندانی در مقادیر ضرایب لیفت و دراگ مشاهده نمی شود. اما با افزایش زاویه حمله و با اضافه کردن مگنت افزایش چشمگیری در ضرایب دیده شده است. اما تفاوت اصلی در این حالت نسبت به حالت دوبعدی این است که زاویه استال در هندسه دو بعدی مشاهده نشد، اما همانطور که در شکل 16 نشان داده شده است زاویه استال برای حالت بدون MHD اتفاق افتاده، که با اضافه کردن مگنت در این هندسه و این شرایط حل این زاویه به تاخیر افتاده است. جداول 11 و 12 مقادیر عددی نمودارها را نشان می دهد.
همچنین برای هندسه سه بعدی، شرایط حل B مورد بررسی قرار داده شده است. شرایط حل برای این حالت در جدول 13 مشاهده می شود. همانطور که از شکل 4-9 مشخص است با افزایش فشار و افزایش ماخ زاویه استال به تاخیر می افتد. همچنین در این ارتفاع و ماخ میزان افزایش ضریب لیفت با اضافه گردن مگنت، نسبت به حالت بدون مگنت افزایش بسیار بیشتری نسبت به شرایط حل A دارد. همچنین با مشاهده تغییرات ضریب دراگ در شکل 17 مشاهده میشود که سیر صعودی ضریب دراگ با افزایش زاویه حمله نسبت به شرایط حل A کمتر میباشد.
جدول 10. شرایط حل A برای هندسه سه بعدی.
| ارتفاع | ماخ |
|
0-0.7 | 9000 m | 6 | 40-0 |
شکل 14. اثر زاویه حمله و اضافه کردن مگنت بر ضریب لیفت در هندسه سه بعدی و شرایط A.
جدول 11. مقادیر عددی اثر زاویه حمله و اضافه کردن مگنت بر ضریب لیفت. | |||||||||||||||||||||||||||||||||||||
NON MHD | With MHD | ||||||||||||||||||||||||||||||||||||
|
|
شکل 15. اثر زاویه حمله و اضافه کردن مگنت بر ضریب دراگ در هندسه سه بعدی و شرایط A.
جدول 12. مقادیر عددی اثر زاویه حمله و اضافه کردن مگنت بر ضریب دراگ.
NON MHD | With MHD | ||||
0.006 | 0 | 0.0067 | 0 | ||
0.012 | 5 | 0.013 | 5 | ||
0.031 | 10 | 0.037 | 10 | ||
0.0732 | 15 | 0.087 | 15 | ||
0.133 | 20 | 0.175 | 20 | ||
0.268 | 25 | 0.281 | 25 | ||
0.398 | 30 | 0.427 | 30 | ||
0.581 | 35 | 0.598 | 35 | ||
0.668 | 40 | 0.757 | 40 |
جدول 13. شرایط حل B برای هندسه سه بعدی.
| ارتفاع | ماخ |
| |||||||||||||||||||||||||||||||||||||
0-0.7 | 50000 m | 8 | 40-0 | |||||||||||||||||||||||||||||||||||||
| ||||||||||||||||||||||||||||||||||||||||
شکل 16. اثر زاویه حمله و اضافه کردن مگنت بر ضریب لیفت در هندسه سه بعدی و شرایط B. | ||||||||||||||||||||||||||||||||||||||||
جدول 14. مقادیر عددی اثر زاویه حمله و اضافه کردن مگنت بر ضریب لیفت. | ||||||||||||||||||||||||||||||||||||||||
NON MHD | With MHD | |||||||||||||||||||||||||||||||||||||||
|
| |||||||||||||||||||||||||||||||||||||||
| ||||||||||||||||||||||||||||||||||||||||
شکل 17. اثر زاویه حمله و اضافه کردن مگنت بر ضریب دراگ در هندسه سه بعدی و شرایط B. |
جدول 15. مقادیر عددی اثر زاویه حمله و اضافه کردن مگنت بر ضریب دراگ | |||||||||||||||||||||||||||||||||||||
NON MHD |
With MHD | ||||||||||||||||||||||||||||||||||||
|
|
|
شکل 18. کانتور ماخ برای زاویه حمله 10 درجه بدون MHD برای شرایط حل A و هندسه سه بعدی. |
|
شکل 19. کانتور ماخ برای زاویه حمله 10 درجه با MHD برای شرایط حل A و هندسه سه بعدی. |
در شکل های 18 و 19 کانتور ماخ برای شرایط حل B به صورت دو بعدی نشان داده شده است. همانطور که مشاهده می شود جریان سکون در پشت هندسه در حالت جریان MHD کاهش پیدا می کند. قابل ذکر است که این کانتور در مقطع وسط هندسه رسم شده است.
5- نتيجهگيري
در این پژوهش به بررسی اثر جریان MHD و اثر زاویه حمله بر مقادیر ضریب دراگ و ضریب لیفت هندسه دو بعدی پرداخته و مشاهده شد در زوایای حمله پایین تفاوت چندانی در مقادیر ضرایب لیفت و دراگ مشاهده نمی شود. اما با افزایش زاویه حمله و با اضافه کردن مگنت افزایش چشمگیری در ضرایب مشاهده شد. همچنین نشان داده شد در هیچکدام از حالت های بدون مگنت و با مگنت زاویه استال رخ نمی دهد که این امر ناشی از شرایط حل در فشار و ماخ بالا می باشد. خلاصه ای از نتایج به شرح زیر است:
· در کانتورهای ماخ حول هندسه دوبعدی نشان داده شد که در انتهای هندسه با جریان MHD جریان سکون کمتری به وجود می آید.
· در قسمت بررسی اثر جریان MHD و اثر زاویه حمله بر مقادیر ضریب دراگ و ضریب لیفت هندسه سه بعدی، دو شرایط حل بررسی شد یکی شرایط A ، فشار 9000 و ماخ 6 و دیگری شرایط B، فشار 50000 و ماخ 8. برای حالت A در زوایای حمله پایین تفاوت چندانی در مقادیر ضرایب لیفت و دراگ مشاهده نشد. اما با افزایش زاویه حمله و با اضافه کردن مگنت افزایش چشمگیری در ضرایب نشان داده شد.
· برای شرایط حل B نیز مشاهده شد که با افزایش فشار و افزایش ماخ زاویه استال به تاخیر می افتد. همچنین در این ارتفاع و ماخ میزان افزایش ضریب لیفت با اضافه گردن مگنت، نسبت به حالت بدون مگنت افزایش بسیار بیشتری نسبت به شرایط حل A داشت.
· با مشاهده تغییرات ضریب دراگ مشاهده شد که سیر صعودی ضریب دراگ با افزایش زاویه حمله نسبت به شرایط حل A کمتر می باشد.
· در کانتورهای ماخ برای شرایط حل B برای هندسه سه بعدی جریان سکون در پشت هندسه در حالت جریان MHD کاهش پیدا کرد.
· با مقایسه دو هندسه سه بعدی و دو بعدی مشاهده شد که زاویه استال در هندسه دو بعدی در حالت ارتفاع 9000 و ماخ 6 رخ نداده اما در هندسه سه بعدی و همین شرایط حل (شرایط حل A) زاویه استال برای حالت بدون MHD اتفاق افتاده، که علت این امر را می توان در خط جریانی بودن هندسه دو بعدی دانست. البته در ادامه مشاهده شد که با اضافه کردن مگنت در هندسه سه بعدی و شرایط حل مذکور این زاویه به تاخیر افتاده است.
· بیشترین درصد افزایش ضریب لیفت برای هندسه سه بعدی و شرایط ماخ 8 و ارتفاع 50000 متر رخ داده است که 77.5% افزایش در ضریب لیفت ناشی از میدان مغناطیسی مشاهده شده است.
فهرست مراجع
[1] A. Salmasi, A. Shadaram, M. Mirzaei, A. Sh. Taleghani, Numerical andexperimental investigation on the effect of a plasma actuator on NLF0414 airfoil’s efficiency after the stall, Modares Mechanical Engineerin, Vol. 12, No. 6, pp. 104-116, 2013.
[2] Seyed Saeed Hoseininezhad, Nima Amanifard*, Hamed Mohaddes Deylami, Farid Dolati, Numerical study of low characteristics around a NACA 4412 asymmetric airfoil under the influence of electric field, Modares Mechanical Engineerin,2014.
[3] Xisto, Carlos M., et al. "Implementation of a 3D compressible MHD solver able to model transonic flows." Proc. ECCOMAS CFD 2010-V European Conference on Computational Fluid Dynamics. 2010.
[4] L. Zhao, K. Adamiak, Numerical simulation of the electrohydrodynamic flow in a single wire plate electrostatic precipitator, IEEE Transactions on Industry Applications, Vol. 44, pp. 683-691, 2013.
[5] Y. Tada, A. Takimoto, D. Ueda, Y. Hayashi, Heat Transfer Enhancement a Convective Field with Corona Discharge, Transactions of the Japan Society of Mechanical Engineers, Part B, pp. 217-228, 2003.
[6] K. T. Hyun, C. H. Chun, The wake flow control behind a circular cylinder using ion wind, Experiments in Fluids, Vol. 35, pp. 541-552, 2014.
[7] L. Leger, E. Moreau, G. Artana, G. Touchard, Influence of a DC corona discharge on the airflow along an inclined fat plate, Journal of Electrostatics, Vol. 51-52, pp. 300-306, 2010.
[8] S. Roberto, A. Guillermo, Steady control of laminar separation over airfoils with plasma sheet actuators, Electrostatics, Vol. 64, pp. 604-610, 2016.
[9] E. Esmaeilzadeh, M. Aghazeinali, Numerical investigation of flow around cirular cylinder affacted by an electrohydrodynamic flush mounted wireplate actuator, Journal of Faculty of Eng, Vol. 32, No. 2, 2006.
[10] Gh. Tathiri, E. Esmaeilzadeh, The effect of EHD phenomenon on the control of sheet boundary layer over a cylinder using two electrode wires and two attached tape electrodes to the surface of the cylinder, International Journal of Engineering, Vol. 19, No. 3, pp. 51-60, 2015. (in persian)
[11] Yi-Nan, Zhang, Cao Hui-Jing, and Zhang Ming-Ming. "A calculation method for modeling the flow characteristics of the wind turbine airfoil with leading-edge protuberances." Journal of Wind Engineering and Industrial Aerodynamics 212 (2021): 104613.
[12] Zhu, Linyang, Weiwei Zhang, Xuxiang Sun, Yilang Liu, and Xianxu Yuan. "Turbulence closure for high Reynolds number airfoil flows by deep neural networks." Aerospace Science and Technology 110 (2021): 106452.
[13] Zhi, Haolin, Zhenhao Zhu, Yujin Lu, Shuanghou Deng, and Tianhang Xiao. "Aerodynamic performance enhancement of co-flow jet airfoil with simple high-lift device." Chinese Journal of Aeronautics (2021).
[14] Lee, Young Mo, Jae Hwa Lee, Lawrence Prince Raj, Je Hyun Jo, and Rho Shin Myong. "Large-eddy simulations of complex aerodynamic flows over multi-element iced airfoils." Aerospace Science and Technology 109 (2021): 106417.
[15] Albers, Marian, and Wolfgang Schröder. "Lower drag and higher lift for turbulent airfoil flow by moving surfaces." International Journal of Heat and Fluid Flow 88 (2021): 108770.
[16] Holman, Jiří, and Jiří Fürst. "Numerical simulation of separation induced laminar to turbulent transition over an airfoil." Journal of Computational and Applied Mathematics (2021): 113530.
[17] Ryakhovskiy, Aleksey Igorevich, and A. A. Schmidt. "MHD supersonic flow control: OpenFOAM simulation." Труды Института системного программирования РАН 28, no. 1 (2016)