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

ЕСТЕСТВЕННЫЕ И ТОЧНЫЕ НАУКИ

УДК 621.94.084

РЕШЕНИЕ НЕОСЕСИММЕТРИЧНОЙ ЗАДАЧИ ТЕРМОУПРУГОСТИ ДЛЯ НЕРАВНОМЕРНО НАГРЕТОГО ДЛИННОГО ЦИЛИНДРА

В УСЛОВИЯХ ПОЛЗУЧЕСТИ

Докт. физ.-мат. наук, проф. КУЛИКОВ И. С., магистрант ШИРВЕЛЬ П. И.

Белорусский национальный технический университет

Одним из главных факторов, определяющих работоспособность ядерного реактора АЭС, являются температурные условия работы его тепловыделяющих элементов (ТВЭЛ), имеющих, как правило, форму длинных цилиндров. Жесткость условий работы ТВЭЛов (большие тепловые и радиационные нагрузки, высокое внешнее давление теплоносителя) предполагает повышенные требования к ним. Это требует от механики деформируемого твердого тела развития методов определения напряженно-деформированного состояния (НДС) тел цилиндрической геометрии, учитывающих влияние нелинейных деформаций.

Одной из причин появления напряжений в сплошном теле является именно неравномерный нагрев. Заметим, что ТВЭЛы активной зоны реактора работают при высоких температурах. Поэтому расчетно-теоретическое исследование НДС ТВЭЛов и включает в себя в первую очередь определение температурных полей и возникающих термонапряжений. При этом предполагаем, что теплофизические и механические характеристики материала остаются неизменными (или меняются в соответствии с полем температуры). Это относится и к прочностным свойствам материала. Абсолютное значение и характер распределения поля температуры определяются мощностью внутреннего тепловыделения, теплофизическими свойствами материала и условиями теплосъема с поверхности цилиндра. Заметим, что ТВЭЛы чаще всего представляют собой длинные ци-

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

Исследуем неосесимметричное НДС однородного бесконечно длинного, сплошного цилиндра (рис. 1), находящегося в неравномерном температурном поле T(r, 0) и подвергающегося действию приложенной нагрузки P (внешнее давление) с учетом деформации тепловой ползучести. В реальности это модель НДС топливного цилиндрического сердечника ТВЭЛа активной зоны реактора АЭС. Уравнение равновесия для такой модели имеет вид:

дог + 1 Эа ге

r de

, 1Э°8 r эе

r до

r

= 0;

r

+ 2-

r

^ = 0. r

(1)

Компоненты тензора деформаций связаны с компонентами тензора напряжений зависимостями (2):

r

Е

Е

(2)

Е

1 + V

-I

£

'ге

Р

Рис. 7. Модель топливного сердечника

Предполагая условие плоского деформирования, принимаем: £_ = 0; егв = = ег2 = = е2Г =0. Тогда из (2) имеем ог, о0, ог, огв. Так как деформация ползучести происходит при постоянном объеме материала, компоненты ес связаны условием несжимаемости: £сг + £сн +

+ есг —0. Использование условия несжимаемости упрощает запись выражений для главных напряжений. В результате компоненты тензора напряжений имеют вид:

Е

=-

' (1 + у)(1-2У)

х

х((1 - у)ег + У£е - < (1 + V) - есг (1 - 2у));

Е

(1 + у)(1-2У) х((1 - У)£е + уег - < (1 + V) - (1 - 2у)); Е

(1-2у)(1 + У) х у(£г+£е)-<(1 + у) + «+£се)(1-2у) ;

(3)

аг0 = °вг =

Е

1 + V

еге £ге •

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

В (3) компоненты тензора деформации связаны с компонентами вектора перемещений геометрическими уравнениями (соотношения Коши):

_ Эи _ и 1 ЭтЭ-

£'" + _За '

дг г г аЬ

~~ ~~ 2

\ди_ + Эд г Э6 дг

г

(4)

У

где и, д - компоненты вектора перемещения в радиальном и окружном направлениях.

Объемные (температурные) деформации в цилиндре определим следующим образом: £гг=егв=е1=£ггв=ег =аТ(г,В), где Т{г, 9) -заданная функция в [1, 2]. Соответственно для

цилиндра Л>,е) = 75(е) + -^(Я2-г2), где Я _

4Х

радиус цилиндра; X - коэффициент теплопроводности материала цилиндра; qv - тепловыделение в единице объема материала; а - коэффициент линейного расширения материала цилиндра; У\(9) - температура поверхности цилиндра, определяется из задачи теплообме-

г

г

на между теплоносителем и поверхностью цилиндра.

Зададим механические (необратимые) деформации температурной ползучести. При известном законе изменения интенсивности скорости ползучести гси = /(ом, /. Т, Ф) деформации ползучести на каждом временном шаге определим из выражений (5):

Рс - Рс -4- А Рс • 1) 1 111] ■

■ с _ 3 ¿и(п-1) . Ч/(и-1) 1 с у ' аи(и-1) (5)

£си=Дои,^ Т, Ф);

/з

V 2 '

где = о,; - 5,;а; о - : 5(; - символ Кроне-

кера; /, / = г, 0: Ф = ф/ - интегральный нейтронный поток; ф- плотность нейтронного потока; ^ - время. Отметим, что в случае одноосного напряженного состояния согласно (5) интенсивность скорости ползучести равна скорости ползучести материала: гси = = = /(ом, /, Т, Ф) = Ао„"', что значительно упрощает процедуру расчета деформаций ползучести.

Учитывая соотношения Коши (4), перейдем от деформаций к перемещениям в полученных выражениях (3) и, подставляя эти соотношения в уравнения равновесия (1), получаем систему уравнений равновесия, выраженную через перемещения (6):

Э 2и 1 ди Эг2 г Эг

(

-гм +-:

г2 2г(1-у)

1 - 2у д и д2$ 1 , 0 ЭдЛ --7 +-+ - 4У-3 —

г эе2 дгдв г эе

V У

1 + удг1 1-2у

1—V Э г

(1-у)г

эе 1

Эе'о Эе"

эе э г

д2Ь 1 эа

Эг2 г Эг

^ 2(1 — V) Э2д Э2и 3 - 4у Эм Л +-+

(6)

г( 1-2У)

х

г

„Т

эе2 ЭгЭе

эе

2 1 + у дет

~ дг г 1-2у эе

+ £г0. г

Определение НДС поставленной задачи сводится к решению системы уравнений равновесия (6) со следующими граничными условиями:

и = 0 при г = 0.

(7)

Ввиду симметрии напряженно-деформированного состояния относительно радиальной плоскости, так как температурное поле цилиндра симметрично (0 < 9 < 9°), справедливы также условия:

с)м

ф=о, — = 0 при е = о ие = е°, (8)

где 9° может принимать различные значения в зависимости от закона распределения температуры по периметру (период функции распределения температуры).

Принимая, что в начале координат £г ~ £(). и учитывая условие (8), получаем:

д = 0 при г = 0;

(9)

1

ке ег

г

4

ог - -Р при г = 0. (10)

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

ний известно лишь для единичных случаев. Поэтому решение нелинейных задач сводится к использованию численных методов, т. е. к тем результатам, которые получены после приближения-замены исходного уравнения в области существования решения. Главное преимущество этих методов заключается в том, что они позволяют получить решение задачи в любом случае (и тогда, когда неизвестно аналитическое решение). Заметим, что частные случаи решения уравнения теплопроводности для различных тел как с источниками теплоты, так и без них, рассмотрены во многих работах. Однако эти решения могут быть получены лишь в том случае, если функция распределения источников теплоты допускает интегрирование уравнения. Примером может служить вывод формул для термических напряжений в полом цилиндре [3]. Недостаток таких решений температурных задач состоит еще и в том, что эти формулы получены из интегральных соотношений, выведенных применительно к телам, не испытывающим внешних нагрузок и влияния температурной ползучести. В действительности же топливный сердечник испытывает напряжение не только вследствие нагрева, но и из-за действия внешнего давления со стороны газового зазора и ползучести материала. Поэтому аналитические интегральные формулы решения нельзя принять для проводимого расчета.

Для интегрирования уравнений (6) используем метод конечных разностей. В частности, в [4] предлагался численный метод решения подобной задачи. Суть разработанного метода состоит в следующем. На первом этапе разрешаем неосесимметричную задачу конечно-

разностным методом при значени"

il=0

равных

l,m l,m-1 — О

-Ц.-1 =0

и° =0, fpj =0

78 = д.? = о

нулю. Для этого в области изменения переменных 0 < 9 < 9°, 0<r <R строим сетку из линий 9 = const и /• = const (рис.2).

7 —> Г

/->е

Разобьем otj/cjur 0, R равноотстоящими

Л

точками: г0 = 0, rn=R, /• = ihr, h. =—.

/ =

= -1, 0,1,...,«. Аналогично отрезок 9 равноотстоящими точками: 90 = 0, 0/н = 9°; 9; = , 9"

\ =—, 7 = 0, 1, ..., т, где Ъг и Ае- шаги сет-т

ки в радиальном и окружном направлениях. Выпишем неявную разностную схему с весовыми коэффициентами [5]. Это значит, что для вычисления производных по 9 в точке г, разностной схемы используются значения сеточных функций на (п + 1)-м и (п - 1)-м слоях по

радиусу. Так,

Ъ\

{г, 9) = А0-2Л2(аи'+1 +ри'-1).

с)9 .......' ''

причем а + (3 = 1. Производные по радиусу и угловой координате аппроксимируем конечно-разностными выражениями. Для первого уравнения системы (6) имеем:

2V-1 Э2-

Э2-

13« 1 1

---гИ+-

V V

г Эг г2 2r(l-v)

X

X

(Э2д 1 Ш

-+ — (4v - 3)—

ЭгЭ9 г Э9

л\

-С1(ея)

/у

(12)

Рис. 2. Схема численного определения НДС

п

c

2У-1

2г2(1-У)

/ге 2Л2 (а и /+1 + 0и ) + ^ (м, А) =

и'/1 -2и)+и1'1

(14)

и1;1 - 2 и) + и'71 = ^ 1 X

1 3 3 2г (1-у)

х(аи 1+1 + рм/'-1) +

(15)

XI =

2у-1 /г2 2г2(1-У) \

2У-1 /г2 2г2(1-У) \

-Р. (16)

Соответствующая система, отвечающая первому уравнению равновесия (6), может быть записана в следующем виде:

= -2«; +«;.-1(1 -зс2л2) - ^ (17)

Аналогично можно представить второе уравнение равновесия системы (6):

= -2Ц +Ъ1~1(\-ЪА^-Р^и,(18)

$1 =

2(у -1) /г2 г2(1-2У) ^

■а; =

г г

\\

1ЭА 1 а

дг г

2(У~1)

г\ 1-2У) % 1

■Р; (19)

г(1-2У)

-х

эгэе г эе

\Л

-С2(ея)

/у

(20)

С2(гн) - 2а—+

Э Т 2 1 + V

Э г гм 1-2у

-х

дТ 4 ха— +—(аТ + е^). ив г^.

(21)

Слагаемые !•]'. (и. 1'}) нелинейные относительно и';, г)';, поэтому для решения уравнений

непосредственно неприменимы методы, развитые для линейных систем. Однако решение может быть построено по методу итераций [6]. На каждом шаге последовательных приближений величины Р1]1 (и, тЭ) вычисляются по зна-

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

Граничные условия (7)-(10) в разностной форме примут вид:

М;°=0, А°=0, т^о = 0, ^=0,

и\-иг0= 0, и'„ =0;

п п-1 и1 ~и1

(1—у) + у

'и) Л

Я

ма

(22)

- аТ" (1 + у) + ес" (1 - 2у) -

Р(1 + у)(1-2У) Е '

Начальная система уравнений (6) стала эквивалентной системам (17) и (18) с граничными условиями (22). Системы уравнений (17), (18) имеют матрицы трехдиагональной структуры и могут быть решены по методу прогонки [7], где , %2г, , £,2г зависят только от текущего радиуса г, а правые части уравнений (17), (18) - функции от перемещений на предыдущих, уже известных, слоях радиуса (/, г — 1) , т. е. не содержат независимые переменные (и- 1. т)'; 1) для рассматриваемого (г +1) -го слоя

по радиусу и являются экстраполированными числовыми значениями.

Иллюстрация схемы численного решения представлена на рис. 2. Расчет начинается с 1 — 0. Учитывая граничные условия (22)

и и ^ = 0; х)/1 ~ 0, вычисляем сначала (и, г)). а потом - правую часть (17). Затем разрешаем (17) методом прогонки вдоль кривой 0.. Найдя значения » , переходим к решению

(18): вычисляем /ч', (и, 1')). правую часть (18), и решаем это уравнение вдоль кривой 6 ■ относительно Фу. После первого слоя переходим ко

второму - устанавливаем / = 1 и проделываем действия, аналогичные описанным выше. Подобным образом, двигаясь по радиусу г от

2

Л

Г

Г

Г

к

Г

г

предыдущих слоев к текущему и «пробегая» вдоль кривой 0,, методом прогонки от / = О

до у=/и осуществляем вычисления и'- 1 по

первому уравнению равновесия и г'}',11 - по

второму, до г = п — 1. Таким образом, находим численное неосесимметричное решение задачи. Зная вектор перемещения точек цилиндра и'- 1

и г)',11. по (4) и (3) можем построить тензоры

деформаций и напряжений в любой точке по периметру цилиндра с дальнейшим пересчетом характеристик НДС через промежуток А/ для учета температурной ползучести по (5).

В Ы В О Д

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

неполным без учета облучения. Отсюда вытекает необходимость моделирования расчета НДС длинного сплошного цилиндра, в котором будет рассматриваться поведение материалов ТВЭЛов при облучении. А полученные в работе выражения и предложенная численная схема решения задачи определения неосесимметрич-ного НДС длинного цилиндра, подверженного неравномерному нагреву в условиях температурной ползучести, может послужить основой для разработки и проведения подобных практических расчетов на прочность.

Л И Т Е Р А Т У Р А

1. Куликов, И. С. Прочность элементов конструкций при облучении / И. С. Куликов, Нестеренко, Б. Е. Твер-ковкин. - Минск: Наука и техника, 1990. - 144 с.

2. Куликов, И. С. Прочность тепловыделяющих элементов быстрых газоохлаждаемых реакторов / И. С. Куликов, Б. Е. Тверковкин. - Минск: Наука и техника, 1984. - 104 с.

3. Тимошенко, С. П. Теория упругости / С. П. Тимошенко, Дж. Гудьер. - М., 1979. - 551 с.

4. Ширвель, П. И. О неосесимметричном НДС неравномерно нагретого длинного сплошного цилиндра, подверженного нейтронному облучению / П. И. Ширвель, И. С. Куликов // Машиностроение: респ. межвед. сб. -Минск, 2008. - Вып. 24. - Т. 1. - С. 185-191.

5. Яненко, Н. Н. Метод дробных шагов для решения многомерных задач математической физики / Н. Н. Яненко. - Новосибирск: Наука, 1967. - 195 с.

6. Вольмир, А. С. Оболочки в потоке жидкости и газа /

A. С. Вольмир. - М., 1976. - 416 с.

7. Годунов, С. К. Разностные схемы / С. К. Годунов,

B. С. Рябенький. - М.: Наука, 1973. - 440 с.

Поступила 03.03.2009

УДК 517.4

ПРЕОБРАЗОВАНИЕ МЕЛЛИНА КАК ЧАСТНАЯ РЕАЛИЗАЦИЯ ОБЩЕЙ СХЕМЫ ПОСТРОЕНИЯ ИНТЕГРАЛЬНЫХ ПРЕОБРАЗОВАНИЙ

Канд. физ.-мат. наук ГАХОВИЧА. С.

Белорусский национальный технический университет