CC BY
54
3. Коссович Л.Ю., Кушеккалиев А. Н. Анализ приближений в задаче Лэмба для бесконечного упругого слоя // Известия вузов. Северо-Кавказский регион. Естеств. науки. 2003. № 5. C. 10-22.
4. Ewing W.M., Jardetzky W.S., Press F. Elastic waves in layered media. N.Y., 1957.
5. Miklowitz J. Elastic wave propagation // In Applied mechanics surveys (Ed. H. N. Abramson, H. Liebowitz, J. M. Crowley, S. Juhasz). Washington D.C., 1966.
6. Каплунов Ю.Д., Коссович Л.Ю. Асимптотическая модель для вычисления дальнего поля волны Рэлея в случае упругой полуплоскости // Докл. АН. 2004. Т. 395, № 4. C. 482-484.
7. Коссович Л. Ю, Кушеккалиев А. Н. Поле Рэлея в бесконечном упругом слое // Математика. Механика : сб. науч. тр. Саратов, 2003. Вып. 5. С. 159-161.
8. Ковалев В. А., Коссович Л.Ю., Таранов О. Г. Дальнее поле волны Рэлея для упругой полуполосы при действии торцевой нагрузки // Изв. РАН. МТТ. 2005. № 5. С. 89-96.
9. Коссович Л.Ю, Ковалев В. А., Таранов О. Г. Поле Рэлая в задаче Лэмба для цилиндрической оболочки // Известия вузов. Северо-Кавказский регион. Естеств. науки. Спецвыпуск. 2004. C. 52-54.
10. Ковалев В. А., Таранов О. Г. Расчленение нестационарного НДС цилиндрических оболочек при ударных торцевых воздействиях нормального типа // Смешанные задачи механики деформируемого твердого тела : материалы V Рос. конф. с междунар. участием / под ред. акад. Н.Ф. Морозова. Саратов, 2005. С. 78-82.
11. Коул Дж. Методы возмущений в прикладной математике. М., 1974. 274 с.
12. Гринченко В. Т., Мелешко В. В. Гармонические колебания и волны в упругих телах. Киев, 1981. 284 с.
13. Лаврентьев М.А., Шабат Б. В. Методы теории функций комплексного переменного. М., 1987. 688 с.
14. Kaplunov J.D., Kossovich L.Yu., Nolde E.V. Dynamics of Thin Walled Elastic Bodies. San Diego, 1998. 226 p.
УДК 539.374
СВЯЗАННЫЕ ТЕРМОУПРУГИЕ ВОЛНЫ ТРЕТЬЕГО ТИПА ЗАДАННОГО АЗИМУТАЛЬНОГО ПОРЯДКА В ВОЛНОВОДЕ С ПРОНИЦАЕМОЙ ДЛЯ ТЕПЛА СТЕНКОЙ
В. А. Ковалев, Ю. Н. Радаев*
Московский городской университет управления Правительства Москвы, кафедра прикладной математики;
* Институт проблем механики им. А. Ю. Ишлинского РАН, Москва E-mail: vlad_koval@mail.ru, y.radayev@gmail.com
Работа посвящена изучению распространения обобщенных связанных термоупругих волн заданного азимутального порядка в длинном цилиндрическом волноводе кругового поперечного сечения. При этом предполагается, что стенка волновода свободна от нагрузок и является проницаемой для тепла. Исследование проводится в рамках теории связанной обобщенной термоупругости третьего типа (GNIII), согласующейся с основными принципами термомеханики.
Данная теория сочетает оба известных типа распространения тепла в твердых деформируемых телах: термодиффузионный и волновой. Предельными случаями обобщенной термоупругости типа III являются классическая термоупругость (GNI/CTE) и гиперболическая термоупругость (GNII), которые могут быть сформулированы в терминах классической теории поля. Дифференциальные уравнения поля в этом случае принадлежат гиперболическому аналитическому типу. Методом разделения переменных в связанных уравнениях линейной термоупругости третьего типа получено их замкнутое аналитическое решение, которое удовлетворяет необходимым краевым условиям на боковой стенке волновода, в том числе условию конвективного теплообмена с окружающей средой. Установлено, что краевые условия на поверхности волновода выполняются отдельно для
Cross-Coupled Type-III Thermoelastic Waves of a Given Azimuthal Number in a Waveguide under Sidewall Heat Interchanging
V. A. Kovalev, Yu. N. Radayev*
Moscow City Government University of Management, Chair of Applied Mathematics * Institute for Problems in Mechanics RAS, Moscow E-mail: vlad_koval@mail.ru, y.radayev@gmail.com
The paper is devoted to a study of cross-coupled type-III generalized thermoelastic waves of a given azimuthal order propagating via a long cylindrical waveguide with circular cross-section. Sidewall of the waveguide is assumed free from tractions and permeable to heat. The study is carried out in the framework of coupled generalized theory of type-III thermoelasticity (GNIII) consistent with the fundamental principles of continuum thermomechanics. The type-III theory combines the both possible mechanisms of heat transfer: thermodiffusion and wave. Type-III generalized thermoelasticity includes classical thermoelasticity (GNI/CTE) and the theory of hyperbolic thermoelasticity (GNII) as limiting cases.
каждой из волн фиксированного азимутального порядка, поэтому волны различного азимутального порядка распространяются в волноводе независимо друг от друга. Для термоупругой волны заданного азимута построен частотный детерминант. Выполнен численный анализ частотного уравнения на предмет поиска его комплексных корней. При этом в частотном уравнении произведено выделение всех возможных однозначных ветвей квадратных радикалов. Детально описана схема локализации корней частотного уравнения и определены волновые числа связанных термоупругих волн, в частности, первого и седьмого азимутального порядков. Приведены результаты численного анализа в случае связанной волны азимутального порядка 70. Обсуждаются различные аспекты численной реализации предлагаемого подхода.
Ключевые слова: термомеханика, термоупругость, GNШ, частотное уравнение, волновод, волновое число, форма волны, азимутальный порядок, теплообмен.
The GNII-theory can be formulated as a field theory and differential field equations are of hyperbolic analytical type. Closed solution of the coupled linear GNIII-thermoelasticity partial differential equations satisfying the required boundary conditions on the surface of waveguide including convective heat interchanging condition is obtained by the separation of variables technique. For a given azimuthal number the frequency equation is derived. A numerical analysis of frequency equation is carried out by Mathematica. A scheme of frequency equation roots localization is proposed and wavenumbers of the coupled type-III thermoelastic waves of the first and seventh azimuthal numbers are computed. A numerical study of the coupled thermoelastic waves of the 70th azimutal number is also presented. Some aspects of numerical realization of the proposed approach are discussed.
Key words: thermomechanics, thermoelasticity, type-III thermo-elasticity, frequency equation, waveguide, wavenumber, wave mode, azimuthal number, heat interchanging.
1. ВВОДНЫЕ ЗАМЕЧАНИЯ И ОСНОВНЫЕ УРАВНЕНИЯ
В настоящее время термин «термоупругость» охватывает достаточно широкий спектр таких явлений, как теплопроводность, термические напряжения, связанные термоупругие деформации, затухание тепловых и упругих импульсов в твердых телах и, конечно же, тепловые волны «второго звука» в деформируемых твердых телах. С чисто теоретической точки зрения термоупругость следует рассматривать как одну из важнейших составляющих термомеханики и теории поля. В частности, термоупругость второго типа GNII допускает формулировку в терминах классической теории поля [1-3].
Первые исследования температурных напряжений и деформаций в деформируемых твердых телах в линейном приближении восходят к J.M.C. Duhamel (1837) и F.Neumann (1885) [4]. Одна из первых работ на русском языке, посвященная температурным напряжениям и содержащая все известные к тому времени достаточно законченные формы решений прикладных задач термоупругости, — монография Н.Н. Лебедева [5].
Классическая теория термоупругости (CTE, conventional thermoelasticity) основывается на законе теплопроводности Фурье, содержание которого, как известно, сводится к следующему: как только в заданном месте пространства возникает градиент температуры V9, мгновенно (т. е. без всякой задержки по времени) устанавливается противоположно направленный поток тепла:
h = -A, V9,
где h — вектор потока тепла (heat flux), Л, — коэффициент теплопроводности (thermal conductivity), 9 — температура, V — пространственный оператор Гамильтона. Таким образом, закон теплопроводности Фурье устанавливает коллинеарность вектора потока тепла и пространственного антиградиента температуры. Известно, что соответствующее закону Фурье уравнение теплопроводности принадлежит параболическому аналитическому типу; уравнения такого типа допускают бесконечные скорости распространения возмущений, а в том случае, когда решение имеет «волновой» характер, тепловые волны будут иметь затухающие амплитуды. Совершенно ясно, что возможность распространения теплового сигнала с бесконечно большой скоростью явно не соответствует реальным физическим наблюдениям и нарушает принцип причинности.
По существу, парадокс бесконечной скорости распространения температурных возмущений впервые обсуждался J.C. Maxwell [6]. M.A. Biot [7] впервые методами термодинамики необратимых процессов корректно построил теорию связанной термоупругости; связанные соотношения этой теории, классифицируемые в настоящее время как CTE, включают векторное уравнение движения, принадлежащее гиперболическому типу, и параболическое уравнение теплопроводности. Очевидно, что теория Biot в полной мере обладает указанным выше недостатком, поскольку предсказывает конечную скорость распространения упругой волны и — бесконечную для теплового импульса. «Полевые»
уравнения CTE для связанных (сопряженных) перемещений и температуры в линейном приближении состоят из одного векторного и сопряженного ему скалярного связанных дифференциальных уравнений в частных производных:
J (Ди + (Л + ^)VV ■ u - aV9 - pU = 0,
( Л*Д9 - к9 - a9oV ■ иi = 0. ( )
В этих уравнениях Л— упругие постоянные Ламе (для изотермической деформации); Д = V ■ V — оператор Лапласа; и — вектор перемещения; а — термомеханическая постоянная (а = 1/3(3Л+2()в*); в* — коэффициент объемного теплового расширения; 9 — абсолютная температура, точнее, ее превышение (инкремент) над отсчетной (референциальной) температурой 90; 90 — отсчетная (равновесная) температура (при температуре, равной 90, отсутствуют деформации и напряжения); к — теплоемкость (на единицу объема) при нулевой постоянной деформации. Связанность уравнений CTE-теории обусловливается взаимным влиянием нагрева (охлаждения) и объемной деформации тела, проявляется через определяющую постоянную a и является следствием основных принципов термодинамики (см.
[7]).
Второе из уравнений системы (1), сформулированное для недеформируемого проводника тепла, представляет собой одно из наиболее важных линейных параболических дифференциальных уравнений математической физики — уравнение диффузии:
Д9 - Л*9 = 0. (2)
За последние десятилетия было предложено достаточно много различных обобщений теории связанной термоупругости Biot с тем, чтобы устранить парадокс бесконечной скорости распространения тепла и допустить возможность волнового характера транспорта тепла в форме незатухающих волн (undamped heat waves). Последнее из указанных явлений называется «вторым звуком». Обзор работ по проблеме «второго звука», включая весьма интересный исторический аспект, дан в статьях [8]. Отметим, что волны и волновые процессы очень широко распространены в природе, однако дать точное и ясное определение того, что же считать волной, не совсем просто, поскольку грань между волновым и неволновым поведением часто оказывается размытой. Традиционный отличительный критерий волнового процесса — перенос энергии в пространстве, не сопровождающийся переносом массы.
Эксперименты, проведенные в разные годы различными исследователями, подтверждают вывод о том, что тепло при определенных условиях может распространяться как незатухающая волна. Первоначально (1946 г.) термические волны «второго звука» были экспериментально обнаружены в жидком гелии и позднее (1966 г.) — в твердом гелии. Феномен второго звука наблюдается при весьма низких температурах в кристаллах высокой чистоты. Например, в работе [9] указывается на обнаружение «второго звука» в кристалле NaF. «Второй звук» в NaF наблюдается при температуре около 15K [10]. В работах [11] и [12] также подтверждается существование волн второго звука в NaF; согласно [13] скорость волн «второго звука» составляет 1953.1 м/с. В публикации [14] приводятся результаты измерений волн «второго звука» в висмуте (Bi) при температуре около 3K. Скорость волн «второго звука» в висмуте оказывается равной 78 ± 5 м/с.
Остановимся на нескольких наиболее часто встречающихся в прикладных задачах термомеханики обобщениях основополагающих уравнений CTE. Одно из них было дано H.Lord и Y. Shulman (LS-theory) [15]. Ими было получено гиперболическое уравнение распространения тепла, основанное на новом законе теплопроводности, заменившем классический закон теплопроводности Фурье. В законе LS-теплопроводности введена дополнительная определяющая постоянная — время релаксации Trei в соответствии с уравнением C. Cattaneo:
h + Treldh = -Л* V9. (3)
dt
Время релаксации представляет собой (всегда весьма малое) время запаздывания после возникновения температурного градиента, необходимое для того, чтобы в элементе тела сформировался установившийся поток тепла. Для металлов величина Trei заключена в пределах от 10-12 с до 10-9 с.
Независимо от C. Cattaneo аналогичное уравнение было получено в работе [17]. Уравнение (3) часто выводится на основании уравнения «запаздывания» (the lag equation)
h(x, t + Trei) = -A* V0(x, t), (4)
которое формулируется для заданного места в пространстве x и, как уже отмечалось, устанавливает, что поток тепла запаздывает по отношению к пространственному градиенту температуры. Принимая во внимание, что
d h
h(x, t + Trei) = h(x, t) + Trei (x, t) + . . . ,
где частное дифференцирование по времени выполняется при фиксированном положении x, и ограничиваясь лишь первыми двумя слагаемыми в приведенном разложении, приходим к уравнению (3).
«Уравнения поля» LS-теории для связанных перемещений и температуры в линейном приближении состоят из одного векторного (совпадающего с соответствующим уравнением CTE-теории) и сопряженного ему скалярного гиперболического связанных дифференциальных уравнений в частных производных:
J ^Au + (Л + ^)VV ■ u - aV0 - pU = 0,
\ Л* Ав - к(в + Trei0) - a^a(V ■ U + TrelV ■ U) = 0. ( )
Второе из этих уравнений, сформулированное для недеформируемого тела, оказывается одним из важнейших гиперболических дифференциальных уравнений математической физики — телеграфным уравнением (telegraphist's equation):
Ав - ^в - -^0 = 0. (6)
Л* Л*
Тем самым, в LS-теории как уравнение движения, так и уравнение транспорта тепла, принадлежат гиперболическому типу, что и обеспечивает конечную скорость распространения и тепловых, и упругих волн деформации. В отечественной научной литературе первое упоминание о времени релаксации теплового потока и конечной скорости распространения теплового возмущения было сделано А. В. Лыковым в 1967 г. [18].
Вторым хорошо известным обобщением связанной теории термоупругости Biot, допускающим «второй звук», является теория, предложенная A.E. Green, K.A. Lindsay (GL-theory) [19]. GL-теория характеризуется тем, что вектор потока тепла в термоупругом теле h зависит от скорости изменения абсолютной температуры в и градиента температуры; закон GL-теплопроводности имеет форму
h = -ьв - л* ■ ve,
где Ь — антисимметричный вектор, Л* — тензор теплопроводности. Тензор второго ранга Л* в GL-теории симметричен: Л* = Л^. Заметим, что в случае центральной материальной симметрии выполняется равенство Ь = 0. Уравнение транспорта тепла оказывается гиперболическим. GL-теория в случае Ь = 0 допускает «второй звук», не нарушая при этом классического закона теплопроводности Фурье.
Определяющее уравнение GL-теории связанной термоупругости имеет вид
о = 2^s + Atr £ - а + т* Jt) (в - во)1, (7)
где о — тензор напряжений, £ — тензор малых деформаций, I — единичный тензор.
Замкнутая система дифференциальных уравнений GL-теории состоит из уравнения движения
^Au + (A + ^)VV ■ u - a ( 1 + т*Ve - pU = 0 (8)
dt
и уравнения распространения тепла
д
Л* Ав - аво V ■ U - к( 1+ т* —J 0 = 0. (9)
Определяющие постоянные ОЬ-теории т* и т*, имеющие смысл времен релаксации, подчиняются неравенствам т* > т* > 0, гарантирующим неотрицательность внутреннего производства энтропии при распространении волн деформаций и температуры. Если т* = т* = 0, то ОЬ-теория сводится к
Заметим, что теории ЬБ и ОЬ принципиально отличаются друг от друга и не получаются одна из другой как частные случаи.
Уравнение (9), сформулированное для недеформируемого теплопроводящего тела, представляет собой телеграфное уравнение вида
В начале 90-х гг. XX в. A.E.Green и P.M.Naghdi [20,21] была развита теория обобщенной термоупругости (GN-theory). Предложенная теория сочетает в себе как свойства классической термоупругости CTE, выстроенной согласно закону теплопроводности Фурье, так и свойства недиссипативной термоупругости, предполагающей отсутствие производства энтропии и волновой характер распространения теплового сигнала. GN-теория была сформулирована в трех различных термодинамически корректных вариантах: GNI, GNII и GNIII. В линейном приближении первый вариант (GNI) приводит к закону теплопроводности Фурье и параболическому уравнению теплопроводности; второй (GNII)
— предлагает считать распространение тепла как волновой недиссипативный процесс (dissipationless thermoelasticity), не сопровождающийся внутренним производством энтропии, и приводит к гиперболическому уравнению транспорта тепла; третий (GNIII) — наиболее общий из рассматриваемых
— включает GNI и GNII в качестве предельных случаев. Именно по этой причине GN-теория в состоянии моделировать значительно более широкий круг явлений, связанных с переносом тепла, по сравнению с теорией Фурье.
Классическая связанная теория термоупругости (CTE) широко применяется в различных прикладных задачах термомеханики в тех ситуациях, когда речь идет о быстропротекающих переходных процессах или когда происходит интенсивный нагрев тела (например, с помощью импульсного лазерного излучения). В таких ситуациях температурное поле связано с упругим полем и связанная форма уравнений термоупругости должна выглядеть наиболее адекватной. Систематическое изложение динамической теории связанной (сопряженной) CTE-термоупругости дано в известной монографии В. Новацкого [22].
Отличительной чертой теории GNII гиперболической термоупругости является то, что она полностью согласуется с основными принципами термодинамики (и в этом смысле она термодинамически корректна) и предсказывает нулевое внутреннее производство энтропии, т. е. отсутствие рассеяния энергии, при распространении термической волны «второго звука». Отсутствие диссипации энергии в GNII-термоупругом теле позволяет дать вариационную (а следовательно, и полевую) формулировку теории с целью дальнейшего поиска вариационных симметрий связанной системы дифференциальных уравнений в частных производных термоупругости. Лагранжева и гамильтонова полевые формы нелинейной теории термоупругости типа GNII рассматривались в работах G.A. Maugin, V.K. Kalpakides
Анализ литературных источников показывает, что к настоящему времени опубликовано сравнительно немного работ, выполненных с привлечением связанных уравнений ОЫШ-теории. Например, распространение плоских гармонических связанных термоупругих волн в рамках ОЫШ-теории изучалось лишь сравнительно недавно в работах [27] и [28], и нельзя сказать, что эта проблема полностью разрешена в плане определения нормальных волновых чисел указанного типа волн. В статье [29] в рамках теории обобщенной термоупругости ОЫШ исследуется термоупругие волны в неограниченном однородном изотропном теле, вызванные линейным источником тепла. Полученные в указанной работе результаты показывают, что теория ОЫШ в целом прогнозирует диффузионный механизм распространения тепла и только в отдельных случаях допускает волновой механизм теплопроводности. В [1,30] получено решение задачи о рапространении связанной ОЫШ-термоупругой волны вдоль теплоизолированного цилиндрического волновода. В работе [31] исследуется частотное уравнение указанных волн в случае достаточно высоких азимутальных чисел.
Целью настоящей работы является изучение распространения обобщенных связанных ОЫШ-тер-
CTE.
Д0
(10)
[23-26].
моупругих волн вдоль цилиндрического волновода со свободной теплопроницаемой стенкой. Одной из задач исследования является определение влияния теплообмена через стенку цилиндрического волновода на величину волновых чисел распространяющихся в волноводе нормальных волн и форму волн деформаций и температуры.
2. НЕЛИНЕЙНЫЕ УРАВНЕНИЯ ТЕРМОУПРУГОСТИ ТИПА III И УРАВНЕНИЯ ТЕРМОДИНАМИЧЕСКОЙ ОРТОГОНАЛЬНОСТИ
Для представления деформации термоупругого тела воспользуемся классическим лагранжевым отсчетным описанием: x = x(X, t), x — пространственное положение точки, которая занимала место X в отсчетном положении. В дальнейшем изложении систематически будет использоваться отсчетный оператор Гамильтона Vr, ассоциированный с лагранжевой переменной X.
Система основных нелинейных диффференциальных уравнений связанной термоупругости типа GNIII состоит из:
— уравнения баланса массы
% =о; (ii)
dt X
— уравнения баланса импульса
PRx = Vr ■ S, (12)
где S = JF-T ■ T — первый тензор напряжений Пиола—Кирхгофа, T — тензор напряжений Коши, F = VR ® x — градиент деформации, J = detF — якобиан деформации;
— уравнения баланса внутренней энергии
e = -Vr ■ hR + tr(S ■ FT) + e, (13)
где e — плотность (в расчете на единицу объема в отсчетном состоянии) внутренней энергии, hR — вектор потока тепла (в единицу времени через единицу площади в отсчетном состоянии), e — плотность объемных источников тепла (лучистое тепло);
— уравнения баланса энтропии
S = -Vr ■ jR + a + С, (14)
где s — плотность энтропии (в расчете на единицу объема в отсчетном состоянии), jR — вектор потока энтропии (в единицу времени через единицу площади в отсчетном состоянии), a — внешнее производство энтропии, С > 0 — внутреннее производство энтропии.
В качестве термодинамических переменных состояния в термоупругости GNIII выбираются переменные (tf — температурное смещение (thermal displacement))
tf, tf, Vr^, Vr^, F. (15)
Далее полагается, что допустимы лишь такие процессы, когда выполняется уравнение баланса энтропии (14) при ограничениях (11)—(13). Возможная двойственная формулировка: выполняется уравнение баланса энергии (13) при ограничениях (11), (12), (14). Исходная и двойственная формулировка выступают как проявление принципа непротиворечивости: баланс энтропии (энергии) не должен нарушать баланса массы, импульса, энергии (энтропии).
Ограничения в форме дифференциальных уравнений (11)—(13) учитываются с помощью множителей Лагранжа. В результате можно получить соотношения (0 — абсолютная температура, обратное значение которой («холодность») выступает в качестве множителя Лагранжа)
¿jR = hR, 0а = e. (16)
На основании уравнений (11)—(14) и соотношений (16) выводится уравнение
-(ф + s0) + tr (S ■ FT) - jR ■ Vr0 = ¿С, (17)
где ф = e - s0 — плотность (в расчете на единицу объема в отсчетном состоянии) свободной энергии Гельмгольца.
Дифференциальные уравнения термоупругости должны быть дополнены определяющими уравнениями. Определяющие уравнения связанного термоупругого поля обязаны иметь форму
^ = ^(tf,tf, Vrtf, Vrtf, F),
hR = hR(tf,tf, VRtf, VRtf, F), S = S(tf,tf, VRtf, VRtf, F),
s = s(tf,tf, Vr tf, Vr tf, F), ( )
С = C(tf,tf, VRtf, VRtf, F),
9 = 9(tf).
Дальнейшие ограничения на форму определяющих уравнений (18) могут быть получены следующим образом. Выполнив подстановку определяющих зависимостей (18) в уравнение (17) и замечая, что коэффициенты при tf, Vr tf, F должны обращаться в нуль, заключаем, что свободная энергия в действительности не может зависеть от переменной состояния VRtf:
^ = ^(tf,tf, VRtf, F);
кроме того, должны выполняться следующие равенства:
s = - ж—, (19)
Sdtf
S = i- (20)
Для внутреннего производства энтропии, следовательно, остается выражение
9С = - (jR + dVRtf) ■ VRtf - dtftf- (2,)
Предполагая, что определяющие уравнения (18) не содержат явных вхождений температурного смещения tf, находим ( )
9С = - (jR + JRi) ■ VRtf- (22)
Вводя затем потенциал рассеяния (точнее говоря, сопряженный потенциал рассеяния) согласно 9С = D = D(tf, VRtf, F; VRtf) и рассматривая его как функцию от термодинамической силы -VRtf, приходим к определяющему уравнению в форме соотношения термодинамической ортогональности1
в пространстве термодинамических сил термодинамического потока jr + —-- поверхности уров-
д VRtf
ня потенциала рассеяния D(tf, VRtf, F; VRtf) = const. В гладких точках поверхности уровня имеем
следующее уравнение градиентальности:
. + = Л* д D (23)
Jr + q =-Л •, (23) д VRtf д VRtf
или
. = Л* д D (24)
Jr = q - Л • - (24)
д VRtf д VRtf
Множитель Л* в уравнении (23) вычисляется на основании определяющего потенциала рассеяния D равенства
9С = D = Л* (VRtf) ,
дVR tf
'Оригинальный подход был развит H. Ziegler в ряде публикаций начиная с 1958 г. как обобщение линейной теории L. Onsager, восходящей к 1931 г. См.: Циглер Г. Экстремальные принципы термодинамики необратимых процессов и механика сплошной среды. М., 1966. 134 с.; Ziegler H. Proof of an orthogonality principle in irreversible thermodynamics // Zeitschrift fur Angewandte Mathematik und Physik (ZAMP). 1970. Vol. 21, № 6. P. 853-863; Ziegler H. Discussion of some objections to thermomechanical orthogonality // Archive of Applied Mechanics. 1981. Vol. 50, № 3. P. 149-164.
откуда сразу же следует, что
Л* =-D. (25)
(Vrtf) ■ ^ } d Vr tf
Таким образом, подставляя (25) в уравнение (24), закон теплопроводности в теории термоупругости типа GNIII с учетом принципа термодинамической ортогональности будет иметь следующий вид:
. =_ J±___D д D (26)
jR = - дVrtf - (vrtf) дVrtf. (26)
^ dVRtftf
Выполнение неравенства С > 0 обеспечивается выпуклостью поверхностей уровня D(tf, VRtf, F; Vri?) = const в пространстве термодинамических сил.
Полученное выше уравнение (26) устанавливает определяющий закон теплопроводности в GNIII теории, основанный на принципе ортогональности термодинамических потоков и сил. Очевидно, что термодинамически корректен и такой вариант, когда процесс термоупругого деформирования не сопровождается внутренним производством энтропии. В этом случае D = 0, и закон теплопроводности приобретает форму
jR = - (27)
д VRtf
характерную для теории GNII. Это определяющее уравнение для вектора потока энтропии — одно из самых замечательных в термомеханике континуума, поскольку поток энтропии, а следовательно, и поток тепла определяются согласно (27) заданием свободной энергии.
Теория связанной термоупругости GNI/CTE основывается на допущении о независимости свободной энергии ф от температурного смещения tf и градиента температурного смещения Vtf. Поэтому следует полагать, что
ф = ф^, F).
В результате внутреннее производство энтропии вычисляется как
0С = -jR ■ Vrtf.
Уравнение термодинамической ортогональности, обобщающее закон теплопроводности Фурье, имеет форму
i = D дР (28) jR =--^^ ^. (28)
(Vrtf) дVRtf
^ R J д^^Rttf
Итак, в теории GNIII, дополненной принципом термодинамической ортогональности, вектор потока энтропии определяется, как следует из уравнения (26), через два потенциала: свободную энергию ф и потенциал рассеяния D. Он естественным образом распадается на две составляющие: обратимую
дф
д Vr $
и необратимую
D д D
(VR^^ dVR$'
V ^ dVR$
3. ЛИНЕЙНЫЕ СВЯЗАННЫЕ УРАВНЕНИЯ ТЕРМОУПРУГОСТИ ТИПА III
На практике обычно применяются линеаризованные уравнения термоупругости. Не вдаваясь в подробности, приведем замкнутую систему уравнений линейной теории термоупругости третьего типа (GNIII).
В линейном приближении GNIII термоупругое поведение тела подчиняется определяющему закону Дюгамеля - Неймана:
& = 2де + (Atr е - а(0 - 0q)) I, (29)
где a — тензор напряжений, г — тензор малых деформаций, I — единичный тензор. Мы по-прежнему используем следующие обозначения: Л, д — упругие постоянные Ламе; a — термомеханическая постоянная (a = 1/3(3Л + * ); в * — коэффициент объемного теплового расширения (coefficient of thermal expansion); в — абсолютная температура; в0 — отсчетная (равновесная) температура. Отметим, что при температуре, равной в0, отсутствуют деформации и напряжения.
Линейная теория GNIII-термоупругости характеризуется тем, что вектор потока тепла h линейно зависит как от градиента температуры, так и от градиента температурного смещения:
h = -Л* ve - AVtf, (30)
где Л — характерная скорость теплопроводности (thermal conductivity rate), V — трехмерный оператор Гамильтона, tf (tf = в) — температурное смещение.
В рамках линейной теории следует предполагать линейную зависимость между термодинамическим потоком:
д-0
V7 + д^
и термодинамической силой
равно как и квадратичную зависимость свободной энергии ф от градиента температурного смещения В результате определяющее уравнение для теплового потока будет иметь вид (30).
Таким образом, уравнение (30) получается, как следствие (26), при условии, что свободная энергия квадратично зависит от градиента температурного смещения а потенциал рассеяния — от градиента температуры и, следовательно, не противоречит принципу термодинамической ортогональности.
Заметим, что внутреннее производство энтропии в модели ОЫП-термоупругости исчезает £ = 0 в силу
дф
j =
д V
Классическая связанная термоупругость ОЫ1/СТЕ основывается на допущении о независимости свободной энергии ф от температурного смещения $ и градиента температурного смещения в силу чего внутреннее производство энтропии вычисляется как
0£ = -7 ■
Кроме того, полная система соотношений ОЫШ-теории включает уравнения движения:
ё1у а — ри = 0, (31)
где р — плотность, и — вектор перемещений; уравнение баланса энтропии:
5 + V ■ 7 = а + £, (32)
где 5 — плотность (на единицу объема) энтропии, 7 — вектор потока энтропии, а — внешнее производство энтропии, £ > 0 — внутреннее производство энтропии, или
05 + V ■ (07) = 7 ■ ^ + 0(а + £); (33)
и уравнение баланса энергии, которое мы примем в форме
— (ф + 50) + Ъ(а ■ ¿) — Ь ■ V0 = 0£, (34)
где ф — плотность (на единицу объема) свободной энергии Гельмгольца. Здесь необходимо учесть, что вектор потока тепла и вектор потока энтропии связаны уравнением
Ь = 07, (35)
а для внутреннего производства энтропии справедливо соотношение
0С=- (j+sVy ■V • (36)
Приведем также соотношения Коши, связывающие тензор малых деформаций е и градиент вектора перемещения V 0 u:
2е = V 0 u +(V 0 u)T . (37)
Условие конвективного теплообмена с окружающей средой через поверхность с единичной нормалью n в линейном приближении имеет вид
П ■ h = а (0 - 0env), (38)
где а — коэффициент теплообмена, 0env — температура окружающей среды. В дальнейшем будем полагать, что температура окружающей среды совпадает с отсчетной 0env = 0О.
Для сокращения записи уравнений в дальнейшем через 0 будем обозначать превышение температуры над отсчетной (равновесной) температурой 0О, т.е. символ 0 в последующем изложении следует понимать как разность 0 - 0О.
Приведем линейные связанные уравнения движения и теплопроводности GNIII-термоупругости [1,2] в той форме, в которой они используются в настоящей работе:
J ^Au + (Л + д) VV ■ u - aV0 - pU = 0,
\ ЛА0 + Л* A0 - кв - авоV ■ u = 0. ( )
Здесь A — трехмерный оператор Лапласа, к — (specific heat of the unit volume) теплоемкость (на единицу объема) при отсутствии деформации2.
Если теплопроводящая среда недеформируема и определяющий параметр Л* ^ 0, то термическое поле согласно уравнениям GNIII удовлетворяет классическому волновому уравнению:
A0 - Кв = 0. (40)
Л
Разделим второе уравнение системы (39) на Л и 0О. В дальнейшем постоянные Л, Л* и к будут считаться отнесенными к отсчетной температуре 0О. Таким образом, можно минимизировать число постоянных, необходимых для формулировки связанных уравнений GNIII-теории. В итоге приходим к системе уравнений
^Au + (Л + ^)VV ■ u - aV0 - pu = 0,
Л* • к •• а_ .. (41)
A0 + — A0 - -0 - — V ■ u = 0. Л Л Л
Как упоминалось ранее, классическая связанная термоупругость (GNI/CTE) и гиперболическая термоупругость (GNII) являются предельными случаями GNIII-теории. Заметим, что система (41) исключительна удобна при переходе к уравнениям гиперболической GNII-теории: для этого необходимо положить Л* = 0. Ясно также, что переход Л ^ 0 к классической термоупругости (GNI/CTE) в уравнениях (41) требует соблюдения ряда мер предосторожности.
Заключая этот раздел работы, заметим, что из приведенных выше прямых тензорных соотношений (29), (37) находятся физические компоненты тензора напряжений и деформаций в цилиндрических координатах r, (, z
л ^ диг > f 1 ur \
arr = Ле + 2д—--ав, а^^ = Ле + 2д - —^ +--- ав,
дг \r д( r )
а„ = Ле + 2д^ - ав, а„, = д - ^ + ^ , (42)
дг V дг r r <У( )
(ди2 д^ \ /1 ди2 дп,п
arz = Д + Ь = Д +
\ дr дг ) \r д( дг
2Поскольку рассматривается линейная теория, к не зависит от температуры для ее значений в окрестности референциальной температуры.
(43)
£ £rr + + £zz ;
dur 1 du^ ur duz
£rr — j — - + , £zz — j
dr r r dz
1 duz du^ duz dur du^ u^ 1 dur
2£œz — Ô + ""^ j 2£rz — ""^ + ^ j 2еги — ""^ + т; •
r dz dr dz dr r r
4. РАЗДЕЛЕНИЕ ПРОСТРАНСТВЕННЫХ ПЕРЕМЕННЫХ В УРАВНЕНИЯХ GNIII-ТЕОРИИ В ЦИЛИНДРИЧЕСКИХ КООРДИНАТАХ
Будем рассматривать гармоническую зависимость перемещений и температуры от времени: u — Ue-iWt, 9 — ©e-iWt, где u — циклическая частота; U, © — комплексные амплитуды. Представим вектор комплексной амплитуды U в виде разложения Гельмгольца:
U — VФ + V х Ф, (44)
где Ф — скалярный потенциал, Ф — векторный потенциал. При этом необходимо учесть условие калибровки векторного потенциала
V ■ Ф — 0. (45)
Исходная система уравнений (41) будет удовлетворяться, если потенциалы Ф, Ф и амплитуда © будут удовлетворять следующей системе уравнений:
' ^2рФ + дДФ = 0,
(А + 2д)ДФ - а© + w2рФ = 0, (46)
л ^л • Л* л ^ aw2
Д© - Д© + —:-© + —— ДФ = 0. Л Л Л
Потенциал Ф удовлетворяет независимому уравнению, векторному уравнению Гельмгольца:
ДФ + k\ Ф = 0, (47)
где — квадрат волнового числа чисто упругой поперечной волны.
В цилиндрической системе координат в векторном уравнении (47) удается разделить пространственные переменные. Для определения векторного поля Ф будем искать его компонеты в виде
Фг (r,p,z) = Ar (r)i sin ^ !e±ifcz, Ф^ (r,p,z) = Á9 (r){ COS ^ 1 e±ifcz, I cos n^ I I — sin n^ I
Фz (r,p,z) = Az (r)i Sin ^ je^, I cos n^ I
где n — азимутальный порядок волны, k — волновое число распространяющейся в волноводе обобщенной термоупругой волны заданного азимутального порядка n. Физические компоненты Ar (r), A^(r), Az (r) находятся в форме линейных комбинаций модифицированных функций Бесселя (см. [1]). Окончательно, компоненты Фг, Ф^, Фz принимают вид
Фг(r,^,z) = (CsIn-1 (q2r)+ C4/n+i(q2r))i sin ^ 1 e±ikz,
I cos n^ I
ф„(г,^) = (сз1п-1(q2r) — C4In+1 (q2r)) i ^^ 1 e±ikz, (48)
I — sin n^ I
Фz (r,^,z) = C5In (q2r)( sin ^ 1 e±ikz, I cos n^ I
где In(C) — модифицированная функция Бесселя целого порядка n; = k2 — k2; C3, C4, — произвольные постоянные.
Два последних уравнения в системе (46) связывают скалярный потенциал Ф и комплексную амплитуду ©. Следуя работе [1], определим Ф и © через новый скалярный потенциал О согласно
Ф = аО, © = ЬО;
при этом постоянные а и Ь могут принимать в точности два различных значения
а з = р
,
К 2
Ьз = Л Л-
Здесь мы вводим следующие обозначения:
р2 = к2 - 72, £2 = к2 - ЛЛ2Л-2^,
Л = Л21 + ¿Л3
1 + Л3
(0 = 1, 2).
Л2 = :±
Л2 = 12 ,
(49)
(50)
Л3 =
Л-^
"Г"'
7-2 = К
1 = Л'
где постоянные 7з- (0 = 1, 2) имеют смысл волновых чисел и определяются ниже, кЦ — квадрат волнового числа чисто упругой продольной волны. Потенциал О удовлетворяет скалярному уравнению Гельмгольца:
ДО + 720 = 0 (7 = 71,72). (51)
Остановимся на определении постоянной 7. Подставляя равенства (49) во второе и третье уравнения (46), приходим к линейной алгебраической системе уравнений относительно а и Ь (см. [1,28])
(72 (Л + 2д) - ра>2) а + аЬ = 0,
' X
22
г^ л а + i7 -1 л
7
Ь = 0.
(52)
Для того чтобы существовало нетривиальное решение относительно а и Ь этой системы линейных алгебраических уравнений, необходимо, чтобы определитель системы обращался в нуль. Вычисляя определитель системы, приходим к биквадратному уравнению, из которого находятся волновые чис-
ла 7:
1
Ь № -1) + ё (Л? - ¿Л§) - Л2 = 0,
к4
к2
л? = 1+л2+а-.
1 2 рЛ
(53)
Квадраты волновых чисел, следовательно, вычисляются в виде (см. [1,21])
272 _ ¿Л§ - Л2 + л/(Л? - ¿Л2)2 +4Л2(гЛ2 - 1)
к2
¿л3 -1
Извлекая квадратный корень, приходим к выражению
2т^ = ¿Л - Л? + а?, 2 + ¿Ь? 2 2 к2 =
¿Л3 -1
(54)
где а? 2 и Ь? 2 вычисляются согласно
а? ,2 = ±
- л4 - 4Л2 + ^/(Л4 - лз - 4Л2)2 + 4Л4 (2Л2 - л?) У
Выражения для самих волновых чисел 7 будут следующими:
2 (Л3 + 1) 7 = а?, 2; 3, 4 + ¿Ь?, 2; 3, 4,
Ь? 2 =
Л2 (2Л2 - Л2) (55)
а?, 2
(56)
где введены обозначения [1]
П, 2;3, 4
=
£2 + Ь?, 2Л2 - а?, 2 + л/(^2 + Ь?, 2Л2 - а?, 2)2 + (Т2 - а?, 2^2 - Ь?, 2)2
Ь?, 2; 3, 4 =
Т2 - а?, 2Л3 - Ь?, 2
2а
3
?, 2; 3, 4
(57)
2
2
2
при этом S2 и T2 выражаются в соответствии с формулами
s 2 = + h2, T2 = h2(h2 -1).
(58)
Следует отметить, что 72 2 есть квадраты волновых чисел плоской связанной термоупругой волны третьего типа (ОЫШ) [28].
Вернемся к определению потенциала О. В цилиндрической системе координат дифференциальное уравнение в частных производных (51) принимает вид
д2 О I дО 2^-0
dr2 r2 д^2 дг2 r дг Y
а соответствующие решения для О находятся в виде
(j -1, 2),
(59)
0|Y=Yj - in(Pjr)
cos n^
e
(j -1, 2).
sin n^ I
Окончательно скалярный потенциал Ф и комплексная амплитуда © определяются в форме
(60)
Ф - (С (р2 - g2) In(pir) + C2 (р2 - g2) in(P2r)) cos ^ L e±ikz
sin n^
© - h^ (Ci7l2in(pir) + C2Y22in(P2r) ) ( COs^ } e± A* — sin n^
(61)
где C1 и C2 — произвольные постоянные.
5. ГРАНИЧНЫЕ УСЛОВИЯ НА БОКОВОЙ ПОВЕРХНОСТИ ВОЛНОВОДА И ЧАСТОТНОЕ УРАВНЕНИЕ ТЕРМОУПРУГОЙ ВОЛНЫ ЗАДАННОГО АЗИМУТА
В рассматриваемой задаче предполагается, что боковая поверхность волновода свободна от нагрузок, т. е. должны удовлетворяться граничные условия
Отт 1г = Д - 0' |Г = Д --|г = Л - °
(62)
где игг, — физические компоненты тензора напряжений.
С окружающей средой через боковую стенку волновода происходит конвективный теплообмен по линейному закону (38). Положив 9епу = 90, разделим уравнение (38) на 90, затем разность 9 — 90 обозначим, как было оговорено ранее, через 9 и продифференцируем обе части полученного уравнения по времени. В результате, приходим к условию
Л— + Л* — + о9 дr дr
- 0,
(63)
r=R
и — коэффициент теплообмена, отнесенный к отсчетной температуре 90.
На основании (44) физические компоненты перемещений в цилиндрической системе координат выражаются (с точностью до гармонической экспоненты е-г^) через потенциалы Ф и Ф
дФ 1 дФ,
ит - ——I—
дr r д^
дф
дг
- 1 дФ д Фт
u,n — — Ь -г r д^ дг
дФ, u - + + 1 дФт
д^ z дг r дr r д^
(64)
Принимая во внимание соотношения (48), (61), с помощью (64) получаем представления для перемещений и температуры в гармонической связанной ОЫШ-термоупругой волне:
Ci (p2 - g2) (nin(pir) + piin+i(pir)) + C2 (p2 - g2) (^in^r) + p2in+i(p2r)) +
n 2n
+^Сбin(q2r) + C3—In(q2r) + C3in+i(q2r) - C4in+i(?2r)
COs I ±ifcz-¿Wt sin n^ I
—
(±гк) (С3 —/п(92г) + Сз/п+1(д2г) + С^/п + ^г) | - С5 (П/п(92г) + 92/п+1 (92г)) -V 92 г у г
-- {б1 (р2 - А2) 1п(Р1 г) + б2 (р2 - д2) 1п (Р2 г)}
в1П П^ I соэ п^ (
иг —
(±*к) {б1 (р? - А2) /п(р1 г) + б2 (р2 - А2) 1п (Р2 г)} +
+ (бз - 64)921п(92Г)
СОЭ 1±гк*-г^
Э1п п^
в — ^ (С1 (к2 - р1)/п(Р1 г) + 62(к2 - р2)1п(Р2г))
А:
СОЭ 1±гк*-г^
э1п П^ I
Далее, удовлетворяя граничным условиям (62), условию конвективного теплообмена (63), а также калибровочному условию (45), приходим к линейной однородной алгебраической системе уравнений относительно пяти постоянных Сз:
Аз С — 0 (м — 1, 5).
(65)
Приравнивая нулю определитель системы (65), который мы в дальнейшем будем обозначать через получаем частотное уравнение Б — 0. В случае теплоизолированного цилинидрического волновода элементы частотного определителя приведены в [1,30]. В рассматриваемой задаче с конвективным теплообменом через боковую поверхность волновода элементы частотного определителя Б41 и Б42 будут отличаться от приведенных в указанных работах. Для удобства ниже выписываются все элементы частотного определителя Б (п — 1, 2,...):
Бц —
1
(Р2 - 52) ((П - П + Р1 )/п (Р1) - Р1/п + 1(Р1^ +
к2 - 92
+ (292к"- ^^ ^ ) ((П - П2 - к2 )1п(Р1) + Р11п + 1 (Р1 )) - ^2 (к2 - Р1)/п (Р1 ),
Б12 — 7 2 1 2 (Р2 - 52) ((п2 - П + Р2Кп(Р2) - Р2/п + 1(Р2)) + к 91
+ (2^ ^^^) ((П - П2 - к2)1п(Р2) + Р2/п+1 (Р2)) - ^2(к2 - р2)1п(Р2),
Б13 — ^(¿к)—^ рп2 -2П + 92 (/п-1(92) - /п+1 (92)) + (П - 1)/п+1ы1 ,
к2 - 92
2П
Б14 — т(гк) к2 - 92 (-921'п(92) + (п + 1)^+1(92)), £15 — - ((п2 - п)/п(92) + 92п/п+1 (92)) , Б21 — 2(р2 - 52) ((п - п2)/п(р1) - ПР1/П+1(Р1^ , Б22 — 2(р2 - 52) ((п - п2)/п(Р2) - ПР21п+1 (Р2)) , Б23 — ±(^к) ("4п2 -4П + 92 (/'-1 (92) - 1п+1(92)) + (2п - 2)/п+1 (92)) ,
2п
Б24 — ±(гк) (92/'(92) - (2п + 2)/п+1 (92)), £25 — (2п - 2п2 - 9^)/п(92) + 292^+1(92),
£31 — ±(2гк)(р1 - д ) (п/„(р1) + Р1/п+1(Р1)), £32 — ±(2гк)(р2 - д ) (п/„(р2) + Р2/п+1 (Р2))
2^2 + 92
Дзз —-^ (/п-1 (92) - /п+1 (92)) + (92 + к2)/п+1(92),
2
Бз4 — -92п/п(92) - (92 + к )/п+1 (92), Бз5 — ±(гк)п/п(92),
Б41 — (к2 -Р2) (п/п(р1) + Р1 /п+1 (Р1)) + (к2 -Р2) /п(Р1), Б42 — (к2 - Р2) (п/п(р2) + Р2/п+1 (Р2)) + ^7 (к2 - Р2) /п(Р2), £43 — 0, £44 — 0, £45 — 0, £51 — 0, £52 — 0, £53 — 92, £54 — 92, £55 — ±(гк),
*
где принято обозначение д2 — к кц, 5* — —— и введены следующие безразмерные величины:
рЛ*
^ = , = ^, ^ = 1 + гщ
волны, указывающий на отсутствие физической размерности, для краткости опущен. Заметим, что
Л-4 2 а2 Д = , ^4 =
^5
= —^, к = кД, ку = куД, р, = р, Д, д.,- = д,Д, и символ
—
к
,2 = Л* С = ЛЯ
рЛ* Сг
. В дальнейшем будет использоваться еще одна определяющая постоянная
6. ВИЗУАЛЬНАЯ ЛОКАЛИЗАЦИЯ КОРНЕЙ ЧАСТОТНОГО УРАВНЕНИЯ НА КОМПЛЕКСНОЙ ПЛОСКОСТИ (Ие к, 1т к)
Переходя к анализу частотного уравнения, отметим, что в рассматриваемой задаче всего имеется пять независимых безразмерных материальных постоянных, не зависящих от частоты: , ,2, , , к_|_/кц = С /с^ . Изменение частоты будет определяться варьированием еще одного (шестого) безразмерного параметра ку.
Для отыскания решений частотного уравнения и анализа полученных данных применялась система символьных вычислений МаШетаИса 6.0. С ее помощью численно исследовалось уравнение Б = 0 в случае различных азимутальных чисел п и осуществлялся поиск безразмерных волновых чисел » = кД при заданных значениях пяти независимых безразмерных определяющих постоянных и различных значениях волнового числа чисто упругой продольной волны ку; проводилась визуальная локализация на комплексной плоскости (Ие 1т») корней частотного уравнения Б = 0 и были построены соответствующие графические образы, иллюстрирующие процесс локализации волновых чисел.
Предваряя описание процедуры визуальной локализации корней частотного уравнения, заметим, что используемые в формулах для элементов частотного детерминанта величины , д2, р!, р2, д в силу своего определения многозначны. Последнее обстоятельство существенно осложняет вычисление частотного детерминанта в системах типа МаШетаИса 6.0. В частности, для получения полного набора линий нулевого уровня значений Ие Б(Ие 1т ») и 1т Б(Ие 1т ») требуется исследование всех 25 вариантов выбора значений двузначных величин д!, д2, р!, р2, д. С целью контроля их значений произведем явное вычисление соответствующих квадратных корней. Для этого воспользуемся тем, что извлечение квадратного корня в выражении
д
= л/к2—С2
(66)
приводит к двум следующим значениям:
Ие д = ±1
'р + 7 р 2 + 4я2
1т д =
Я
Ие д'
(67)
где введены обозначения Р = (Иек)2 — (1тк)2 — (ИеС)2 + (1тС)2, Я = (Иек)(1тк) — (ИеС)(1тС).
Ясно, что приведенные только что формулы вычисления значений у/к2 — С2 пригодны только в том случае, когда Ие л/к2 — С2 = 0. Это условие не выполняется только тогда, когда разность к2 — С2 есть отрицательное вещественное число. Поэтому в том случае, когда разность С2 — к2 вещественна и положительна, получаем д = л/к2 — С2 = ±2\/С2 — к2.
На основании (1), (67) для величин
д! = А /к2 — к2,
к2 =
получим
Ие д! = ±\
д2 = к2 - к 2 ,
Ру + А/Р|2 + 4ЯУ
к2 = к1 = „2
"С2
1т д! =
Яу
Ие д1
2
2
2
С
Р|| — (Ивк)2 - (1тк)2 - к2, О" — (Ивк)(1тк);
не,2 — ± Рг+ УРГ+Ж
1т 92 —
Ог
Ив 92'
Рг — (Ивк)2 - (1тк)2 - кГ, Ог — (Ивк)(1тк).
Квадратный корень в выражении
5
— л/к2-^, С2 — да
2 — к2^2 1 + ¿/2 21 + /4,
поскольку
Ив С — ±^2ку ^ извлекается в следующем виде:
Ив 5 — ±
1 + уТТ/3 2(1 + /3) ■
1Р +у/ Р2 + 402
1т С —
к|| /2/2
2(1 + /4 )Ив С'
(68)
1т д —
О
Ив д'
Р — (Ивк)2 - (1тк)2 - (ИвС)2 + (1тС)2, О — (Ивк)(1тк) - (ИвС)(1тС).
В формуле (68) необходимо выбрать один из знаков (все равно какой, так как дальнейшие вычисления не зависят от конкретизации знака). По смыслу величина С есть комплексное волновое число плоской чисто термической (т.е. не сопровождающейся деформацией) волны. Поскольку термические волны в чистом виде (т. е. отдельно от волн деформации) в термоупругой среде не существуют, то значение С имеет указанный выше смысл только с этой поправкой.
Остается рассмотреть квадратные корни в
где
Р1 — V к2 - ТГ, Р2 — \/к2 - 72,
271, 2 ¿/2 - /1 + «1, 2 + ¿&1, 2
к2
¿/3 - 1
/1 - /4 - 4/2 + уда - /3 - 4/2)2 + 4/3(2/2 - /Ц?
ш2 — ±\/-~-
. ¿3 (2/2
01, 2 —
/2(2/2 - /2)
«1, 2
Поскольку
У2
^1 + /4 Ив 71,2
—
/4 + /301, 2 + /2 - «1, 2 + л/(/4 + /2&1, 2 + /2 - «1, 2)2 + (/1/2 - «1, 2/3 - /3 - 01, 2)
2
"к"2 \/1 + /31т 71,2 —
/2/3 - «1, 2/2 - /2 - 01, 2
22
2
У2
^1 + /4 Ив 71,2
(69)
значения Р1, Р2 вычисляются как
ИвР1 — ±
+ у/ рР2 + 4Я2
1т р1 —
О1
Ив Р1
Р1 — (Ивк)2 - (1тк)2 - (Ив71)2 + (1т71)2, 01 — (Ивк)(1тк) - (Ив71)(1т71);
Ив р2 — ±1
Р2 + л/ Р22 + 402
2
1т р2 —
О2 Ив Р2'
Р2 — (Ивк)2 - (1тк)2 - (Ив72)2 + (1т72)2, О2 — (Ивк)(1тк) - (Ив72)(!т72).
к
2
к
В формуле (69) также необходимо указать один из знаков (все равно какой, так как интересующие нас вычисления не зависят от того, какой именно знак выбран). Заметим, что постоянные 7!, 72 — скорости распространения плоских связанных ОЫШ-термоупругих волн.
Учитывая изложенное, все программы в кодах МаШетаИса 6.0 были написаны так, чтобы обеспечить различные варианты выбора знаков перед знаками радикалов в данных выше формулах для вещественных и мнимых частей величин д!, д2, р!, Р2, д. В качестве примера на рис. 1 для связанной ОЫШ-термоупругой волны азимутального порядка п = 3 приводится картина линий нулевого уровня значений Ие Б(Ие 1т ») и 1т Б(Ие 1т »), полученная в результате объединения всех 25 вариантов возможного выбора значений двузначных величин д!, д2, р!, р2, д. Значения безразмерных определяющих постоянных приняты равными =0.1, ,2 = 100.0, ,4 = 0.01, ,6 = 21, о /С* = 1.9; ку =0.1. Построение выполнено с помощью процедуры ContourPlot, PlotPoints = 450. Осуществлен учет 25 = 32 возможных вариантов выбора однозначных ветвей квадратных радикалов в выражении для Б.
В табл. 1 сгруппированы некоторые наиболее близкие к нулю безразмерные волновые числа кк термоупругой волны первого азимутального порядка, найденные с помощью процедуры FindRoot при следующих значениях безразмерных определяющих постоянных: = л/Ло/Ля) = 0.1, = сг/1 = 100.0, = «д/ДДрЛ*сг) = 0.01, = л/<гсг/Л = 1.1, сг/с^ = 1.9. В выражении для параметра следует учитывать, что Л* относится к отсчетной температуре 00.
Рис. 1
Таблица 1
к\\ Волновые числа к = кй к\\ Волновые числа к = кй
0.1 -0.011814+ 6.085313г -0.19 - 3.39887 х 10-18г 1.98107 х 10-13 + 2.81428г 1.71312 х 10-14 + 0.456983г 0.2 -0.2 + 4.98733 х 10-18г 0.38 - 1.3892 х 10-22г 1.33648 х 10-13 + 2.79345г 6.96231 х 10-15 + 0.605541г
0.3 0.57 - 1.90096 х 10-18г 0.3 - 1.0262 х 10-17г 2.20287 х 10-15 + 0.689327г 1.74377 х 10-15 + 6.6713г 0.4 0.4 - 1.5378 х 10-17г 0.76 + 4.57084 х 10-18г 2.35725 х 10-15 + 0.731949г -4.41931 х 10-18 + 2.70864г
0.5 2.967 х 10-15 + 0.741351г -3.51351 х 10-15 - 0.741351г 0.95 - 4.0535 х 10-19г 0.5 - 3.07992 х 10-17г 0.6 2.12995 х 10-15 + 0.719426г 2.67966 х 10-15 + 0.719426г 1.14 - 3.57361 х 10-18г 0.6 - 4.43641 х 10-17г
0.7 1.69503 х 10-15 + 0.663357г -1.60072 х 10-15 +0.663357г 1.33 + 7.8288 х 10-20г 0.7 - 6.05126 х 10-17г 0.8 1.21172 х 10-16 + 0.562791г 1.44482 х 10-15 + 0.562791г 0.8 - 7.91034 х 10-17г 1.52 + 1.94567 х 10-24г
Далее остановимся на процедуре визуальной локализации корней частотного уравнения. Для определенности выберем первый из корней частотного уравнения, включенных в табл. 1: » = = —0.011814314799580788+ 6.085313119171459г. Сначала построим кривые нулевого уровня значений Ие Б(Ие 1т ») и 1т Б(Ие 1т »). Линии уровня Ие Б(Ие 1т к) = 0 и 1т Б(Ие 1т ») = 0 выстраиваются совместно на одном и том же рисунке с помощью стандартной процедуры МаШетаИса 6.0 ContourPlot с заданным значением PlotPoints = 450.
На рис. 2-4 с различной степенью детализации совместно построены нулевые линии уровня Ке Б(Ке к, 1т к) = 0 и 1т Б(Ке к, 1т к) = 0 в окрестности волнового числа к' = -0.011814314799580788+6.085313119171459«. По-прежнему рассматривается термоупругая волна первого азимутального порядка при следующих значениях безразмерных определяющих постоянных: = ^Л*С1 /(АД) = 0.1, Н2 = с/ = 100.0, НА = «уя/(рЛ*сг) = 0.01, Н6 = л//Л = 1.1, с/с* = 1.9. Для безразмерной постоянной кц принято значение 0.1. В выражениях для вещественных частей дь выбраны отрицательные знаки.
6.8
-а
6.6
6.4 6.2
6.0
5.8
-0.3 -0.2 -0.1 0.0 0.1 0.2 0.3
Рис. 2
Рис. 4 дает наиболее детальную картину. Точки, в которых пересекаются указанные линии уровня, и являются корнями частотного уравнения Б = 0. После чего, собственно говоря, процесс визуальной локализации волнового числа к = к' можно считать завершенным.
Подставляя найденное волновое число к = к' в частотный определитель Б, получаем значение Б = 5.4090618305084216-10-11 -1.1823313693210316х х10-1°весьма близкое к нулю, что свидетельствует о достаточно точном нахождении корня частотного уравнения Б = 0.
Остановимся также на локализации волновых чисел связанной термоупругой волны седьмого азимутального порядка. Численно можно определить следующие два корня частотного уравнения:
_..........................................
-
"Л -
Кек
6.8 | -
Ь
6.6 - -
6.4 ------^
6.2 -
6.0 - -
5.8 Яек
-0.10 -0.05 0.00 0.05 0.10
Рис. 3
-а 6.8 1 -
6.6 - -
6.4 - -
6.2 -
6.0 -
-0.02 -0.01 0.00 0.01 0.02
Рис. 4
к = 0.35813631327764267+0.13995345934802972«, к = -0.3300000112936643 + 0.14999999521083596«.
По-прежнему задаются следующие значения безразмерных определяющих постоянных: = у/Л.сг/(ЛД) = 0.1, /2 = сг/1 = 100.0, /4 = аУЯ/(рЛ*сг) = 0.01, /е = л/^/Л = 1.1, сг/с* = 1.9. Для безразмерной постоянной кц принято значение 0.1. В выражениях для вещественных частей д1, д2 выбраны отрицательные знаки. Совместное изображение нулевых линий уровня вещественной и мнимой части частотного определителя Б в зоне локализации указанных выше корней частотного уравнения представлено на рис. 5, а. На рис. 5, б приводится более детальное изображение этой же зоны. Построения выполнены с заданным значением Р1с^Ро1^з = 450.
Для контроля точности нахождения волновых чисел необходимо вычислить значения частотного детерминанта Б для определенных численно корней частотного уравнения. Для безразмерного
волнового числа к = 0.35813631327764267 + 0.13995345934802972« имеем весьма малое значение
Б = 3.6385463758821164 * 10-23 + 8.860204458854764 * 10-22 для к = -0.3300000112936643 + 0.14999999521083596« также получается чрезвычайно малое значение Б = 2.733596645 5513133 * 10 -22 - 2.6113 56 1 644 592067 * 10- 23г.
е о
-2
-2
0
Кек
-0.2
-0.4
0.4 Кек
е о
-2
Рис. 5
С ростом азимутального числа п связанной термоупругой волны процесс локализации волновых чисел существенно усложняется. Это обстоятельство прекрасно иллюстрируется картиной (рис. 6) совместно построенных изолиний Ке Б(Ке к, 1тк) =0 и 1т Б(Ке к, 1тк) = 0 для волны достаточно высокого азимутального порядка п = 70. Указанные изолинии в выбранном масштабе представления графических данных практически сплошь заполняют комплексную плоскость (Кек, 1тк), не проникая, однако, внутрь четко обозначенной на рис. 6 «окружности». Графические построения проводились при следующих значениях безразмерных определяющих постоянных:
= у/Л,с/(ЛД) = 0.1, /2 = сг/1 = 100.0, /4 = = «уд/(рА*сг) = 0.01, /е = л/^г/Л 2.0, Рис. 6 сг /с* = 1.9; для безразмерной постоянной кц при-
нято значение 0.1; задано значение PlotPoints = 450. В выражениях для вещественных частей р1, д1, д были выбраны отрицательные знаки.
Полученные численные и графические данные позволяют также заключить, что при фиксированной частоте в волноводе, по-видимому, не могут распространяться термоупругие волны третьего типа, которым соответствуют волновые числа, имеющие модуль, меньший некоторого порогового значения Л*: |к| < Л*. Это значение Л* = 2.934, оцениваемое численно, представляет собой характерный
-2
о
Кек
4
2
4
1
2
б
а
радиус «окружности» (см. рис. 6), внутри которой не отображаются линии нулевого уровня вещественной и мнимой части частотного определителя. Геометрия самой «окружности», как показывают более точные вычисления, оказывается достаточно сложной (рис. 7-9).
2.05
2.04
2.03
2.02
2.01
2.00
2.00 2.02 2.04
2.10 2.12
2.00 2.05
2.0150 -к
е ;
2.0149 2.0148 2.0147
ш,
2.0146
2.06
2.07 2.08
Рис. 8
2.09 Ыей 2.10
.....
2.06 2.08
Рис. 7
На рис. 7 изображена структура границы «окружности», внутрь которой не проникают регулярные линии уровня Ие Б(Ие1т») = 0 и 1т Б(Ие 1т ») = 0. Графические данные построены с заданным значением PlotPoints = 750. На рис. 8 более детально изображен участок границы «окружности», приведенной на рис. 7. Наряду с регулярными линиями уровня присутствуют также линии с хаотичной геометрией. Графические данные построены с заданным значением PlotPoints = 650. Дальнейшее исследование осложняется тем, что в окрестности границы «окружности» линии уровня характеризуются достаточно хаотичной геометрией, детализация которой, как показывает рис. 9 (значение PlotPoints = 750),из-за сопутствующего «серого фона» практически невозможна.
В табл. 2 приводятся результаты численного анализа влияния безразмерной определяющей постоянной (волна азимутального порядка п = 70) на характерный радиальный размер «окружности».
Таблица 2
ЙШ
ЙИ ЙЙ1
шщщ йЩ
ШЙЩЩ
ШШ
шй
НИ
2.0145 ;
2.07555 2.07560 2.07565
Рис. 9
¡¡¡¡¡Я
2.07570
2.07575 ЫеА;
Значения Характерный Значения Характерный Значения Характерный
постоянной линейный постоянной линейный постоянной линейный
Л,6 размер Л* Л,6 размер Л* Л,6 размер Л*
2 2.9342 16000 2.738 120000 2.695
2000 2.817 32000 2.737 240000 2.679
4000 2.793 48000 2.7265 480000 2.659
6000 2.782 60000 2.718 960000 2.643
8000 2.775 96000 2.708 1920000 2.631
Следует заметить, что значения частотного определителя внутри указанной «окружности» весьма малы. Так, при = 2.0 (п = 70) согласно данным вычислений имеем:
к = 0.0950440 + 0.000084г значение Б = —2.044685 х 10- 620 + 1.34 0236 х 10-620г; к = 1.783787 + 1.782571г значение Б = 0.0 + 0.0г; к = 0.792573 + 0.693172г значение Б = 0.0 + 0.0г.
10
По этой причине локализация корней частотного уравнения на фоне чрезвычайно малых значений самого частотного детерминанта оказывается невозможной. Вопрос о построении нулевых изолиний ИеБ(Ие1т») = 0 и 1тБ(Ие1т») = 0 внутри «окружности» будет рассматриваться отдельно.
Картина изолиний, подобная рис. 6, характерна лишь для сравнительно больших азимутальных чисел п, поскольку, например, при п = 20, как свидетельствуют проведенные вычисления, она не наблюдается (рис. 10). На рис. 10 изображены изолинии вЫШ-термоупругой волны азимутального порядка п = 20. Значения безразмерных определяющих постоянных приняты равными =0.1, ,2 = 100.0, ,4 = 0.01, ,6 = 2.0, о/С* = 1.9; к = 0.1 (значение PlotPoints = 450). При переходе к гиперболической вЫП-термоупругой волне3 и малых значениях безразмерной постоянной «окружность» вытягивается вдоль оси Ие» так, как показано на рис. 11, 12 (значение PlotPoints = 450). На рис. 12 изображены изолинии гиперболической вЫП-термоупругой волны азимутального порядка п = 47; значения безразмерных определяющих постоянных приняты равными =0, = 100.0, 5 = 0.01, = 0.001, сг /с* = 1.9; к» = 0.1. Менее детальное по сравнению с рис. 11 изображение приводится на рис. 12.
В 0[
10
10
0
Кек
Рис. 10
10
е о
Кек
Рис. 11
о
Кек
Рис. 12
5
5
4
4
4
2
2
4
7. ПОСТРОЕНИЕ НУЛЕВЫХ ИЗОЛИНИЙ НА КОМПЛЕКСНОЙ ПЛОСКОСТИ (Ие к, 1т к) В ОБЛАСТИ МАЛЫХ |Б|
Как было отмечено, локализация корней частотного уравнения на фоне чрезвычайно малых значений самого частотного детерминанта сопряжена с дополнительными трудностями. На рис. 13-15 на комплексной плоскости представлены нулевые линии уровня вещественной и мнимой частей частотного детерминанта Б в случае достаточно больших азимутальных чисел п, а именно п = 70, построенные при заданных значениях определяющих постоянных = 0.1, = 100.0, = 0.01,
3Соотношения гиперболической термоупругости вЫП получаются из соотношений вЫШ-теории по следующей схеме: Но ^ 0, Нз ^ 0, Н ^ 0, Н4 ^ ж, Н0Н4 = в, Л.2 = 1 + ^2 + в2. Безразмерная постоянная в* = Н^/^/ку сама по себе в соотношения вЫП-теории не входит; Н(в* )2 ^ в2.
ЩШшШШШШу
Ш
/е = 2, сг/с* = 1.9; кц = 0.1. Использовано значение PlotPoints = 150. В выражениях для вещественных частей квадратных радикалов р1, д1, д выбраны отрицательные знаки.
Из рис. 13 видно, что при больших значениях азимутального числа п, нулевые линии уровня вещественной и мнимой части величины Б в области, близкой к началу координат на комплексной плоскости (Кек, 1тк), не отображаются. Это обстоятельство обусловлено тем, что в указанной области, ограниченной «окружностью», величина Б принимает очень малые значения порядка 10-295.
На рис. 14 представлены нулевые линии уровня вещественной и мнимой части величины Б, построенные в той же области, что и на рис. 13. Заданы значения = 0.1, /2 = 100.0, /4 = 0.01, /е = 2, сг/с* = 1.9; кц = 0.1. Вычисления с нормированным значением Б, PlotPoints = 150. При построении изолиний значения Ке Б(Ке к, 1тк) и 1т Б(Ке к, 1т к) нормировались множителем 103°°. Тем самым показано, что внутри «окружности» нулевые линии уровня на самом деле существуют и даже могут быть до определенной степени выстроены численно.
Нулевые линии уровня вещественной и мнимой части частотного детерминанта Б внутри «окружности» в заданном квадрате наиболее детально изображены на рис. 15. Видно, что указанные линии достаточно плотно заполняют рассматриваемую область и характеризуются хаотической геометрией.
В 0
-2
-2
0
Rek
Рис. 13
3FT 1
2-1
-1
-2-1
-3
-3 -2
* S
-1 0 1
Рис. 14
2 Refc 3
0.80
-se
В
0.75 0.70 0.65
0.60
/ г
У , У / // /V
У'////у / ' 'V V Л'
/ / / / О А £, / < у / t « . у ,
/ / 1 , / / у У / /*-*■ у УУУУ, У
^///////''/ Ж;
Ï-ÏA Ъ
0.60
0.75 Refc 0.80
Рис. 15
Работа выполнена при частичной финансовой поддержке РФФИ (проект 10-01-00184-a).
Библиографический список
1. Ковалев В. А., Радаев Ю. Н. Волновые задачи теории поля и термомеханика. Саратов, 2010. 328 с.
2. Ковалев В. А., Радаев Ю.Н. Элементы теории поля : вариационные симметрии и геометрические инварианты. М., 2009. 156 с.
3. Ковалев В. А., Радаев Ю. Н. Волновые задачи теории поля и термомеханика // Математическая физика и ее
приложения : материалы второй междунар. конф. (под ред. чл.-корр. РАН И. В. Воловича и проф. Ю. Н. Рада-ева). Самара, 2010. С. 165-166.
4. Duhamel J. Second Mémoire sur les Phenomenes Thermo-Mecanique // J. de L'Ecole Polytech. 1837. Vol. 15. P. 1-57; Duhamel J. Memoire sur le Calcul des Actions Moleculaires Developpees par les Changements de Temperature dans les Corps Solides // Memoirs
par Divers Savants. A l'Acad. Roy. des Sci. de l'Inst. de France. 1838. Vol. 5. P. 440-498; Neumann F. Vorlesungen über die Theorie der Elasticitat der festen Körper und des Lichtathers. Breslau, 1885.
5. Лебедев Н. Н. Температурные напряжения в теории упругости. М.; Л., 1937. 110 c.
6. Maxwell J. C. On the Dynamical Theory of Gases // Phil. Trans. Royal Soc. Lond. 1867. Vol. 157. P. 49-88.
7. Biot M.A. Thermoelasticity and irreversible thermodynamics // J. Appl. Phys. 1956. Vol. 27(3). P. 240-253.
8. Joseph D. D., Preziozi L. Heat waves // Rev. Modern Physics. 1989. Vol. 61, № 1. P. 41-73; Joseph D.D., Preziozi L. Addendum to the paper «Heat waves» // Rev. Modern Physics. 1990. Vol. 62, № 2. P. 375-391.
9. McNelly T.F., Rogers S.J., Channin D.J., Rollef-son R. J., Goubau W. M., Schmidt G. E., Krumhansl J. A., Pohl R. O. Heat pulses in NaF: Onset of second sound // Phys. Rev. 1970. Vol. 24(3). P. 100-102.
10. Jackson H. E, Walker C. T., McNelly T. F. Second sound in NaF // Phys. Rev. Letters. 1970. Vol. 25(1). P. 26-28.
11. Rogers S. J. Transport of heat and approach to second sound in some isotopically pure Alkali-Halide crystals // Phys. Rev. B. 1971. Vol. 3(4). P. 1440-1457.
12. Pohl D. W., Irniger V. Observation of second sound in NaF by means of light scattering // Phys. Rev. Letters. 1976. Vol. 36(9). P. 480-483.
13. Hardy R.J., Jaswal S.S. Velocity of second sound in NaF // Phys. Rev. B. 1971. Vol. 3(12). P. 4385-4387.
14. Narayanamurti V., Dynes R. C. Observation of second sound in Bismuth // Phys. Rev. Letters. 1972. Vol. 28. P. 1461-1464.
15. Lord H, Shulman Y. A generalized dynamical theory of thermoelasticity // J. Mech. Phys. Solid. 1967. Vol. 15. P. 299-309.
16. Cattaneo C. Sur une forme de l'equation de la chaleur eliminant le paradoxe d'une propagation instantanee // J. of Comptes-Rendus Hebdomadaires des Seances de l'Academie des Sciences. 1958. Vol. 247. P. 431-433.
17. Vernotte P. Les paradoxes de la theorie continue de l'equation de la chaleur // J. of Comptes-Rendus Hebdomadaires des Seances de l'Academie des Sciences. 1958. Vol. 246. P. 3154-3155.
18. Лыков А. В. Теория теплопроводности. М., 1967. 600 c.
19. Green A.E., Lindsay K.A. Thermoelasticity // J. Elasticity. 1972. Vol. 2. P. 1-7.
20. Green A.E., Naghdi P.M. On undamped heat waves in an elastic solid // J. Thermal Stresses. 1992. Vol. 15. P. 253-264.
21. Green A.E., Naghdi P.M. Thermoelasticity without energy dissipation // J. Elasticity. 1993. Vol. 31. P. 189208.
22. Новацкий В. Динамические задачи термоупругости. М., 1970. 256 c.
23. Maugin G.A. Towards an analytical mechanics of dissipative materials // Rend. Sem. Mat. Univ. Pol. Torino. 2000. Vol. 58, № 2. Geom., Cont. and Micros., II. P. 171-180.
24. Maugin G.A., Kalpakides V.K. The slow march towards an analytical mechanics of dissipative materials // Technische Mechanik. 2002. B. 22, H. 2. S. 98-103.
25. Maugin G.A., Kalpakides V.K. A Hamiltonian formulation for elasticity and thermoelasticity // J. Phys. A: Math. Gen. 2002. Vol. 35. P. 10775-10788.
26. Kalpakides V.K., Maugin G.A. Canonical formulation and conservation laws of thermoelasticity // Reports in Mathematical Physics. 2004. Vol. 53. P. 371-391.
27. Puri P., Jordan P.M. On the propagation of plane waves in type-III thermoelastic media // Proc. Royal Soc. Lond. A. 2004. Vol. 460. P. 3203-3221.
28. Ковалев В. А., Радаев Ю.Н. Волновые числа плоских GNIII-термоупругих волн и неравенства, обеспечивающие их нормальность // Изв. Сарат. ун-та. Нов. сер. 2010. Т. 10. Сер. Математика. Механика. Информатика, вып. 3. С. 46-53.
29. Dhaliwal R. S., Majumdar S. R., Wang J. Thermoelastic waves in an infinite solid caused by a line heat source // Intern. J. Math. & Math. Sci. 1997. Vol. 20, № 2. P. 323-334.
30. Ковалев В. А., Радаев Ю.Н. Распространение связанных гармонических GNIII-термоупругих волн в длинном цилиндрическом волноводе // Вестн. Чувашского гос. пед. ун-та им. И.Я. Яковлева. Сер. Механика предельного состояния. 2010. № 2(8), ч. 2. С. 207-255.
31. Ковалев В. А., Радаев Ю.Н., Романов А.Е. Прохождение теплового GNIII-волнового сигнала с высокой окружной гармоникой через цилиндрический волновод // Актуальные проблемы прикладной математики, информатики и механики : сб. тр. междунар. конф., посвящ. 80-летию д-ра физ.-мат. наук, проф. Д. Д. Ивлева. Воронеж, 2010. С. 173-180.
