CC BY
65
ВЕСТНИК ТГГПУ. 2010. №4(22)
ФИЗИКО-МАТЕМАТИЧЕСКИЕ НАУКИ
УДК 532.5.296:534.12
РАСЧЕТ СИЛОВОГО ВОЗДЕЙСТВИЯ КАВИТАЦИОННОГО ПУЗЫРЬКА НА УПРУГОЕ ТЕЛО
© А.А.Аганин, В.Г.Малахов, Т.Ф.Халитова, Н.А.Хисматуллина
В статье рассматривается задача об осесимметричном силовом воздействии на упругое тело высокоскоростной струи, образующейся на поверхности кавитационного пузырька в процессе его схлопывания вблизи тела. Тело моделируется изотропным упругим полупространством, струя -цилиндрическим столбом жидкости. Волны в жидкости описываются по модели линейной акустики, волны в полупространстве - по модели линейно-упругого тела. Решение находится методом С.К.Годунова первого порядка точности. Проведено исследование возможностей сокращения расчетной области при использовании в качестве краевых условий на внешних искусственных границах в области упругого полупространства условий свободного края, неподвижной абсолютно твердой стенки и неотражающих условий. Установлено, что неотражающие условия более предпочтительны.
Ключевые слова: кавитационное разрушение, кавитационная эрозия, численное моделирование гидроупругого взаимодействия, ударное воздействие струи.
1. Введение
Термин кавитация был введен в 1894 году британским инженером Р.Фрудом для описания процесса образования газовых пузырьков в жидкости. Как правило, кавитационные пузырьки появляются при понижении давления жидкости ниже уровня давления ее насыщенных паров. Если в последующем давление жидкости становится выше давления насыщенных паров, то кавитационные пузырьки с силой схлопываются. Схлопывание кавитационных пузырьков создает шум, вызывает вибрацию, неблагоприятно отражается на функционировании различных работающих с жидкостями машин и механизмов (гидротурбин, жидкостных насосов, гребных винтов кораблей, подводных звукоизлучателей, жидкостных систем высотных самолётов и т.д.), снижая их коэффициент полезного действия, вызывая в них повреждения или даже разрушения.
Для понимания природы кавитационного разрушения необходимо изучение особенностей воздействия отдельного пузырька на тело. Этому вопросу посвящено много работ (например, [16]). Несмотря на это, режим ударного воздействия отдельных пузырьков, при котором в финальной стадии сжатия пузырька в жидкости возникают ударные волны, остается малоизученным. Растущие возможности экспериментального оборудования и вычислительной техники позволяют приступить к его детальному исследованию.
В настоящее время считается, что одним из основных механизмов силового воздействия пузырьков на тело, приводящих к кавитационному
разрушению, является высокоскоростная струя [2], образующаяся на поверхности пузырька в процессе его схлопывания вблизи стенки. Образование струи сопровождается сильными деформациями поверхности пузырька. Для изучения динамики пузырьков с большими деформациями его поверхности применяют методы конечных разностей и конечных элементов. Эти методы позволяют достаточно полно учесть эффекты вязкости жидкости, ее теплообмена с пузырьком и т. д.
В настоящей работе предлагается методика численного исследования осесимметричного воздействия струи жидкости, возникающей при схлопывании кавитационного пузырька, на линейно-упругое изотропное тело. Струя моделируется цилиндрическим столбом жидкости. Волны в жидкости описываются по модели линейной акустики, волны в полупространстве - по модели линейно-упругого тела. Получающаяся при таком подходе система линейных уравнений первого порядка в частных производных решается методом С.К.Годунова первого порядка точности. Приводятся результаты исследования возможностей сокращения расчетной области при использовании в качестве краевых условий на внешних искусственных границах в области упругого полупространства условий свободного края, неподвижной абсолютно твердой стенки и неотражающих условий.
2. Математическая постановка задачи
Рассматривается ударное воздействие на тело высокоскоростной струи, образующейся на поверхности пузырька при его схлопывании у
стенки тела. Задача полагается осесимметричной, ось симметрии г направлена перпендикулярно плоской поверхности тела. Область, занимаемая струей, принимается в виде прямого кругового цилиндра. Такую форму струи можно получить в результате упрощения формы пузырька по схеме, приведенной на рис.1.
Линейные размеры кавитационных пузырьков обычно намного меньше размеров тела. С учетом этого рассматриваемое тело считается упругим полупространством. Область, занимаемая струей и телом, представлена на рис.2а, расчетная область для численного решения задачи -на рис.2 Ъ.
а))
Рис.1. Упрощение формы струи (слева направо)
Динамика упругого тела при осевой симметрии описывается уравнениями [7]
ди 1
дґ р
да.
даг дтг
дг дг
/о , ч ду I ди и
-(Л + 2м)~ + Ц~дГ + ~
дґ дг \дг г
ду _ 1 I даг дґ р\ дг
дтгг
дг
(1)
дґ
ди ду
даш (ди ду\ „ ч и
—-_ЛI — + — \ + (Л + 2и)~ ;
дґ {дг дг \ v ’ г
а динамика жидкости - уравнениями
ди 1 др ду 1 др
— +---— _ 0, — +---------— _ 0,
дґ р0 дг дґ р0 дг
1 др ди и ду
------2---------1-----1-----1----
р0с0 дґ дг г дг
(2)
_ 0.
Здесь г, г - цилиндрические координаты точки среды; t - время; и, V - компоненты скорости по осям г и г соответственно; аг, с2, а(р, тГ2 - компоненты напряжений в упругом теле; X = р (с2 - 2с2 ), /и = рс\ - параметры Ламе в упругом теле, где с\ и с2 - продольная и сдвиговая скорости звука в упругом теле; р - невозмущенная плотность упругого тела; р - давление жидкости; р0, с0 - невозмущенные значения плотности и скорости звука в жидкости.
Ь)
£
/?
тело
+ — а г _-р
у + = у “
і Л,
Рис.2. Область, занимаемая струей и телом (а), и расчетная область (б).
За начальный момент времени принимается момент касания струи поверхности тела.
Условия на границе раздела фаз жид-кость-твердое тело (г = 0, 0 < г < Я) имеют вид
5 + = 0 (3)
Здесь знак плюс указывает, что соответствующий параметр относится к жидкости, а минус - к телу.
На оси симметрии г = 0 ставятся условия и = 0, тгг = 0 при -да < г < 0 (область тела) и = 0, при 0 < г < да (область жидкости).
На внешней границе струи (г = Я, 0 < г < да) имеем р = р0, где р0 - невозмущенное давление в жидкости.
На свободной поверхности тела (Я < г < да, г
= 0) 02 = - р0, Тгг = 0.
На бесконечно удаленных границах тело и жидкость считаются невозмущенными.
В начальный момент времени t = 0 давление в струе однородно и равно р0, компоненты скорости равны V = V) и и = 0; в упругом цилиндре деформации в радиальном и окружном направлениях отсутствуют, при этом реализуется статическое напряженное состояние
.0 „ _0____0 _ 0 , о .Л-1 „ 0
с =-P0, с = ар =-х(х + 2и) Po, т ^
и0 = V0 = 0.
Для решения задачи (1)-(5) применяется разностная схема Годунова первого порядка точности [8].
С целью оценки напряженного состояния в упругом теле используется величина интенсивности напряжений ог,
с, =^2^ (С 2 -Сг )2 +(С 2 -Ср)2 +(Сг -Ср)2 + 6X22 .
При исследовании напряженно-деформированного состояния металлических упругих тел экстремальные значения интенсивности напряжений аг и их положение в теле представляют обычно наибольший интерес, поскольку величина ог- применяется во многих критериях разрушения.
Для анализа максимальных по времени возмущений интенсивности напряжений используется безразмерная величина
>-тах / \
С i (г, г ) =
сг -с0
тах ' а* -а0
где oi - начальное значение интенсивности напряжений в теле, ai - максимальное значение этого параметра в одномерной задаче об ударе струи бесконечного радиуса по упругому полупространству.
3. Результаты расчетов При расчетах были взяты следующие данные: р = 8-103 кг/м3, р0 = 103 кг/м3, с0 = 1500 м/с, С = 5742.9 м/с, с2 = 3069.7 м/с, р0 = 0,1 МПа, Я = 20 мкм. Они соответствуют удару водяной струи по стальному полупространству. Расчетная область представлена на рис. 2Ь. Струя занимает цилиндрическую область 0 < г < Я, 0 < * < hj , а
упругое тело - область 0 < г < г8, К < г < 0. Внешние границы расчетной области г = hj,
г = г и г = К являются "искусственными". Они выбираются столь далекими, чтобы отраженные от них волны не оказывали заметного влияния на величину экстремальных напряжений в теле, которые возникают в окрестности оси симметрии вблизи контактной поверхности. Удаленность искусственных границ зависит от условий, принятых на этих границах.
Шаг по времени определяется из условия устойчивости [9]
5 Аг • Аг
1 =-----1----Г- ,
с1 Аг + Аг
где Аг, Аг - минимальные шаги сетки по радиальной и осевой координате, 5 = 0,95.
Расчеты показывают, что если положение искусственных границ задать равенствами г8 = ~К = 5Я, К = Я, то решение рассматриваемой задачи в области, где достигаются максимальные значения сжимающих напряжений в теле, не будет зависеть от условий на этих границах. В настоящей работе при изучении возможностей сокращения размеров расчетной области данное численное решение используется в качестве эталонного.
Распределение максимальных по времени
возмущений интенсивности напряжений ^тах в
эталонном решении показано на рис.3. Видно, что наибольшие значения возмущения интенсивности напряжений достигаются в сравнительно небольшой области на оси симметрии вблизи поверхности контакта (г/Я ~ 0,75).
Рис.3. Изолинии Стах - максимальных по времени возмущений интенсивности напряжений. Приведено 15 изолиний (для Стах =0,2^2,8)
Для получения результата удовлетворительной точности с использованием схемы первого порядка требуются довольно мелкие сетки. Расчет на таких сетках связан с большим объемом оперативной памяти и процессорного времени. Численный эксперимент показал, что применение неравномерных сеток (с увеличением размеров ячеек в направлении внешних границ) несколько повышает эффективность алгоритма. Однако более существенного сокращения процессорного времени удается добиться за счет уменьшения размеров расчетной области при соответствующем выборе граничных условий на искусственных границах. При этом наибольший эффект достигается при уменьшении расчетной области в той ее части, которая относится к телу, поскольку в области тела приходится применять более мелкую, чем в области жидкости, сетку, что обусловлено различием в скоростях распространения возмущений в теле и жидкости. В силу этого в настоящей работе изучаются возможности сокращения расчетной области только за счет приближения ее внешних искусственных границ г = г , г = К , ограничивающих область тела. В качестве краевых условий на этих границах рассматриваются следующие: и = 0 при г = г ; V = 0 при г = К ; (4)
аг = а°, Тг = 0 при г = ГS ; (5)
0 (5)
гг = аг , тгг = 0 при г = К рс1и + аг = рс1и° + г ,
рC2V + тгг = Pc2v 0 + Т ^и г = ^
р^ + а2 = рс° + а г,
(6)
(7)
рсги + Тгг = рсги + ТГг при г = К .
Условия (4) соответствуют жесткой стенке, условия (5) - свободному краю. Условия (6), (7) представляют собой один из простейших вариан-
тов неотражающих условий для двумерных задач Положение искусственных границ варьиру-
[8]. ется следующим образом:
а)
В
і-.Г
О
-3
гт
-і
о -з
г/Л
-і
г„ = -К = ЪЯ,2Я-, 1,75/г.
Л Ь) л с) ЗА -/ '/Л ! \ I \
3 / V, 2У у\ / \ V. Ф 2 \£у ч- ^ У/ \
1 1 1 1 1 Г т~ -Г- - Г-—т < Т | Г , Г
о
-3
2ІК
-1
о
Рис.4. Распределение максимальных по времени возмущений интенсивности напряжений на оси симмет-
рии. Результаты получены на расчетных областях с искусственными границами г = - кх = 5Я (кривые 0, эталонное решение); 3Я (кривые 1); 2Я (кривые 2); 1,75Я (кривые 3) при использовании на этих границах (а) неотражающих условий (6), (7), (Ь) условий свободного края (5) и (с) жесткой стенки (4)
Время расчета при таком сокращении расчетной области по сравнению с затратами на эталонное решение уменьшается соответственно в 4,3; 6,2 и 10,1 раз. Естественно, что эти числа практически не зависят от того, какое из трех перечисленных в (4)-(7) условий применяется, поскольку в общем объеме вычислений их расчет составляет лишь ничтожно малую долю. Результаты, полученные на большой области при г8 = —к3 =5Я и на областях с г = - к = 3Я; 2Я; 1,75Я для трех рассматриваемых видов условий на искусственных границах, показаны на рис.4-7.
На рис.4-5 приведено распределение максимальных по времени значений возмущений интенсивности напряжений на оси симметрии. В интервале —3Я < г < 0 графики функции ^тах(о, г), рассчитанные на большой области с кj = 2,5Я , г3 = —к3 = 5Я (кривые 0), для всех
трех типов краевых условий совпадают (поэтому данное решение и принято за эталонное). Максимальное значение функции ^™ах (0, г) в эталонном решении составляет 2,97 и достигается при г/Я ~ 0,75.
Рис.5. Распределение максимальных по времени возмущений интенсивности напряжений ^Гах на оси симметрии при использовании неотражающих условий: кривая 0 - эталонное решение, кривая 4 соответствует гх = - Их = 1,25Я, а кривая 5 - гх = - ^ = 1,05Я
Из рис.4 следует, что наиболее близкие к эталону результаты получаются с неотражающими условиями (а). Влияние неотражающих условий при дальнейшем уменьшении расчетной области иллюстрирует рис.5. Согласно рис.4, 5, поведение функции <^гтах (0, г) в окрестности максимума (-1 < г/Я < - 0,5) остается неизменным не только при Гц = - кц = 3Я; 2Я; 1,75Я, но и с г = -И = 1,25Я (при этом время расчета уменьшается по сравнению с затратами на эталонное решение в 20 раз). Лишь при гх = - к = 1,05Я (рис.5) значение максимума несколько снижается (на 2%).
Рис.6. Изолинии функции (г,2) на фрагменте расчетной области: (а) - эталонное решение, (Ь)-(й) - решения, полученные на области с гя = - кя = 2Я при использовании неотражающих условий (Ь), условий свободного
края (с) и жесткой стенки (й)
При использовании условий свободного края (рис.4Ь) и жесткой стенки (рис.4с) в окрестности искусственной границы возникают довольно большие отклонения от эталона. Вместе с тем в относительно небольшой окрестности максимума функции ^™ах (0, г) расхождения с эталоном остаются незначительными вплоть до области с Г = - к = 1,75Я. При этом в случае свободной границы в области с гя = - кя = 1,75Я появляются отклонения от эталона вблизи контактной поверхности, а при условии жесткой стенки возникают дополнительные локальные максимумы, так что в области с г я = -кя = 1,75Я локальный максимум, вызванный отражением от искусственных границ, становится больше эталонного. В
определенных случаях такие погрешности могут оказаться неприемлемыми.
На рис.6 и 7 представлены изолинии функции СГаХ(г, г), рассчитанные при г = -кя = 2Я (рис.6) и Гц = - кц = 1,75Я (рис.7). На этих рисунках показан наиболее интересный фрагмент расчетной области тела, демонстрирующий поле интенсивности напряжений в окрестности его максимума. Рис.6, 7 позволяют более наглядно оценить влияние граничных условий при уменьшении размеров области. Видно, что результаты расчета с неотражающими границами (рис.6Ь, 7Ь) практически не отличаются от эталонного (рис.6а, 7а).
Рис.7. Изолинии функции СГХ (г, z) на фрагменте расчетной области: (а) - эталонное решение, (Ь)-(й) - решения, полученные на области с г. = - hs = 1,75Л при использовании неотражающих условий (Л), условий свободного края (с) и жесткой стенки (й)
Условия свободного края при достаточном удалении границ г. > 2Я (рис.бс) также по-
зволяют получить удовлетворительный результат. С уменьшением размеров области до г.. = -Ь = 1,75Я картина изолиний в окрестности искусственных границ при использовании этих условий несколько искажается (рис.7с). Наибольшие отклонения изолиний от эталона получаются при условии жесткой стенки (рис.бй и 7й). Здесь качественные различия с эталоном в окрестности нижней искусственной границы наблюдаются уже при г.! = - hs = 2Я.
4. Заключение
Проведено численное моделирование ударного воздействия на упругое тело высокоскоростной струи жидкости, имеющей форму столбика с плоским концом. Входные данные являются характерными для водяной струи, возникающей на поверхности маленького кавитационного пузырька при его схлопывании около поверхности большого стального тела (скорость струи 100 м/с, радиус 20 мкм; давление в полости пузырька
1 бар).
Показано, что при воздействии водяной струи на стальное тело наиболее напряженная область возникает в окрестности оси симметрии на удалении примерно 0,75 радиуса струи от по-
верхности тела. В этой области максимальные значения интенсивности напряжений оказываются в два раза выше, чем на поверхности тела.
Выполнено численное исследование возможности сокращения размеров расчетной области и повышения за счет этого экономичности вычислений при использовании на искусственных границах неотражающих граничных условий, условий свободного края и жесткой стенки. Наибольшего сокращения линейных размеров расчетной области (более чем в два раза) удалось достичь при использовании неотражающих граничных условий.
Работа выполнена в рамках программы РАН и при поддержке РФФИ (номер проекта 08-01-00215).
1. Rayleigh Lord. On the pressure developed in a liquid on a collapse of a spherical cavity // Phylos Mag. -1917. - Vol.34. - №200. - P.94-97.
2. KornfeldM. and Suvorov N. On the destructive action of cavitation // App. Phys. - 1944. - Vol.15. - P.495-506.
3. Plesset M.S., Chapman R.B. Collapse of an initially spherical vapour cavity in the neighbourhood of a solid boundary // J. Fluid Mech. - 1971. - Vol.47. -P.283-290.
4. Lauterborn W., Bolle H. Experimental investigations of cavitation-bubble collapse in the neighbourhood of a solid boundary // J. Fluid Mech. - 1975. - Vol.72. -P.391-399.
5. Аганин А.А., Малахов В.Г., Халитова Т.Ф., Хис-матуллина Н.А. Силовое воздействие кавитационного пузырька на упругое тело// Проблемы нелинейной механики деформируемого твердого тела: Тр. Второй междунар. конф. (Казань, 8-11 декабря 2009 г.). - Казань: Казан. гос. ун-т, 2009.
- С.20-24.
6. Mettin R., Luther S., et. al. Investigations of cavitation bubble dynamics by means of fast Cinematography // The Intern. Conf. on multiphase systems.
Gilem Publisher, Pol Publisher. - Ufa, 2000. - P.279-287.
7. Чебан В.Г., Навал И.К., Сабодаш П.Ф., Чередниченко Р.А. Численные методы решения задач динамической теории упругости. - Кишинев: Изд-во Штиинца, 197б. - 22б с.
8. Ильгамов М.А., Гильманов А.Н. Неотражающие условия на границах расчетной области. - М.: ФИЗМАТЛИТ, 2003. - 240 с.
9. Годунов С.К., Забродин А.В., Иванов М.Я., Крайко А.Н., Прокопов Г.П. Численное решение многомерных задач газовой динамики. - М.: Наука, 197б. - 400 с.
COMPUTATION OF THE CAVITATION BUBBLE IMPACT ON THE ELASTIC BODY
A.A.Aganin, V.G.Malakhov, T.F.Khalitiva, N.A.Khismatullina
A problem of axially symmetric impact of a high-speed jet on an elastic body is considered in the article.
The jet arises on the surface of the bubble during its collapse near the body. The body is modeled as an elastic semi-space, the jet as a cylindrical column of the liquid. The waves in the liquid are described by linear acoustics, the waves in the semi-space by the linearly elastic body. The problem is computed by Godunov’s method of first order of accuracy. Possibilities of computational domain reduction using the conditions of the free boundary, the non-moving absolutely rigid wall and the non-reflecting conditions as boundary conditions for the artificial boundaries in the area of the elastic semi-space have been investigated. The non-reflecting conditions are found to be more preferable.
Key words: cavitation destruction, cavitation erosion, numerical simulation of hydro-elastic interaction, jet impact.
Аганин Александр Алексеевич - доктор физико-математических наук, профессор кафедры вычислительной математики и информатики Татарского государственного гуманитарнопедагогического университета.
E-mail: aganin@kfti.knc.ru
Малахов Владимир Георгиевич - кандидат физико-математических наук, старший научный сотрудник Института механики и машиностроения Казанского научного центра РАН.
E-mail: malahov@kfti.knc.ru
Халитова Талия Фаритовна - кандидат физико-математических наук, научный сотрудник Института механики и машиностроения Казанского научного центра РАН.
E-mail: taliny@mail.ru
Хисматуллина Наиля Абдулхаевна - кандидат физико-математических наук, доцент кафедры вычислительной математики и информатики Татарского государственного гуманитарно-педагогического университета.
E-mail: nailya_hism@mail.ru
Поступила в редакцию 01.09.2010
