Диссертация (781854), страница 21
Текст из файла (страница 21)
Проведено моделирование зон тепловыделяющего слоя. В частности, плавление частиц стали, а затем топлива учтено путем моделирования стоков тепла в тепловыделяющем слое.137Поставлена и решена задача движения пузыря переменной массы в жидкости. Полученное решение используется для определения источников теплапри конденсации паров натрия.Получена формула для стоков тепла в зоне с теплообменниками.Поставлена и решена задача определения напряженного состояния вверхней и нижней плитах напорной камеры.2.
Разработанная математическая модель реализована в программе БРУТ.3. Выполнена верификация программы БРУТ путем сопоставления результатов расчетов с экспериментальными данными и данными аналитическихтестов. Разработан тест для верификации модуля «Расчет проплавления внутриреакторных конструкций» программы БРУТ. Проведенные верификационныерасчеты показали адекватность моделирования процессов и явлений при исследовании тяжелой аварии.4.
Разработана математическая модель в одномерном приближении. Наоснове этой модели создана программа БРУТ – О.5. Выполнен расчет аварии UTOP для вариантов активной зоны реактораБН большой мощности с нитридным топливом и с MOX-топливом с помощьюпрограммы БРУТ. Показано, что во всех вариантах расплав удерживается вкорпусе реактора.6.
Проведен расчет аварии ULOF, при которой происходит разрушение 18ТВС первого и второго рядов активной зоны в реакторе МБИР. В данном случае частично разрушенная активная зона удерживается в корпусе реактора.138Глава 3Математическое моделирование процессов разгона реакторана мгновенных нейтронахПри обосновании безопасности реактора типа БН необходимо определитьколичество энергии, выделяющейся при неконтролируемом аварийном ростереактивности. Важной прочностной задачей является также расчет напряженнодеформированного состояния корпусов реакторов в условиях запроектнойаварии.
При запроектных авариях динамическое нагружение корпуса можетпроисходить при возникновении вторичной критичности в разрушеннойактивной зоне и соответственно при большом импульсе энерговыделения,приводящем к плавлению и испарению части топлива.Анализ аварии с разгоном реактора на мгновенных нейтронах требует совместного решения уравнений термодинамики, гидродинамики, теплофизики инейтронной кинетики. У истоков решения рассматриваемой задачи находитсяшироко известный метод Бете - Тайта, имеющий различные модификации[250]. Первые расчеты, проведенные в аналитической форме Бете и Тайтом,имели целью ориентировочно оценить максимально возможное количествоэнергии, выделяемое при разрушении активной зоны.
Несмотря на множестводопущений, принятых в расчете, полученные результаты стали основой для детальных исследований аварийных процессов.Развитие вычислительной техники позволило перейти к численному решению задачи. Удалось отказаться от некоторых упрощающих предположений,снизить степень консервативности результатов расчета, что является вполнеоправданным по отношению к расчетному анализу тяжелых запроектных аварий.Отечественные программы, предназначенные для расчета энерговыделения в реакторе при разгоне его на мгновенных нейтронах, были разработаныболее тридцати лет назад и использовались для анализа тяжелых аварий в про-139цессе проектирования отечественных реакторов на быстрых нейтронах.
Описание соответствующих математических моделей содержится в основном в работах ГНЦ РФ-ФЭИ и РНЦ “КИ”. Все отечественные программы были одномерными и позволяли рассчитывать энерговыделение в реакторе при вводе избыточной реактивности с заданной скоростью.Решению рассматриваемой задачи посвящены зарубежные программыMAX, ODEXCO, VENUS, MARS, POOL, KADIS [189]. При создании математической модели программы, предназначенной для расчета разгона реактора намгновенных нейтронах, был учтен опыт разработки перечисленных программ.Для расчета энерговыделения в реакторе при разгоне его на мгновенныхнейтронах был создан код ANPEX [74] (сокращение от английских слов analysisof power excursions), в математической модели которого движение материаловреактора описывается в цилиндрической системе координат в осесимметричном приближении. Математическая модель программы описана ниже.3.1.
Математическая модель разгона реактора на мгновенных нейтронах3.1.1. ГидродинамикаКак известно [160], при лагранжевом методе описания движений жидкости за основу берется движение фиксированных «жидких частиц», прослеживаемое начиная от некоторого начального момента времени t=t0. Под «жидкимичастицами» при этом понимаются объемы жидкости, линейные размеры которых очень велики по сравнению со средним расстоянием между молекулами,но все же настолько малы, что скорость и давление внутри частицы можно считать практически постоянными, и что в течение рассматриваемых промежутковвремени эти частицы можно считать перемещающимися «как одно целое».«Жидкие частицы» фактически представляют собой математические точки,плывущие вместе с жидкостью.
При эйлеровом методе описания движенийжидкости рассматривается поток среды, пересекающий границы конечноразностных элементов; при этом форма элементов остается неизменной [189]. Всистеме отсчета Эйлера элементы считаются неподвижными, а поток среды140проходит сквозь них. В дальнейшем используется лагранжев подход для описания движения материалов реактора.Уравнение сохранения массы записывается как 0 V0,V(3.1)где – плотность, V – объем фиксированной массы материала, 0 и V0 –значения и V в момент времени t = 0. Таким образом, задав метод расчетаобъема лагранжевой массы, изменяющегося вследствие перемещения ее границ, можно использовать уравнение (3.1) для определения плотности .В координатной системе Лагранжа уравнение сохранения количествадвижения имеет видDvρ X 1 v -p ρg X 2 ,Dt(3.2)где X 1 , X 2 ˗ дополнительные члены, выражающие перенос массы и количествадвижения от одной фазы к другой.Инерционными силами и силами тяжести в уравнении (3.2) можно пренебречь по сравнению с ростом давления, которое может составлять десяткиатмосфер и более при продолжительности аварийного процесса 5÷20 мс.
Крометого, X 1 X 2 0. Тогда из (3.2) следуют уравнения для определения составляющих скоростиu r 1 P; r1 Pυ z , zгде P = p + q, p – давление, рассчитанное из уравнения состояния, q - псевдовязкое давление.Существует вероятность возникновения больших градиентов давления впроцессе аварии, что затрудняет решение уравнений гидродинамики в конечноразностной форме. В данной модели используется методика [271], согласно ко-141торой в уравнения вводится новый параметр: псевдовязкое давление, определяемое следующим образом:1,44 A 02 dv dv<0, , еслиdtv dt 2qq 0 , еслиdv 0,dtгде А – площадь ячейки.Если в рассматриваемом процессе происходит сжатие жидкости, в уравнениях энергии и количества движения к действительному значению p добавляется q.
В результате градиент давления распространяется на несколько элементов жидкости (эффект размазывания ударной волны), что приводит к повышению устойчивости численного решения.Граничные условия, используемые в модели кода ANPEX:а) Материал на оси симметрии движется только в аксиальном направлении;б) Разработано и применяется в расчетах граничное условие, которое основанона втором законе Ньютона, используемом для определения ускорений в точкахна боковой границеaPS,mгде a – ускорение на границе активная зона – боковой экран, P – давление вприграничной ячейке активной зоны, S – площадь взаимодействия приграничной ячейки с фиктивной ячейкой, m – масса фиктивной ячейки.Под фиктивной ячейкой понимается ячейка вне активной зоны, в которойсосредоточена масса части бокового экрана на шаге z , с шириной r , равнойначальной ширине приграничной ячейки.Получение конечно-разностных уравнений.
Конечно-разностная схема,использованная при аппроксимации уравнений гидродинамики, является версией метода Кольски [238]. Данная схема получила название “метод средней точки”.142Нанесем на область решения задачи, которой является активная зона, конечную координатную сетку, определяющую границы лагранжевых частиц(рис. 3.1).
Введем следующие обозначения: RI начальная координата R после Iго смещения R , Z J начальная координата Z после J-го смещения Z . Значениякоординат r и z в момент времени t для точки ( R I , Z J ) обозначим через ri , j и z i , j .Положения, скорости и ускорения рассматриваются в вершинах конечноразностной сетки, плотности и давления – как средние по ячейке и обозначают1 1I ,J 2 2и т.д.zся rРис.3.1. Координатная сеткаПри получении конечно-разностных уравнений значительно удобнее использовать обозначения на рис. 3.1.
Затем уравнения можно преобразовать вуравнения с индексами.Конечно-разностное представление уравнения сохранения массы получим для ячейки AHDO, показанной на рис. 3.1. Текущую площадь AHDO можно выразить через координаты точек в видеA1( z H zO )(rD rA ) (rH rO )(z A z D );2(3.3)143начальная площадь Ao в момент времени t=0 определяется по формуле (3.3) заменой r и z на R и Z соответственно.
Для умеренных искажений радиус центраячейки AHDO представляется следующим образом:1r (rH rD rO rA );4начальное значение радиуса ro рассчитывается аналогично.Таким образом, уравнение сохранения массы в конечно-разностном видезаписывается как (см. формулу (3.1))1 1I ,J 2 2Ar1 11 11 1o, I , J o, I , J o, I , J 2 22 22 2Ar1 1 1 1I ,J I ,J 2 2 2 2Давления известны в дискретных точках 1,2,3 и 4, окружающих точку O(см.
рис. 3.1). Градиенты давления в точке O необходимо оценить.Точки 5,6,7 и 8 - средние точки OA,OB,OC и OD соответственно, принимается, что соответствующие давления будутP5 1( P1 P2 ),2P6 1( P2 P3 ),2P7 1( P3 P4 ),2P8 1( P4 P1 ).2Координаты точек 5,6,7 и 8 выражаются какr5 1(rA rO ),2z5 1( z A z O ),2r6 1(rB rO ),2z6 1( z B z O ),2r7 1(rC rO ),2z7 1( z C z O ),2r8 1(rD rO ),2z8 1( z D z O ).2Разложения в ряд Тейлора разностей давления между точками O и 5, O и6, O и 7, O и 8 можно записать в видеP5 PО (r5 rO )PP ( z5 zO ) ...,rz144P6 PO (r6 rO )PP ( z6 zO ) ... ,rzP7 PО (r7 rO )P8 PO (r8 rO )PP ( z7 zO ) ...
,rzPP ( z8 z O ) ....rzПренебрегая членами второго и более высокого порядков малости и выполняя алгебрагические преобразования, получим системуPP P5 P7 (r5 r7 ) r ( z 5 z 7 ) z ; P P (r r ) P ( z z ) P ,86868 6rzрешения которой сутьP ( P5 P7 )( z6 z8 ) ( P6 P8 )( z5 z7 )r( z6 z8 )( r5 r7 ) ( z5 z7 )( r6 r8 )P ( P6 P8 )( r5 r7 ) ( P5 P7 )( r6 r8 )z ( z6 z8 )( r5 r7 ) ( z5 z7 )( r6 r8 )Используя выражения для координат точек 5,6,7 и 8, получимP ( P1 P3 P2 P4 )( z B z D ) ( P2 P4 P3 P1 )( z A zC )r( z B z D )( rA rC ) ( z A zC )( rB rD )P ( P2 P4 P3 P1 )( rA rC ) ( P1 P3 P2 P4 )( z B z D )z( z B z D )( rA rC ) ( z A zC )( rB rD )Преобразовав числители дробей, будем иметьP ( P1 P3 )( z A zC z B z D ) ( P2 P4 )( z A zC z B z D )r( z B z D )( rA rC ) ( z A zC )( rB rD )P ( P1 P3 )( rA rC rB rD ) ( P2 P4 )( rA rC rB rD )z( z A zC )( rB rD ) ( z B z D )( rA rC )Соответственно, конечно-разностная аппроксимация уравнений движения имеет вид ( r ) I , J 4 / 1 1 1 1 1 1 1 1 ( P 1 1 P 1 1 )( z I , J 1 z I , J 1 I ,J I ,J I ,J I ,J I ,J I ,J 22222222 2222 z I 1,J z I 1,J ) ( P11I ,J 22P1 1I ,J 2 2)(z I ,J 1 z I ,J 1 z I 1,J z I 1,J )/( z I 1,J z I 1,J )(rI ,J 1 rI ,J 1 ) ( z I , J 1 z I , J 1 )( rI 1, J rI 1, J )(3.4)145 ( z ) I , J 4 / 1 1 1 1 1 1 1 1 ( P 1 1 P 1 1 )( rI , J 1 rI , J 1 rI 1, J I ,J I ,J I ,J I ,J I ,J I ,J 222222 2222 2 2 rI 1,J ) ( P11I ,J 22P11I ,J 22)(rI ,J 1 rI ,J 1 rI 1,J rI 1,J )/(rI 1,J rI 1,J )(z I ,J 1 z I ,J 1 ) (3.5) ( rI , J 1 rI , J 1 )( z I 1, J z I 1, J )Расчет ускорений осуществляется по формулам (3.4),(3.5).Для расчета скоростей используются следующие соотношения:r n 1/ 2 r n 1/ 2 rn tz n 1/ 2 z n 1/ 2 zn t ,11nn1где t ( t 2 t 2 ) .2Значения координат r и z при t t n1 рассчитываются по формуламr n1 r n r n1/ 2 t n1/ 2z n 1 z n z n1/ 2 t n 1/ 2Таким образом, перемещения определяются в результате двойного интегрирования ускорений на шаге по времени.3.1.2.