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

УДК 539.3, 534.1

ОБ ОДНОЙ ЗАДАЧЕ МАТЕМАТИЧЕСКОЙ ФИЗИКИ С ПОДВИЖНОЙ ГРАНИЦЕЙ © 2008 г. Л.О. Кечеджиян, Н.А. Пинчук, А.М. Столяр

Method of asymptotic integration has been developed for solution of initial boundary value problem with a moving boundary. A comparison of asymptotic results with results of numerical integration of initial problem with moving boundary has been carried.

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

Рассматриваемая задача моделирует продольные колебания троса, который разматывается с заданной скоростью, и подвешенного на него груза. Анализ задач, связанных с движением тел на тросе переменной длины, имеет давнюю историю. Первые работы об исследовании напряжений в тросах переменной длины появились более 60 лет назад [1] и были вызваны потребностями угольной отрасли в проектировании надежных подъемных устройств. В работах А.Ю. Иш-линского [2, 3] получено уравнение в частных производных, описывающее колебания груза на тросе переменной длины, которое сведено при упрощающих предположениях к интегро-дифференциальному уравнению относительно функции времени. В трудах Г.Н. Савина и О.А. Горошко [4, 5] многочисленные задачи динамики объектов переменной длины сводятся к решению дифференциальных и интегро-дифференциаль-ных уравнений с изменяющимися во времени пределами интегрирования и ядрами. Теории волновых процессов в упругих системах с изменяющимися во времени параметрами посвящена монография А.И. Вес-ницкого [6], в которой, кроме того, содержится обширная библиография по изучаемым проблемам.

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

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

Используемые в настоящей работе численные алгоритмы впервые предложены в [9] для расчета продольных колебаний троса переменной длины, а в [10]

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

Постановка задачи

Рассмотрим поведение системы, схематически изображенной на рис. 1. Это трос, намотанный на катушку и разматывающийся с неё с заданной скоростью, и висящий на нём груз. Начало подвижной системы координат Ох расположим в точке крепления груза к тросу и направим ось Ох вертикально вверх. Будем изучать продольные колебания троса. Поведение такой системы опишем следующей начально-краевой задачей:

Рис. 1.

pF

/д 2u d 2Л

dt2 dt2 ,

J

d _ - EF du

dt2 dx

,д 2u

= EF-- -pFg. дх 2

+ mg,

х = 0

4(t) = i(t) + u(i (t), t),

u( x, o| x = 0 = 0 • u(x, t)| t = 0 = ф (x) ,

du( x, t)

(1)

dt

= Ф2 (x) ■

t = 0

l(t)\ t = 0 = l0 •

du(l, t)

dl_ dt

1 +

dl

= ey (t) •

(2)

(3)

Здесь и( х, /) - смещение сечения х троса в момент времени t; р, ¥, Е - плотность, площадь попереч-

m

ного сечения и модуль Юнга для троса; g - ускорение силы тяжести; m - масса груза; l(t) - длина троса в момент времени t в недеформированном состоянии;

- фактическое расстояние от точки схода троса с катушки до груза. Соотношения (1) для моделирования системы, изображенной на рис. 1, впервые, по-видимому, даны в работе А.Ю. Ишлинского [2]. Функции <(x) и (piix) описывают начальные смещение и скорость сечения x троса; еу(t) - скорость разматывания троса, е - малый числовой параметр, о котором подробнее речь идет ниже - в конце текущего раздела. Соотношение (3) выражает равенство между скоростью сечения троса x = l(t) и линейной скоростью точки обода колеса в месте схода троса с

колеса. Отметим также, что слагаемое

dl du (l, t)

в (3)

dt дI

появляется как раз вследствие того, что длина троса -величина переменная [2]. Кроме того, предполагается, что выполняется условие согласования начального и граничного условий: ^ (0) = 0 .

Сделаем замену переменных и (х, ^ = £,($) — х—ы(х^). После соответствующих преобразований задачу (1), (2) перепишем в виде

d 2U dt2 d 2U

2 d2U

= a+ g ,

dt

U

x = 0

dx 2 _ EF m

2 E a =—.

P

dU

x = l (t)

= 0,

dx

dU

+ 1

x=0

+ g,

U|t = о =$i(x)

= Ф2(Х) ,

(4)

(5)

(6)

t = 0

где ф^(х) = ¡0 +Ф\(1о) — х — ф\(х), Фг(х) = —<Рг(х).

Примем, что скорость разматывания троса, т.е. величина d¡ / dt, мала по сравнению со скоростью распространения волны в тросе (реальная скорость распространения волны может быть величиной порядка 4000 - 5000 м/с). Малый параметр е введём по формуле е = V# / У00, где У* - максимальная (за время рассматриваемого процесса) скорость возрастания длины троса; У00 =^[Щр - скорость распространения волны в тросе. Кроме этого, в равенстве (3) можно

положить — = , т.е. пренебречь величиной де-

dt

формации — по сравнению с 1. Как показывают чис-

дх

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

Асимптотическое интегрирование начально-краевой задачи с подвижной границей

где W (x,l (t)) = 1 ^g 12(t) - x2 ]+ ^ mg + ^ (t) - x] • (7)

Функция W при l(t) = const (случай троса постоянной длины) зависит только от переменной x и является решением статической задачи, соответствующей данной динамической. Полагаем, что в случае троса переменной длины функция W зависит еще и от l(t), которую можно считать «медленной переменной», имея в виду, что в силу (3) скорость изменения l(t) пропорциональна малому параметру е •

Подставляя (7) в (4), (6), с учетом (3) получаем уравнение, начальные и граничное условия относительно функции w(x, t):

d 2 w dt2 d2w

d2w

+ F (t)s+ F (t )e = a 2 —:

dx

2

dt

w

x = 0 t = 0 =фll(x)

9 EF dw + Fi(t )e + F1(t)e2 = ef

m dx

(8)

x=0

fit = 0 =ф22( x) , где Flit) ==if,0 +mg + 1

d^ dt '

F2(t) =

PPg_ E

v2{t) + i(t) ^ dt

ll (t) = J V(T)d? ,

0

Ф11 =Ф - W(x,l0),

Ф22(x) = Ф2(x)+ ^ +1 \y{0)e .

(9)

Отдельно остановимся на преобразовании граничного условия (5). Подставляя (7) в (5), находим

> = 0.

(10)

w х = Щ) ''

где для I(Г) имеем, исходя из (3), выражение

¡(0 = ¡0 + ¡^)е. (11)

Подставим (11) в (10), учтем, что параметр е мал, и разложим функцию w(¡о + ¡\^)е,1) в ряд Тейлора в окрестности точки х = ¡0 . Из равенства (10) имеем

w

x=ln

dw dx

x=ln

x l1 (t)e + —

1 d2 w

2 dx2

(12)

x=lr

x l12(t)e2 +... = 0.

Задача (8), (12) содержит малый параметр е в уравнении и граничных условиях (она является регулярно возмущенной - по терминологии [11]). Применяя метод асимптотического интегрирования, её решение построим в виде ряда по степеням малого параметра е .

Будем искать решение задачи (3-6) в виде U (x, t) = w(x, t) + W(x, l(t)),

w = 2 wt (x,t)e'. i=0

2

2

x

Подставим (13) в (8), (12) и сгруппируем коэффициенты при одинаковых степенях е . Приравнивая коэффициенты при одинаковых степенях е, для определения членов ряда (13) получаем начально-краевые задачи уже на постоянном отрезке [0, /д ] изменения переменной х.

0 12

Приравнивая коэффициенты при е , е , е ,

ек (к > 2), находим

при е

о .

d 2w0 dt2 - a2 ^ dx2 ' (14)

d w0 EF dw0

dt 2 x=0 m dx x=0 '

w0 x= = 0, l0

wo\t=0 =^1<X>>

W = ф2( x). dt l'=° V2 )

при ei :

2 2 d_wL + m = a2 d-Wl

dt2 d wi I

dt Wi

dx

2

(15)

2 I x=°

X=lo

+ Fi(t) =

dw° dx

EF dwi m dx

, li(t) ,

w = 0, 1 t=0

dwi I

_drlt=° = "

при e2:

f+mg+iW°);

d 2 w2

dt d w2

^ ^ ч 2 d2w2 2 + F2(t) = a2

dt w2

2

dx'

+ F2(t) = EF dw2

x=0 m dx

(16)

x=0

x=ln

=-li(t)

dwi

dx

x=lr

i l 2(t) d2w0 -2li (t)l^

x=ln

w = 0 .

2 t=0

dwo

dt

t=0

k,

= 0.

при ek (k > 2):

2 2

d wk = a2 d wk

dt2

d wk

dx

2

dt wk

2

_ EF dwk x=0 m dx

x = 0

k-i [li(t)]k-i dk-iw1

x=lo

¿=0 (k - i)

dx'

k-i

x=L

wk = 0. k t=0

dwi.

dt

t=0

= 0.

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

Задачи (14) - (17) допускают аналитическое решение методами математической физики (см., например, [12]). Однако с ростом номера задачи (и номера члена ряда (13)) её решение приводит ко все более громоздким математическим выкладкам. Ниже будет проведено сравнение результатов численного решения задачи с применением ряда (13) и результатов непосредственного численного интегрирования исходной задачи с подвижной границей.

Замечание 1. Поведение главного члена разложения (13) определяется начальными условиями задачи (14), т.е. функциями фц(х) и . Поэтому в случае фц(х) = 0, ф2(х) = 0, очевидно, имеем Wo(x,г) = 0 .

Замечание 2. В уравнении и граничном условии задачи (15), определяющей первый член разложения -функцию w1(x, г), в качестве «внешнего воздействия» на систему выступает функция F1 (?), которая описывает влияние на колебательный процесс масс троса и груза и зависит от скорости разматывания троса: чем больше масса троса и груза, чем больше ускорение разматывания, тем это влияние существеннее.

Основные факторы, влияющие на второй член разложения, - это скорость и ускорение разматывания троса (функция F2(t), определяемая по (9)).

Численное интегрирование начально-краевой задачи с подвижной границей

Для численного решения задачи (1) - (3) применялись алгоритмы, основанные на методах конечных разностей и Рунге-Кутта, которые впервые были предло-, жены в [9]. Ниже будем называть их МКР-П и РК-П -

соответственно метод конечных разностей и метод Рунге-Кутта для системы с подвижной границей. Они зарекомендовали себя в качестве эффективных методов решения задач с подвижными границами - как модельных, так и реальных инженерных [7, 8]. Коротко опишем применение МКР-П.

В каждый момент времени г у разобьем трос дли-

(17) ной / (г) на N равных частей узлами х1 (г у) = ¡к у , где г у = ¡Ы; Ы - шаг интегрирования по времени; у -

соответствующий временной слой сетки, /=0,1, ..., Ы; к у = /(г у)/ N - шаг сетки по пространственной переменной в момент времени г у . Производные, входящие в задачу (1) - (3), аппроксимируем центрально-

x=0

0

разностными формулами. На рис. 2 показан шаблон сетки, которая получается при описанном способе дискретизации задачи.

Крестиками отмечены точки сетки, в которых вычисляются значения неизвестной функции при конечно-

Рис. 2. Шаблон метода МКР-П

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

со временем, узлы сетки на предыдущих слоях по времени, участвующие в вычислениях, являются гипотетическими, фиктивными (они и отмечены крестиками), и значения функции в них определяются при помощи формул интерполяции. Некоторые расчетные формулы алгоритмов МКР-П и РК-П приведены в работах [9, 10].

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

Для проверки эффективности асимптотики проведено сравнение результатов численного интегрирования исходной задачи (1) - (3) с результатами асимптотического интегрирования с использованием ряда (13). Расчет проводился для следующих значений параметров и функций задачи:

F = 4,41х10"6я- i 2;

E = 1011 1/i

*"=E if'2 + m-

; f = 5197 ea/i 3 ;

= 0

l0 = 50 i dl

dt

= Kt . (18)

Вид функций p(x) и <P2(x) соответствует следующей постановке: в начальный момент времени трос растянут весом груза и собственным весом, и его начальная скорость равна нулю. Из этого положения начинается разматывание троса с заданной скоростью

—; коэффициент пропорциональности K - это уско-

dt

рение, с которым происходит разматывание. Нетрудно показать, что в указанном случае w0 = 0. Это вполне предсказуемый и объясняемый с механической точки зрения результат. Дело в том, что при указанных начальных условиях решение задачи в случае троса постоянной длины должно совпадать с решением соответствующей статической задачи, которое определяет функция W в (7) при / (t) = /о = const. В случае троса переменной длины динамический процесс

генерируется только разматыванием троса и учитывается последующими членами разложения (13) -функциями х,/) , ^(х,г) , м>у.(х,г)(к > 2) .

В [9], [10] при расчете задачи (1) - (3) и в [10] при расчете задачи (4) - (6) показано, что оба численных алгоритма - МКР-П и РК-П дают весьма близкие результаты на большом отрезке времени (равном примерно 3,5 периодам колебаний), а на отрезке времени, равном периоду колебаний, результаты практически совпадают.

Отметим, что совпадают также результаты расчетов задач (1) - (3) и (4) - (6) между собой и результаты расчетов задачи (1) - (3) для различных вариантов начальных условий, задаваемых функциями х) и <Р2(х), а не только для выражений (18). Проводились расчеты для следующих функций, задающих начальные условия: , с 2 '

<М -) =

fgx

+, *2(x) = 0;

2E FE

p(x) = 0, q>2(x) = k(x-10) ; ^(x) = 0, ^(x) = 0;

x) = fj2, %(x) = 0.

(19)

В частности, на рис. 3 дано сравнение результатов расчетов исходной задачи для начальных условий, которые задаются функциями (19). Условия (19) означают, что в начальный момент времени трос не деформирован, находится в покое, на него подвешивается груз и начинается разматывание с заданной скоростью —. Кривая 1 получена с применением МКР-

П, а кривая 2 - РК-П.

и(Щ I), ц _

U(l(t), t)

Рис. 3. Сравнение методов МКР-П и РК-П

Ввиду данного успеха взаимотестирования методов МКР-П и РК-П, которые можно считать альтернативными, сравнение с результатами асимптотики проводилось уже только с применением одного метода, а именно РК-П. На рис. 4 приводятся графики зависимости от времени функции ы(1 (/),/) - удлинения троса в точке его схода с колеса. Нетрудно показать, что

и(1(1), 1) = и(0, 1) - ¡(1:). Обозначим её через и. Номером 1 отмечена кривая, полученная численным интегрированием задачи (4) - (6) при помощи метода РК-П; 2 и 3 - кривые, отвечающие функции VI, которая вычисляется по формуле VI (1) = ^(0,1) + W(0, ¡(1)) - ¡(1), где w(х, 1) определяется рядом (13). Как отмечалось выше, для начальных условий, которые задаются по формулам (18), главный член разложения (13) - функция Wo(x, 1) = 0 . Функции w■í(x, 1) и W2(x, 1) опреде-

И1 м

лялись численным решением задач (15) и (16). Рис. 4 отвечает значениям К = 2 и К = 20 соответственно (коэффициент К в формуле (18) равен ускорению, с которым возрастает скорость разматывания троса). Эти значения К отвечают значениям е» 0,001 и е» 0,01. При вычислении кривой 2 на рис. 4 учтен только первый член ряда (13) Wl(х, 1)е; в кривой 3 учтены два

члена разложения (13), т. е. сумма w1(x, 1)ш2(х, 1)е2.

Ui м

16 13

Рис. 4. Сравнение асимптотических и численных результатов при а - К=2, б - К=20

Как мы видим, на рис. 4а при среднем ускорении разматывания троса (другими словами, возрастания области интегрирования) все три кривые, как численная, так и асимптотические, практически совпадают. Кривые на рис. 4а представляются физически достоверными, так как функция и(1,(1),1) остается положительной во всё время движения системы. При высоком ускорении разматывания троса (К = 20 на рис. 4б) кривые с течением времени начинают расходиться, функция VI колеблется около нуля, становясь периодически отрицательной. Это означает (с механической точки зрения), что в соответствующие промежутки времени трос стремится сжаться в продольном направле*"™ В реальности это соответствует образованию ненапряженного участка. Физически этот результат также достоверен: при К=20 ускорение разматывания троса превышает ускорение силы тяжести, т.е. колесо стремится опередить трос, что и приводит к его провисанию.

При «сверхвысоком» ускорении разматывания троса (К=200) описанный эффект наступает раньше.

Анализ данных расчетов позволяет сделать вывод, что во всякий фиксированный момент времени 1 с уменьшением е численные и асимптотические результаты сближаются. Метод асимптотического интегрирования в сочетании с численным интегрированием позволяет эффективно решать широкий класс задач о колебаниях троса переменной длины. В этот класс входят задачи с различным набором начальных условий при ускорении разматывания троса, не превышающем 10 м/с2 (что соответствует е«0,005).

Авторы считают своим приятным долгом выразить благодарность И.М. Бермусу [7], который внес решающий вклад в разработку численных алгоритмов, использованных для расчетов в настоящей статье.

Авторы признательны В.Б. Левенштаму за внимание к работе.

Работа выполнена при поддержке РФФИ, грант № 06-01-00287.

Литература

Неронов Н.П. // Определение напряжений в подъемных канатах: Тр. совещания по шахтным подъемным канатам. М., 1944.

Ишлинский А.Ю. // Укр. мат. журн. АН УССР. 1953. Т. 5. № 4. С. 370-374.

Ишлинский А.Ю. // Докл. АН СССР. 1954. Т. 95. № 5.

Савин Г.Н., Горошко О.А. Динамика нити переменной длины. Киев, 1962.

Горошко О.А., Савин Г.Н. Введение в механику деформируемых одномерных тел переменной длины. Киев, 1971.

Весницкий А.И. Волны в системах с движущимися границами и нагрузками. М., 2001.

Бермус И.М. и др. // Докл. VII Всесоюз. съезда по тео-рет. и прикл. механике. М., 1991. С. 47.

Прошкин С.Г., Столяр А.М., Бернштейн Л.А. // Современные проблемы математического моделирования: Докл. IX Всерос. школы-семинара. Абрау-Дюрсо, 2001. С. 287-293.

6.

t, c

t, c

б

а

9. Бермус И.М. и др. // Докл. расширенных заседаний се-

минара по прикладной математике им. И.Н. Векуа. Тбилиси, 1989. Т. 4. № 3. С. 29-32.

10. Пинчук Н.А., Столяр А.М. // Современные проблемы

механики сплошной среды: Докл. X Междунар. конф. Ростов н/Д, 2006. С. 217-221.

11. Васильев А.Б., Бутузова В.Ф. Асимптотические методы

в теории сингулярных возмущений. М., 1990.

12. Тимошенко С.П., Янг Д.Х., Уивер У. Колебания в инже-

нерном деле. М., 1985.

Южный федеральный университет_20 марта 2007 г.