CC BY
18
мума на нижнем основании указанного параллелепипеда или его боковых границах, где в силу условия (16)-(17), будет выполняться
H^k){x,y,t)> c0 = const> 0, к, к = 1,2,...,N. Заключение
Авторами предложена нелинейная математическая пространственно-двумерная модель транспорта наносов. Исследование указанной модели проводится путем ее линеаризации. Получены условия положительности решений линеаризованной начально-краевой задачи, представленные в виде соответствующей теоремы. Кроме того, получена априорная оценка положительного решения.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Sukhinov, A.I., Chistyakov, A.E., Alekseenko. E.V. Numerical realization of the three-dimensional model of hydrodynamics for shallow water basins on a high-performance system // Mathematical Models and Computer Simulations, October 2011, Volume 3, Issue 5, pp 562-574. doi:10.1134/S2070048211050115.
2. Alekseenko, E., Roux, B., Sukhinov, A., Kotarba, R., Fougere, D. 2013. Coastal hydrodynamics in a windy lagoon. Computers and Fluids, 2013, 77, pp. 24-35 doi:10.1016/j.compfluid.2013.02.003.
3. Леонтьев, И.О. Прибрежная динамика: волны, течения потоки наносов. - М.: ГЕОС, 2001. - 272 с.
4. Сидорякина В.В., Сухинов А.И. Существование и единственность решения линеаризованной двумерной задачи транспорта наносов // В книге: Теория операторов, комплексный анализ и математическое моделирование. Тезисы докладов XIII Международной научной конференции. 2016. - С. 184-185.
5. Сухинов А.И., Сидорякина В.В. О единственности решения линеаризованной двумерной начально-краевой задачи транспорта наносов // Вестник Таганрогского института имени А.П. Чехова. 2016. № 2.- С. 270-274.
6. Сухинов А.И., Чистяков А.Е., Проценко Е.А. Математическое моделирование транспорта наносов в прибрежной зоне мелководных водоемов//Математическое моделирование. 2013. N. 25, № 12, С.65-82.
7. Сухинов А.И., Чистяков А.Е., Проценко Е.А. Математическое моделирование транспорта наносов в прибрежных водных системах на многопроцессорной вычислительной системе // Выч. мет. Программирование. 2014. 15:4, 610-620.
8. Sukhinov, A.I., Chistyakov, A.E. Adaptive modified alternating triangular iterative method for solving grid equations with a non-self-adjoint operator// Mathematical Models and Computer Simulations, July 2012, Volume 4, Issue 4, pp 398-409. doi:10.1134/S2070048212040084.
9. Сухинов А.И., Проценко Е.А., Чистяков А.Е., Шретер С.А. Сравнение вычислительных эффективностей явной и неявной схем для задачи транспорта наносов в прибрежных водных системах // Вычислительные методы и программирование: новые вычислительные технологии. 2015. Т. 16. № 3. - С. 328-338.
10. Годунов С. К. Уравнения математической физики. — М., 1979. - 392 с.
И.В. Яковенко
ОПТИМИЗАЦИЯ ПАРАМЕТРОВ ДИСКРЕТНОЙ МАТЕМАТИЧЕСКОЙ МОДЕЛИ ТРАНСПОРТА ТЕПЛА В ГАЗОТУРБИННОЙ УСТАНОВКЕ
Аннотация. Целью работы является расчет параметров математической модели транспорта тепла в паровой турбине на основе вычислительных методов, позволяющих значительно повысить точность расчетов. Практическая значимость состоит в разработке и программной реализации модели тепловых процессов на корпусе газотурбинной установки.
Ключевые слова: математическая модель, термодинамика, паровая турбина, вычислительные эксперименты.
I.V. Yakovenko
OPTIMIZATION OF THE PARAMETERS OF THE DISCRETE MATHEMATICAL MODEL OF TRANSPORT OF HEAT IN THE GAS TURBINE INSTALLATION
Abstract. The aim of this work is to calculate the parameters of the mathematical model for heat transport in PA-world turbine-based computational methods to significantly improve the accuracy of calculations. Practical value consists in development and software implementation of models of thermal processes on the housing of the gas turbine installation.
Key words: mathematical model, thermodynamics, steam turbine, computational experiments.
Теория теплопроводности широко используется для решения различных производственных задач, связанных с определением температурных полей и тепловых поток. Как правило, используются методы математической физики, рассматривающие протекание наблюдаемого процесса в произвольно выбранном элементарном объеме в течение бесконечно малого промежутка времени. Такой подход позволяет учитывать только существенно вляющие величины и значительно упростить вычислительный процесс.
Первоначально тело принимается однородным, изотропным, с постоянными физическими характеристиками, такими как X, с (теплоемкость) и р (плотность), и равномерно распределенными в нем внутренними источниками тепла. Источники тепла характеризуются плотностью теплового потока - величиной qv.
Одним из основных уравнений, которое используется при описании теплового процесса, является уравнение теплопроводности. Это уравнение выведено в соответствии с законом сохранения энергии: вся теплота, выделенная внутренними источниками dQeн и внесенная извне в элементарный объем путем теплопроводности dQm за время dt, меняет внутреннюю энергию вещества, содержащегося в этом объеме:
dQm + dQm = dU . (1)
Произвольным образом выделим элементарный объем тела в виде параллелепипеда с ребрами dx, dy, dz, параллельными осям системы координат Оxyz. Тогда количество теплоты, проходящее внутрь выделенного объема в направлении оси Ох через элементарную площадку dy•dz за время dt, равно:
дТ
dQх1 = Чх • dy ■ dz ■ dt = —Х---Су . (2)
дх
При этом температура на противоположной грани параллелепипеда получит приращение
--dx и будет равна Т +---dx. Количество тепла, отведенного через эту грань:
дх дх
д ( дТ Л
dQx2 = —X---\ Т н---Сх у dy ■ dz ■ dt. (3)
дх ^ дх )
Таким образом, разность количества теплоты, подведенного к элементарному параллелепипеду и отведенного от него, равна количеству теплоты, внесенного путем теплопроводности в направлении оси Ох:
д^т
дх2
Аналогично:
д2Т
Полное количество теплоты, внесенное в выделенный параллелепипед, равно:
'д 2Т д 2Т д Т ~дхГ + ~ду2 + 1к2
Количество теплоты, выделенное в элементарном объеме за счет внутренних источников, определяется следующим образом:
dQвн = qv с ■ dt, (7)
где dv=dx•dy•dz - объем элементарного параллелепипеда.
Приращение внутренней энергии определяется через массу параллелепипеда рС\, теплоемкость с, приращение температуры: дТ
Си = с ■ р■ С\---а. (8)
дt
Подставляя (6) - (8) в уравнение (1), получим уравнение теплопроводности:
дТ _ (д2Т д2Т д2Т
-= Х-\ -- +-- +-:
дt { дх2 ду2 дz
Сумма вторых производных функции Т есть оператор Лапласа, который обозначается:
д 2Т д2Т д 2Т =^Т дх2 ду2 дz2 " .
Величина —■ называется коэффициентом температуропроводности и в математической Р-с
литературе обозначается буквой а. Тогда уравнение (9) примет вид:
— = а ■У2Т . (10)
дt с ■ р
Уравнение (10) называется дифференциальным уравнением теплопроводности (или уравнением Фурье). Оно является одним из основных уравнений математической теории теплопроводно-
^х = — ^ х 2 = ■ Сх ■ СУ ■ ^ ■ А . (4)
dQy = X--- ■ Сх ■ Су ■ dz ■ dt. (5)
dQm = dQx +dQy + dQz = X] — + —г + —г [ Сх ■ Су ■ dz ■ dt. (6)
с ■ Р-— = X-\— + —г + —T\ + Ч\. (9)
сти. Коэффициент температуропроводности a является физическим параметром вещества: как видно из уравнения (10) для любой точки тела изменение температуры во времени пропорционально величине a.
Уравнение Фурье применимо для всех без исключения задач теплопроводности. Но при решении конкретной задачи это уравнение необходимо рассматривать совместно с математическим выражением особенностей процесса. Эти данные называются краевыми условиями задачи, или условиями однозначности.
Существуют различные краевые условия.
• Геометрические: характеризуют форму и размеры тела, в котором протекает процесс теплопроводности.
• Физические: характеризуют физические свойства тела.
• Временные: характеризуют распределение температуры тела в начальный момент времени.
• Граничные: характеризуют взаимодействие тела с окружающей средой. В свою очередь граничные условия различают трех родов:
1) первого рода, задается распределение температуры на поверхности тела в функции времени;
2) второго рода, задается плотность теплового потока для всей поверхности тела в функции времени;
3) третьего рода, задаются температура окружающей среды ж и закон теплоотдачи между поверхностью тела и окружающей средой — закон Ньютона—Рихмана:
= а\Тс -Тж)• dF■ Л. (11)
где Тс — температура поверхности тела, К; а — коэффициент пропорциональности, называется коэффициентом теплоотдачи, Вт/(м2К). Численно коэффициент теплоотдачи равен количеству теплоты, отдаваемому (или воспринимаемому) единицей поверхности в единицу времени при разности температур между поверхностью тела и окружающей средой в один градус. Коэффициент учитывает все особенности явлении теплообмена, которые происходят между поверхностью тела и окружающей средой. Тогда плотность теплового потока, передаваемого от поверхности тела в окружающую среду, равна
q =-Л2% = а-(Тс-Тж) . (12)
dF■dt
Учитывая закон сохранения энергии, теплота, передаваемая от поверхности тела в окружающую среду, будет равна теплоте, подводимой к поверхности изнутри тела путем теплопроводности:
а <т--^--"Ш' <13)
Преобразуем последнее равенство:
(£ 1 ~ ат - тж). <">
Равенство (14) является математическим представлением граничных условий третьего рода. Таким образом, в результате решения уравнения теплопроводности с учетом граничных условий можно найти температурное поле, а на основании закона Фурье — и соответствующие ему тепловые потоки. Рассмотрим однородную плоскую стенку толщиной д. Известно, что на внешних поверхностях стенки поддерживаются постоянные температуры Тс1 и Тс2, коэффициент теплопроводности стенки постоянен и равен X. При стационарном режиме = 01 и отсутствии внутренних источников теплоты (ду = 0) дифференциальное уравнение теплопроводности (9) примет вид:
д2Т д2Т д2Т
-Т +-Г + —2
дх2 ду2 дz
При заданных условиях температура будет изменяться только в направлении, перпендикулярном плоскости стенки (вдоль оси Ох). Тогда:
дТ=дТ=о. (16)
ду дz
Граничные условия для уравнения (10) имеют вид:
Т (0) = ТС1, Т (8) = Тс 2. (17)
Интегрируя уравнение (15), находим:
.У+ТГ + ТТ = (15)
f = ^ dx
После второго интегрирования получим:
T = C ■ x + C2. (19)
Из граничных условий (17) определим постоянные С и С2:
ТС1 = с ■ 0 + C2, Тс 2 = C1 ■ S + C2, отсюда
T - T
f _ c1_c2_ f _
C1 _ ~ , c2 ~ Т c1 .
О
Подставляя значения C¡ и С2 в уравнение (19), получим уравнение распределения температуры по толщине стенки:
T, -T2
T = T . —^-— ■ x .
1
О
Для определения плотности теплового потока, проходящего через стенку в направлении оси
i dT
ox, используем закон Фурье: q = -А--.
dx
dT TA - Tc 2
Учитывая: — = ц =--, получим:
dx 1 О
а
q=-(Tc 1 - Tc2). (20)
О
Таким образом, общее количество теплоты, которое передается через поверхность стенки F за время t:
А
Q, =-{TA -Tc2)■ F■ t. (21)
О
а „ О
Отношение — называется тепловой проводимостью стенки, а обратная ей величина RA = —
О А
- термическим сопротивлением теплопроводности. Величина X зависит от температуры. Следовательно, в уравнения (20) - (21) необходимо подставить коэффициент теплопроводности Xc, взятый при средней температуре стенки.
При решении некоторых задач транспорта тепла целесообразно уравнение (10) записать в цилиндрической системе координат. Уравнения связи декартовой системы координат и цилиндрической системы координат имеют вид:
x = x r = ^1 y2 + z2
y = r cosm » z
. m = arctg —
z = r sinm y
Тогда функция T = T(x, y, z) в цилиндрической системе координат T = T(x, r,m). Получим выражение оператора Лапласа в цилиндрической системе координат.
Переходя к цилиндрическим координатам функция T'y примет вид:
dT dT dx dT dr dT dm dT dT dT sinm
— =----1-----1----, — =--cosm----.
dy dx dy dr dy dm dy dy dr dm r
Аналогично T'z :
dT dT . dT cosm
— = — sin m +----.
dz dr dm r
Запишем производные второго порядка. По переменной y :
d2T (dT dT sinm
dy2 \dr dm r
d (dT dT sinm I d ( dT dT sin®! sinm =—i— cosm----I cosm—I— cosm---
dr \dr dm r ) dm\dr dm r ) r
d2T 2 dT sin2m d2T sin2m dT sin2 m d2T sin2 m —- ■ cos2 m+---r~----+---+—г —т^'
dr2 dm r2 dmdr r dr r dm2 r2
По переменной z :
д2Т (dT . dT cos®
pc— = 2 —- + 2 —- +---+ —--- + f. (22)
г я..2 я.. „2 J v 7
2 = I--sin® +-
dz удг д® r
д (дТ . дТ cos® | . д (дТ . дТ cos® | cos®
= —I--sin®+----I sin®+--1--sin®+----I--=
дг \дг д® r ) дф\ дг д® r ) r
д2Т . 2 дТ sin 2® д2Т sin 2® дТ cos2 ® д2Т cos2 ®
= —- ■ sin2 ®----+---- +---- +—---.
дr д® r дфдr r дr r д® r
Запишем сумму вторых производных по координатам y и z:
д2Т д2Т = д2Т +1 Т +_1_ д2Т ду2 дz2 dr2 r дr r2 д®2
Тогда рассматриваемая задача в цилиндрической системе координат представится следующим уравнением:
dT „ д2Т „ д2Т 2 дТ 2 д2Т
— = 2—- + 2—- +---+ —---
dt дх дr r дr r д®
В случае осевой симметрии уравнение (22) запишется в виде:
pcr^ = r(2Т;)'Х +(2гТ;1 + rf . (23)
где X - коэффициент турбулентного обмена; f - функция, описывающая интенсивность и распределение источников
Стоит отметить, что при численном моделировании многих физических процессов, как правило, приходится решать задачу диффузии-конвекции-реакции. При этом для построения разностных схем часто используется интегроинтерполяционный метод.
Поставим задачу оптимизации параметров дискретной математической модели транспорта тепла в паровой турбинной установке. Ниже представлена аппроксимация операторов диффузионного и конвективного переноса. Вернемся к уравнению (23):
pcrjt(Т )=r )' х +2гт:)' r + rf.
п _ с1Т
Преобразуем оператор —: dt
pcd_ {Т)=д1с1 +dTdx +cTdr = + иТ, + (24)
dt dt dt дх dt Сг dt
где {u,v} - компоненты вектора скорости вдоль осей х и r соответственно. Преобразуем уравнение (24):
pcrr;+тТ' + rvVr = r (2Т')'х +(2гТr)r + rf . (25)
Будем рассматривать уравнение (25) с граничными условиями третьего рода вида: Т'(х, у, t) = апТ + /Зп . (26)
Расчетная область вписана в прямоугольник. Для численной реализации дискретной математической модели задачи вводится равномерная сетка [1]:
wh ={tn = nz,х, = ihx,rj = jhr; n = Щ,i = O,,j = ^Ñ; Nz = lt,Nxhx = lx,Nrhr = lr}, где z - шаг по времени, hx, hr - шаги по пространству; Nt - верхняя граница по времени; Nx, Nr - границы по пространству.
Аппроксимации уравнения (25) по временной переменной выполняется на основе схем с
весами:
Ф_т
pcr-+ тТ' + rvVy = r )^ + (2rf;)^ + rf, (27)
где Т = аТ +(1 -а)Т , сте[0,1] - вес схемы.
В результате получаем дискретный аналог уравнения диффузии-конвекции-реакции (25), записанного в цилиндрической системе координат с граничными условиями третьего рода (26):
253
Т. — т т — т т — т
9> \ j Р, А Г + (),,' jrj + (92 j —"2,Г 'j 2к ~1'j
т — т. т — т т — т.
+ ( 9з ), ' V, '+1/2 Г+ ( 9 4 ^ j V, j—1/2 Г'—1/2 Л ' 2к ''j—1 = ( 91 X, ' ^+1/2, Л
— (92 ),' Л—1/2,'Г' — 1(91 ^' — (92 ^' I Л'Г' ^К + + (^' Л'+1/2Г+1/2 Т''^'
г
т т сс т + ((
— (94),,' Л,'—1/2Г'—1/2 '' —1(93),,' — (94),,' | Л,'Г' Сг'/; 3 + (90 ),' ОХ,' , (28)
где 91,' = 0,4 - коэффициенты заполненностей контрольных областей [2].
Дискретные аналоги операторов конвективного и диффузионного переноса, записанные в цилиндрической системе координат, в случае частичной заполненности ячеек:
т — т т — т
'+1,' ',' , /"„ \ ,, „ ',' '—1,'
(9о)','игт □ (91)',' и'+1/2,л '1 + (92)',' и—1/2,2Ь
гр гр гр гр гр Д
(9о)',' □ (91)',' Л+1/2,— (92)',' Л—— |( 91)',' —(92)' Л,' а'и 3
X
т —т т —т
(90 ' (ЛгтГ)'г П (9з )',' Л,'+1/2Г+1/2 " 1+1к2 ''1 — (94 )',' Л,'—1/2' 1/2 '' —
— |(9з),' —(94),'Л,Г . (29)
Для внутреннего узла (91 = 92 = 1) дискретные аналоги операторов конвективного и диффузионного переноса имеют следующий вид:
т — т т — т т — т т — т
Т^' П ' + 1,' '',' , '',' ' — 1,' /ОТ^Л' Г! 0 '+1,' '',' 0 '',' ' — 1,'
игтх □ и'+1/2, А—'ъ-+ и—1/2, /! - ' (Лтх ) X П Л+1/2,' -^2--Л—1/2,' -Й2- '
т — т т — т.
(Лгт;)'г □ Л,'+1/2Г'+1/2 "'+'2 и] — Л,'—1/2Г'—1/2 "' ^^'—1. (30)
К К
Данные конечно-разностные аналоги обладают вторым порядком погрешности аппроксимации по пространственной координате. Для того чтобы показать, что аппроксимации (27) обладают вторым порядком погрешности, нужно доопределить задачу путем расширения расчетной области путем ввода фиктивных узлов.
Поставим задачу нахождения оптимального параметра схемы с весами для трехмерного уравнения диффузии [3]. При помощи этого уравнения описывается широкий класс задач математической физики. Одними из наиболее распространенных методов решения подобных задач являются сеточные методы. При использовании метода сеток возникает необходимость оценки погрешности аппроксимации уравнения в зависимости от шага по временной переменной [4-6]. В случае большого шага точность решения сильно падает. В случае малого шага - возрастает трудоемкость вычислений. Поэтому в любом случае возникает необходимость выбора оптимального шага.
Рассмотрим трехмерное уравнение диффузии [7-9]
т; = d'v (agrad (т )) + / (31)
с начальным условием
т (0, X, у, 2) = ф( X, у, 2). (32)
где источник / - задан кусочно-постоянной функцией.
Аппроксимируем уравнение (31) по пространственной переменной, при этом используем равномерную прямоугольную сетку. В результате получим:
т'=(*т )х +(*т \ +(*т I + /, (33)
т —т т —т
' = а+1/2,1',k-2--а'—1/2,1',k-2- - дискретный аналог оператора диффу-
(лМ) Кх Кх
зии.
Выражение (33) можно записать в матричной форме:
254
Т'= -Лт+г, t , tи+1 ]
(34)
где Л - самосопряженный, положительно определенный оператор (Л = Л* > 0). Найдем погрешность аппроксимации выражения (34) по временной переменной. Для решения задачи используем разностные схемы с весами.
Оценим погрешность для уравнения (34). Для этого разложим вектора с и f по ортонорми-рованной системе, составленной из собственных векторов оператора Л . Данная система является полной и образует базис [7].
Т =^а1Х1 , Г РЧХ, , (35)
ния:
где X. - собственный вектор оператора Л, для которого справедливы следующие выраже-
1, ' = }
(36)
ЛХ, = "X,, (X.,X. ) = 8.. ,8.. = Для аппроксимации уравнения (34) по временной переменной используем схемы с весами:
К
.=-" (а+(1 -^)ап)+рп
(37)
Для оператора Л в уравнении (37) имеют место неравенства:
Л > 0, 2D -Л > 0, (38)
где О - диагональная часть оператора Л . Оценим оператор О максимальным значением:
(„ , „ „ , „ „ , „ Л
Е. (39)
О < тах
а.+М2, ¡,к + а.-М2, ¡,к + а',.+1/2,к + а', .-1/2,к + а.,. ,к+1/2 + а., ¡,к-1/2
К2 Н2 И2
\ "х "у "z
В результате получим оценку максимального собственного числа оператора Л :
(„ , „ „ , „ „ , „ Л
"тах < 2таХ
а.+1/2, .,к + а.-1/2, .,к + а', .+1/2, к + аг, .-1/2,к + аг,. ,к+1/2 + а',. ,к-1/2
К К
к2
Обозначим t = "тахг . Тогда а[ примет вид:
Следовательно:
г V " Я" (а) =—— а +
(40)
(41)
(42)
Таким образом, решение уравнения (31) сводится к решению задачи (42). Точное решение задачи (42) имеет вид:
а ('-1 )=Ь ^)-Кт\ •" К+ Р
". I "
I У I
Аппроксимация этого уравнения с учетом схем с весами в тех же обозначениях:
= -" (*аГ +(1 -а)а" ) +
К
"
К "
(43)
(44)
Стоит отметить, что явная схема (а = 0) монотонна при К < 1, устойчива при К < 2 . Аппроксимация вида (44) может быть записана и в другой форме:
" "Я"
а"+1 -а" =--'-К аа" -(1 -а)—!-к а" + К —
' ' " ' " ' "
или
1 —
1 + " К а
Я"
К —
1 + " К а
(45)
где а"1 - приближенное значение функции а(х) на текущем временном слое, а" на предыдущем.
а -а"
а -а"
Значение погрешности на n -ом временном слое через точное at (г") и приближенное значения функции а" поля можно выразить через функцию:
у" = а" -а, (г"). (46)
Погрешность аппроксимации по временной переменной на n -ом слое будет меньше либо
равна max \у"\. / I I
Обозначим:
*=Т к
max
Перепишем (45) в обозначениях (47):
( R"-1 \
*
(1 + *&)
-а"-1 + Я
т
Подставим (46) в выражение (48):
( \
\ =
i—
*
(1+*&), С учетом (43) получим:
у"-1 + а (г"-1 )|i-
*
(1+*&),
-а (г")
" \ + *
(1+*&)
(р"-1 \ v ~Х J
У" =
1-
*
(1+*&),
у"-1 +
-л Я
а (г"-1)-
<-л
т
1-
*
(1+*&)
- e
Погрешность на " -ой итерации будет иметь вид:
У" =|1-
*
= 11 --
*
1 + *&
(
у"-1 +
-л Я
а (г"-1 )-
<-л (
т
1-
*
(1+*&)
-e
а г
1+*& -л Я
1-
*
1 + *&
У"-2 +
( П"-2\ (
а [г"-2)- Я
[г"-1 У
т
+...+| 1 -
*
I / (
1 -
*
1+*&
(1+*&) г-)-Я
- e
т
1 --
1-
*
(1 + *&)
-e
; J
*
1 + *&
Я Л
а (о)-Т
v Т J
я
J
1-
*
\
(1+*&)
Получим оценку погрешности численного решения поставленной задачи:
\у"\ < max
I I к=0"-1
)-т
max
жФА ]
1 —
*
1 + *&
-1
к=0 V 1
*
< max
к=0,"-1
= max
k =0,"-1
а гк)-Т
а
г)-—
max *ФА]
max
*фл ]
1 —
*
1 -
1 + *& *
1 + *&
+ *&
k
*
У1 1 —
к = 0 V 1 + *&
1 + *&
*
Полученная оценка имеет место в случае
1 —
1 + К &
(47)
(48)
(49)
(50)
(51)
< 1. Отсюда следует ограничение Кг <
(52) 2
1 - 2&
Вернемся к переменной К . Относительно данной переменной погрешность имеет вид:
< max
к=0,"-1
а, (гк)-Т
Т
max
*Що,к ]
1-
*
Л
-e
1 + *&
1 + ^) ) X Нетрудно убедиться, что полученная оценка достигается. Найдем параметр а*, при котором относительная погрешность минимальна:
(С .. \ V
(53)
& = arg min
max
*е(0Д]
1 --
*
-e
1+*&
*
(54)
+
+
+
х
+
k
-*
- e
-*
-e
-*
— e
В таблице 1 приведены значения шага и оптимального веса, зависящие от величин относительных погрешностей.
Таблица 1.
Значения оптимальных соотношений погрешности <р,
шага сетки hT и параметра а
р К а Р К а
0.0001 0.08468 0.5058 0.006 0.7203 0.5478
0.0002 0.1204 0.5083 0.007 0.7856 0.5519
0.0003 0.148 0.5102 0.008 0.8476 0.5558
0.0004 0.1715 0.5117 0.009 0.907 0.5595
0.0005 0.1924 0.5132 0.01 0.9642 0.563
0.0006 0.2113 0.5144 0.02 1.47 0.5923
0.0007 0.2288 0.5156 0.03 1.925 0.6162
0.0008 0.2452 0.5167 0.04 2.374 0.6373
0.0009 0.2607 0.5178 0.05 2.837 0.6567
0.001 0.2754 0.5188 0.06 3.333 0.6748
0.002 0.3964 0.5268 0.07 3.882 0.692
0.003 0.4925 0.5332 0.08 4.508 0.7084
0.004 0.5757 0.5386 0.09 5.25 0.7243
0.005 0.6508 0.5434 0.1 6.166 0.7397
На рис. 1 представлены функции зависимости относительной погрешности аппроксимации ф от шага ^ в случаях: а = 0,5 и оптимального веса.
Рис. 1. Зависимость ошибки аппроксимации р от шага сетки К в случаях симметричной схемы а = 0,5 (сверху), оптимального веса (снизу).
Из рис. 1 видно, что кривая, которая соответствует оптимальному весу, проходит в 5 раз ниже кривой, которая соответствует а = 0.5. В случае, если задана минимально допустимая погрешность, например, 1%, при использовании схем с оптимальным весом временной шаг можно выбирать в 2,5 раза больше по сравнению симметричной схемой, а это в свою очередь приводит к уменьшению времени счета.
Разработанное в ходе решения поставленной задачи программное обеспечение на базе ЭВМ предназначено для математического моделирования транспорта тепла.
В программе задаются следующие параметры:
- шаги по пространственным координатам (double hx=0.0164; double hr=0.0164) и по времени (double ht=1);
- размеры расчетной сетки (int Nx=819; int Nr=249);
- временной интервал (double lt=3600);
- погрешность вычисления сеточных уравнений (double eps=1e-3);
- плотность метала (double po=7800);
- коэффициент теплоемкости металла (double c=0.462);
- вес схемы (double sigma=0.563).
В программе выделяется память под следующие массивы:
• одномерный вещественный массив для описания поля концентрации C[N];
• одномерные вещественные массивы для функции, описывающей интенсивность и распределение источников f[N];
• одномерные вещественные массивы для поля коэффициента теплопроводности landa[N];
• одномерные вещественные массивы для компонентов вектора скорости u[N], v[N];
• одномерные вещественные массивы для функции, описывающей заполненность ячеек O[N];
• одномерные вещественные массивы для коэффициентов сеточных уравнений на текущем временном слое A[N], B1[N], B2[N], B3[N], B4[N];
• одномерные вещественные массивы для коэффициентов сеточных уравнений на предыдущем временном слое B5[N], B6[N], B7[N], B8[N], B9[N];
• одномерный вещественный массив для правых частей сеточных уравнений F[N]. В программе формируются следующие переменные:
• целочисленная переменная, в которой хранится размерность массивов N=Nx*Nr;
• целочисленные переменные под счетчики цикла i,j;
• вещественные переменные под значение времени t;
• целочисленные переменные под номера узлов m0,m1,m2,m3,m4,m24;
• вещественные переменные, характеризующие заполненности контрольных областей
q1,q2,q3,q4,q0.
Для работы программы также необходима функция расчета сеточных уравнений SLAY(A,B1,B2,B3,B4,F,T,eps) [2-6], входными параметрами которой являются: коэффициенты A,B1,B2,B3,B4 и правая часть F сеточных уравнений, погрешность вычисления сеточных уравнений eps. Выходным параметром является поле температур Т.
Таким образом, получена оценка погрешности аппроксимации трехмерного уравнения диффузии по временной переменной. Получены зависимости оптимального веса и шага схемы от величины относительной погрешности.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Самарский, А.А. Теория разностных схем / А.А. Самарский — М.: Наука, 1989. — 616 с.
2. Сухинов А.И., Тимофеева Е.Ф. Чистяков А.Е. Построение и исследование дискретной математической модели расчета прибрежных волновых процессов // Известия ЮФУ. Технические науки. 2011. № 8 (121). С 22-32.
3. Белоцерковский О. М., Гущин В. А., Щенников В. В. Метод расщепления в применении к решению задач динамики вязкой несжимаемой жидкости // Журнал вычислительной метаматематики и математической физики. 1975. 15:1. С. 197-207.
4. Сухинов А.И., Чистяков А.Е., Проценко Е.А. Математическое моделирование транспорта наносов в прибрежной зоне мелководных водоемов // Математическое моделирование. 2013. Т. 25, № 12. С. 65-82.
5. Сухинов А.И., Чистяков А.Е., Проценко Е.А. Математическое моделирование транспорта наносов в прибрежных водных системах на многопроцессорной вычислительной системе // Вычислительные методы и программирование. 2014. Т. 15. С. 610-620.
6. Сухинов А.И., Никитина А.В., Чистяков А.Е., Семенов И.С. Математическое моделирование условий формирования заморов в мелководных водоемах на многопроцессорной вычислительной системе // Вычислительные методы и программирование. 2013. Т. 14. С. 103-112.
7. Сухинов, А.И. Параллельная реализация задач транспорта веществ и восстановления донной поверхности на основе схем повышенного порядка точности / А.И. Сухинов, А.Е. Чистяков, А.А. Семенякина, А.В. Никитина // Вычислительные методы и программирование: новые вычислительные технологии. — 2015. — Т. 16. — C. 256-267.
8. Четверушкин Б.Н. Пределы детализации и формулировка моделей уравнений сплошных сред // Математематическое моделирование. 2012. Т. 24. № 11. С. 33-52.
9. Сухинов А.И., Проценко Е.А., Чистяков А.Е., Шретер С.А. Сравнение вычислительных эффективностей явной и неявной схем для задачи транспорта наносов в прибрежных водных системах // Вычислительные методы и программирование: Новые вычислительные технологии. - 2015. - Т 16, №3. - С. 328-338.
10. Сухинов, А.И. Адаптивный модифицированный попеременно-треугольный итерационный метод для решения сеточных уравнений с несамосопряженным оператором / А.И. Сухинов, А.Е. Чистяков // Математическое моделирование. — 2012. — Т. 24, №1. — С. 3 - 20.
11. Чистяков, А.Е. Теоретические оценки ускорения и эффективности параллельной реализации ПТМ скорейшего спуска / А.Е. Чистяков // Известия ЮФУ. Технические науки. — 2010. — №6(107). — С 237-249.
12. Чистяков, А.Е. Библиотека параллельных итерационных методов решателей СЛАУ для задачи конвекции-диффузии на основе декомпозиции по одному пространственному направлению / А.Е. Чистяков, Д.С. Хачунц, А.В. Никитина, Е.А. Проценко, И.Ю. Кузнецова // Современные проблемы науки и образования. — 2015. — № 1. URL: http://www.science-education.ru /121-19510 (дата обращения 4.06.2015).
13. Самарский, А.А. Методы решения сеточных уравнений / А.А. Самарский, Е.С. Николаев — М.: Наука, 1978. — 592 с.
14. Коновалов, А.Н. К теории попеременно-треугольного итерационного метода / А.Н. Коновалов // Сибирский математический журнал. — 2002. — 43:3. — С. 552-572.
15. Сухинов, А.И. Оценка погрешности решения уравнения диффузии на основе схем с весами / А.И. Сухинов, А.Е. Чистяков, Ю.С. Бондаренко // Известия ЮФУ. Технические науки. — 2011. — №8 (121). — С. 6-13.
