Спектральные разложения в динамических задачах вязкоупругости Текст научной статьи по специальности «Математика»

Лычева Т.Н., Лычев С.A. Спектральные разложения в динамических задачах вязкоупругости // Вестник Пермского национального исследовательского политехнического университета. Механика. - 2016. - № 4. - С. 120-150. DOI: 10.15593/perm.mech/2016.4.08

Lycheva T.N., Lychev S.A. Spectral decompositions in dynamical viscoelastic problems. PNRPUMechanics Bulletin. 2016. No. 4. Рр. 120-150. DOI: 10.15593/perm.mech/2016.4.08

ВЕСТНИК ПНИПУ. МЕХАНИКА

№ 4,2016 PNRPU MECHANICS BULLETIN

http://vestnik.pstu.ru/mechanics/about/inf/

001 10.15593/регш.шесЬ/2016.4.08 УДК 539.3

СПЕКТРАЛЬНЫЕ РАЗЛОЖЕНИЯ В ДИНАМИЧЕСКИХ ЗАДАЧАХ ВЯЗКОУПРУГОСТИ

Т.Н. Лычева1, С.А. Лычев1,2,3

1 Институт проблем механики им. А.Ю. Ишлинского РАН, Москва, Россия

2Московский государственный технический университет им. Н.Э. Баумана, Москва, Россия

Национальный исследовательский ядерный университет «МИФИ», Москва, Россия

О СТАТЬЕ

АННОТАЦИЯ

Получена: 05 августа 2016 г. Принята: 04 декабря 2016 г. Опубликована: 30 декабря 2016 г.

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

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

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

© ПНИПУ

© Лычева Татьяна Николаевна - инженер-математик,е-таИ: taniko@mail.ru

Лычев Сергей Александрович - доктор физико-математических наук, ведущий научный сотрудник,е-таН: lychevsa@mail.ru Tatiana N. Lycheva - Applied Mathematician, e-mail: taniko@mail.ru

Sergeii A. Lychev - Doctor of Physical and Mathematical Sciences, Leading Researcher, e-mail: lychevsa@mail.ru

SPECTRAL DECOMPOSITIONS IN DYNAMICAL VISCOELASTIC PROBLEMS T.N. Lycheva1, S.A. Lychev1,2,3

11nstitute for Problems in Mechanics RAS, Moscow, Russian Federation Bauman Moscow State Technical University, Moscow, Russian Federation 3National Research Nuclear University MEPhI, Moscow, Russian Federation

ABSTRACT

Theoretical relations obtained from solutions of dynamic problems of viscoelasticity represent an effective framework for experimental identification of dynamic rheological properties of materials. For the construction of such relations, closed solutions of boundary value problems (i.e. written in the form of convergent series or integrals) are preferred, because they(unlike solutions obtained by numerical methods) allow strict error estimates. However, the construction of analytical solutions is associated with the following difficulties. 1. As usual, the hypothesis of proportionality is accepted for relaxation operators corresponding to the first and second Lamé moduli, which is equivalent to the hypothesis of a constant Poisson's ratio. This significantly reduces the generality of consideration. 2. Representation of solutions for three-dimensional problems in the form of expansions in eigenfunctions makes it necessary to taking into account the large eigenvalues which in the vast majority of problems can be found only numerically, as the roots of transcendental equations, thus, it is likely to skip closely spaced and multiple roots. 3. Constructed series converge slowly.

In this paper we suggest ways to overcome these difficulties. Solutions of initial boundary value problems are presented in the form of spectral expansions, but in contrast to the classical method of Fourier decomposition they are expanded over biorthogonal system of eigenfunctions of mutually conjugate pencils of differential operators. This pencils define generalized Sturm-Liouville problem with a polynomial spectral parameter. This eliminates the hypothesis of proportionality relaxation operators. Effective relations for the terms of spectral(in particular normalization factors) coordinate functions and asymptotic formulas for the initial approximations of eigenvalues excluding their omission in calculations are obtained. Power related ranking of elements of the spectral decomposition is proposed which allows achievingthe required accuracy of calculations on the partial sums of a low order.

© PNRPU

Введение

Математическому моделированию нестационарных колебаний линейновязкоупругих тел посвящены многочисленные исследования зарубежных [1, 2, 3, 4, 5] и отечественных авторов [6, 7, 8, 9]. Современное состояние теории вязкоупругости в рамках общей теории сред с затухающей памятью представлено в [10, 11]. Методы построения точных решений уравнений вязкоупругости приведены в работах [12, 13, 14]. Способы идентификации динамических моделей описаны в [15].

Современные методы описания феноменов вязкости и упругости твердых тел основаны на идеях, сформулированных в первой половине XIX века Томсоном, Максвеллом и др. (исторические аспекты и оригинальные работы, см. [16-19]). Строгая математическая формулировка проблемы для малых деформаций была дана Фойгтом в 1892 г. В его работе[20] сформулирован линейный закон состояния, выражающий напряжения о через деформации е и скорости деформаций е :

де

о = Е: е + ЕЕ: е, е = —. (1)

д/

ARTICLE INFO

Received: 05 August 2016 Accepted: 04 December 2016 Published: 30 December 2016

Keywords:

linear viscoelasticity, rate-type model, dynamics, closed solutions, spectral decompositions, biorthogonality, asymptotic representations for the eigenvalues.

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

е = и = 2 ( + Уит),

(2)

где и = 1 к ик - перемещения из натуральной формы; У = 1к —— - пространственный опе-

• Л

дхк

ратор Гамильтона; 1 к, к = 1,2,3 - орты декартовой координатной системы: 1 к • 1 = 5кр (точка (•) - скалярное произведение); хк - декартовы координаты; / - время; Е - тензор

модулей упругости; Е - тензор модулей вязкости.

Замечание 1. Дальнейшие построения удобно вести в терминах общих криволинейных координат дк = дк (х1,..., х3). Символами ек = 13 будем обозначать локальные

базисы криволинейной координатной системы, а символами е = 1 .дд / дх3 - взаимные к ним базисы (ек • е3 = Ь3к ); = ei • е., gil = ег • е1. Оператор Гамильтона определяется выражением: У = екдк, дк = д / ддк. Двоеточие в записи Е : е обозначает двойное внутреннее произведение. Если определены тензоры Е = Е'зк1е1 0 е3 0 ек 0 е1 и е = вптеп 0 ет, то их двойное внутреннее произведение определяется соотношением Е: е = Е3вnmglmglmei 0 е3. Далее мы придерживаемся нотации Гиббса [19]. Для векторного поля V = вур имеем следующее определение дивергенции:

V = У-V = (ек дк) • (е рУр ) = д^к + ур укрк = дк (ик ^),

где урх - символы Кристоффеля, соответствующие выбранной системе координат, g = ёе! (3.). Дивергенция тензорного поля о определяется соотношением

(V-о )п = У-( о • п ),

которое полагается справедливым для любого постоянного векторного поля п [20].

Если полагать, что закон состояния (1) задает изотропное отношение, то для определения тензоров Е, Е , в силу теоремы о представлении изотропных симметрических тензорных функций [21], потребуются четыре константы (реологических модуля А, д, А, д ):

Е = 2д1 + А10 I, ЕЕ = 2д 1 + А10 I,

где А, д - модули упругости Ламе; А, д - модули вязкости; 1 - тензорная единица, т.е. тензор четвертого ранга, отображающий любой тензор второго ранга А в себя; 1: А = А ; I - единичный тензор второго ранга, отображающий любой вектор V в себя, I • V = V . В полиадных разложениях эти операторы имеют вид I = ек 0 ек, 1 = ек 0 е 0 ек 0 ер .

Подстановка закона состояния (1) в уравнение баланса импульса

V-Т + рК-рй = 0,

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

V•(E : е)-1^-(Ё: е) + рК-ри = 0. (3)

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

Прежде чем переходить к общей задаче, поясним проблематику на примере уравнения (3). Полагаем, что тело обладает натуральной формой. Тогда с учетом соотношений Коши (2) для изотропного материала уравнение (3) преобразуется к виду ( А = V • V )

Г ~ дЛд

ц + ц— |Ди +

V дг)

А + ц + (Х

д_ дг

VV•u + рК-рй = 0. (4)

Уравнение (4) можно записать в более удобной форме, если ввести обозначения для

_ д ~ д

дифференциальных операторов по времени В = ц + ц—, ВА = А + А—, а именно:

дг дг

ВцАи + (Вц + и + рК-ри = 0. (5)

Левая часть уравнения (5) представляет собой сумму произведений (композиций) операторов по пространственным переменным и по времени и в общем случае не допускает разделения переменных в классическом смысле (по Фурье). Такое разложение возможно, если только имеется пропорциональность ц = ац, А = аА. Тогда уравнение (5) преобразуется к виду

{ д Л

1 + а—)[цАи + (А + ц)VV•й] = рй-рК. (6)

В левой части уравнения (6) оказалась композиция дифференциального оператора по времени (оператора релаксации) и пространственного оператора Ламе V :

Ти = цАи + (А + ц^-и. (7)

Такая факторизация - результат пропорциональности операторов: Вц = ц/ АВц, или,

эквивалентно, постоянства коэффициента Пуассона V. В самом деле, если по аналогии

X 1

с классической формулой V = —-- ввести оператор В, = 1/ 2ВХ (Вх + Вц ) , то при этих

2 (X + ц)

условиях он окажется постоянным: Вч, =—-т——т I, где I - тождественный оператор. То-

2 (А + ц)

гда решение начально-краевой задачи, определяемой уравнением (6), получается классическим методом Фурье. Действительно, пусть О - конечная область с кусочно-гладкой

границей 5^ ; на 5^ заданы некоторые краевые условия, например, соответствующие жесткому закреплению (т.е. и |эп = 0 ). Тогда выражение (7) совместно с краевыми условиями

определяют регулярный дифференциальный оператор V с областью определения, плотной в пространстве функций, интегрируемых с квадратом по области ^ - в пространстве Ь2 . Этот оператор самосопряженный и отрицательно определен, следовательно, его собственные значения отрицательные и действительные [24]. Пусть & - одно из собственных значений, а и0 - соответствующая собственная функция, т.е.

Vu0 = &и0, &< 0.

Полагаем К = 0, а начальные условия такими, что и = и0е^ - решение начально-краевой задачи, определяемой уравнением (6) и соответствующими краевыми условиями. При этом параметры & и должны быть связаны алгебраическим соотношением

(1 + а£)&-р£2 = 0. Разрешая это уравнение относительно приходим к выражению

с & ^ &2 2 &

£ = —аа± Л—а2 +-.

2р \4р

Если а <

V

4Д , то решение определяет затухающие колебания, в частности

1&1

• & I & 2 & и = и 0е ^тю?, у =--а, /ю = . —- а +—.

2р V 4р2 р

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

В общем случае решение представляется разложением по счетному множеству собственных функций и оператора Ламе V :

ад

и = &к ()ик (х) ,

к=1

где фк () - координатные функции, которые находятся из последовательности независимых задач Коши

п (л 52 _ . 0 5 . 1

&к I1 + аа у Фк-р5Г7 Фк = Ф к, Фк ^=0 = Ф°, Фк ^ = Фк.

Здесь Фк =<-рК, и к >, ф0 =< У0, и к >, фк =< У1, и к >; У0, У1 - начальные смещения и скорости, (•, •> - скалярное произведение в Ь2 (. Каждая такая задача определяет независимый ангармонический осциллятор, а отклик всего вязкоупругого тела на динамическое воздействие, задаваемое ненулевыми начальными данными и динамическим полем массовых сил К , представляется суперпозицией независимых ангармонических колебаний.

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

Перейдем теперь к общему случаю. Для изотропной вязкоупругой трехмерной среды уравнения состояния темпового типа могут быть записаны в виде [23]

Л с = 2ВЦ г + Шх Не, (8)

где Л, Вц, В, - дифференциальные операторы,

п т т ^

к=0 ш к=0 ш к=0 ш

ак, , Хк - константы, определяющие реологические свойства материала ( ап = 1). Используя очевидные преобразования, выразим в явной форме напряжения с:

с = 2Ме + Юге,

где М , С - интегродифференциальные операторы

м = л, с=л -1в,,

и деформацию е:

е=2м(с-Юго), М=ВГЧ С=В,(2ВЦ-3В,)-1.

Предполагая существование натуральной конфигурации, определим поле деформаций е кинематическими соотношениями Коши:

е = 2(и + (Уи)т), е = 8^ ®еу, е/ = 2( + дпикЕ^п),

где и = ие/ - поле перемещений относительно указанной натуральной конфигурации. Тогда напряжения с могут быть определены с помощью следующего оператора:

Е : и ^ с, Е (и) = 2МёеГи + 1СУ- и. Оператор Е далее будем называть оператором напряжений. 1. Формулировка начально-краевой задачи

Пусть В - исследуемое вязкоупругое деформируемое тело, отклик которого требуется определить. Как принято в теориях о малых деформациях, тело будем отождествлять с множеством геометрических мест его материальных точек, полагая, что деформации тела приводят к пренебрежимо малым изменениям их мест. Предполагая использование теоремы Стокса, условимся, что граница тела Г = дВ регулярная в смысле Келлога [25].

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

отображения замкнутого интервала). Регулярная граница может содержать конечное число ребер и конических точек [25, 26]. Это ограничение не слишком обременительное и соответствует интуитивным представлениям о «сплошном теле».

Кроме того, полагаем, что отсчетная форма свободна от напряжений. Тогда из баланса импульса

V- о — + рК = 0,

где р - плотность массы; К - плотность массовых сил, и закона состояния (8) вытекает уравнение движения

,^2и + £VV• и — ри + рК = 0, (9)

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

В^ 2и + ( + В})VV• и — рЛи + рЛК = 0. (10)

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

п - 2(и) |г = 2 - и |г, (11)

где Е - тензорное поле второго ранга, заданное на Г и определяющее упругие характеристики закрепления тела на его границе; п - внешняя единичная нормаль к поверхности Г . Далее рассматривается несколько упрощенная постановка начально-краевой задачи: будем полагать, что Е определяется постоянными матрицами Еп и Е12 такими, что

п-2(и)-Е^ = 0, и -= 0, (12)

где Гг - попарно непересекающиеся части границы (I = 1,...,Ь), образующие в совокупности ее разбиение:

Г=и г, ((* рН Г П Г р

I=1

а для матриц Еп, Е12 выполняются условия

гапкЕг 1 + гаикЕг 2 = 3.

Эти условия обеспечивают регулярность краевой задачи [27]. Начальные условия могут быть записаны следующим образом:

| У ^и | у Эки | у )

и \г=0 = У0, =0 = У1,"', |г=0 = Ук, (13)

Э от

где к = тах(т,п + 2) ; У0,У^...,Ук - начальные смещения, скорости, ускорения и т.д.

Уравнения (10), краевые условия (11) и начальные условия (13) определяют формулировку исследуемой начально-краевой задачи. Для дальнейших построений систему уравнений (10) удобно представить в виде

Му 2и + М2 УУ- и + М3и = $. (14)

Здесь М, Мг, М - дифференциальные операторы,

т (к т (к п (к+2

М, М=Е(^к+^к )—г, мз=-рТак-п+г,

к=0 ( к=0 ( к=0 (

а $ - заданная вектор-функция, которая вычисляется следующим образом:

п (к

$ = -рЕак(-т К.

к=0

Уравнения (14) совместно с соответствующими начальными и краевыми условиями определяют начально-краевую задачу, которая может быть представлена как задача Коши с операторными коэффициентами:

г дк

УЦк дТ и = $,

где г = тах (т, п + 2), а Ц - дифференциальные операторы по пространственным переменным, подобные оператору Ламе (7):

Ци = ДкУ2и + (дк к )УУ и - ра^и. (15)

Здесь {д к } {А, к }к=0 - реологические модули, упорядоченные на одном и том же множестве индексов (0,..., г) с добавлением нулей там, где это необходимо:

^ к =

Vk, к < m - ГАk, k < m

Л

0, k > m

А k =

v 0, к > m

k

0, k (0 л ^и'

Уравнения движения упрощаются, если операторы С и M пропорциональны:

£ = kM, k = const. (16)

Это условие классифицируется в [28, с. 58] как условие разделения переменных. В большинстве работ, рассматривающих замкнутые решения для вязкоупругих тел, принимается (16).

При условии (16) уравнение (9) может быть представлено в виде

MVu -pii + pK = 0, (17)

где V - классический оператор Ламе (7), область определения которого задается краевыми условиями (12). Заметим также, что в квазистатическом приближении (без учета сил инерции) уравнения (17) оказываются эквивалентными уравнениям упругости

Vu + pK = 0,

где K = M-1K, т.е. u = -pV-1K.

2. Спектральное представление решения

Пусть }2 (В) - гильбертово пространство комплекснозначных векторов-функций, определенных в области В со скалярным произведением:

Уи, V е }2 (В) <ц, у> = ]и • уйУ,

в

где черта означает комплексное сопряжение. Будем отыскивать решения и (х, /) краевой задачи (12)-(14) в форме разложения по функциям пространственных переменных х, образующих базис {ек (хв пространстве } (В), т.е.

u (xt )=&( К (x )•

(18)

k=0

Функции ф(/), выступающие в роли коэффициентов разложения, параметрически

зависящих от I, будем называть координатными функциями. Далее для сокращения письма аргументы базисных и координатных функций опускаем. Идеология используемого представления подобна методу Фурье, однако, в отличие от него, непосредственная подстановка (18) в уравнения (12)-(14) не позволяет разделить переменные.

Разделение переменных можно осуществить, если преобразовать эти уравнения к системе уравнений первого порядка по переменной I, а именно:

V v 0 + V v 0 + Р2 v 1 + Рз v 2 +... + Vm v m_i = f,

vi " v0 = 0 v2 " v 1 = 0,

v - v = 0

v m-1 v m-2

(19)

где функция v0 совпадает с искомым решением (т.е. u = v 0). При этом краевые условия могут быть записаны в виде

((n • def vk + В,n V- vk ) • H,1 |Гк = ЛНka • vk Irk, k = 0,1,., m -1. (20)

Начальные данные задаются соотношениями

v I = V v I = V

v0 It =0 T0> V1 It=0 M>

v I = V

T m-1 It=0 ' m-1'

(21)

Для разделения переменных (отделения переменной I) рассмотрим обобщенную задачу Штурма-Лиувилля, которая определяется системой уравнений

Г Р0е0 + ^ + е, + ^ +... + уТ^ = 0,

e1 -ve0 = 0, e2 - ve1 = 0,

e , -ve 2 = 0

m-1 m-2

и краевыми условиями

(( п• ^е + 4пУ-е)) |гк = 4 • ег ^, / = 0,1,т -1,

4 = ^к Vк, ^=1^ V*, 4=!^

V

к=0

к=0

к=0

Здесь v - комплекснозначный спектральный параметр. Уравнения удобно представить в матричной форме

Г V ' 0 0 0 . . 0 Л Г V V; Vз . V ' т-1 VVЛ

0 I 0 . . 0 -I 0 0 0 0

0 0 I . . 0 + V 0 -I 0 0 0

V 0 0 0 . 1У V 0 0 0 .. -I 0 У

из этих уравнений исключить , ет-1 , выразив их через

(е \

, ^

0.

(24)

V ет-1 У

ференциальное уравнение, в которое параметр v входит в различных целых степенях. Дифференциальный оператор, порождаемый этим уравнением в области из I2 (ВВ), определяемой краевыми условиями (23), будем называть полиномиальным пучком Vv,

4= V ^ V2V V. (25)

Заметим, что спектральный параметр V входит в выражение ^ полиномиально.

В рамках рассматриваемой задачи о движении вязкоупругого тела выражение для ^ упрощается. Поскольку все V, г = 0,...,т, имеют вид (15), то уравнение = 0 может быть записано в виде

2Х ^)У2е + Z2 ^)УУ • е-р23 (V) е = 0 или, полагая Z1 (V) ^ 0, в виде

У2е + К (V) УУ- е - Р (V) е = 0, (26)

Z (V) Z (V)

где Z, (V), Z? (V), Z3 (V) - многочлены от v, а К (V) = Р (V) = 3, ч .

^ (V) 4 ' Z1 (V)

Замечание 3. Изучение спектра пучка операторов теории упругости началось в работах Эжен и Франсуа Коссера. В 1898-1901 гг. они опубликовали серию статей, в которых исследовались собственные числа и собственные векторы пучка операторов У2 + vУёгу (v -спектральный параметр) и рассматривались их приложения к решению основных задач теории упругости. Работы Коссера по спектру пучка операторов теории упругости были продолжены в 1967 г. С.Г. Михлиным и В. Г. Мазья (более подробно см. [24]).

Для определения дифференциальных операторов, порождаемых расширенной системой (22), введем расширенное пространство I2 (В), элементы которого т компонентные

кортежи векторов (и0,их,...,ит-1) со скалярным произведением

т-1

(((, и1,..., ии-1), ( VI,..., Vт-1 ))) = • V.

Вк=0

Дифференциальное выражение для Н (24) и краевые условия (20) определяют в Е2 (В) дифференциальный оператор, который будем обозначать тем же символом Н . Оператор НД, сопряженный оператору Н , может быть найден из соотношения

((НДи о, u1,., и т-1 ) (Уо, y1,., У т-1 ))) = (((и о, U1,., и т-1 ), (У 0, ^^^^^^^ У т-1 )))

Используя теорему о дивергенции, приходим к следующему виду для сопряженного оператора Н*:

Н* =

д

V * ' 0 0 0

0 I 0

0 0 I

0 0 0

+ д

V -I 0 V* 0 -I

V* 0

V V 0

V т

0

0

*

е0

*

е1

Н* * е2

У * V ет-1 У

определенного в области, задаваемой краевыми условиями Вк*, т.е.

п • Вк*ек |г = 0, к = 0,1,..., т -1.

Здесь V* - дифференциальное выражение, сопряженное с Vk, а В*, - сопряженные операторы краевых условий. Последние находятся из условия обращения в ноль интегралов по границе тела, получаемых при переходе к сопряженному выражению [29, 30, 31]. Сопряженному оператору соответствует сопряженный полиномиальный пучок

Р0*е* + дЯ*е* + д2 "Р2*е* +... + дт ^е* = 0. (27)

Замечание 4. Вопрос о полноте и базисности собственных функций (а в случае кратных собственных значений - также и присоединенных) несамосопряженных операторов является довольно тонким. В общем случае для несамосопряженных операторов не имеет места утверждение, аналогичное спектральной теореме для самосопряженных операторов [27]. Однако для специальных классов несамосопряженных операторов, в частности для полиномиальных пучков регулярных дифференциальных операторов, рассматриваемых в настоящей работе, полнота является следствием асимптотического поведения их резольвент [32] (см. также [29]).

Пусть X ид - два различных комплексных собственных числа и (е0, е1,., ет-1 ),-

отвечающие им собственные функции прямого Нх и сопряженного НД пучков соответственно. Имеют место соотношения биортогональности:

К V V V V V V е ^ (е* ^

' ' 1 '2 ' 3 ' т-1 ' т С0 с0

-I 0 0 ... 0 0 0 -I 0 ... 0 0

0 0 0

У Vет-1 У

V ет-1 у

V -I 0

V* 0 -I

0

0

0 ¥ е0 0 0

0 IУ1

V * 0 0 .

0 I 0 .

0 0 I .

IV ет-1 У V 0

00

(29)

I е

V Vет-1 у

Замечание 5. Соотношения биортогональности не исчерпываются формулами (28), (29). Общий вид таких соотношений может быть получен следующим образом. Пусть ех,

еД - собственные функции пучков (25) и (27), отвечающие различным собственным

значениям X , д . Тогда имеют место равенства

{Цех, еД) + А(7Х, еД) + А2^ех, еД) + ... + Ат [Ттех, еД) = 0, (30)

(ея, + ^ех, + Д2^, Г^е^ + ... + Дт(ех, ^е^) = 0. (31)

Преобразуем второе уравнение (31)

, е*) + д( Рlex, <) + Д2( Р2ех, е*) + ... + дт( Ттех, еД) = 0. (32)

Уравнения (30) и (32) можно рассматривать как одно и то же алгебраическое уравнение относительно параметра х (или д ). В силу теоремы Виета имеем

('^р , ( V1 + V2 + . + Vт ) = (^, , ^ ), (Р0ev, , eV^Sk (V1, V2 , . , Vт ) = (Рkev, , < ) ,

(^р , 4) V V 2 "^т = (^4 , <) .

Здесь р,д = 1,.,т, v1,V2,...,Vт - корни алгебраического уравнения (30), к числу которых принадлежат х , Д . Линейные комбинации этих равенств определяют различные соотношения биортогональности (см. также [29], теорема 4).

Перейдем теперь к построению решения краевой задачи (19)-(21). Решение и = и (х, г) будем отыскивать в виде разложения

ад

и (х, г ) = 5>( )Ек (х ),

к=0

где Ек = (е0, ек, ..., ект_1) - решение задачи (24), соответствующее значению спектрального

параметра V = Vк (для ограниченного тела В такие решения образуют счетное множество). Подстановка этого представления в уравнения (19) и в начальные данные (21) приводит к следующим соотношениям:

( ц 0 0 0 . . 0 1 Г ц ц ц . ц -рл т (Г1

ад 0 I 0 . . 0 I 0 0 0 0 0

£ к=0 Фк 0 0 I . . 0 Ек+ф к 0 I 0 0 0 Ек> = 0

V 0 0 0 . I) V 0 0 0 .. I 0) V0 ,

£фЕк |г=о =((,У„...,).

к=0

Поскольку функции Ек удовлетворяют уравнению (24), то

(ц о о ... о 1

0 I о 0 0 I

чо о о ... I,

Ек =^к

ц т

0 0

0 0

_1 0 0 0 _1 о

о о о ... I о

Е,

Таким образом, имеем

£

к=0

г ц ц2 ц ц т (Г 1

I 0 о . 0 0 0

Vк Фк) 0 I 0 0 0 Ек . = 0

V 0 0 0 . I 0) V0 ,

(33)

Действуя операторами проектирования на левую и правую части уравнения (33), приходим к последовательности независимых задач Коши относительно координатных функций фк,

Фк ^кФк = Фк, Фк 1|=0 =ф0, Фк =(((,0,..,0)Х)) = <Г,е*>,

(

(ио, и!,., и ),

ц0* 0 0 ... о

0 I о ... о

0 0 I ... о

0 0 0 ... I

= ( V0U0 , > + <U1> ,1 ) + ^ + ,m-1 >•

В случае постоянного коэффициента Пуассона разложение удобно осуществлять по собственным функциям оператора Ламе V . Оператор V как оператор, действующий в L2 (B), самосопряжен и отрицательно определен, следовательно, система его собственных функций образует ортогональный базис в L2 (B) [24]. Собственные функции находятся из задачи Штурма-Лиувилля

Vek = v4, n■2(ek)■ £Jrk = 0, ek ■ S2lrt = 0, <ek,ek> = 1. (З4)

Координатные функции находятся как решения последовательности независимых интегро-дифференциальных задач Коши

dp

^Ф^2-РФ* = fk, fk =-<Uk,pK>, dtpU = Ф^, Ф°,р =<Uk,VP>• (35)

Подействуем оператором А на левую и правую части уравнения (35). Приходим к неоднородному обыкновенному дифференциальному уравнению

Вфкчк -рАфк = f, I = Afk. (36)

Характеристическое уравнение, соответствующее уравнению (36), имеет вид

Zk (k) = v2B(k)-pk2A(k) = 0, degZk (k) = r = max(m,n + 2). (37)

Пусть .p - корни уравнения (37). Будем полагать, что они все различны. Тогда решение уравнения (36) может быть представлено соотношением

Фk = I П ( р . {¿(-1)'SP-' (^^У )+F'(T)drl,

р=1 Ц ,=1 (.k -.k )l '= 0 J

, * P

где sp (kj,...,kk) - симметричный многочлен порядкаp от кортежа (kx,...,kk),

Sp (V- kk)= I Пk„, 1<11 < h <■■■< jp <k ,=1

а sp (,■.., kn) - симметрический многочлен от кортежа Aj,„.,k ,■..,An, в котором опуЩен элемент kp , т.е. Sp (( ■ , kn ) = Sk (k1, ■ , kp-1, kp ■ , kn ) •

Заметим, что выражение для фк понимается в классическом смысле только при условии и-кратной дифференцируемости функций fk. В противном случае все соотношения

следует понимать в обобщенном смысле и использовать методы дифференцирования и интегрирования обобщенных функций. При этом выражения для координатных функций могут содержать дельта-функцию и ее обобщенные производные. Однако, если порядок дифференциального оператора А не превосходит порядка операторов Вц, Bk, т.е. n < m ,

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

3. Колебания вязкоупругого цилиндра

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

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

В настоящей статье рассматриваются конечные круговые цилиндры. Динамическая задача для упругого конечного цилиндра имеет более чем вековую историю: впервые подобная начально-краевая задача исследовалась для бесконечного цилиндра Похгаммером [33] и, независимо, Кри [34], а для конечного цилиндра и специальных классов краевых условий - Файлоном [35] (см. также [36, 37]).

Далее символы г, ф, г будут обозначать цилиндрические координаты, причем ось г совпадает с осью тела, а плоскость г = 0 содержит его нижнее основание. Символами и, V, w обозначим проекции собственных функций на нормированные векторы локального базиса ег, еф, ег.

Пусть тело представляет собой конечный цилиндр, радиусом Я и высотой Н, основания которого закреплены на гладких жестких плитах, а боковая поверхность свободна от нагрузок, т.е. Г = Г1 иГ2 иГ3, где Г2Г1 - основания, Г3 - боковая поверхность. Тогда

Основными элементами решения (18) являются собственные функции, определенные из уравнения (26) и соответствующих краевых условий. По своей структуре уравнение

(34) в предположении, что для решения уравнения (26) следует в полученных соотношениях осуществить указанную подстановку.

Задача Штурма-Лиувилля (34) формулируется следующим образом. Система дифференциальных уравнений (26) принимает вид

w 1Гж = 0 w 1г2 = 0 ° гг 1г3 = 0

(26) подобно более простому уравнению (34) с точностью до подстановки К (у) ^ к, этой причине все дальнейшие построения производятся для уравнения

Ли

г

и

,2

4- у,ф+к [ -(ги), г+-w

г ^ г г

1

1

,г

Лу

г

V

2

Краевые условия на цилиндрической (боковой) поверхности Г3 определяются соотношениями

2иг + (к-1)1 иг + 1уф + и + м,г I \г=К = О,

1 1 ,

~игГ + ^--У \г = Д = 0,

Г г

+ и,г \г =Д = 0

(39)

а краевые условия на основаниях Г1, Г2 задаются равенствами

+ и 2 \г=0,Н = 0

V,г + "\*=0,Н = 0 г

™ \г=0,Н = 0

(40)

Здесь А - оператор Лапласа в цилиндрических координатах,

А/ = / + - / +1 и + и .

Г > г

Запятая используется для обозначения частных производных по указанным после нее переменным.

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

1° Осевая симметрия. Полагаем V = 0 . Система трех дифференциальных уравнений (39) сокращается до двух уравнений:

Аи - и + к | — (ги )Г г V Г

Ам + к | — (ги) г + м,

V г , "

■ м.

= -V 2и,

2

(38')

= -V ж

Краевые условия на цилиндрической поверхности Гз записываются в виде

2и,Г +(к - ^и,Г + и + ™ 2 ^ \г =Д = 0,

(39')

+ и, 2 \г=Д = 0.

а краевые условия на основаниях Г1, Г2 задаются двумя равенствами

Кг + и,2 \г=0,Н = 0

к \г=0,Н = 0

Здесь А - оператор Лапласа в полярных координатах:

(40')

А/- = / +1 / + и .

и и ,ГГ «у ,Г ,22

Г

2° Осевая симметрия и независимость от вертикальной координаты 2 (это условие соответствует плоской деформации). Полагаем V = 0, w = 0. В системе (39) остается лишь одно нетривиальное уравнение

(k+i)í1 (u ),r - u J

V r r2 J

= -v 2u.

(38'')

Краевые условия на основаниях Г1, Г2 удовлетворяются автоматически, а краевое условие на боковой поверхности Г3 задается одним равенством

(k + l)ur +(k-1)

u

Ir=R

= 0.

(39'')

Известны различные формы представления решений системы уравнений (39) [14, 36, 38]. Для дальнейшего анализа удобна следующая форма, представляющая ограниченные на оси r = 0 и периодические по ф решения [30]:

U = ( u, v, w) = Yc Cc,

где Cc =(Cic

, C 2 , C

) - арифметические векторы постоянных интегрирования, точка (.) обозначает вектор-матричное произведение,

т„ ( muz muz . muz Л ,• / • ч

Y = diag cos-,cos-,sin-I .diag (sin nф,cos /Kp,sin жр).У ,

V H H H J

T„ .. ( muz muz . muz Л ,• / • \

Y = diag cos-,cos-,sin-I .diag (sin ^,coswp,sin жр).У ,

V H H H J

f

Уs =

-aJn+i (ra)+ nJn (ra) -nJn (rb)

-mJn+i (rb)

Л

rJn (ra)

-mJn (ra)

bJn+i (rb)--Jn (rb ) mJn+l (rb) r

0 bJn(rb)

Уc =

-aJn+i (ra) + -Jn (ra) -Jn (rb )

r r

-mJn+i (rb)

J Л

--Jn (ra)

r

-mJn (ra)

bJn+l (rb)- -Jn (rb) -mJn+i (rb) r

0 bJn{rb)

(4i)

a =

2 2

v2 m2n2 / 2 m"n b = J v2 -

k +1 H2

H

2

(42)

Здесь diаg (...) - диагональная матрица; Зп (г) - функции Бесселя первого рода порядка п,

] (2) = !

(-1) '

к=0

к !Г (к + п +1)) 2

2к+п

1 ™

— 1с08 (пт- 281Пт) йТ,

которые для Ы > п2 -1/4 имеют асимптотические представления [39]:

V К

2 4 / е

0081 2 - —-- 1 + е^Ю(\ 2 \-1)

Аналогично строятся решения сокращенных систем (38') и (39''). Решения уравнений (38'), удовлетворяющие условиям (40'), могут быть представлены в виде У8ут.С8ут, С8ут = (ут, С2ут ), где

Уут = ^ I 0088 ^ , I 3ут ,

[-а] (га) -т] (гЬ)

(44)

38 т =

8ут |-т]0 (га) Ь]0 (гЬ) а решения уравнений (38'') имеют вид

Уутр8 = ]1 ( Га) .

4. Нахождение собственных значений

у

(45)

Для нахождения нетривиальных решений задач Штурма-Лиувилля (38)-(40), (38')-(40') и (38'')-(39'') следует подставить решения (41), (44), (45) в соответствующие условия на боковой поверхности цилиндра, т.е. в соотношения (38), (38'), (38''). Получаемые при этом однородные системы линейных уравнений относительно элементов векторов С имеют нетривиальные решения при обращении в нуль определителей их коэффициентов. Соответствующие трансцендентные уравнения относительно v будем называть спектральными уравнениями. Известную сложность представляет вычисление достаточного количества их корней. Ниже предлагается методика, позволяющая осуществить этот поиск эффективно в вычислительном смысле.

Начнем с самого простого случая - спектрального уравнения, порождаемого решением (45) и краевым условием (39''). В результате подстановки (45) в (39'') приходим к трансцендентному уравнению:

к (к +1)]0 (к )-2] (к) = 0, к = аД

(46)

Поскольку для действительных аргументов к , согласно (43), функции Бесселя допускают асимптотическое представление

]П (кУ

с

-еоэ

кк

v

пк К

к----

2 4,

то и для больших к уравнению (46) соответствует уравнение

к

(к+1

00Б

к-К I-2э1п

V 4 у

к-К 1 = 0.

(47)

(48)

Следовательно, для больших значений х второй член в левой части уравнения (48) -малая величина в сравнении с первым, и им можно пренебречь. Тогда приближенные значения корней уравнения (46) можно определить как нули '0 р функции J0 (к)•(р е М). Значения '0 табулированы, в том числе в известных системах символьных вычислений,

и могут быть легко найдены. Заметим также, что из (48) можно получить приближения. С помощью уравнения (47) можно построить более точные приближения. Для этого воспользуемся кусочно-линейной интерполяцией для тригонометрических функций:

sina «(-1)

a— l\f % 2

1

a

% I %

2

cosa

(—1)

1-

a

% I %

(49)

где [...] - целая часть, а {...} - дробная часть числа. Это уравнение может быть приближено уравнением

к

(k+о

i—2'

%

%

к--I mod %

v 4 J

г

%

к--I mod %,

v 4 J

корни которого образуют последовательность {кР° }^=0 с членами

(о)

кР =

8 (k +1)

%(k +1)(4 p + 3) + ^%(k + l)(%(k +1))4 p + 3)2 + 64 p —16) + 64 — 8

Используя эти значения как начальные приближения и раскладывая в их окрестностях левую часть уравнений (46) в ряд Тейлора, получаем следующее приближение:

Г (к° ) Г (К ) = Ж 2_{к + ^ ^ 2Ж/ 2 )_.

К 2 (к +1)-2)/ ( к) + к (1 - к)/0 (к)

Более точные приближения можно получить, используя оператор Г рекурсивно:

г \

ку =

(F ° F°...F)( к Р0))

(50)

Этим приближениям соответствуют собственные значения:

>/кТТ

v p = lim к p

R

Практические вычисления показывают, что «разумная точность» достигается при ^ = 3 . В качестве примера вычислим собственные значения при Н = % , Я = 1 методом последовательного шагового поиска и с помощью приближений (50). На рис. 1 показано распределение собственных значений в зависимости от их номера, а на рис. 2 - логарифмы абсолютных значений погрешностей асимптотических приближений различных порядков.

Спектральное уравнение для осесимметричных колебаний получается при подстановке решений (44) в краевые условия (39'). В результате приходим к трансцендентному уравнению:

J

1 (к)[ J0 (к)R (b2H2 - %2m2 )(a2H2 (k +1) + %2m2 (k -1)) - 2J1 (к) aH2 (b2H2 + m2))

+ 4%2т2Н2аЬЮх (к)/0 (К) = 0. Здесь к = аЯ , К = ЬЯ . Замечая, что с учетом соотношений (42) справедливы равенства

Log •

Рис. 1. Распределение собственных значений Fig. 1. The distribution of eigenvalues

- I:

_____ r 10 -•-1

5 0 ......... -........•......... ►-------*-------1 ........»........■

A .............♦........... I T..... .......................^ .....♦....... *

\

í i -----... "■*--_____

Рис. 2. Погрешности асимптотических

приближений различных порядков: • 5 = 1; ■ 5 = 2; ♦ 5 = 3; А 5 = 4; p - номер

собственного значения Fig. 2. Errors of asymptotic approximations of several orders

a2H2 (k +1) + n2m2 (k -1) = b2H2 - n2m2 = v2m2 - 2m2n2, преобразуем уравнение (51) к более простому виду:

J1 (л)[J0 (л)R(b2 - С2)2 -2J1 (л)a(b2 + С2)] + 4С2abRJ1 (л) J0 (л) = 0, (52) где СН = mn.

Для построения асимптотических приближений собственных значений рассмотрим два случая

1°. Пусть v > Vk +1 С , тогда Im a = Im b = 0, а асимптотические выражения для функций Бесселя определяются тригонометрическими функциями (47), что позволяет записать асимптотическое представление уравнения (52) в виде

sin(л-n /4)[cos(л-n/4)R(b2-С2)2 -2sin(л- n/4)a(b2 + С2)] + + 4С2abRsin ( л) cos ( л - n /4) = 0.

(53)

Повторяя рассуждения, приведенные выше, получим уравнения для грубых приближений корней уравнения (51)

J1 (Л(0)) J0 (л(0>) = о,

(54)

>(0)_

n

1

p - 41, Л Р0) = Á,

т.е. ^p = 71,p ражения для собственных значений

n

p + 4 I. Этим приближениям соответствуют вы-

,(0) =

(

= \ (к +1)

(pn-n / 4)2 + m2n2 | v(0)= I(pn + n / 4)2 + m2n2

R2

H2

y

R2

H2

Более точные приближения могут быть получены следующим образом. Полагая а , Ь большими величинами, представляя левую часть уравнения (53) в форме разложения по их степеням и удерживая четвертую и третью степени, приходим к уравнению

sin (bR - n / 4)[cos(aR -n /4)Rb4 - 2sin(aR -n /4)ab2] = 0. (55)

Это уравнение эквивалентно двум уравнениям

sin(bR-n/4) = 0, cos(aR-n/4)Rb2 -2sin(aR-n/4)a = 0. (56)

Для того чтобы сделать эти уравнения независимыми, следует выразить b через a, т.е.

[mn J

\ + l)a2 + l)a2, и подставить результат во второе уравнение. Имеем

cos (aR-n/4)R (k + l)a-2sin (aR-n/4) = 0.

Первое уравнение (56) дает те же приближения, что и (54). Для нахождения приближений из второго уравнения воспользуемся кусочно-линейной интерполяцией (49). Полагая a = aR- n /4 = x + pn, p e N, приходим к соотношению

'l - — V + l)(x + pn + n/4)-— = 0. V n ) n

Положительные корни этого уравнения имеют вид

n(-4kp + s¡(k +1)2(4p + 3)2 + 16n(k +1)(4p -1) + 64n2 + k -4p - 8n +1

xp =■

>(k +1)

Соответственно, приближения для собственных значений могут быть вычислены по формулам

v(0)=£ p R

p 3 8nH2 -J 64n2 H4 +16nH2 (k +1) (4 p -1) + (k +1)2(4 p + 3)2

2 + 8 8 (k +1)

2°. Пусть теперь

£ < v < V k +.

При таких значениях v 1тЬ = 0, Яе а = 0 и функции Бесселя согласно (43) могут быть приближены выражениями

J0 (Л )

а

-cos

nx

V

x-fj, J1(x):

2 . -sin

nx

V

_ n

x--

4 ,

При подстановке этих соотношений уравнение (51) приводится к виду

. ( _ n

sin I x--

V 4)

2H (b2 H2 +n2 m2)

n2km2 - b2H2

k +1

R (b2 H2-n2 m2 )2

. ^ 2„ f~ nJ n2km2 -b2H2

4n bHm RcosI x--IJ-= 0.

V 4 )V k +1

Используя кусочно-линейную интерполяцию (49) для тригонометрических функций, полагая а = х-ж/4 = х + рж, р е N и заменяя значениеЬ вне аргументов

тригонометрических функций приближенным выражением (рж + ж/4)/Я, приходим

к уравнению относительно х :

4НЯ((2Нр + Н)2 + 4т2Я2))^ + ж ((2Нр + Н)2 -4т2Я2)2

ътт г,п „ ч ¡4кт2Я2 -(2Нр + Н)2 л -ж3 Нт2 (2 р + 1))ж-2 х)-(к +\)яр ) = 0

Его решения могут быть записаны следующим образом:

8ж2 т 2дЯ3 7

8R3

хр =

4HRZ(4m2R2 + Q2) + n(Q2 - 4m2R2)2 + 16nm2QR3Z '

где Qp = H(2p +1); Zp деляются по формуле

4km2 - Q2/ R2

к +1

. Приближения для собственных значений опре-

(0)

mH

n

2 Л

+ pn + n /4 л R

Эффективность полученных формул проиллюстрирована на рис. 3, на котором приведены графики левых частей трансцендентных (спектральных) уравнений, стрелками показаны приближения, полученные из асимптотических формул. Вычисления производились для Я = 1, Н = ж , к = 5 • Левому графику соответствует т = 1, правому - т = 10 .

Рис. 3. Графики нормированных левых частей спектральных уравнений и асимптотические приближения для их нулей (показаны стрелками) Fig. 3. Graphs of the normalized LHS of spectral equations and asymptotic approximations for their zeros (indicated by arrows)

Имея начальные приближения, более точные значения корней можно определить численно из уравнения (51).

Спектральное уравнение для трехмерного случая получается в результате подстановки решений (41) в (38) и может быть представлено в виде

Ъ + Ъ + Ъ = о

Ъ = 4д ( (л)- л1п+1 (л)){[4п (п2 -1)- 3,?2 ] Jn+1 (л?(л?)--2(п2 - 1)л/п+1 (л?)2 + л?(2 -3(п -1)п)Jn(л)2},

Ъ2 = (( -л2 )) (л) + пМп (л))х x{Jn (л) Jn (л?)(-(к-1)?(2(п- 1)п-%2) + л2 ((Л + 1)л2 -2(п- 1)п)-2(к + 1)(п- 1)пл2)-

-2л?^п+1 (л?) ((к -1) ? + (к +1) л2 - 2п (п2 -1) +

+ 2л^п+1 (л) [(2п (п2 -1) - л?2) ) (л) - 2 (п2 -1) )п+1)]},

Ъ = nqJn (л?) {^+2 (л?) (Л (л) ((к - 1)д + (к +1) л2 - 2(п -1) п) - 2 л/^ (л)) -

-4п((п-1) Jn (л)-л1п+1 (л)))(п + l)Jn+1 (л?)-&п (л?))},

где ? = (СД)2. Этим выражениям соответствуют асимптотические представления в терминах тригонометрических функций:

81пу 008 у

■2(п2 -1)

л? 81П у +

Ъ = 4? (п 008 х - л 81п X){ 4п (п2 -1)- 3 л?

+ л? [л?2 - 3 (п - 1)п] 0082у},

Р2 = ( пл? 008 у + ( - л?2) 81п у ){оо8 х [008 у (-(к -1)? ( 2 (п - 1)п - л? 2) + + л?2 ((к +1) л2 - 2 (п -1) п) - 2 (к +1) (п -1) пл2) -- 2л? 81п у ((к -1) + (к +1) л2 - 2п (п2 -1)) + + 2л 81п X ((2п (п2 -1)-л?2 ) 008 у - 2 (п2 - 1) 81п у)},

Ъ3 = Щ 008 у {л? 008 у [2л 81п X - ((к - 1) ? + (к +1) л2 - 2 (п -1) п) 008 X -4п ((п - 1) 008 X - л 81п X) ) +1) 81п у - л? 008 у)}.

Главный член асимптотики определяет уравнение

С

J„

Л (

Д V2--

Н2

J.

Л

Д V2 --

Н2

2 2 2 V т к

К\к+1 Н2

V

= 0.

Этому уравнению согласно (43) соответствует асимптотическое уравнение в терминах тригонометрических функций

(

81п

-2 2

2Д V2 --

Н2

- кп - -

(

008

V

„ / V2 к2т2 кп к - Д/---^ + —+—

к +1 Н2 2 4

= 0.

Таким образом, грубые приближения корней могут быть найдены из (57).

5. Спектральные разложения

Начнем с наиболее простого случая - осесимметричного, не зависящего от вертикальной координаты г. Собственные значения найдены, и все, что остается сделать -подставить их в решение (40'), и полученную функцию пронормировать, поделив ее на корень из ее скалярного квадрата. Имеем

ик =■

■А! ^

'я X2(■0(кк)2 + ■ (кк)2)-2{кк)■0 {кк) В осесимметричном случае собственные функции определяются соотношением

(58)

иГ =

^(и )• сг, Ч (СГ)Т • ыг • сг,

11 Л

ыг=2ж||г;т • У^ ,

0 0

где сг = ((т , с2Ук") - алгебраический вектор постоянных интегрирования, компоненты

которого определяются нетривиальными решениями системы однородных алгебраических уравнений

( 2 2 Л Л

2а/ (к)-Я (к + 2)а2 + к^-^г /0 (к) 2тГ/ (к)-к/0 (к)]

I Н У

2ат/х (к)

( „„2_2

\

Н2

- ь2

/1 (к)

(с;г Л (0 Л

/"»Бут

Vс2,к У

V 0 У

Так как при V = Vк определитель этой системы равен нулю, то ее решения имеют вид

^Бут _

С1,к =

ь2 -

т2ж2 Л Н2

2 2 т ж

/1 (к), С2ут = 2а 'Нг /1 (к).

Матрица ЫГ может быть вычислена в замкнутом виде:

ыг =

М М12 ^

УМ12 М22;

м„ =

ж

11

м„ =

2Н жЯ

22

2ЬН

(к2 Н2 + т2 Я2 ж2 ))(к)2 + /0 (к)2)- 2кН/0 (к )/1 (к) к? (к?2 + т2 ж2 ))(к)2 + /0 (к )2) - 2ж2 т2/0 (к)/ (к?)

М12 = -тж2 Я/0 (к) /1 (к).

Таким образом, осесимметричные нормированные собственные функции могут быть представлены следующими выражениями:

Ьуеквуа Т.Ы., Ьуекву БЛ. /РМЯРПМвекст^ БиНвПп 4 (2016) 120-150 1

шпг

и,;

-008-

Н

а

Г шХ н2

81П-

шп

Н

шп

{ 2 2 ш п

Н

Н

- Ь2^/ (л)/ (аг)- 2аН J1 (л)/ (Ьг) -- Ь2^/ (л)/0 (аг)- 2а Н/ (л)/ (Ьг)

ег +

В общем трехмерном случае собственные функции разделяются на два класса: чет-

ные ыск и нечетные п\ по переменной ф :

ик =

1 (( ) Ск> ик 1 ( |у=ук ) >

где Ск = (С18к, С2к, С38к ) - нетривиальные решения однородной системы уравнений

ш11,к ш12,к ш13,к

Л

ш21,к ш22,к ш23,к

^ш31,к ш32,к шз3,к У

^ Л

1,к

С8

2,к

С8

Vй3,к У

0 Л

V 0 У

а Ск = (С1кк, С2кк, С3кк ) - нетривиальные решения системы

ш11,к -ш12,к ш13,к ^ ГС к Л 1,к Г 0 Л

-ш21,к ш22,к -ш23,к С0 ^ 2,к = 0

ш31,к -ш32,к ш33,к у С0 Vй3,к У V0 У

где ш^к определяются следующими выражениями:

ш

1,к =-2лк/п-1 (Лк ) +

2 2 ш п

Л,,

(к + 2) + кЕ2 - 2п (п +1)

ш12,к = 2П [(П + 1) (■Лк ) - ЛкЛ-1 ( Лк )] ,

н

2шп

ш13,к =■

н

[(п + 1)/ )-Л*к/п (и )],

2п ( (-1 ( Лк ) - (п + 1) ¿п ( Лк )) ш22,к = 2л<кЦп-1 ( лк ) + Цп (Лк ) (( - 2п (п + 1)) ,

, шп

ш.

23, к

= -Л?к^"Н'/п+2 (Лк ) ш31,к = 2ш (( (л)к -ЛкЕ/п-1 ("к ))

ш

32,к

= пЕ/п (Л к ),

Л

2

Лк2 ¿п-1 (Лк )/ Е 2 -

Г ~2,и2 2 П2 Л

лV Е - 2—

V

Н2

¿п+1 (Лк )

у

Нормирующие множители вычисляются по формулам

2 п НЕ

-к , = | Я (г)Т • г^^ к

0 0 0 2 п НЕ

— = >/(С)Т • N • С,, N = Ш(7е)Т • Гг^ф^

=ук

0 0 0

к

Элементы матриц N, N могут быть вычислена в замкнутом виде:

'м„,к М12,к мв,кЛ 'мП,к -М12,к М,3,к '

N М13,к М 22,к М 23,к N М12,к М 22,к -М 23,к

VМ12,к М13,к М33,к у Ч-М13,к М 23,к М33,к у

Здесь

М1и = -Н {Н

4п

2Р31 п,п + 2";п +1,п + 1,2п +1;-жк

2п (п!)2

-ж

2Р31 п +1, п + 3; п + 2, п + 2,2п + 2; -ж-к2

2п ((п +1)!)

+ ^ (ж,)[жкН2У + 1(жк)-Н2Л2т2У-1 К)]-У (Ж,)[жIН2У + 2(жк)-Н2Л2т2у (жк)]

м = 4-п-1 ^2 п

^ 22,к

(п + 1)(п - 1)!(п +1)!

' 1 2

2 (п +1 )2 2Р3 ( п, п + -2; п +1, п +1,2п +1; -ж

-ж2 „К I п + —;п +1,п +1,п,п + 2,п + 2,2п +1;-ж,2

к 2 3

к

М33,к =

жЯ2 4Н

{Жк [Л (Жк )2 - Л-1 ( Жк К+1 ( Жк )] + 2^2 [Л+1(Жк )2 - Л (Жк К+2 (ж

]

М2з,к [ Л+^к )2 - Л (-к )Л+2 (Жк )] ,

м1з,к = (^к) Л+2 (жк), м12,к = (n, п+1; жк , жк) )kZ (п -1 п; жк, жк )],

где

(от )у (¿г )<* = £

2 (п, т; ж, ж) = 2 (п, т; аг, Ьг) =

(-1)р Ьта2(р+и)2р (-п - р, р; т +1; Ь2 / а2)

)И+т+2 р

р=0

т! р !(п + р )!(2 р + т + п +1)

Я2 р

+п+ т+1

а р Р, - гипергеометрическая функция Гаусса

Т7 ( и и \ ^ (а1)п (ар )п 2

п=о (Ь>)п-(Ь,)п п!

(а)0 = 1, (а)п = а(а +1)(а + 2)„.(а + п-1).

д е

Рис. 4. Ранжирование элементов спектрального разложения Fig. 4. The ranking of elements of spectral decomposition

Таким образом, все элементы спектрального разложения (18) найдены: собственные функции находятся по формулам подразд. 4, а координатные - подразд. 1. Вместе с тем следует отметить, что для трехмерных по пространственным переменным задач (т.е. в отсутствие каких-либо симметрий) число слагаемых в разложении (18) оказывается большим (согласно грубой оценке n3, где n - число слагаемых в одномерной задаче при заданной точности). Для оптимизации вычислений предлагается использовать энергетическое ранжирование элементов спектрального разложения, выполняемое для заданных массовых сил K и начальных данных V0, V1,..., Vm1 по весам:

sk =<Фк , 2 [фк ]>.

Члены суммы (18) упорядочиваются по мере увеличения усредненных по времени абсолютных значений sk. На рис. 4 приведены результаты вычислительного эксперимента с ранжированием. В качестве исходных данных вновь использовались значения R = 1, H = % , к = 5 , а разложение строилось для осесимметричной ступенчатой функции, соответствующей вертикальной компоненте поля массовых сил. На рисунке показаны: a -распределение весов в зависимости от спектрального параметра v (точками обозначены собственные значения, горизонтальная ось - значения спектрального параметра, вертикальная ось - вес); б - результат ранжирования (tol - порядок частичной суммы по переставленному ряду); в - графическое представление тестовой осесимметричной функции в координатах z, r; г - график частичной суммы 70-го порядка при обычном порядке слагаемых (по возрастанию спектрального параметра); д - график такой же суммы утроенного, 210-го порядка; е - график частичной суммы 70-го порядка при ранжировании слагаемых. На этом простом тестовом вычислении видно, что суммирование с изменением порядка позволяет существенно сократить число членов разложения для достижения заданной точности.

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

Работа выполнена при поддержке РНФ 14-29-00198. Библиографический список

1. Рейнер М. Реология. - М.: Наука, 1958. - 224 с.

2. Фрейденталь А., Гейрингер Ф. Математические теории неупругой сплошной среды. — М.: Гос. изд-во физ.-мат. лит., 1962. — 432 с.

3. Бленд Д.Р. Теория линейной вязкоупругости. — М.: Мир, 1965. — 390 с.

4. Reiner M., Lewis H.K. Advancedrheology. — London, 1971. — 374 p.

5. Alber H.-D. Materials with Memory: Initial-Boundary Value Problems for Constitutive Equations with Internal Variables. - Springer, 1988. - 170 p.

6. Работнов Ю.Н. Ползучесть элементов конструкций. — М.: Наука, 1966. — 752 с.

7. Ржаницын А.Р. Теория ползучести. — М.: Изд-во лит. по строит., 1968. — 416 с.

8. Ильюшин А.А., Победря Б.Е. Основы математической теории термовязкоупругости.— М.: Наука, 1970. — 280 с.

9. Пшеничнов С.Г. Динамические задачи линейной вязкоупругости для кусочно-однородных тел // Изв. РАН. МТТ. — 2016. — № 1. — С. 79-89.

10. Amendola G., Fabrizio M., Golden J.M. Thermodynamics of materials with memory. — Springer, 2012. — 574 p.

11. Дэй У.А. Термодинамикапростыхсредспамятью. - М.: Мир 1974. - 190 с.

12. Полянин А.Д., Лычев С.А. Различные представления решений систем уравнений механики сплошных сред // Докл. Акад. наук. - 2014. - Т. 455, № 2. - С. 162-166.

13. Полянин А.Д., Лычев С.А. Различные способы декомпозиции линейных уравнений механики сплошных сред // Докл. Акад. наук. - 2014. - Т. 458, № 6. - С. 663-666.

14. Polyanin A.D., Lychev S.A. Decomposition methods for coupled 3d equations of applied mathematics and continuum mechanics: Partial survey, classification, new results, and generalizations // Applied Mathematical Modelling. - 2016. - Vol. 40. - No. 4. - P. 3298-3324.

15. Nairn J.A. Measurement of polymer viscoelastic response during an impact experiment // Polym. Eng. Sci. - 1989. - Vol. 29. - No. 10. - P. 654-661.

16. Baily F. On the correction of a pendulum for the reduction to a vacuum: Together with remarks on some anomalies observed in pendulum experiments // Philosophical Transactions of the Royal Society of London. - 1832. - Vol. 122. - P. 399-492.

17. Stokes G.G. On the theories of the internal friction of fluids in motion, and of the equilibrium and motion of elastic solids // Mathematical and Physical Papers. Cambridge University Press (CUP). -1866. - Vol. 1. - P. 75-129.

18. Maxwell J.C. The bakerian lecture: On the viscosity or internal friction of air and other gases // Philosophical Transactions of the Royal Society of London. - 1866. - Vol. 156. - P. 249-268.

19. Thomson W. On the elasticity and viscosity of metals // Proceedings of the Royal Society of London. - 1865. - Vol. 14. - P. 289-297.

20. Voigt W. Ueber innere reibung fester körper, insbesondere der metalle // Ann. Phys. Chem.-1892. - Vol. 283. - No. 12. - P. 671-693.

21. Gibbs J. W. Vector Analysis, a text-book for the use of students of mathematics and physics, founded upon the lectures of J. Willard Gibbs. - New Haven Yale University Press, 1901. - 480 p.

22. Gurtin M.E. The linear theory of elasticity // Linear Theories of Elasticity and Thermoelasticity. - Springer, 1973. - P. 1-295.

23. Truesdell C., Noll W. The non-linear field theories of mechanics. - Springer, 2004. - 602 p.

24. Михлин С.Г. Курс математической физики. - М.: Наука, 1968. - 576 с.

25. Kellogg O.D. Foundations of potential theory // Die Grundlehren der Mathematischen Wissenschaften. Vol. 31. - Springer, 1967. - 386 p.

26. Трусделл К. Первоначальный курс рациональной механики сплошной среды. - М.: Мир, 1975. - 592 с.

27. Наймарк М. А. Линейные дифференциальные операторы. - М.: Наука, 1969. - 526 с.

28. Кристенсен Р. Введение в теорию вязкоупругости. - М.: Мир, 1974. - 338 с.

29. Лычев С. А., Сеницкий Ю. Э. Несимметричные интегральные преобразования и их приложения к задачам вязкоупругости // Вестн. Самар. гос. ун-та. Естественно-научная серия. -2002. - Спец. вып. - C. 16-38.

30. Льтев С. А. Связанная динамическая задача термовязкоупругости // Изв. РАН. МТТ. -2008. - № 5. - С. 95-113.

31. Лычев С. А., Манжиров А.В., Юбер С. В. Замкнутые решения краевых задач связанной термоупругости // Изв. РАН. МТТ. - 2010. - № 4. - С. 138-154.

32. Келдыш М. В. О полште собственных функций некоторых классов несамосопряженных линейных операторов // УМН. - 1971. - T. 26,4(160). - С. 15-41.

33. Pochhammer L. Ueber die fortpflanzungsgeschwindigkeiten kleiner schwingungen in einem unbegrenzten isotropen kreiscylinder // Journal für die reine und angewandte Mathematik. - 1876. -Vol. 81. - P. 324-336.

34. Chree C. The equations of an isotropic elastic solid in polar and cylindrical co-ordinates their solution and application // Transactions of the Cambridge Philosophical Society. - 1889. - Vol. 14. -P. 250.

35. Filon L.N.G. On the elastic equilibrium of circular cylinders under certain practical systems of load // Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences. - 1902. - Vol. 198. - No. 300-311. - P. 147-233.

36. Eringen A.C., Suhubi E S. Elastodynamics. Vol. 1; Finite Motions, 1974. - Vol. 2: Linear Theory. - New York-London, Academic Press, 1975. - 1018 p.

37. Meleshko V. Equilibrium of an elastic finite cylinder: Filon's problem revisited //Journal of Engineering Mathematics. - 2003. - Vol. 46. - No. 3/4. - P. 355-376.

38. Ляв А. Математическая теория упругости. - М.-Л.: ОНТИ, 1935. - 674 с.

39. Ватсон Дж.Н. Теория бесселевых функций. - М.: Изд-во иностр. лит., 1949. - 798 с.

References

1. Reiner M. Rheology. Elasticity and Plasticity. Encyclopedia of Physics. Springer, 1958, vol. 3, 6, pp. 434-550.

2. Freudenthal A.M., Geiringer H. The mathematical theories of the inelastic continuum. Elasticity and Plasticity. Encyclopedia of Physics. Springer, 1958, vol. 3,6, pp. 229-433.

3. Bland D.R. The theory of linear viscoelasticity. Oxford, New York: Pergamon Press, 1960. 125 p.

4. Reiner M., Lewis H.K. Advanced rheology. London, 1971, 374 p.

5. Alber H.-D. Materials with memory: initial-boundary value problems for constitutive equations with internal variables. Springer, 1988, 170 p.

6. Rabotnov Iu.N. Polzuchest' elementov konstruktsii [Creep of structural elements]. Moskow: Nauka, 1966, 752 p.

7. Rzhanitsyn A.R. Teoriia polzuchesti [Theory of creep]. Moskow: Izdatel'stvo literatury po stroitel'stvu, 1968, 416 p.

8. Il'iushin A.A., Pobedria B.E. Osnovy matematicheskoi teorii termoviazko-uprugosti [Foundations of the mathematical theory of thermoviscoplastic]. Moskow: Nauka, 1970, 280 p.

9. Pshenichnov S.G. Dinamicheskie zadachi lineinoi viazkouprugosti dlia kusochno-odnorodnykh tel [Dynamic linear viscoelasticity problems for piecewise homogeneous bodies]. Mechanics of Solids, 2016, vol. 51, no. 1, pp. 65-74.

10. Amendola G., Fabrizio M., Golden J. M. Thermodynamics of materials with memory. Springer, 2012, 574 p.

11. Day W.A. The thermodynamics of simple materials with fading memory. Springer Science & Business Media, 2013, 134 p.

12. Polyanin A., Lychev S. Razlichnye predstavleniia reshenii sistem uravnenii mekhaniki sploshnykh sred [Various representations of the solutions of systems of equations of continuum mechanics]. Doklady Physics, 2014, vol. 59, pp. 148-152.

13. Polyanin A., Lychev S. Razlichnye sposoby dekompozitsii lineinykh uravnenii mekhaniki sploshnykh sred [Various decomposition techniques for linear equations of continuum mechanics]. Doklady Physics, 2014, vol. 59, pp. 487-490.

14. Polyanin A.D., Lychev S.A. Decomposition methods for coupled 3d equations of applied mathematics and continuum mechanics: Partial survey, classification, new results, and generalizations. Applied Mathematical Modelling, 2016, vol. 40, no. 4, pp. 3298-3324.

15. Nairn J.A. Measurement of polymer viscoelastic response during an impact experiment. Polym. Eng. Sci, 1989, vol. 29, no. 10, pp. 654-661.

16. Baily F. On the correction of a pendulum for the reduction to a vacuum: Together with remarks on some anomalies observed in pendulum experiments. Philosophical Transactions of the Royal Society of London, 1832, vol. 122, pp. 399-492.

17. Stokes G.G. On the theories of the internal friction of fluids in motion, and of the equilibrium and motion of elastic solids. Mathematical and Physical Papers, vol. 1, Cambridge University Press (CUP), 1866, pp. 75-129.

18. Maxwell J.C. The bakerian lecture: On the viscosity or internal friction of air and other gases. Philosophical Transactions of the Royal Society of London, 1866, vol. 156, pp. 249-268.

19. Thomson W. On the elasticity and viscosity of metals. Proceedings of the Royal Society of London, 1865, vol. 14, pp. 289-297.

20. Voigt W. Ueber innere reibung fester korper, insbesondere der metalle. Ann. Phys. Chem., 1892, vol. 283, no. 12, pp. 671-693.

21. Gibbs J.W. Vector Analysis, a text-book for the use of students of mathematics and physics, founded upon the lectures of J. Willard Gibbs. New Haven Yale University Press, 1901, 480 p.

22. Gurtin M.E. The linear theory of elasticity. Linear Theories of Elasticity and Thermoelasticity. Springer, 1973, pp. 1-295.

23. Truesdell C., Noll W. The non-linear field theories of mechanics. Springer, 2004, 602 p.

24. Mihlin S.G. Mathematical physics, an advanced course. In: North-Holland series in applied mathematics and mechanics 11. Amsterdam : North-Holland publishing co., 1970, 561 p.

25. Kellogg O.D. Foundations of potential theory. In: Die Grundlehren der Mathematischen Wissenschaften Vol. 31. Springer, 1967, 386 p.

26. Truesdell C.A.A First Course in Rational Continuum Mechanics. Academic Press, 1976, 352 p.

27. Naimark M.A. Lineinye differentsial'nye operatory [Linear differential operators]. Moskva: Nauka, 1969, 526 p.

28. Christensen P.M. Theory of Viscoelasticity: An Introduction, 2nd Edition. Academic Press, 1982, 364 p.

29. Lychev S.A., Senitskii Iu.E. Nesimmetrichnye integral'nye preobrazovaniia i ikh prilozheniia k zadacham viazkouprugosti [The nonsymmetric integral transform and their applications in viscoelasticity]. Vestnik Samarskogo gos. universiteta. Estestvennonauchnaia seriia. 2002, Spetsial'nyi vypusk, 16-38 pp.

30. Lychev S. Coupled dynamic thermoviscoelasticity problem. Mechanics of solids, 2008, Vol. 43, no. 5, pp. 769-784.

31. Lychev S.A., Manzhirov A.V., Joubert S.V. Closed Solutions of Boundary-Value Problems of Coupled Thermoelasticity. Mech. Solids, 2010, vol. 45, no 4, pp. 610-623.

32. Keldysh M.V. O polnote sobstvennykh funktsii nekotorykh klassov nesamosopriazhennykh lineinykh operatorov [To completeness of eigenfunctions of some classes of nonselfadjoint linear operators]. UMN, 1971, vol. 26, no. 4(160), pp. 15-41.

33. Pochhammer L. Ueber die fortpflanzungsgeschwindigkeiten kleiner schwingungen in einem unbegrenzten isotropen kreiscylinder. Journal fur die reine und angewandte Mathematik, 1876, vol. 81, pp. 324-336.

34. Chree C. The equations of an isotropic elastic solid in polar and cylindrical co-ordinates their solution and application. Transactions of the Cambridge Philosophical Society, 1889, vol. 14, pp. 250.

35. Filon L.N.G. On the elastic equilibrium of circular cylinders under certain practical systems of load. Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences, 1902, vol. 198, no. 300-311, pp. 147-233.

36. Eringen A.C., Suhubi E. S. Elastodynamics. Vol. 1; Finite Motions. 1974. Vol. 2: Linear Theory. New York-London, Academic Press, 1975, 1018 p.

37. Meleshko V. Equilibrium of an elastic finite cylinder: Filon's problem revisited. Journal of Engineering Mathematics, 2003, vol. 46, no. 3/4, pp. 355-376.

38. Love A.E.H. A treatise on the mathematical theory of elasticity. Cambridge University Press, 1906, 583 p.

39. Watson G.N. A treatise on the theory of Bessel functions. Cambridge university press, 1966, 804 p.