-۱ مقدمه

راکتورهای تحقیقاتی بدلیل اسـتفادههای وسیع در زمینههای علمی، صنعتی و کاربردهای پزشکی، دارای نقش مهمی در حوضـههای علوم هسـتهای و تکنولوژی، میباشند. هدف راکتورهای تحقیقاتی داشتن شار نوترونی ماکزیمم در محل تابش دهی نوترون میباشـد. این راکتورها علاوه بر شـرط شار ماکزیمم، بایستی پارامترهای ایمنی را نیز برآورده کنند. اسـتفاده از کدهایی که برای آنالیز ایمنی راکتورهای قدرت بکار برده میشود، بدلیل پیچیدگیهای آنها، در راکتورهای تحقیقاتی توصـیه نمیشود .[۱] بکارگیری تعداد زیادی از این کدها، نیاز به تلاش و مهارت زیاد برای تهیه ورودی و همچنین تحلیل خروجی کدها میباشــد. در صــورتی که میتوانیم با تقریبهای مناسـبی معادلات را برای این گونه از راکتورها به سادگی حل کنیم و در نتیجه کد مخصوص این راکتورها را تولید کنیم. اگر چه برخی از این کدها برای بکار گیری در راکتورهای تحقیقاتی سـازگار شدهاند ( [۲] و([۳] ولی بـا این وجود توســـعـه روشهـای ســـاده برای آنـالیز رفتارهای گذرای راکتور بدون نیاز به شـبیهسازی کامل آن بسیار مطلوب است. در سالهای اخیر تلاشهای زیادی در این راستا انجام شده است( [۴] ( [۵], هـدف اصـــلی این تحقیق، تهیه یک مدل عددی برای آنالیز رفتار کاهش جریان خنک کنندگی راکتور تهران میباشد.

۷۰۹

.۲ روش کار

مدل عددی پیشــنهادی مانند هر کد دینامیکی بکار رفته برای تحلیل رفتارهای گذرای یک راکتور هســتهای از سـه بخش نوترونیک، ترموهیدرولیک و محاسـبه فیدبکها تشـکیل شده است . در بخش نوترونیک معادلات مربوط به تولید توان حل شده و با استفاده از توان بدست آمده دماهای مربوط به سوخت، غلاف و خنککننده در بخش ترموهیدرولیک محاســبه میشــود. افزایش یا کاهش دما باعث تغییر نرخ واکنشهای مربوط به تولید توان میشــود. فیدبکهای دمایی ســهم هریک از عناصــر موجود در قلب راکتور (ســوخت، خنککننده) را محاسبه و به بخش نوترونیک اعمال میکند.

.۱,۲ مدل سینتیک نوترونی

محاسـبه دقیق توان تولید شـده در قلب هر راکتور هستهای مستلزم حل معادلات نوترونیک وابسته به مکان و زمان میباشــد. از آنجاییکه حل همزمان وابســتگی زمانی و مکانی بســیار پیچیده اســت، از تقریب ســینتیک نقطهای برای حل معادلات دینامیک اسـتفاده میشـود. تقریب سـینتیک نقطهای راه حل مناسـب و سریع برای محاســبه توان راکتورهایی اســت که قلب کوچکی دارند.[۶] معادلات ســینتیک نقطهای به روشهای متفاوتی قابل حل است و در اکثر مواقع از روش رانجو گوتا مرتبه ۴ استفاده میشود .[۷] اما در این کد از مدل ریاضی که در سال ۲۰۰۴ توسط ماتیو کینارد۱ برای حل معادلات دینامیک راکتور پیشنهاد شد، استفاده شده است .[۸] الگوریتم حـل مـدل PCM2 بـدینصـــورت اســـت که معادلات دینامیک نوترون، در یک بازه کوچک h بین زمانهای ti,ti+1 مسـتقل از زمان فرض میشـوند و مقادیر متوسط آنها در زمان ti+h/2 محاسبه شده و بعنوان مقدار ثابت در این بازه در نظر گرفته میشـود. معادلات به روش ویژه مقداری حل میشوند که منجر به ریشه یابی معادله inhour شـده و با اسـتفاده از این ریشـهها، بردارهای ویژه سـاخته و در نهایت جوابهای مسئله بدسـت میآیند. شـرایط اولیه حل معادلات سینتیک در این مدل بصورت تعادلی در نظر گرفته شده است، که در آن تغیرات زمانی دانسـیته مولدهای نوترونهای تاخیری صفر در نظر گرفته میشود. با استفاده از نرم افزار Matlab این روش در حالت شــش گروهی پیاده ســازی شــد. اگر جملات مربوط به چشــمه تولید نوترون و راکتیویتـه تابع آرامی از زمان باشـــد، روش PCM جوابهای بســـیار دقیقی میدهد .([۸]) برای اطمینان از جوابهای بخش نوترونیک این مدل، نتایج مربوط به راکتیویته نوســانی ذکر شــده در کتاب هتریک به عنوان مرجع، برای راکتوری با مشـــخصـــات = ۰٫۰۰۷۹ ،= ۰٫۰۷۷ ، =۰,۰۰۱s ، T=50 با راکتیویته نوســـانیρ( ) = ۰٫۰۰۵۳۳sin⁡(۵۰ ) بدست آورده شده است.[۶] تغییرات چگالی نوترون در حالت دوگروهی در مدت ۱۰۰ثانیه با بازه زمانی یک ثانیه محاسبه شده است. معیار این اعتبار سنجی، سنجش حداکثر چگالی نوترون در این بازه زمانی اســت. جواب دقیق( حداکثر چگالی نوترون) برای این مســاله در حدود ۴۹,۲۵ اســت که روش PCM مقدار آن را با بازه نســبتا بلند ۱ ثانیه ۵۰ محاســبه میکند که هم خوانی خوبی با جواب مرجع دارد. با انتخاب

۷۱۰

h در حدود میلی ثانیه مدل PCM نیز همان جواب را خواهد داد. بدیهی است هر چقدر گام زمانی را کوچکتر انتخاب کنبم جوابها از دقت بالایی برخوردار خواهند شــد. درصــد خطای محاســبه توان راکتور در روش PCM متناسب با h2 است.[۸]

.۲,۲ مدل ترموهیدرولیک

بخش ترموهیدرولیک این کد به روش مدل توده ای۳ تهیه شــده اســت. از این مدل عموما برای مطالعه رفتار کیفی دینامیک راکتور اسـتفاده میشـود. در این مدل از قانون سـرمایش نیوتن برای انتقال حرارت از غلاف به خنککننده اسـتفاده میشـود .[۹] در بخش ترموهیدرولیک، شبیهسازیها بر اساس مدل دوکانالی در نظر گرفته شــده اســت. کانال داغ نماینده داغ ترین صــفحه ســوخت به همراه کانال خنککننده مجاورش اســت و کانال متوســـط نماینده کل قلب راکتور میباشـــد. در این مدل اثرات دما، تنها برای یک دمای موثر برای هریک از نواحی ســوخت، غلاف و خنککننده بکار میرود. از دماهای محاســبه شــده در کانال متوســط برای محاســبه فیدبکهای دمایی اسـتفاده میشود. برای محاسبه حداکثر دماها در کانال داغ، بایستی ضرایب مربوط به فاکتور مهندسـی۴ و قله توان شـعاعی ۵ را در توان متوسـط ضـرب کرد. وابسـتگی دما به مکان بصورت استاتیک و معمولا بصـورت کسینوسی در نظر گرفته میشود در این مرحله فاکتور قله توان محوری۶ نیز اعمال میگردد. در این مدل راکتور تحقیقاتی MTR تهران که توسط آب سبک خنک میشود برای مطالعه انتخاب شده است. در این مدل از رابطه دیتوس باتلر۷ به عنوان ضریب انتقال حرارتی غلاف به خنککننده در جریانهای آشوبی اســتفاده شــده اســت. دمای اســتخر یا دمای ورودی کانال خنککننده به عنوان یک پارامتر ثابت در نظر گرفته میشـــود. ولی این پارامتر نیز میتواند متغیر باشـــد. در این مدل تغییرات چگالی خنککننده بر حســـب دما بصــورت یک رابطه خطی در نظر گرفته شــده اســت. مدل عددی ارایه شــده تنها به حل معادله بالانس انرژی میپردازد. معادلات وابســته به زمان و مکان مربوط به نواحی ســوخت غلاف و خنککننده راکتورهای MTR در مرجع شـماره ۴ آورده شده است. توان موضعی دارای شکل کسینوسی در راستای محور میباشد که شامل فاکتورهای قله توان محوری و شـــعاعی نیز میباشـــد .[۴] گرما از ســـطح غلاف از طریق روش همرفتی به خنککننده منتقل میشــود. برای ضــریب همرفتی گرمایی برای خنککنندگی اجباری در رژیم زیر ســرد۸ از رابطـه ســـیـدر- تیـت۹ یا دیتوس باتلراســـتفاده میشـــود .[۱۰] در مواقعی مانند خاموش بودن پمپ مدار خنککنندگی ثانویه ، دمای ورودی قلب دیگر پارامتر ثابتی نبوده و با زمان متغیر اســـت. در این موارد معادله موازنه انرژی بین توان تولیدی راکتور و افزایش دمای آب استخر به معادلات فوق افزوده میشود.

.۳,۲ مدل راکتیویته

مدل راکتیویته کل شامل راکتیویته خارجی تزریق شده به همراه راکتیویتههای ناشی از فیدبکهای دمایی میباشد. راکتیویته خارجی میتواند بصورت پلهای یا شیب دار با زمان باشد. همانگونه که میدانیم ضرایب راکتیویته

۷۱۱

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