Об одном подходе к восстановлению неоднородных свойств термоупругих тел Текст научной статьи по специальности «Математика»

УДК 539.3

ОБ ОДНОМ ПОДХОДЕ К ВОССТАНОВЛЕНИЮ НЕОДНОРОДНЫХ СВОЙСТВ ТЕРМОУПРУГИХ ТЕЛ

© 2012 г. А. О. Ватульян, С.А. Нестеров

Ватульян Александр Ованесович - доктор физико-математических наук, профессор, заведующий кафедрой теории упругости, факультет математики, механики и компьютерных наук, Южный федеральный университет, ул. Мильчако-ва, 8а, г. Ростов н/Д, 344090, e-mail: [email protected].

Нестеров Сергей Анатольевич - аспирант, кафедра теории упругости, факультет математики, механики и компьютерных наук, Южный федеральный университет, ул. Мильчакова, 8а, г. Ростов н/Д, 344090, e-mail: [email protected].

Vatulyan Aleksandr Ovanesovich - Doctor of Physical and Mathematical Science, Professor, Head of Department of the Elasticity Theory, Faculty of Mathematics, Mechanics and Computer Sciences, Southern Federal University, Milchakov St., 8a, Rostov-on-Don, 344090, e-mail: [email protected].

Nesterov Sergey Anatolyevich - Post-Graduate Student, Department of the Elasticity Theory, Faculty of Mathematics, Mechanics and Computer Sciences, Southern Federal University, Milchakov St., 8a, Rostov-on-Don, 344090, e-mail: [email protected].

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

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

In general terms, the inverse problem of determining the inhomogeneous characteristics of thermoelastic body on some additional information is formulated. The weak formulation of the problem of vibrations of thermoelastic bodies for the general case of loading is presented. Based on this formulation the iteration method of constructing of the inverse problem solution is provided. As a concrete example the several one dimensional problems of reconstruction of the pair modules of thermoelasticity is considered. The solution of the direct problem is reduced to consequent solving a Friedholm integral equation of the 2nd kind. Results of numerical experiments are given.

Keywords: thermoelasticity, inverse problem, inhomogeneous body, integral equation.

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

алами [8]. Одним из основных подходов к решению КОЗ является способ, основанный на нахождении минимума неквадратичного функционала невязки. При этом для минимизации функционала используют, как правило, градиентные методы [1-3] или генетические алгоритмы. При решении прямых задач авторы опираются в основном на метод конечных элементов [3, 4] или конечных разностей [1, 2]. Обзор результатов по обратным задачам, возникающим при изучении общих задач механики деформируемого твердого тела, дан в [9]. В случае неоднородных тел главная трудность состоит в сложной процедуре построения операторных соотношений, связывающих искомые и измеряемые в эксперименте функции. Отметим, что в последние годы в ряде работ развит альтернативный подход к решению КОЗ теории упругости [9], теплопроводности [6, 10] и термоупругости [11, 12]. В этих работах решение нелинейных обратных задач сводится к итерационной процедуре, на каждом шаге которой решается линейная задача. При этом линеаризация осуществляется либо на основе обобщенной теоремы взаимности [11], либо на основе условия ортогональности [12]. В

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

Слабая постановка задачи термоупругости

Рассмотрим колебания термоупругого тела, имеющего объем V и кусочно-гладкую границу = ^ = ^ & . Общие уравнения движения,

определяющие соотношения и граничные условия, имеют вид [14, 15]:

" I CijkMißk ,idV - Р2+

= суыек1 -Уув , (1)

(ку(М)ву(МЛ), -св(М)в- Ту^ = 0,

в\8т = 0 , = д , (2)

^к = 0 , аупу к = Л •

Здесь ау (М, Л), и, (М, Л) - компоненты тензора напряжений и вектора перемещения; в(М, {) = Т(М, {) - Т - приращение температуры от естественного состояния с температурой Т0; сук1(М) -компоненты тензора упругих постоянных; р - плотность; се - удельная теплоемкость; к^ (М), у (М), п ■, р1 - компоненты тензора теплопроводности, тензора температурных напряжений, единичного вектора внешней нормали к 5 , вектора активной нагрузки, приложенной к телу; д - тепловой поток.

Применим к уравнениям (1) и граничным условиям (2) преобразование Лапласа по времени, полагая тепловые и механические начальные условия нулевыми:

ст.. f = p2ppi, Sit = Cj,ßk, -Yifi

(3)

(kj(M)0 j (M,p))i - pes (M)d - pT0yv~. j = 0, (4)

в Ist = 0' -kueini k = p ' p 1& = 0 ' Sun> k = pi ■

ij J j

ij j

Для получения слабой постановки задачи введем в рассмотрение гладкие пробные функции р и ~, удовлетворяющие главным граничным условиям р ^ = 0, т\8т = 0, и воспользуемся методом множителей

Лагранжа. Для этого умножим уравнение (3) на р , (4) - на ах р , сложим полученные выражения и проинтегрируем по объему V. Аналогично поступим с граничными условиями: условие на части 5 умножим на а2р, на части - на а3р, складываем и

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

укг^1,у к ,1"-

V V

+ I у и у ~ - Рт0а~увР)Ж¥ - а1 I ~ -

V V

- а | рсев тЖУ + | ((1 + а2 )<Ууп. - а2 р )р ЖБ +

V Ба

+ |((а + а3)кувП - а3Р)рЖ8 = 0. (5)

Полагая в (5) ах =-, а2 = -1, а3=--, получим

РТ0 РТ0

I СукРуРк+ Р 2 I Р^^^ + -^Т I Р ^ +

V V рт0 V

+± I сев^ +у - = I I д ж.

Т0 V V Ба рТ0 Бд

В этом случае слабая формулировка приобретает вид А(а,и,= Ъ(у) [13], где А(а,и,~) - трилинейная форма (линейная по каждому аргументу) переменных а, и , ~ ; Ъ(р) - линейная форма. Для рассматриваемой в настоящей работе модели термоупругости эти формы имеют вид

V) = 12ь(щ, 'V,, в, ~, ст, Ру, ку, с^ р)йУ,

V

Ъ(и) = IIр иж,

б а рТ0 Бд

, Р , C¡jkl, Р) = СукР Лу +-1- квУ Р + рТ0

Р Р Р ~ 1 ~Р 2__

+ уЦ(иии -сзви +рр РР • (6)

Пусть а - начальное приближение для нахождения искомых коэффициентов, которое обычно ищется в классе простых функций (линейных, кусочно -постоянных); иР - соответствующее поле, удовлетворяющее слабой постановке А(а0, и0, = Ъ(~).

Тогда последовательность задач для определения иР формируется следующим образом: А(а„, ип, = Ъ(~), а поправка а - ап находится из операторного уравнения 1 -го рода с компактным оператором А(а - ап, ип, ип) = ъ(~ - ~).

Одномерная коэффициентная обратная задача термоупругости

В качестве конкретного примера рассмотрим закрепленный на торце х = 0 термоупругий стержень длиной I, в котором колебания возбуждаются приложением к торцу х = I механической нагрузки р1. Считаем, что модуль Юнга, удельная теплоемкость, плотность, коэффициент теплопроводности, коэффициент температурного напряжения - произвольные положительные функции координаты х. Целью решения обратной задачи является реконструкция двух коэффициентов дифференциальных операторов уравнений термоупругости по дополнительной (входной) информации на торце стержня и(1,Л) = /(Л). Эта задача относится к одномерным КОЗ. Слабая постановка представленной задачи может быть легко получена на основе общего соотношения (6):

JL(ce,р,к,Е,и,v,ß,r)dx = F ,

(7)

где L = EWP' + у(вР' -tu') + — кв'т' + — ceв P +

PTo

T

1 n

JE-i dr-f dx + p2 JP-Uld + JSkni (dx +

1 i ddn —J ¿W в

pT0 о dx

1

i

+ —JSCn-ien-1dx = Pl(UT - un-l)\x=l •

T0 0

(8)

Здесь ит - измеренные торцевые значения трансформанты смещения; и виЧ - найденные значения торцевых смещений и температур из решения прямой задачи на (п -1) -й итерации. Для случая нагрузки р1 = —Н(I) или в трансформантах

р = перейдем в (8) к безразмерным параметрам Р

и функциям, полагая г = — ; кп (г) = п (— ;

I кй

~сп ^ = -^ ; Рп (г) = ; ЁЫ = Еп

Ро =

Ро

соРо12Р . - ч _an (x) .

кс\

ип (z) =■

Wn (z, p) =

Ео

ena0 к0

соРо12

p и к„ p <т к„ n к2

u~n(z,P) = 1uL^; fin(z,p) = -<J4i; = 0

l соро E0C0p0l

C0p0 E0l

,dW„

+ sJSkn-1(-n-1)2 dz + p0 Jöcn-1W 2 n-idz =

о dz о

= ^1(Ut - Un-1)\z==1 •

(9)

Решение прямой задачи для неоднородного термоупругого стержня

Для нахождения безразмерных трансформант и„_г(г,р0) и р0) необходимо решить прямую

задачу [11]:

dfi , 9 p ^ _ dU , p n- = pa2pn-l(z)Un-l, fin-1 = En-1(z)U1 -Y(z)Wn-1,

dz

dz

sd (K-^z) dWn--) = pö^n-1 (z)Wn-1 +SpоУ(z) ^ dz dz dz

(10)

+ р риу , Г = ру .

Обратная задача может быть решена на основе итерационного процесса. Для этого, используя обычную технику, проведем линеаризацию (7), положив

-п = Сп-1 + 3-п-1 , Рп = Рп-1 + ¿Рп-1 , Еп = Еп-1 + 3Еп-1 ,

кп = кп-1 + 3кп-1 , 7п = Гп-1 + % п-1 , ~ = ~п-1 + 3ип-1 ,

в п = в п-1 + 3вп-1, 7 = ип-1, т = в п-1, и сохраняя в уравнении (7) линейные слагаемые, с учетом дополнительной информации получим интегральное уравнение Фредгольма 1-го рода для нахождения поправок восстанавливаемых коэффициентов:

К „ -и ,

! ЯЦ ( п-1

^п-1 ( ,

ах

~ ~ № , ~ 2 Р-^Ро)=КА(0,Р0) = 0, -^(Хр») = 0, й^а р0) = —.

аг ро

Большинство материалов характеризуется слабой связанностью упругих и тепловых полей, и поэтому прямую задачу удобно исследовать методом малого параметра, в качестве которого был выбран параметр ¿1. В этом случае решение прямой задачи на (я-1)-й итерации находят в виде разложения по малому параметру: ^п-! = и™ + ¿¿7™ , = + З^Г™ , = й^] + ¿й^ . Подставив эти суммы в уравнения и граничные условия (10) и приравнивая коэффициенты при одинаковых степенях параметра д1, получают цепочку задач [11] для последовательного

нахождения Ж7^), йч , ^^Х, й^- . Решение каждой задачи находят в виде интегрального уравнения

(1)

VX)

Фредгольма 2-го рода: W(fl = о, &(пи_\ (z, p0) =

№

1

Л

= -pо2 Jw^ JPn-1 p0)dv-

о En-1(V)minxz,V} pо

pо ^

Wn-\( z, p0) J (Cn-1(n) +

£ _ аоТоЕ0 . 2 = —0 -оро Ео

Здесь а0, C0, р0, E0 - некоторые харакгерн^1е коэффициенты теплопроводности, теплового расширения, теплоемкость, плотность и модуль Юнга соответственно.

Для безразмерных параметров уравнение (8) примет вид

1 — -и 1 ~

¿1 ро ¡гЕп-А—т1)2 -г +¿1 р30 \5рп-1и п-1-г + о —г о

£ о

+ б1а2(Л)Еп_1(л)) | ро)—Л-

о кп-1(6)

-31 '¡^^¡Е^о^, р0——л ,

£ о кп-1(Л) V и 111 ~ йп- (г, ро) = -р2 ¡-¿Г— ¡Рп-1 (4)—$й(1\ (V, ро)—Л -

о Еп-1 (71) шт{гЛ}

1 V и

- ро ¡Рп-1 ШЕ(4)ЖП1 (4, рд—£—1.

г о

Точность построения вычислительной схемы проверена путем сравнения численного и аналитического решений исходной задачи для однородного стержня. Относительная погрешность сравниваемых решений составила не более 1 % при равномерном разбиении отрезка интегрирования на 20 частей.

Получение дополнительного интегрального уравнения

Для одновременного восстановления нескольких коэффициентов уравнений термоупругости одного уравнения вида (8) недостаточно. Нужно получить дополнительные уравнения путем изменения вида нагру-жения торца стержня. В данной работе для одновременного восстановления двух коэффициентов рассматривается 2-я задача о колебаниях неоднородного термоупругого стержня, возбужденных тепловой нагрузкой q. Входной информацией здесь служит измеренная температура на торце стержня в(1, ^ = /2 (0. В этом случае, выполняя действия, аналогичные решению рассмотренной выше задачи, построено второе интегральное уравнение для нахождения поправок реконструируемых функций, которое для случая нагрузки q = до Н) после обезразмеривания примет вид

о

а

о

i _ dv i ~

$ po №n-i dz +öxe2 p0 \5pn-Vn-\dz +

o dz o

<€-1,2^. p2 , 1

(11)

+ \5кп_х (ЖТп-1)2 + р0\SOn-Jhdz = —(Тт - Тп-1) и • 0 dz 0 р0

Здесь УпА, Тп_х - безразмерные трансформанты смещения и температуры, вычисленные на (п-1)-й итерации для случая тепловой нагрузки.

Решение интегральных уравнений вида (9), (11) является некорректной задачей и требует регуляризации, которая осуществлялась на основе метода регуляризации А.Н. Тихонова.

Результаты вычислительных экспериментов

В работе натурный эксперимент заменен вычислительным. Коэффициенты дифференциальных операторов уравнений термоупругости восстанавливались в 2 этапа. На 1 -м определялось начальное приближение в классе положительных ограниченных линейных функций на основе минимизации функционала невязки:

Jn =|(Ut(Xpo) -Un-i(1,po))2dpo +

0

да

+ I (Tt (1 po) - Tn-i(1, po))2 dpo ■

(12)

На 2-м этапе определялись поправки реконструируемых функций путем решения системы интегральных уравнений вида (9), (1 1). После нахождения поправок строилось новое приближение и осуществлялся итерационный процесс по схеме ап ^) = ап_х ^) + <5ап_г (z) . Выход из итерационного процесса осуществлялся по достижении некоторого порогового значения функционала невязки (12).

В 1-й серии вычислительных экспериментов одновременно восстанавливались две безразмерные функции - модуль Юнга Е (z) и удельная теплоемкость с (z), остальные считались известными. Были получены системы интегральных уравнений вида (9), (11) для нахождения поправок для следующих типов нагружения: 1) р1 = а0Н(Л), д = д0Н(Л); 2) р1 = ст0£(О,

д = д0ё(Л); 3) р1 = а0Ле- , д = д0Ле- . Так, в случае 1-й нагрузки интегральные уравнения примут вид

1 _ жи 1 ~

ёр ^^(ж^--)2 dz +р{0 ^^-Ж =

0 dz 0

= Лё1(ит - иТп-1)\ ,

1 — 1 Т

I ёЕп-1 (—г-1)2 dz +р01 ёСп-Т-Ж =

0 dz 0

= —(Тт -Тп_1)\г=.

р0

Расчеты проводились при следующих безразмерных параметрах: е = 10-5, ¿1 = 103. Исследовалось влияние шумов, характера нагружения, монотонности функций на результат реконструкции. Монотонные функции восстанавливались лучше немонотонных. Наибольшая погрешность реконструкции возникала на торцах стержня. Задание торцевых значений восстанавливаемых функций значительно снижало погрешность реконструкции и в других точках. В случае, если коэффициенты имели разные законы неоднородности,

восстановление происходило хуже и за большее число итераций, в случае отсутствия зашумления - примерно одинаково при любых типах нагружения; погрешность реконструкции не превышала 15 %. Однако наложение даже 1%-го шума делало процесс реконструкции удельной теплоемкости невозможным.

Во 2-й серии вычислительных экспериментов одновременно восстанавливались 2 безразмерные функции - модуль Юнга Е (z) и коэффициент теплопроводности к (¿), остальные считались известными. Выяснено, что процедура реконструкции в этом случае устойчива даже при наложении 4%-го аддитивного шума. Погрешность реконструкции без шума не превышала 11 %, при 4%-м уровне шума - 16 %.

На рис. 1, 2 представлены результаты одновремен-

— ж

ной реконструкции двух функций к (z) = 0,5 + 51п(— z)

и Е (¿) = 1,1 + cQs(лz) при торцевой нагрузке р = а0Н(Л), д = д0Н(Л) и заданных значениях восстанавливаемых функций на торцах стержня. При этом сплошной линией изображен точный закон, точками - восстановленный. Начальное приближение для

функции к (х) выбрано в виде 0,4 +, для Е (z) -2,2 -1,8z . Погрешность реконструкции на 7-й итерации не превысила для коэффициента теплопроводности 2 %, а для модуля Юнга - 12, параметр регуляризации а = 7,2 -Ю-8. к

1,5 1,3 1,1 -0,9 -0,7 -0,5

0

0,2

0,4

0,6

0,8

Рис. 1. Реконструкция функции k (z) = 0,5 + sin(~ z)

E

2,1 1,6 -1,1 -0,6 0,1

0

0,2

0,4

0,6

0,8

Рис. 2. Реконструкция функции Е (z) = 1,1 + оо$(ж)

На рис. 3, 4 представлены результаты одновременной реконструкции 2 функций Е (т.) = 1 + 0,5е~z и - 3 2

с(z) = — z - z +1 при торцевой нагрузке р1 = ст0ё(0,

д = д0 ё(Л). При этом сплошной линией изображен точный закон, точками - восстановленный. Начальное приближение для Е (z) выбрано в виде 1,4 - 0,3z , с (¿) -0,8 + 04. Погрешность реконструкции на 6-й итерации не превысила для модуля Юнга 2 %, а для удельной теплоемкости - 6, параметр регуляризации а = 4,4 -10 8.

СО

0

с

1,5 -1,4 -1,3 -1,2 -1,1 -1

0,9 -0,8

0

0,2

0,4

0,6

0,8

Рис. 3. Реконструкция функции c (z) = — z2 - z +1 E

1,5 1,4 -1,3 -1,2

1,1

1 z

0 0,2 0,4 0,6 0,8

Рис. 4. Реконструкция функции E(z) = 1+0,5e z

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

Работа выполнена при частичной поддержке Российского фонда фундаментальных исследований (№ 10-01-00194-а), ФЦП «Научные и научно-педагогические кадры инновационной России» на 2009 -2013 годы (госконтракт П596).

Литература

1. Алифанов ОМ., Артюхин ЕА., Румянцев С.В. Экстремальные методы решения некорректных задач. М., 1988. 288 с.

2. Кабанихин С.И., Гасанов А., Пененко А.В. Метод градиентного спуска для решения коэффициентной обратной задачи теплопроводности // Сиб. журн. выч. мат. 2008. Т. 11, № 1. С. 41 - 54.

3. Jarny Y., Ozisik N., Bardon J.P. A general optimization method using adjoint equation for solving multidimensional inverse heat conduction // International J. of Heat and Mass Transfer. 1991. Vol. 34, № 11. P. 2911 - 2919.

4. Kravaris C., Seinfeld J.H. Identification of spatially-varying conductivity from point observation as an in-verse Sturm-Liouville problem // SIAM J. Control & Optim. 1986. P. 522 - 542.

5. Dimitriau G. Parameter identification in a two-dimensional parabolic equation using ADI based solver // Computer science. 2001. Vol. 79. P. 479 - 486.

6. Xu M.H., Cheng J.C., Chang S.Y. Reconstruction theory of the thermal conductivity depth profiles by the modulated photo reflectance technique // J. of Applied Physics. 2004. Vol. 84, № 2. P. 675 - 682.

7. Cardenas-Garcia J.F., Shabana YM., Medina R.A. Thermal loading and material property characterization of a functionally graded plate with a hole using an inverse problem methodology // J. of Thermal Stresses. 2006. Vol. 29, issue 1. P. 1 - 20.

8. Ломазов В.А. Задачи диагностики неоднородных термоупругих сред. Орел, 2002. 168 с.

9. Ватульян А.О. Обратные задачи в механике деформируемого твердого тела. М., 2007. 223 с.

10. Ватульян А.О., Нестеров С.А. Об одном подходе к восстановлению коэффициентов переноса // Изв. вузов. Сев.-Кавк. регион. Естеств. науки. 2009. № 3. С. 39 - 43.

11. Ватульян А.О., Нестеров С.А. Коэффициентные обратные задачи термоупругости для неоднородных тел // Эколог. вестн. науч. центров ЧЭС. 2009. № 3. С. 24 - 30.

12. Ватульян А.О., Нестеров С.А. Об особенностях идентификации неоднородных свойств термоупругих тел // Эколог. вестн. науч. центров ЧЭС. 2011. № 1. С. 29 - 36.

13. Ватульян А.О. К теории обратных коэффициентных задач в линейной механике деформируемого тела // ПММ. 2010. № 6. С. 911 - 918.

14. Новацкий В. Динамические задачи термоупругости. М., 1970. 256 с.

15. Подстригач Я.С., Ломакин ВА., Коляно Ю.М. Термоупругость тел неоднородной структуры. М., 1984. 386 с.

Поступила в редакцию

19 марта 2012 г.

z