Численное моделирование задач термоупругости для конструкции с внутренним источником Текст научной статьи по специальности «Физика»

Математические заметки СВФУ Июль—сентябрь, 2017. Том 24, № 3

УДК 519.63

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ЗАДАЧ

ТЕРМОУПРУГОСТИ ДЛЯ КОНСТРУКЦИИ

С ВНУТРЕННИМ ИСТОЧНИКОМ

М. В. Васильева, П. Е. Захаров, П. В. Сивцев, Д. А. Спиридонов

Аннотация. Рассматривается численное моделирование термомеханического состояния конструкции, состоящего из источника тепла, газового зазора и оболочки. Математическая модель описывается нелинейной системой уравнений для температуры и перемещений. Выделение тепла происходит в подобласти источника. Возникающие перемещения за счет градиента температур рассчитываются в области источника и отдельно в оболочке и могут быть описаны как моделями линейной упругости, так и нелинейными моделями пластичности. Численная реализация построена на основе метода конечных элементов. Представлены результаты численного моделирования нелинейной модельной задачи в двумерной и трехмерной областях.

DOI 10.25587/SVFU.2018.3.10889

Ключевые слова: задачи термоупругости, тепловое расширение, теплоперенос, задача линейной упругости, модели пластичности, нелинейные задачи, метод конечных элементов, математическое моделирование.

Введение

Математическое моделирование мультифизичных задач широко используется для прогнозирования и понимания процессов, происходящих в конструкциях с внутренним источником тепла, например, в тепловыделяющих элементах [1—5]. Тепловыделяющий элемент является главным конструкционным элементом активной зоны ядерного реактора и состоит из топливных таблеток, газового зазора и оболочки. Процессы, происходящие в них, описываются, в частности, связанной задачей термоупругости [6-9]. Для моделирования таких мультифизичных процессов для нужд атомной энергетики разработаны различные программы, например, такие как BISON, FALCON и SVECHA [10-14]. Такое программное обеспечение, как правило, разрабатывается в национальных лабораториях. Следует отметить, что при построении математической модели необходимо также учитывать и различные мультимасштабные процессы, такие

Работа выполнена при поддержке Мегагранта Правительства Российской Федерации №14.Y26.31.0013.

© 2017 Васильева М. В., Захаров П. Е., Сивцев П. В., Спиридонов Д. А.

как выделение продуктов (осколков) деления, которое приводит к значительному изменению теплопроводности материала за счет накапливания слабопрово-дящих газов, а также появлению разрывов и трещин [15].

В данной работе рассматривается численное моделирование распределения температуры в конструкции с источником тепла и возникающего напряженно-деформированного состояния. Математическая модель термоупругого состояния описывается связанной системой уравнений для температуры и вектора перемещений. Для решения таких задач авторами были разработаны методы расщепления, позволяющие проводить последовательный расчет температуры и возникающих тепловых расширений [16]. Для аппроксимации напряженно-деформированного состояния используется метод конечных элементов, позволяющий строить дискретную модель для сложных геометрических областей. Геометрическая конфигурация включает в себя источник, оболочку и зазор, заполненный газом. Расчет температурного поля проводится во всей области. В каждой подобласти задаются характерные теплофизические параметры. Последующий расчет напряженно-деформированного состояния за счет термических расширений происходит в подобласти источника и отдельно в подобласти оболочки в предположении, что при расширении источник и оболочка не будут касаться друг друга.

Работа состоит из четырех разделов. В разд. 1 приводится постановка задачи. В разд. 2 представлен вычислительный алгоритм и конечно-элементная аппроксимация. В разд. 3 приводятся результаты численного решения задачи для осесимметричной постановки для области, содержащей один источник. Результаты моделирования для несимметричного источника приводятся в разд. 4 для трехмерной задачи.

Для построения геометрии и генерации сетки используется свободно распространяемая программа СМ8И [17]. Программная реализация численного расчета построена с использованием свободно распространяемой библиотеки ЕЕшСБ [18].

1. Математическая модель

Рассмотрим математическую модель термомеханического состояния конструкции с внутренними источниками тепла, которая описывается нелинейной системой уравнений для температуры и перемещений [2-4, 6, 7].

Пусть конструкция состоит из источников, газового зазора и оболочки О = О1 и О2 и О3, где О — подобласть источников, О2 — подобласть газового зазора и О3 — оболочка. Термомеханическое состояние описывается следующей связанной системой уравнений [19-22]: дТ д div и

= хеп, - а + 7 grad(T - Тг^) = 0, х е О, где С = ср — объемный коэффициент теплоемкости, с = с(х) — удельный коэффициент теплоемкости, р = р(х) —плотность, Т — температура, 7 = 7т(3А+2^),

Рис. 1. Схема области

Тт(ж) — температурный коэффициент линейного расширения, Л = Л(х), ц = ^(х) — параметры Ламе, и — вектор перемещения, к = к(х) — коэффициент теплопроводности, Q = Q(x) — мощность внутреннего источника тепла, и = с(и) — тензор наряжения и Tref — постоянная абсолютная температура, при которой тело находится в начальном равновесном состоянии.

Коэффициенты и внутренний источник для системы уравнений (1) зависят от подобласти и определяются следующим образом:

С(х) = Сг, А(ж) = А г, ^(ж) = к(х) = кг, Q(x) = Qi, х € , г = 1,2,3. (2)

В общем случае теплофизические свойства подобластей зависят и от других параметров, но в данной работе мы ограничимся только зависимостью от температуры [10,14, 20].

Система уравнений дополняется начальным условием

T(x, 0) = T0, x е tt, и граничными условиями для температуры

dT

-к^Г = a(T - Tout) dn

-к— dn

0,

x е гД,

x е rL и гт и Гв,

(3)

(4)

где a — коэффициент конвективного теплообмена, Tout — внешняя температура.

Для расчета напряженно-деформированного состояния будем использовать следующие граничные условия:

ui = 0, i = 1,..., d, x е гт,

an = g(x), x е rL U Гв U Гд U Г/ U Г0,

(5)

где й — размерность задачи.

Далее будем предполагать, что нагрев происходит в основном за счет внутреннего источника тепла, и будем пренебрегать слагаемым сИу ^ в уравнении (1). Такое предположение развязывает систему уравнений (1) и позволяет последовательно проводить расчет температуры и перемещений [8, 9,16].

Рассмотрим далее две модели для расчета напряженно-деформированного состояния.

Модель линейной упругости. В случае рассмотрения задач линейной упругости для тензора напряжения а и тензора деформаций е имеем

а (и) = И : е = 2 це(и) + ЛсИу(и) У, в (и) = — (grad и + grad ит) ,

где ^ — единичный тензор, А и р — упругие постоянные Ламе,

_ Е Еи

М= 2(1+г/)' (1 + г/)(1 — 2г/)

и Е — модуль упругости, V — коэффициент Пуассона.

Модель пластичности. В качестве примера решения задач нелинейной упругости [23, 24] рассматривается модель идеальной пластичности, для которой тензор напряжения зависит только от упругой части тензора деформаций ее:

а = Б : ее, е = ее + ер, (6)

где ер — тензор пластических деформаций.

Для классической модели идеальной пластичности тензор напряжения должен удовлетворять условию текучести Мизеса:

здесь / — пластический потенциал, ау — предел текучести.

В теории пластического течения скорость изменения пластической деформации задается с помощью введения параметра пластичности А как

дер = -д£ дЬ да'

Здесь член А определяет величину, а ^ задает направление пластической деформации.

С учетом уравнений (6) имеем

а = Б : (е — ер).

Таким образом, в общем случае задача нелинейна и для ее решения применяем метод Ньютона.

2. Конечно-элементная аппроксимация и вычислительный алгоритм

Для аппроксимации системы уравнений (1) будем использовать метод конечных элементов. Для аппроксимации по времени в уравнении для температуры будем использовать устойчивую неявную разностную схему. При аппроксимации нелинейных коэффициентов уравнения теплопроводности будем брать

коэффициенты для температуры с предыдущего временного слоя. Вариационная постановка задачи (1) записывается следующим образом: найти Т € V и и € Ш такие, что

ТП+1 _ ТП

Т п+1 _т п с

С(х, Тп)-г!х+ / к{х,Тп)УТ ■

_ 11(Т"+1 _ Тоиг)д ¿з = I Q(x)q ¿х, (7)

о о

пп+ 1

и

J ст(и) : е(у) ¿х + У 7grad(T _ Tref) • V ¿х = J д • «¿з, (8)

о о г^

где q € "К и V € Ш произвольные. Здесь V = V = Н1 и

Ш = Ш = {и € [Н^ : и = 0, х € },

— размерность пространства (^ = 2, 3), Н1 = Н 1(0) — пространство интегрируемых функций (V и V« конечно интегрируемы в О) [4,5]. Мы будем использовать линейные базисные функции.

При реализации вычислительного алгоритма расчет температуры будем проводить во всей области О, содержащей три подобласти О1, О2 и О3, для которых заданы коэффициенты. Расчет напряженно-деформированного состояния будем проводить только на двух подобластях О1 и О3. Отметим, что в данной математической модели мы не учитываем возможный контакт оболочки (подобласти О3) с источником (подобласть О1). При решении задачи упругости будем предполагать, что в подобласти зазора О2 задано некоторое давление рд и на внешней границе области О течет жидкость с заданными температурой и давлением

Тем самым алгоритм решения поставленной задачи можно представить следующим образом:

• строим геометрическую область О и генерируем расчетную сетку;

• генерируем сетку для подобласти источника О1 из общей сетки для области О;

• генерируем сетку для подобласти оболочки О3 из общей сетки для области О;

• для каждого временного слоя:

(1) вычисляем распределение температуры Т в области О = О1 и О2 и О3 с неоднородными коэффициентами;

(2) проводим расчет напряженно-деформированного состояния источника в подобласти О1 ;

(3) проводим расчет напряженно-деформированного состояния оболочки в подобласти О3;

(4) сохраняем результаты в файл;

г т

(5) переходим на следующий временной слой.

Представим результаты численного моделирования и разработанный программный интерфейс.

3. Численные результаты для осесимметричного случая

Приведем результаты численного моделирования задачи в двухмерной постановке для осесимметричного случая (см. рис. 1). Уравнения (1) запишем в виде

дТ 1 д ( дТ\ д ( дТ\

С(г,г,Т)~---— гк(г,г,Т)— -— к{г,г,Т)— \ =Я(г,г),

у ' ' 7 дЬ гдг\К,,'дг) дг\К,,'дг)^К,> (9)

— div а (и) + 7 grad(T — Т^) = 0

с начальным условием Т = То при Ь = 0 и 7 = в(3А + 2р). Отметим, что в данной модели мы не учитываем влияние перемещений на температуру и рассматриваем случай задачи линейной упругости.

В качестве граничных условий для температуры зададим

дТ дТ

-к(т)— = а(т-тоиг), (г,г)егд, -к{т)— = о, (г,г)егьигтигв,

дп дп

и для перемещений на области источника

иг = 0, аг = 0, (г, г) е , аг = 0, иг = 0, (г, г) е Гд,

—ап = рдп, (г, г) е дО1/(Г2 и Г^),

—ап = рдп, (г, г) е Г|, — ап = р0^п, (г, г) е Г|,

ау = 0, иг = 0, (г, г) е Г| и Г|.

Проведем численное моделирования рассмотренной задачи при следующих параметрах [12-14]:

Р1 = 65.3845 • 109 [Па], А1 = 98,0769 • 109 [Па], в1 = 10-5 [К-1]

Р3 = 27.7272 • 109 [Па], А3 = 40.404 • 109 [Па], вз = 2 • 10-5 [К-1].

Коэффициенты с^ и р^ (г = 1, 2, 3) задавались в виде функций от температуры и представлены на рис 2.

В качестве граничных условий использовались следующие параметры:

Тоиг = 373 [К], а = 7500

Вт

2К

Рд = 2 • 106 [Па], Роиг = 15 • 106 [Па].

Расчет проводился при ^ = 2 • 108 [Дж] с шагом по времени т = 0.1 [с], максимальное время Ьтах = 15 [с], То = 373 [К].

Моделирование проводилось для области со следующими размерами: радиус источника 4.3 • 10-3 м, высота источника 5.5 • 10-3 м, толщина зазора 0.03 •

Рис. 2. Коэффициенты кг, рг, сг, где г = 1, 2, 3 (слева направо)

1 = 2 сек. I = 8 сек. 1 = 15 сек.

Рис. 3. Распределение температуры в различные моменты времени £ = 2, 8 и 15 е

Рис. 4. Распределение перемещений по г и по г и напряжения для таблетки и оболочки

10 3 ми толщина оболочки 0.5 • 10 3 м (рис. 3 слева). Расчеты проводились на треугольной сетке, содержащей 19486 элементов и 9627 узлов (см. рис. 3 слева).

Результаты расчетов представлены на рис. 3,4. Рис. 3 иллюстрирует нагрев подобласти источника в разные моменты времени за счет наличия внутреннего источника тепла и охлаждения за счет конвективного теплообмена с окружающей средой посредством задания граничных условий третьего рода. Возникающие перемещения и напряжение в источнике и в оболочке за счет тепловых расширений проиллюстрированы на рис. 4 и представлены на конечный момент времени £тах = 15 е. При этом представленные перемещения достаточно малы. В случае больших тепловых расширений необходимо будет использовать другие модели для описания упругого состояния, в том числе с учетом возможного кон-

Рис. 5. Геометрия и расчетная область в разрезе. Распределение температуры в разрезе для различных временных слоев Ь = 2 сек, Ь = 8 сек, Ь =15 сек

такта источника с оболочкой. В данном случае мы проводили моделирование перемещений в рамках базовой модели линейной упругости.

4. Численное моделирование трехмерной задачи

Решим трехмерную задачу термоупругости с внутренним источником с несимметричной геометрией. Геометрия и расчетная область в разрезе представлены на рис. 5.

Для численных вычислений использовались следующие значения коэффициентов: для расчета перемещений

Е\ = 170 ■ 109 [Па], VI = 0.3, в = 10-5[К-1],

Е2 = 60 ■ 109 [Па], »2 = 0.32, в2 = 2 ■ 10-5[К-1]

и теплопереноса

(ki ,k2 ,кз) = (7,0.2,23)

Вт

м • К

(ci,c2,c3) = (300,5200, 300)

Дж

кг • К

(Р1 ,Р2,Рз) = (104,0.1,6550)

Моделирование проводилось при Qf = 2 • и г = 0.1 [с] на tr

Дж

= 15 с.

В качестве граничных условий задавалось a = 7500 ) Tout = Tq = 373[А']

для температуры и pg = 2 • 106 [Па], pout = 15 • 106 [Па] для перемещений.

Результаты расчетов распределения температуры для различных временных слоев t = 2 сек, t = 8 сек, t = 15 сек представлены в разрезе на рис. 5.

Рис. 6. Распределения модуля перемещений в срезе для значений предела текучести 1.3 ■ 108, 1.2 ■ 108 и 1.1 ■ 108 Па(слева направо) в момент времени Ь = 15 сек

Рис. 7. Распределения модуля перемещений в срезе решений задач пластичности для значений предела текучести =3 ■ 108 Па, 100 ■ 108 Па и для задачи линейной упругости в момент времени Ь = 15 сек.

Для определения влияния величины предела текучести на нелинейное решение задачи пластичности было проведено исследование распределения модуля перемещения для различных значений предела текучести около предельной величины. Распределения модуля перемещений в срезе для значений предела текучести, равных 1.3 • 108, 1.2 • 108 и 1.1 • 108 Па, представлены на рис. 6 и иллюстрируют увеличение перемещений при уменьшении предела текучести. Из рисунков распределения перемещений видно, что при приближении к некоторому предельному значению предела текучести непропорционально увеличивается модуль перемещения и оболочка сильно разбухает, что показывает

целесообразность использования материалов с большим значением предела текучести. Также из используемой модели пластичности следует, что при увеличении предела текучести решение нелинейной задачи должно стремиться к решению задачи линейной упругости, что и показывают результаты на рис. 7, где представлены сравнения распределения модуля перемещений при значениях предела текучести, равных ay =3 • 108 Па, 100 • 108 Па, и распределение модуля перемещений для линейной задачи.

5. Заключение

В работе представлена модель термоупругого состояния конструкции с внутренним источником. Вычислительная реализация построена с использованием метода конечных элементов. Результаты моделирования представлены для двумерной осесимметричной постановки и для трехмерного случая.

В дальнейшем планируется рассмотреть термоупругое состояние с учетом наличия трещин и течения газа, которые существенно изменяют теплофизиче-ские свойства теплового источника и, следовательно, могут привести к большим градиентам температур в местах особенностей задачи.

ЛИТЕРАТУРА

1. Newman C., Hansen G., Gaston D. Three dimensional coupled simulation of thermomechanics, heat, and oxygen diffusion in UO2 nuclear fuel rods //J. Nuclear Materials. 2009. V. 392, N 1. P. 6-15.

2. Williamson R. L., Hales J. D., Novascone S. R., Tonks M. R., Gaston D. R., Permann C. J., Andrs D., Martineau R. C. Multidimensional multiphysics simulation of nuclear fuel behavior // J. Nuclear Materials. 2012. V. 423, N 1. P. 149-163.

3. Kang C. H., Lee S. U., Yang D. Y., Kim H. C, Yang Y. S. 3D FE simulation of the nuclear fuel rod considering the gap conductance between the pellet and cladding // Proc. KNS Fall Meeting (Kyungju, Rep. Korea, Oct. 23-25, 2013). Daejeon, Rep. Korea: KNS, 2013.

4. Kang C. H., Lee S. U., Yang D. Y., Kim H. C., Yang Y. S. 3D finite element analysis of a nuclear fuel rod with gap elements between the pellet and the cladding // J. Nuclear Sci. Technology. 2016. V. 53, N 2. P. 232-239.

5. Philip B., Berrill M. A., Allu S., Hamilton S. P., Sampath R. S., Clarno K. T., Dilts G. A. A parallel multi-domain solution methodology applied to nonlinear thermal transport problems in nuclear fuel pins // J. Comput. Phys. 2015. V. 286. P. 143-171.

6. Ramirez J. C., Stan M., Cristea P. Simulations of heat and oxygen diffusion in UO2 nuclear fuel rods // J. Nuclear Materials. 2006. V. 359, N 3. P. 174-184.

7. Mihaila B., Stan M., Ramirez J., Cristea P. Simulations of coupled heat transport, oxygen diffusion, and thermal expansion in UO2 nuclear fuel elements // J. Nuclear Materials. 2009. V. 394, N 2. P. 182-189.

8. Brown D. L., Vasilyeva M. A. A generalized multiscale finite element method for poroelasticity problems II: Nonlinear coupling // J. Comput. Appl. Math. 2016. V. 297. P. 132-146.

9. Brown D. L., Vasilyeva M. A. A generalized multiscale finite element method for poroelasticity problems I: Linear problems // J. Comput. Appl. Math. 2016. V. 294. P. 372-388.

10. Hales J. D. et al. BISON theory manual. The equations behind nuclear fuel analysis. Idaho Nat. Lab., 2013.

11. Rashid Y., Dunham R., Montgomery R. Fuel analysis and licensing code: FALCON MOD01. EPRI Rep. EPRI, 2004.

12. Veshchunov M. S. et al. Code Package SVECHA: Modeling of core degradation phenomena at severe accidents // Proc. 7th Int. Topical Meeting on Nuclear Reactor Thermal Hydraulics NURETH-7 (Saratoga Springs, NY, Sept. 10-15, 1995). 1995. P. 1914.

13. Berdyshev A. V., Boldyrev A. V., Palagin A., Shestak V., Veshchunov M. S. SVECHA/QUENCH code for the modeling of reflooding phenomena in severe accidents conditions // Proc. 9th Int. Topical Meeting on Nuclear Reactor Thermal Hydraulics (NURETH-9) (San Francisco, CA). 1999.

14. Hagrman D. L., Reymann G. A. MATPRO-VERSION 11. Handbook of materials properties for use in the analysis of light water reactor fuel rod behavior. Idaho Falls (USA): Idaho Nat. Eng. Lab., 1979.

15. Hales J. D. et al. Asymptotic expansion homogenization for multiscale nuclear fuel analysis // Comput. Materials Sci. 2015. V. 99. P. 290-297.

16. Вабищевич П. Н., Васильева М. В., Колесов А. Е. Схема расщепления для задач поро-упругости и термоупругости // Журн. вычисл. математики и мат. физики. 2014. Т. 54, № 8. C. 1345-1355.

17. Software GMSH. Available at: http://geuz.org/gmsh/ (accessed 06.05.2016).

18. Logg A., Mardal K. A., Wells G. (ed.). Automated solution of differential equations by the finite element method: The FEniCS book. New York: Springer Sci. & Business Media, 2012.

19. Вабищевич П. Н., Самарский А. А. Вычислительная теплопередача. М.: Едиториал УРСС, 2003.

20. Васильева М. В., Стальнов Д. А. Математическое моделирование термомеханического состояния тепловыделяющего элемента // Вестн. СВФУ. 2016. № 1. C. 45-59.

21. Вабищевич П. Н., Васильева М. В. Численное моделирование задач термоупругости // Вестн. СВФУ. 2013. Т. 10, № 3. C. 5-9.

22. Sivtsev P. V., Vabishchevich P. N., Vasilyeva M. V. Numerical simulation of thermoelasticity problems on high performance computing systems // Proc. Int. Conf. Finite Difference Methods. Berlin: Springer, 2014. P. 364-370.

23. Simo J. C., Hughes T. J. R. Computational inelasticity. New York: Springer, 1998. (Interdis-cip. Appl. Math.; V. 7).

24. De Souza Neto E., Peric D., Owen D. R. J. Computational methods for plasticity: Theory and applications. New York: John Wiley & Sons, 2008.

Статья поступила 7 августа 2017 г.

Васильева Мария Васильевна, Захаров Петр Егорович, Сивцев Петр Васильевич, Спиридонов Денис Алексеевич Северо-Восточный федеральный университет им. М. К. Аммосова, Международная научно-исследовательская лаборатория

«Многомасштабное математическое моделирование и компьютерные вычисления», ул. Кулаковского 42, Якутск 677891

[email protected], [email protected], [email protected], [email protected]

Математические заметки СВФУ Июль—сентябрь, 2017. Том 24, № 3

UDC 519.63

NUMERICAL MODELING OF THERMOELASTICITY

PROBLEMS FOR CONSTRUCTIONS

WITH INNER HEAT SOURCE

M. V. Vasilyeva, P. E. Zakharov, P. V. Sivtsev, and D. A. Spiridonov

Abstract: We consider the numerical simulation of the thermomechanical state of a structure consisting of a heat source, a gas gap and a shell. The mathematical model is described by a nonlinear system of equations for temperature and displacements. The heat is released in the subdomain of the heat source. The resulting displacements due to the temperature gradient are calculated in the heat source region and separately in the shell, and can be described by both linear elasticity models and nonlinear plasticity models. The numerical implementation is based on the finite element method. The results of numerical modeling of a nonlinear model problem in two- and three-dimensional domains are presented.

DOI 10.25587/SVFU.2018.3.10889 Keywords: thermoelasticity problems, thermal expansion, heat transfer, linear elasticity problem, plasticity models, nonlinear problems, finite element method, mathematical modeling.

REFERENCES

1. Newman C., Hansen G., and Gaston D., "Three dimensional coupled simulation of thermome-chanics, heat, and oxygen diffusion in UO2 nuclear fuel rods," J. Nuclear Materials, 392, No. 1, 6-15 (2009).

2. Williamson R. L., Hales J. D., Novascone S. R., Tonks M. R., Gaston D. R., Permann C. J., Andrs D., and Martineau R. C., "Multidimensional multiphysics simulation of nuclear fuel behavior," J. Nuclear Materials, 423, No. 1, 149-163 (2012).

3. Kang C. H., Lee S. U., Yang D. Y., Kim H. C, and Yang Y. S., "3D FE simulation of the nuclear fuel rod considering the gap conductance between the pellet and cladding," in: Proc. KNS Fall Meeting (Kyungju, Rep. Korea, Oct. 23-25, 2013). KNS, Daejeon, Rep. Korea (2013).

4. Kang C. H., Lee S. U., Yang D. Y., Kim H. C., and Yang Y. S., "3D finite element analysis of a nuclear fuel rod with gap elements between the pellet and the cladding," J. Nuclear Sci. Technology, 53, No. 2, 232-239 (2016).

5. Philip B., Berrill M. A., Allu S., Hamilton S. P., Sampath R. S., Clarno K. T., and Dilts G. A., "A parallel multi-domain solution methodology applied to nonlinear thermal transport problems in nuclear fuel pins," J. Comput. Phys., 286, 143-171 (2015).

6. Ramirez J. C., Stan M., and Cristea P., "Simulations of heat and oxygen diffusion in UO2 nuclear fuel rods," J. Nuclear Materials, 359, No. 3, 174-184 (2006).

The authors were supported by the by the Government of Russian Federation Megagrant No. 14.Y26.31.0013.

© 2017 M. V. Vasilyeva, P. E. Zakharov, P. V. Sivtsev, D. A. Spiridonov

7. Mihaila B., Stan M., Ramirez J., and Cristea P., "Simulations of coupled heat transport, oxygen diffusion, and thermal expansion in UO2 nuclear fuel elements," J. Nuclear Materials, 394, No. 2, 182-189 (2009).

8. Brown D. L. and Vasilyeva M. A., "A generalized multiscale finite element method for poroelasticity problems II: Nonlinear coupling," J. Comput. Appl. Math., 297, 132-146 (2016).

9. Brown D. L. and Vasilyeva M. A., "A generalized multiscale finite element method for poroelas-ticity problems I: linear problems," J. Comput. Appl. Math., 294, 372-388 (2016).

10. Hales J. D. et al., BISON Theory Manual. The Equations Behind Nuclear Fuel Analysis, Idaho Nat. Lab. (2013).

11. Rashid Y., Dunham R., and Montgomery R., Fuel Analysis and Licensing Code: FALCON M0D01, EPRI Rep., EPRI (2004).

12. Veshchunov M. S. et al., "Code Package SVECHA: Modeling of core degradation phenomena at severe accidents," in: Proc. 7th Int. Topical Meeting on Nuclear Reactor Thermal Hydraulics NURETH-7 (Saratoga Springs, NY, Sept. 10-15, 1995), 1914 (1995).

13. Berdyshev A. V., Boldyrev A. V., Palagin A., Shestak V., and Veshchunov M. S., "SVECHA/ QUENCH code for the modeling of reflooding phenomena in severe accidents conditions," in: Proc. 9th Int. Topical Meeting on Nuclear Reactor Thermal Hydraulics (NURETH-9) (San Francisco, CA) (1999).

14. Hagrman D. L. and Reymann G. A., MATPRO-VERSION 11. Handbook of Materials Properties for Use in the Analysis of Light Water Reactor Fuel Rod Behavior, Idaho Nat. Eng. Lab., Idaho Falls (USA) (1979).

15. Hales J. D. et al., "Asymptotic expansion homogenization for multiscale nuclear fuel analysis," Comput. Materials Sci., 99, 290-297 (2015).

16. Vabishchevich P. N., Vasilyeva M. V., and Kolesov A. E., Splitting scheme for poroelasticity and thermoelasticity problems," Comput. Math. Math. Phys., 54, No. 8, 1305-1315 (2014).

17. Software GMSH, Available at: http://geuz.org/gmsh/ (accessed 06.05.2016).

18. Logg A., Mardal K. A., and Wells G. (eds.), Automated Solution of Differential Equations by the Finite Element Method: The FEniCS Book, Springer Sci. & Business Media, New York (2012).

19. Samarskij A. A. and Vabishhevich P. N., Vychislitel'naia Teploperedacha [in Russian], Editorial URSS, Moscow (2003).

20. Vasilyeva M. V. and Stal'nov D. A., "Mathematical modeling of the thermomechanical state of a heat-inducing element [in Russian]," Vestn. SVFU, No. 1, 45-59 (2016).

21. Vabishhevich P. N. and Vasilyeva M. V., "Numerical modeling for thermoelasticity problems [in Russian]," Vestn. SVFU, No. 3, 5-9 (2013).

22. Sivtsev P. V., Vabishchevich P. N., and Vasilyeva M. V., "Numerical simulation of thermoelasticity problems on high performance computing systems," in: Proc. Int. Conf. Finite Difference Methods, pp. 364-370, Springer, Berlin (2014).

23. Simo J. C. and Hughes T. J. R., Computational Inelasticity, Springer, New York (1998). (Interdiscip. Appl. Math.; V. 7).

24. De Souza Neto E., Peric D., and Owen D. R. J., Computational Methods for Plasticity: Theory and Applications, John Wiley & Sons, New York (2008).

SSubmitted August 7, 2017

Maria V. Vasilyeva, Petr E. Zakharov,

Petr V. Sivtsev, and Denis A. Spiridonov

M. K. Ammosov North-Eastern Federal University,

Institute of Mathematics and Informatics,

42 Kulakovsky Street, Yakutsk 677000, Russia

[email protected], [email protected], [email protected], [email protected]