Научная статья на тему 'Решение обратной задачи теплопроводности по идентификации начального условия краевой задачи'

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

CC BY
48
12
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
ЗАДАЧА ТЕПЛОПРОВОДНОСТИ / HEAT CONDUCTION PROBLEM / АНАЛИТИЧЕСКОЕ РЕШЕНИЕ / ANALYTICAL SOLUTION OF THE DIRECT PROBLEM / ПРЯМАЯ ЗАДАЧА / ОБРАТНАЯ ЗАДАЧА / INVERSE PROBLEM / АППРОКСИМАЦИЯ / APPROXIMATION / ИДЕНТИФИКАЦИЯ / IDENTIFICATION / НАЧАЛЬНОЕ УСЛОВИЕ

Аннотация научной статьи по математике, автор научной работы — Кузнецова Анастасия Эдуардовна, Скворцова Марина Петровна, Стефанюк Екатерина Васильевна

С использованием аналитического решения прямой задачи нестационарной теплопроводности, а также значений величин искомой функции, полученных численным методом, выполнена идентификация начального условия краевой задачи. Аппроксимация полученной из численного решения температуры в одной из точек пространственной переменной выполнялась кубической параболой в некотором диапазоне времени нерегулярного режима. После подстановки аппроксимационной функции в аналитическое решение и интегрирования полученного выражения в диапазоне времени аппроксимации относительно величины начального условия было получено алгебраическое линейное уравнение, из решения которого следует, что точность идентификации составляет 0,018 %.

i Надоели баннеры? Вы всегда можете отключить рекламу.

Похожие темы научных работ по математике , автор научной работы — Кузнецова Анастасия Эдуардовна, Скворцова Марина Петровна, Стефанюк Екатерина Васильевна

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

The inverse problem solution of heat conduction for the initial conditions identification of the boundary value problem

Analytical solution of the direct problem of nonstationary heat conduction, as well as the values of the unknown function obtained by the numerical method is used. Identification of the initial conditions of the boundary value problem is made. The approximation obtained from the numerical solution temperature at one point in the spatial variable was performed by the cubic parabola in a certain irregular time range. After substitution of approximation function in the analytical solution and integrating the obtained expression in the time range of approximation relative to the value of the initial conditions the transcendental equation was obtained. The result is the iteration method. The identification accuracy is 0,018%.

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

ВЕСТН. САМАР. ГОС. ТЕХН. УН-ТА. СЕР. ТЕХНИЧЕСКИЕ НАУКИ. 2014. № 3 (43)

УДК 536. 2 (075)

РЕШЕНИЕ ОБРАТНОЙ ЗАДАЧИ ТЕПЛОПРОВОДНОСТИ ПО ИДЕНТИФИКАЦИИ НАЧАЛЬНОГО УСЛОВИЯ КРАЕВОЙ ЗАДАЧИ*

А.Э. Кузнецова, М.П. Скворцова, Е.В. Стефанюк

Самарский государственный технический университет Россия, 443100, г. Самара, ул. Молодогвардейская, 244

E-mail: totig@yandex.ru

С использованием аналитического решения прямой задачи нестационарной теплопроводности, а также значений величин искомой функции, полученных численным методом, выполнена идентификация начального условия краевой задачи. Аппроксимация полученной из численного решения температуры в одной из точек пространственной переменной выполнялась кубической параболой в некотором диапазоне времени нерегулярного режима. После подстановки аппроксимационной функции в аналитическое решение и интегрирования полученного выражения в диапазоне времени аппроксимации относительно величины начального условия было получено алгебраическое линейное уравнение, из решения которого следует, что точность идентификации составляет 0,018 %.

Ключевые слова: задача теплопроводности, аналитическое решение, прямая задача, обратная задача, аппроксимация, идентификация, начальное условие.

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

* Работа выполнена при финансовой поддержке Министерства образования и науки РФ в рамках базовой части государственного задания ФГБОУ ВПО «СамГТУ» (код проекта: 1273).

Анастасия Эдуардовна Кузнецова, аспирант.

Марина Петровна Скворцова, аспирант.

Екатерина Васильевна Стефанюк (д.т.н.), доцент кафедры «Теоретические основы теплотехники и гидромеханики».

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

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

Обратные задачи подразделяются на ретроспективные, граничные, коэффициентные и геометрические. Ниже будут представлены результаты решения ретроспективной обратной задачи, или задачи с обратным временем, в которой восстанавливается временная предыстория физического процесса. К числу таких задач относится задача по идентификации начального условия краевой задачи [1, 2]. Рассмотрим последовательность восстановления начального условия на примере решения следующей прямой задачи теплопроводности [3, 4]:

—^—=—(ро>0; 01; (1)

©&0)=0; (2) дЭ(0,Ро)/д^-В^[0(0^) + Б] = 0; (3)

д®(1,Ро)/ д^ + Bi2[©(1,Fo) -1] = 0, (4)

где ©(^, Fo) = [Т(^,Fo) - Т ]/(Т2 - 7^) - относительная избыточная температура; Б = (Т0 - Т1)/(Т2 - Т0);

% = х / 5 - безразмерная координата;

5 - толщина пластины;

Fo = ах/52 - число Фурье;

а - коэффициент температуропроводности;

х - время;

В^ = а^/Л', Bi2=«25/Л' - числа Био;

^, «2 - коэффициенты теплоотдачи;

Л' - коэффициент теплопроводности;

Т1, Т2 - температуры сред;

Т) - начальная температура.

Решение задачи (1) - (4) представим в виде суммы двух функций

©(^) = ©с© + ©н №), (5)

где ©с(^) - решение стационарной задачи при неоднородных граничных условиях вида (3), (4). Оно принимается в виде

©с© = + ^ ,

где ^ - неизвестные коэффициенты, определяемые из граничных условий (3), (4). Формулы для них будут (полагая, что Т = Т0 и Б = 0)

^ =-^-; ^ = ^ . (6)

Б11 + Бг2 + БцБ12

Функция ©н (%^о) является решением нестационарной задачи теплопроводности при однородных граничных условиях. Математическая постановка задачи для нее имеет вид

д©н (^о) 52©н (^^о) , ч

' =-(Ро > 0; 0 <%< 1); (7)

дБо д%

©н (%,0)=-©с (%,); (8) д©( (0,Fo)/д% - Bi1 ©( (0,Fo) = 0; (9)

д©( (l,Fo)/д% + Bi2 ©( (1 ^о) = 0. (10)

Следуя методу разделения переменных, решение задачи (7) - (10) принимаем в виде

©( (%^о) = ф^оМ%), (11)

где ф (Fo) - неизвестная функция времени; у(%) - неизвестная функция пространственной переменной.

Подставляя (11) в (7), находим

dф( Fo)/ dFo + Уф^о) = 0; (12)

d:V(%)/й%2 + уу(%) = 0. (13)

где у - коэффициент.

Решение уравнения (12) известно и имеет вид

ф( Fo)= АехрН^о). (14)

Подставляя (11) в (9), (10), находим граничные условия для функции у(%):

0)/ й%- Б^у(0) = 0; (15)

1) / й% + Б12у(1) = 0. (16)

Решение краевой задачи Штурма - Лиувилля (13), (15), (16) принимается в виде

т

у(%)=£Ьк%к , (17)

к=0

где Ьк - неизвестные коэффициенты, определяемые из основных и дополнительных граничных условий (Ь0 = 1); %к (к = 0, т) - координатные функции.

Если ограничиться, например, десятью членами ряда (17), то будем иметь девять неизвестных коэффициентов Ьк (Ь = 1), а граничных условий для их определения только два - (15) и (16). В связи с этим необходимо добавить еще семь дополнительных граничных условий, которые находятся из уравнения (13) путем выполнения самого этого уравнения, а также производных от него различного порядка, в граничной точке % = 0. Такие дополнительные граничные условия в данном случае имеют вид

й 2у( 0)/ й%2 +уу(0) = 0; (18)

^ 0)/ + vdУ-2(0)/d^-2 = 0. (г = 3, 4,5, 6,7, 8)

Подставляя (17) в основные (15), (16) и дополнительные (18) граничные условия, относительно неизвестных коэффициентов Ък (к = 0,9) будем иметь следующую цепочную систему девяти алгебраических линейных уравнений:

¿1- В^Ь0 = 0; 2Ь2 = 0; 6Ь3 +^= 0; 12Ь4 + vЪ2 = 0; 20Ь5 + vЪ3 = 0;

30Ь6 + vЪ4 = 0; 42Ь7 + vЪ5 = 0; 56Ь8 +^6 = 0;

Ь0В12 + ¿1(1 + В12) + ¿2(2 + В12) + ¿3(3 + В12) + ¿4(4 + В12) + ¿5(5 + В12) + + ¿6(6 + В12) + Ьу(7 + В12) + ¿8(8 + В12) + ¿9(9 + В12) = 0. (19)

При известной величине ¿0 = 1 все остальные неизвестные коэффициенты ¿1 (г = 1,9) системы (19) легко находятся из решения лишь одного уравнения (соответствующего данному коэффициенту) с одним неизвестным. Найденные таким путем коэффициенты имеют вид

В^; ¿2 /2; ¿5 =-vBi1/6; ¿4 = V2/24;

¿5 = V2В^/120; ¿6 = -V3/720; ¿7 = v3Bi1/5040; ¿8 = V4/40320; ¿9 = -1 + 9v / 22 - / 440 + 83v3 / 55440 - v4 / 44352.

Коэффициент ¿9 найден из условия В^ = 3; В^ = 2.

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

Для определения собственных чисел составим интеграл взвешенной невязки уравнения (13)

1

{[¿2у(^)/^ 2 + d \ = 0. (20)

0

Подставляя (17) в (20), относительно v к с учетом найденных значений коэффициентов Ък получаем следующее алгебраическое уравнение:

v5 283v4 37v3 27v2 197v

+ _--+ = 0. (21)

1995840 2217600 3850 88 55

Первые два корня уравнения (21) имеют вид: V! = 3,40427% V2 = 17,593928. Ввиду того, что уравнение (13) удовлетворяется лишь при некоторых дискретных значениях vк (собственных значениях), остальные корни этого уравнения отбрасываются как не удовлетворяющие ему.

Для уточнения первых двух собственных чисел составляется невязка уравнения (13) и требуется ортогональность невязки к собственной функции (17):

1

{[¿4(0/+ d \ = 0. (22)

0

После определения интегралов в (22) относительно собственных чисел получаем степенное алгебраическое уравнение, корни которого, удовлетворяющие уравнению (13), следующие: Vl = 3,404290, V2 = 17,689298.

Сравнивая собственные числа с их точными значениями (см. таблицу), можно заметить, что требование ортогональности уравнения (13) к собственной функции (17) позволяет существенно уточнить их значения. Собственные функции по известным собственным значениям находятся из (17):

¥п (%) = £Ьк(у)% (п = 1, 2). (23)

к=0

Число приближений Собственные числа

V 2 V 4 V 5

2 3,404278 17,593927

3,404290 17,689297

3 3,404290 17,683688 48,120431

3,404290 17,683680 48,586633

4 3,404290 17,683680 48,569379 95,497616

3,404290 17,683680 48,569017 98,594028

5 3,404290 17,683680 48,569017 98,357200 166,909671

3,404290 17,683680 48,569017 98,357137 167,648755

Точные значения Vп 3,404025 17,682025 48,566961 98,346889 167,650704

Подставляя (14), (17) в (11), находим частное решение уравнения (12):

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

©н (%, Fo) = Ау(%)ехр(-уБо). (24)

Составляя сумму частных решений, получаем

т

©н (%, Ро) = Е Ап¥п (%)ехр( - у Бо), (25)

п=1

где т - число приближений, равное количеству собственных чисел.

Для нахождения неизвестных коэффициентов А п составляется невязка начального условия (8) и требуется ортогональность невязки к каждой собственной функции:

1 2

1 [ Е АпУп (%) + ^ + ^2%]у; (%)й% = 0 . (26)

0 п=1

Ввиду ортогональности собственных функций получаем следующую формулу для определения коэффициентов А п (п = 1, 2):

1 1 2 Ап =1 (^ + ^2%)Упй%/ 1у2(%)й% . (27)

0 0

Из формулы (27) находим А1 =-0,27528; А2 = 0,13713. Отметим, что коэффициенты А1 и А 2 найдены из уточненных значений собственных чисел.

После нахождения решений нестационарной

©н

и стационарной

©с соотношение (5) принимает вид

Fo) = ^ + ^ \ + Е ©ехр(^п Fo). (28)

П=1

На графиках (рис. 1) представлены результаты расчетов безразмерной температуры по формуле (28) в пятом приближении (т = 5). При этом было использовано 30 членов ряда (17). Значения собственных для различного числа приближений в сравнении с их точными значениями [5] представлены в таблице. Здесь даны значения собственных чисел для случаев, когда интегрируется лишь невязка уравнения (13) (верхняя строка каждого приближения) и когда требуется ортогональность невязки к собственной функции уп (нижняя строка).

Анализ полученных результатов позволяет заключить, что безразмерные температуры, полученные по формуле (28), в пятом приближении для чисел Fo > 0,001 практически совпадают с их точными значениями, которые определялись по формуле (24) из [5]. При этом было использовано 30 членов ряда точного решения.

Анализируя изменение невязки уравнения (1) по координате х (рис. 2) (Fo = 0,02), можно заметить, что на отрезке 0 < х < 0,75 уравнение (1) выполняется точно (8 = 0). Для х > 0,75 невязка возрастает, и в точке х = 1 она принимает максимальное значение. Изменение невязки уравнения (1) в точке х = 1 во времени (рис. 3) показывает, что она уменьшается с увеличением числа Фурье и при Fo > 0,03 становится практически равной нулю.

Анализ невязки начального условия (2) показывает, что максимальная ее величина достигается также в точке х = 1. Невязка начального условия уменьшается с увеличением числа приближений (рис. 4).

0,5

0

0,4

0,6

0,2

0,1

0 1

Рис. 1. Распределение безразмерной температуры:-точное решение [5];

о расчет по формуле (28) (пятое приближение, т = 5)

Преимущества полученного выше аналитического решения вида (28) заключаются в том, что при использовании экспериментальных (или расчетных) данных по изменению во времени температуры в какой-либо одной точке пластины (например ^ = 1) путем решения обратной задачи теплопроводности

могут быть найдены физические характеристики среды, коэффициенты граничных условий, геометрические условия, а также начальные условия краевой задачи. И, в частности, используя численное решение задачи (1) - (4), с помощью соотношения (28) выполним идентификацию начального условия (2) при следующих исходных данных: а,1 = 10 Вт/(м2 / К); а1 = 20 Вт/(м2 / К);

5 = 0,06 м; Х' = 45 Вт/(м / К); а = 10

8

-6 „2

м2 / с; Т0 = Т1 = 100иС; Т2 = 20иС.

0,0033 0

- 0,0033 - 0,01

- 0,0166

0

0,25 0,5 0,75 £ 1

0

8 0

- 0,14 -0,38

- 0,62 -0,86

-1,1

0

0,01 0,02

0,04

Рис. 2. Невязка 8 уравнения (1) при т = 5 (пять приближений) (Fo = 2,02, £ = 1)

Рис. 3. Изменение невязки 8 уравнения (1) во времени при т = 5 (£ = 1)

8

- 0,063

0,038

0,014 0

- 0,011 - 0,035

0 0,25 0,5 0,75 £ 1

Рис. 4. Невязка 8 начального условия при т = 5 (пять приближений) (Fo = 0)

Значения температуры в точке £ = 1 (х = 0,006м), полученные из решения данной задачи численным методом, в диапазоне времени 72 с < т < 288 с аппроксимируем следующей функцией:

3

Т (0,006; т) = Е а т', (29)

г=0

где а (г' = 0,3) - неизвестные коэффициенты.

Записывая соотношение (29) для четырех значений температур (Т = 71,23 °С; Т2 =74,14 °С; Т3 =77,5 °С; Т4=82 °С), наблюдающихся в точке х = 0,006м соответственно в моменты времени т = 72 с ; т2 = 144с ; Т3 = 216 с;

Т4 = 288с, относительно неизвестных коэффициентов аг- (г' = 0,3) будем иметь систему четырех алгебраических линейных уравнений, из решения которой находим а1= 88,3; а2= -0,103819; а3= 0,000243; а4= -3,081061. Проверка результатов аппроксимации позволяет заключить, что значения температур,

получаемых по формуле (29), во всем диапазоне времени 72 с < т < 288 с практически совпадают с численным решением. Запишем формулу (28) в размерном виде:

T(x, т) = T0 + (T2 - To)

/ \k

^ x 5 9 ( x \ ax.

Fi + F2 j + s А Z bk (vn) - I exp(-vn -2)

о n=i k=o I 5 ) 52

(30)

Подставляя (29) с учетом найденных значений коэффициентов й1 (г = 0,3) в левую часть соотношения (30), относительно искомого начального условия 70 будем иметь алгебраическое линейное уравнение 151,2158170 - 15124,366359 =

0. из решения которого находим 70 =100,01841. Таким образом, погрешность определения начального условия из решения обратной задачи теплопроводности составляет 0,018 %.

БИБЛИОГРАФИЧЕСКИЙ СПИСОК

1. Алифанов О.М. Обратные задачи теплообмена. - М.: Машиностроение, 1988. - 297 с.

2. Алифанов О.М. Идентификация процессов теплообмена летательных аппаратов. - М.: Машиностроение, 1979. - 216 с.

3. Кудинов В.А., Дикоп В.В., Назаренко С.А., Стефанюк Е.В. Метод координатных функций в несимметричных задачах теплопроводности // Вестник Самарского государственного технического университета. Сер. Математическая. - Вып. 22. Дифференциальные уравнения и их приложения. № 2. - 2003. - С. 136-141.

4. Кудинов В.А. Метод координатных функций в нестационарных задачах теплопроводности // Известия Российской академии наук. Энергетика. - 2004. - № 3. - С. 82-100.

5. Григорьев Б.А., Маньковский О.Н. Инженерные задачи нестационарного теплообмена. - Л.: Энергия, 1968. - 84 с.

Статья поступила в редакцию 5 января 2014 г.

THE INVERSE PROBLEM SOLUTION OF HEAT CONDUCTION FOR THE INITIAL CONDITIONS IDENTIFICATION OF THE BOUNDARY VALUE PROBLEM

A.E. Kuznetsova, M.P. Skvortsova, E.V. Stefanyk

Samara State Technical University

244, Molodogvardeyskaya st., Samara, 443100, Russian Federation

Analytical solution of the direct problem of nonstationary heat conduction, as well as the values of the unknown function obtained by the numerical method is used. Identification of the initial conditions of the boundary value problem is made. The approximation obtained from the numerical solution temperature at one point in the spatial variable was performed by the cubic parabola in a certain irregular time range. After substitution of approximation function in the analytical solution and integrating the obtained expression in the time range of approximation relative to the value of the initial conditions the transcendental equation was obtained. The result is the iteration method. The identification accuracy is 0,018%.

Keywords: the heat conduction problem, analytical solution of the direct problem, inverse problem, approximation, identification.

Anastasiya E. Kuznetsova, Postgraduate Student. Marina P. Skvortsova, Postgraduate Student. Ekaterina V. Stefanyk (Dr. Sci. (Techn.)), Professor.

i Надоели баннеры? Вы всегда можете отключить рекламу.