Нелинейная динамика. 2011. Т. 7. № 3. С. 635-647. Полнотекстовая версия в свободном доступе http://nd.ics.org.ru
УДК: 532.3,519.635.8
М8С 2010: 76Б05, 35Q30, 65Ш6
Движение тела с переменной геометрией масс в безграничной вязкой жидкости
В статье изложена постановка задачи о движении тела в жидкости, вызванного перемещением внутренней материальной точки. Математическая модель содержит уравнения движения твердого тела и гидродинамические уравнения Навье-Стокса. Задача решается численно с применением конечно-разностной схемы.Проведенные численные исследования выявили определяющее влияние вязких сил и моментов на траекторию движения тела.
Ключевые слова: самопродвижение, уравнения Навье-Стокса, вязкое вихревое движение, численные методы
1. Введение
В работах [1, 2] изучено движение твердых тел в идеальной жидкости за счет изменения геометрии масс. Обнаружено, что закон изменения внутренней геометрии масс можно подобрать таким образом, чтобы осуществить перемещение тела в произвольную точку. При этом, в отличие от подавляющего большинства более ранних и современных работ, внешняя
Получено 1 сентября 2011 года После доработки 24 сентября 2011 года
Работа выполнена в рамках гранта Правительства РФ для государственной поддержки научных исследований, проводимых под руководством ведущих ученых в российских образовательных учреждениях ВПО (дог. №11.034.31.0039); гранта Президента РФ поддержки ведущих научных школ НШ-8784.2010.1 «Динамические системы классической механики»; ФЦП «Научные и научно-педагогические кадры инновационной России» на 2009-2013 годы, мероприятие 1.1. «Научно-образовательный центр "Регулярная и хаотическая динамика"» ГК №02.740.11.0195.
Рамоданов Сергей Михайлович [email protected]
Институт компьютерных исследований Удмуртский государственный университет 426034, Россия, г. Ижевск, ул. Университетская, 1
Тененев Валентин Алексеевич [email protected]
Ижевский государственный технический университет 426069, Россия, г. Ижевск, ул. Студенческая, 7
С. М. Рамоданов, В. А. Тененев
оболочка тела остается неизменной. Управляемость движения обеспечивается эффектом присоединенных масс. В реальной жидкости не действует парадокс Даламбера и наличие вязкости приводит к изменению характера обтекания тел жидкостью. В работе ставится задача о расчете движения тел с изменяемой внутренней геометрией в вязкой жидкости. Проведено тестирование применяемого численного метода на задаче о движении круглого цилиндра в вязкой жидкости. Расчетное исследование показало сильное влияние сил сопротивления, обусловленных вихревым движением жидкости, на траекторию движения тела.
2. Уравнения движения тела и жидкости
Так же, как в работах [1, 2], рассматриваем плоское движение тела с жесткими границами в бесконечной вязкой жидкости. Тело состоит из внешней оболочки и материальной точки конечной массы, которая перемещается в пределах оболочки по заданному закону. Движение всей системы начинается из состояния покоя.
Наряду с неподвижной системой координат (х,у) вводится подвижная система, связанная с телом. Оси этой системы координат (£, ц) ориентированы по главным осям инерции (рис. 1).
Рис. 1. Системы координат.
Внутренняя точка с массой m перемещается по некоторому заданному закону (£(t), n(t)). Для плоской постановки задачи кинетическая энергия системы (тело + жидкость) записывается в виде
Т = ± [Ami + Mvt + Ви2 + m(u2m + v2m)],
где Ai = M + X\, A2 = M + A2, B = J + A3, M, J — масса тела и момент инерции без учета материальной точки, Ai, A2, A3 — коэффициенты присоединенных масс, и — угловая скорость вращения, üb = (щ, vb) — вектор скорости движения точки Oi тела, um = = C(t), vm = n(t) — проекции скорости движения материальной точки в подвижной системе
координат. Тогда уравнения движения тела в жидкости запишем следующим образом:
(1Рх
-т- - иР2 = аЬ
(Р2 , о Т?
-Ж+шР1=Р2> (2.1)
Ш.+иьР2-чР1 = С, Р1(0) = Р2 (0) = к (0) =0,
где Р1 = А1иь + т(иь + {(Ь) — шп), Р2 = А2иь + т(иь + п(Ь) + — компоненты импульса в подвижной системе координат, К = Вш — ш(щ + {(Ь) — шп)п + т(^ь + п(Ь) + -
кинетический момент.
Для вязкой жидкости Е = (^1, — сила, действующая со стороны жидкости на тело, О — момент этой вязкой силы.
Движение вязкой несжимаемой жидкости описывается уравнениями Навье-Стокса:
V-V = 0, (2.2)
^ + (У-У)У = -^+г/У2У, (2.3)
где V = (и, V) — вектор абсолютной скорости жидкости, р, р, V — давление, плотность и кинематический коэффициент вязкости.
Уравнение (2.3) представим в «дивергентной» форме и запишем в подвижной системе координат [3]:
^+V•((V-W)V) = -^ + г/V2V-ГixV (2.4)
с граничными условиями на поверхности тела
V - п = (иь + П х г,) - п, п х V = п х (иь + П х г,)
и на бесконечности V = 0.
Здесь W = Щ + П х Г1 — вектор переносной скорости частицы жидкости в точке Г1; П = (0, 0, ш)Т — вектор угловой скорости. Силы реакции жидкости на тело и момент вязких сил определяются интегралами по поверхности [4]:
ЕС(Ь) = У [—р1 + рv (VV + ТVV)] п(в, (2.5)
Сс(Ь) = У г, х [—р1 + рv ^ + ТVV)] п (в. (2.6)
В соответствии с теоретическим анализом, проведенным в работах [5, 6], суммарная сила воздействия жидкости на движущееся твердое тело состоит из трех частей.
1. Первый вклад определяется действием присоединенных масс для соответствующего потенциального течения.
2. Второй вклад связан с генерацией завихренности в тонком слое около тела при граничных условиях прилипания на поверхности.
3. Третий вклад в суммарную силу определяется вязкими напряжениями трения в вязкой жидкости.
В правой части системы уравнений (2.1) находятся сила и момент, обусловленные вязкостью жидкости, т. е. второй и третий вклады в суммарную силу воздействия жидкости на тело. Первый вклад соответствует слагаемым в левой части уравнений (2.1), содержащих коэффициенты присоединенных масс Х\, Л2, A3. Гидродинамические уравнения (2.2), (2.4) решаются численно с использование конечно-разностной сетки. На каждом шаге по времени At скорости тела U и П считаются постоянными и вклад ускорения тела в сопротивление
учитывается левой частью уравнений (2.1) Сила Fc(i), вычисленная
по формуле (2.5) при значениях p(t), V(t), полученных из численного решения уравнений Навье-Стокса, не равна в точности сумме второго и третьего вкладов. Как показано в работе [4], в этой силе присутствует вклад от интегрирования давления по поверхности тела, связанный с движением и вращением тела. Действительно, запишем уравнения (2.4) в форме Громеки-Ламба для абсолютного движения в подвижной системе координат
9'V + V (- W • V + ^ ) = z/V2V + (V - W) х rot V
dt \ 2 Р^
и рассмотрим случай потенциального течения жидкости с потенциалом Чф = V и граничными условиями непротекания на поверхности тела. Интеграл Коши-Лагранжа в этом случае имеет вид
dt V 2 vvv^ Р '
и
Р = Рос-Р^-Р^~ wvj = 0.
Интегрирование по поверхности тела дает силу f p ■ n
äs = Л^ + Л(П х U5), где Л —
матрица коэффициентов присоединенных масс. Первая часть связана с ускоренным движением тела, вторая — с вращением. Следовательно, при численном решении уравнений гидродинамики (2.4) сила Fc(t) содержит составляющую Л(П х Ub), которая также имеется в левой части уравнений (2.1). С учетом сказанного, уравнения движения тела в жидкости принимают вид:
dPi
- Lü(Mub + т(уъ + ф) + wO) = Fei,
dp
— - ш(Мщ + т{щ + Í(t) - шг])) = Fc2, ^2 ^
^ + иь(Миь + т(иь + Ш - urj)) - Vb{Mvb + т{щ + ф) + = Gc,
Pi(0) = P2(0) = K (0) =0,
где силы и момент вычисляются по формулам (2.4), (2.5). Эти уравнения следует дополнить кинематическими соотношениями, определяющими абсолютные координаты x, y точки Oi тела:
Х = ub cos a — vb sin a, y = ub sin a + vb cos a, a = ш, (2.8)
где a — угол поворота жестко связанной с телом системы координат относительно неподвижных осей. Условие того, что движение началось из состояния покоя, имеет вид
Pi(0)= P2 (0) = K (0)=0. (2.9)
3. Случай идеальной жидкости
Приведем некоторые этапы анализа самопродвижения тела в идеальной жидкости из [1]. Уравнения движения тела имеют вид (2.1), в которых F\ = F2 =0 и G = 0. Они предста-вимы в виде уравнений Кирхгофа
/&ГХ дТ _ /&ГХ дТ _
\д иь) dvb ' \dvbJ диь
(dT X дТ дТ
+ иь о--г'ь = 0,
\дш J dvb диь
где
Т = + A2V% + Вш2) + Ц1Щ + Ц2Щ + + 36
представляет собой кинетическую энергию системы (тело + жидкость). Функции ,№, ж,
описывающие изменение формы тела и перераспределение масс внутри него, считаются известными. В рассматриваемом случае движения внутри тела единственной точки они легко могут быть выражены через ее координаты £(t), n(t) и их производные £(t), n(t).
Поскольку тело начало движение из состояния покоя, импульсы и кинетический момент
дТ дТ дТ
Pi = ?- = 0, Р'2 = 4- = О, К = %- = О (3.1)
диь дvь дш
являются интегралами движения.
Выразив иь, Vb, ш из (3.1) и подставив полученные выражения в (2.8), получим уравнения, описывающие движение точки Oi тела:
№1 , № . Vi . №2 V
х = —— cos а + — sin а, у = —— sin а--— cos а, а = —.
Аг А2 ' у Аг А2 В
Полагая Vk = -Vk/Ak (k = 1, 2) и п = —v/B, получим
z = veia, где z = x + iy, v = v1 + iv2.
Пусть теперь v{t) и rj(t) периодичны с периодом Щ-. Поскольку а = г), то a(t) = Ш + Ф(£), где Q = (п), $(t) — 2п/ш периодична по t. Отсюда z(t) = Z(t)expiQt, где Z(t) периодична по t с периодом 2п/ш. Разлагая эту функцию в ряд Фурье и интегрируя по t, получим координаты точки O1 тела:
z
Z(t)= У const.
v у ^ i(nu + Q)
Значит, при Q = пш, n Е Z, ряд расходится и движение тела не ограничено. Подобные резонансные соотношения реализуются, например, при движении точки по траектории, изображенной на рисунке 6.
4. Численный метод
Для численного решения двумерных уравнений (2.2), (2.4) применяется разностная сетка, координатными линиями которой являются линии тока и эквипотенциальные линии для потенциального обтекания тела. Сетка строится комплексным методом граничных элементов [7]. Для решения уравнений (2.2), (2.4) в криволинейной системе координат применяется проекционный метод [7, 8]:
_ ук
к+-
м
V-V
к+1
+ V- ((V - W)V) 2 = -
0,
2
+ + V*) - (Пх V)" ' 2
к+7
реализуемый в три этапа.
у* _ у к
дГ
1. т г.т + V • ((V - W) У)А'+2 = + г/У2 (V* + УА') - (о х УА'+2
2.
ук+1 = у* _
V-(у* - АЖфк+1^) = 0.
(4.1)
3. //' ' 2 = рк + фк+1 - '.
Разностные уравнения (2) решаются совместно с системой уравнений (2.7), для которых применяется метод Рунге-Кутты.
Применяемый алгоритм тестировался на задаче о движении круглого цилиндра в жидкости. Цилиндр приобретает скорость V = 1 при Ь = 0. Для числа Рейнольдса Ее = 200 изменение во времени сил , показано на рисунке 2.
-1
Р
Рис. 2. Изменение сил сопротивления на поверхности тела во времени.
С момента времени Ь = 40 начинаются существенные колебания продольной и поперечной сил сопротивления. Это связано с нестационарным характером отрыва вихрей от цилиндра, это иллюстрируется рисунками 3 и 4. На рисунке 3 показаны изолинии завихренности
(rot V) около движущегося цилиндра при t = 10. В этот момент времени течение почти симметричное и стационарное. С увеличением времени нарастает неустойчивость течения жидкости и вихри сходят с цилиндра поочередно с каждой стороны (рис. 4).
Рис. 3. Поле завихренности около движущегося цилиндра (t = 10).
Образуется хорошо известная «дорожка Кармана». Несимметричные отрывы вихрей приводят к сильным колебаниям поперечной силы, как это видно на рисунке 2.
Вихри отстают от цилиндра и постепенно диффундируют. Средняя величина коэффициента лобового сопротивления Сд = т^ / д ^ хорошо согласуется с известными
экспериментальными и расчетными данными [8, 9]. При Ее = 200 Си = 1-25, для Ее = = 2000 Си = 1-03- Также согласуется с известными данными и рассчитанное число Стру-хала 5 = 0-18- Со временем режим течения стабилизируется, о чем свидетельствует выход на предельную циклическую фазовую диаграмму представленную на рисунке 5.
Циклическая диаграмма имеет период, соответствующий числу Струхала.
5. Результаты расчетов
Изучалось движение эллиптического цилиндра с полуосями а, Ь- Сначала рассмотрим движение круглого цилиндра (а =1, Ь =1) при движении внутренней точки массой т =
Рис. 5. Фазовая диаграмма движения цилиндра в жидкости.
= 0.5М. Траектория движения точки с постоянной скоростью приведена на рисунке 6 (контур 1).
Рис. 6. Контуры движения внутренней материальной точки.
Радиус окружности гт = 0.56. При нулевых силах трения = ^2 = С = 0) траектория движения тела является замкнутой, как это и показано в работе [2]. При учете вязкого характера движения жидкости траектория движения тела существенно изменяется. На рисунке 7 показана траектория движения цилиндра при Ее = 2000 в сравнении с траекторией в идеальной жидкости. В левой части рисунка 7 приведены траектории тела в идеальной жидкости (и = 0) ив вязкой (Ее = 2000) в окрестности точки старта, а в правой части показана траектория движения круглого цилиндра (А1 = А2) в вязкой жидкости до времени г = 500.
При достаточно большом времени движения траектория тела стремится к некоторой замкнутой кривой. О выходе на предельный цикл также свидетельствует фазовая диаграм-
— Ев =20 ООО
— 0
0.2 0.4
Рис. 7. Траектории движения цилиндра в вязкой и идеальной жидкости.
ма (Р1, Р2), приведенная на рисунке 8. Амплитуды импульсов Р1, Р2 сначала увеличиваются из начального состояния Р1 = 0, Р2 = 0, но затем стабилизируются в предельном цикле. Предельный цикл также устанавливается для составляющих силы сопротивления (рис. 8).
0.4 0.3 -0.2 -0.1
о н
д
-0.1
жДИк
«1
ди
V/ о уу
0.05 -
0.4 0.3 0.2 0.1
0
Р1
-0.1
-0.05 -
-0.1
Рис. 8. Фазовые траектории для импульсов (Р1 ,Р2) и сил сопротивления (^1,^2).
Расчеты для прямоугольного контура движения внутренней точки (контур 2 на рис. 6) также показали переход траектории движения кругового цилиндра к замкнутой кривой (рис. 9).
Для переменных (Р1 ,Р2), ) также устанавливаются предельные циклы. Величи-
на вязкости или числа Рейнольдса изменяют траекторию движения тела. Конфигурации предельных циклов показаны на рисунке 10 при разных значениях числа Рейнольдса. Предельное состояние достигается раньше для больших значений вязкости (низкое число Рейнольдса). С увеличением Ее время достижения предельного значения увеличивается. Как следует из рисунка 10, для Ее = 2000 предельное состояние при Ь = 500, строго говоря, еще не достигнуто, а достигается при времени Ь = 850- В идеальной жидкости траектория движения цилиндра является замкнутой в окрестности точки начала движения (х = 0,у = = 0), в вязкой жидкости тело удаляется от начала координат. Тем не менее, при некотором
У 1 т
О
-1-2
-3
-4-56-1 О 1 2 3 4 5
Рис. 9. Траектории движения кругового цилиндра с прямоугольным контуром.
удалении от точки старта траектория все равно становится замкнутой. Для дальнейшего движения тело с круговой симметрией следует привести в состояние покоя, а затем повторить движение до выхода на предельный цикл.
Р
— Ев = 200
— Ев= 2 ООО
— Ев= 20 ООО
"'"О 0.1 0.2 0.3 0.4 1
Рис. 10. Влияние числа Рейнольдса на движение кругового цилиндра.
Тело при поворотах попадает в возмущенную область жидкости, что также влияет на его траекторию, если характерное время затухания движения в жидкости сравнимо с периодом движения материальной точки. Проведены расчеты для круглого цилиндра, двигавшегося в течение времени от г = 0 до г = 8 по траектории, как на рисунке 9. Это время соответствует одному циклу движения внутренней точки. Затем тело останавливается и жидкость успокаивается. Максимальная скорость движения жидкости в окрестности тела уменьшается в несколько раз за характерное время Т ~ 2. Это означает, что движущееся тело при изменении направления движения находится все время в возмущенном потоке.
При движении круг совершает периодическое вращение. Изменение угловой скорости от времени приведено на рисунке 11. Здесь также представлена угловая скорость для случая идеальной жидкости (зависимость 0). Достаточно большая вязкость (Ее = 200) не оказы-
вает существенного влияния на вращательное движение кругового цилиндра. Тем не менее, наличие вязкости при поступательном и вращательном движении кругового цилиндра приводит к появлению силы, существенно изменяющей траекторию движения тела.
Далее рассмотрим движение эллиптического цилиндра с соотношением осей 5:1. В идеальной жидкости тело совершает движение со смещением, как в работе [2]. Учет вязкости также сильно изменяет траекторию движения (рис. 12).
У
Рис. 12. Траектории движения эллиптического цилиндра.
Я,
-0.2
Рг
Рис. 13. Фазовые диаграммы для движущегося эллипса.
В идеальной жидкости общее направление линейное. В вязкой жидкости общее направление движения также стремится к прямолинейному, но, как видим, эллипс в вязкой жидкости движется в другую сторону. Угловая скорость и углы поворота эллипса невелики,
т. е. он ориентирован примерно по оси х- Предельный цикл для импульсов наступает значительно позднее, чем для плохо обтекаемого кругового цилиндра (рис. 13 слева). Для силы сопротивления предельный цикл наступает быстрее (рис. 13 справа). Движение внутренней материальной точки по прямоугольному контуру 2 сохраняет эти особенности движения, только предельное состояние наступает немного быстрее. Это связано с тем, что сила сопротивления при движении внутренней материальной точки по прямоугольному контуру увеличивается за счет больших ускорений при смене движения.
Заключение
Анализ результатов проведенных численных оценок показывает существенное влияние сил и момента вязкого сопротивления на траекторию движения тел с изменяемой внутренней геометрией масс в вязкой жидкости. В частности, самопродвижение тела с одинаковыми присоединенными массами (круговой цилиндр), невозможное в идеальной жидкости, оказывается (во многом благодаря несимметричному сходу вихрей с поверхности тела) вполне осуществимым в вязкой жидкости. Существенной отличительной чертой движения тела в вязкой жидкости по сравнению с идеальной является то, что скорость частиц жидкости определяется теперь не только мгновенным распределением скоростей на поверхности тела, но и еще всей «предысторией» движения. В результате тело, попадая в возмущенную на ранних этапах движения область жидкости, демонстрирует поведение (выход на предельный цикл), хорошо описываемое теорией конечномерных динамическими систем, что указывает на применимость различных феноменологических моделей (например, задание трения с помощью функции Рэлея) для описания движения тела. Это позволило бы широко применять в изучении данной задачи не только численные, но и аналитические методы. Рассмотренная в статье модель может послужить основой для использования принципиально новых способов самопродвижения при создании подводных роботов.
Авторы благодарят А.В.Борисова за постановку задачи и полезные дискуссии.
Список литературы
[1] Козлов В. В., Рамоданов С.М. О движении в идеальной жидкости тела с жесткой оболочкой и меняющейся геометрией масс // Докл. РАН, 2002, т. 382, №4, с. 478-481.
[2] Козлов В. В., Онищенко Д. А. О движении в идеальной жидкости тела, содержащего внутри себя подвижную сосредоточенную массу // ПММ, 2003, т. 67, №4, с. 620-633.
[3] Кочин Н. Е., Кибель И. А., Розе Н. В. Теоретическая гидромеханика: Ч. 1. М.: Физматлит, 1963.
[4] Mougin G., Magnaudet J. The generalized Kirchhof equations and their application to the interaction between a rigid body and an arbitrary time-dependent viscous flow // IJMF, 2002, vol. 28, pp. 1837-
[5] Howe M. S. On the force and moment on a body in an incompressible fluid, with application to rigid bodies and bubbles at high and low Reynolds numbers // Quart. J. Mech. Appl. Math., 1995, vol. 48, pp.401-426.
[6] Biesheuvel A., Hagmeijer R. On the force on a body moving in a fluid // Fluid Dynam. Res., 2006, vol. 38, pp. 716-742.
[7] Бендерский Б. Я., Тененев В. А. Пространственные дозвуковые течения в областях со сложной геометрией // Матем. моделирование, 2001, т. 13, №8, с. 121-127.
583 с.
1851.
[8] Mittal R., Dong H., Bozkurttas M., Najjar F. M., Vargas A., vonLoebbecke A. A versatile sharp interface immersed boundary method for incompressible flows with complex boundaries //J. Comput. Phys., 2008, vol.227, pp. 4825-4852.
[9] Шлихтинг Г. Теория пограничного слоя. М.: Наука, 1974. 711с.
Motion of a body with variable distribution of mass in a boundless viscous liquid
Sergey M. Ramodanov1, Valentin A.Tenenev2
1 Institute of Computer Science, Udmurt State University Universitetskaya 1, Izhevsk, 426034 Russia
2 Izhevsk State Technical University Studencheskaya 7, Izhevsk, 426069 Russia [email protected], [email protected]
In the paper we consider the motion of a rigid body in a boundless volume of liquid. The body is set in motion by redistribution of internal masses. The mathematical model employs the equations of motion for the rigid body coupled with the hydrodynamic Navier-Stokes equations. The problem is mostly dealt with numerically. Simulations have revealed that the body's trajectory is strongly governed by viscous effects.
Keywords: self-propulsion, Navier-Stokes equations, viscous vortical motion, numerical methods
Received September 1, 2011, accepted September 24, 2011
Citation: Rus. J. Nonlin. Dyn., 2011, vol. 7, no. 3, pp. 635-647 (Russian)
MSC 2010: 76D05, 35Q30, 65N06