УДК 536.25:519.6
ЧИСЛЕННАЯ МОДЕЛЬ КРУТИЛЬНОГО ВИСКОЗИМЕТРА, ЗАПОЛНЕННОГО ньютоновской ЖИДКОСТЬЮ
А.Е. Коренченко., O.A. Головня, ВЛ. Бескачко
Численными методами решается задача о движении крутильного вискозиметра, заполненного ньютоновской жидкостью, за пределами приближений, принятых в стандартных аналитических теориях. Расчеты проведены для случая бесконечного цилиндра и для осесимметричных течений в конечном цилиндре. Проведено сравнение с результатами аналитических вычислений.
Введение. Крутильный вискозиметр является широко распространенным инструментом для изучения вязкости жидкостей, в особенности, агрессивных, например, металлических расплавов. Он представляет собой прямой круговой цилиндр, заполненный исследуемой жидкостью и подвешенный вдоль своей геометрической оси на упругой нити. Непосредственно наблюдаемыми величинами являются период Т и декремент 3 режима установившихся затухающих колебаний. Эти две величины зависят от выбора параметров установки, а также от свойств заполняющей цилиндр жидкости - ее (кинематической) вязкости у и плотности р. Задача вискозиметрической теории заключается в установлении связи между наблюдаемыми параметрами колебаний (Т и 8), параметрами установки и параметрами жидкости (у и р). С формально-математической точки зрения она представляет собой сопряженную задачу о движении твердого тела и заполняющей его жидкости, точное решение которой в конечном аналитическом виде в настоящее время не найдено, да и вряд ли существует вообще. Известен лишь ряд приближенных решений, справедливых при некоторых предположениях о характере течения жидкости, первое из которых было найдено в [1]. Практически все имеющиеся сегодня экспериментальные данные обработаны с использованием методики, предложенной в этой работе или подобных ей и отличающихся друг от друга вычислительными деталями, но не исходными положениями. Накопившиеся, однако, противоречия в вискозиметрических данных делают привлекательной попытку решения той же задачи за пределами принятых приближений, полагаясь только на фундаментальные принципы механики. Сделать это в настоящее время можно только численными методами, что и составляет цель настоящей работы, первые результаты которой были уже представлены в [2].
Математическая модель. Пусть цилиндр с замкнутыми сверху и снизу торцами полностью заполнен исследуемой жидкостью, которая полагается ньютоновской. Уравнение его движения при условии пренебрежения затуханием колебаний, связанным с процессами внутреннего трения в нити и трением о воздух, имеет вид
/ £ ц л
Здесь /ц, / и <р- момент инерции, угловая скорость и угол поворота цилиндра, к - крутильная жесткость нити, Мтр - момент вязких сил, действующий на цилиндр со стороны жидкости:
Мтр=Мх+М2,
обусловлены трением жидкости о боковую поверхность цилиндра и его
= -к<р +Мтр
(1)
где моменты Мх и М2 торцы соответственно:
я
Мх =2м?И2 Í!
oV '
К
Mr
2 nr¡ j
dV„
д
<p_
r y
Rf
dz,
r=R
¿V'
r dr.
rгdr--2л;r|| -г=0 (Л- - J
Здесь Уу - - азимутальная компонента поля скорости V = (УГ,У<Р,У2) (введена цилин-
дрическая система координат), Н и Л - высота и радиус цилиндра. В работе ограничимся изуче-
нием только осесимметричных режимов течения, когда ни одна из интересующих нас переменных не зависит от азимутального угла ср .
Движение цилиндра возбуждает течение заполняющей его жидкости, подчиняющееся системе уравнений Навье - Стокса и неразрывности. Введем следующие безразмерные переменные и
параметры: * = —, у = - компоненты радиус-вектора, ~ безразмерное время, х = (и, у, IV)
Я Я Я
- компоненты вектора скорости, у = V—, р - возмущение давления сверх pgz, отнесенное к
V
С = ^ К - безразмерный коэффициент жесткости нити на кручение, тпг- момент вязких
сил. Тогда движение цилиндра и заполняющей его жидкости описывается следующей системой уравнений:
V ^
ат
^ + + (2)
дт
1 = 0
у ду дх
Граничными условиями для (2) служат условия прилипания жидкости ко всем твердым границам:
w(Ux,r) = fy
w(y, 0, т) = w(y, Я / Л, г) = / -у, (За)
v(y, 0, т) = и(у, 0, r) = v(y, H/R, т) = и(^ H/R, r) = v( 1, г) = и(1, г) = 0. Кроме того, на оси сосуда выполняются условия
—- —= 0, v = 0, w = (36)
ду ду
обусловленные осевой симметрией течения. Жидкость в начальный момент времени покоится, а цилиндр удерживается в состоянии покоя в положении, повернутом на угол <р0 относительно положения равновесия. В соответствии со сказанным выше начальными условиями для уравнений (2) являются:
v(*,jO = 0, ДО)-0, р(0) = (р0. (4)
Задача заключается в отыскании закона движения цилиндра <р(т).
Численное решение. Численное решение системы (2) с граничными условиями (3) и начальными условиями (4) находилось методом конечных разностей. Использовались равномерные пространственные сетки с максимальной размерностью 36x36 в осевом и радиальном направлениях. Дискретизация гидродинамических уравнений производилась по схеме центральных разностей с точностью (Ах)2 по пространственным переменным. Линеаризация получаемых в результате нелинейных разностных уравнений проводилась методом Ньютона, а решение линеаризованных систем уравнений - методом исключения Гаусса. Для контроля консервативных свойств решения проверялось равенство нулю суммарного потока в каждом поперечном сечении цилиндра или замкнутость линий тока в осевом сечении в любой момент времени. Расчеты проводились для цилиндров с размерами, обычно используемыми на практике: внутренний радиус изменялся в пределах от 1 до 3 см, а высота - от 1 до 7 см. Отношение моментов инерции цилиндра и «замороженной» жидкости варьировалось в диапазоне от 2 до 20, а вязкость исследуемой жидкости
выбиралась в интервале (КГ7 м2/с < v <8-10"6м2/с), свойственном легкотекучим жидкостям (вода, металлические расплавы).
Коренченко А £., Головня О.А., Бескачко В П.
Численная модель крутильного вискозиметра, заполненного ньютоновской жидкостью
Результаты и обсуждение. Вернемся к приближениям, используемым в стандартных виско-зиметрических теориях и оценим их справедливость, располагая численным решением рассматриваемой задачи
1 Предположение о малости радиальной и осевой компонент скорости жидкости по сравнению с азимутальной не всегда верно для цилиндра конечной длины В работе [2] обнаружено, что в осевой плоскости вискозиметра образуются течения, длящиеся во все время существования крутильных колебаний Однако вопрос об абсолютных величинах осевой и радиальной компонент в сравнении с азимутальной в [2] не обсуждался На рис 1 показано поле скорости в осевой плоскости для цилиндра с Я^Н =1 см, заполненного жидкостью с вязкостью 5 10~7 м2с 1 Видно образование четырех вихрей, причем наибольшая осевая скорость тбпю^тся вблизи оси цилиндра На рис 2 приведены временные зависимости осевой скорости в жидкости для отношений у=Н/К~1 и у =2, и указано амплитудное значение азимутальной скорости Из рисунка видно, что абсолютное значение осевой скорости увеличивается при уменьшении отношения у, а отношение наибольших значений осевой и азимутальной скоростей для у =7 превышает 0,1 Это не дает достаточных оснований для пренебрежения нелинейными слагаемыми в уравнениях (2), во всяком случае, на начальном этапе колебаний
02п
» I г / /
V, \ \ \ \ Н I / ( / / / ^ ^ » V .
N > М II У\ I \ \ / > , % ч
/ / / / / И \ X М \ N „ ,
/ ✓ / / | \ ч V •ч^,-..____
—- ч \ N ч I ( ( у -
" м | I I I I ^ V_
- у // Г ) I 1 Ч X V -
00 -02 -04 -06 -08 1 0 -12 14 •16
*иа12(6е)
Я/Л-2
«V
, , Чя/ди
30
40
ко
Рис 1 Распределение скоростей в осевой плоскости цилиндра Наибольшая осевая скорость отвечает точке А (0 8Я, 0)
Рис 2 Осевая скорость в области наиболее интенсивного движения, / / /= 20
ц/ ж
'7м2/с
2 Для разделения переменных в уравнениях (2) используется предположение о том, что жидкость внутри цилиндра совершает колебательные движения около положения равновесия с частотой, равной частоте колебаний цилиндра и фазой, определяемой скоростью передачи момента импульса к внутренним слоям Это предположение выполняется не всегда На рис 3 изображены зависимости от времени азимутальной скорости на небольших расстояниях (г=Я/4) от оси вращения для жидкостей с различной вязкостью Условия подобраны так, чтобы частота колебаний была одинаковой для всех трех рассмотренных случаев Как видно из рисунка, предположение о гармоническом характере изменения азимутальной скорости выполняется лишь асимптотически и тем позднее, чем
-0 2-
Рис. 3. Зависимость азимутальной скорости вблизи оси от времени, /ц //ж — 2
меньше вязкость жидкостей. Чем дальше от оси (и ближе к стенкам), тем быстрее устанавливается режим затухающих колебаний.
Сравним, наконец, предсказания аналитических и численных расчетов параметров колебаний вискозиметра. Для этого, выбрав некоторые параметры установки и зафиксировав плотность жидкости, мы выполнили серию расчетов зависимости ср(т) для значений вязкости, типичных для низковязких жидкостей. Далее эти данные были подогнаны под зависимость вида
(5)
0011 л 0 010 -0,009 0,008 0 007-0 006-0,005 -0,004 -0 003-0,002 -
0,001
Численный расчет,
/ц//ж=20; H/R=3
<р(т)~А-е 5тsm[cor + \f/)
Численный расчет,
/ц//ж=20; Я/Д=7
А А А А А
iili*1*
. .........
* д *А . И:,«11"мик::::::
А ,|1Я" /
Аналитический расчет
/ц//ж=20; Я/Я-З
Численный расчет, /ц/^^20; бесконечный цилиндр
0,000000
0,000002
0,000004
0,000006
ГП \ Ч-Ч" I Г 0,000008
0,07-
0,06
0,05 4
0 04-
0,03
0,02-
0,01 0,000000
а)
Численные расчеты /ц/7^2; бесконечный цилиндр
А A *V * * * А А 4 .......
v(m /с)
Численные расчеты /Ц//Ж=2;Я«=7
--------
Аналитические расчеты
I "1"1 | Г "Г I—I VI I | -|—|—р-Т-Т
0,000002 0,000004
1 1 1 I ' ' • 0,000006
1—1—г—1—ГТ—J—Г"П 0,000008
v(M"/C)
б)
с целью определения частоты колебаний со и коэффициента затухания 6. Подгонка осуществлялась методом наименьших квадратов с минимизацией методом Розенброка [3]. На рис. 4а, б изображены зависимости коэффициента затухания от вязкости жидкости. Квадратами во всех случаях помечена аналитическая зависимость, полученная в [1], треугольниками показаны результаты численных расчетов. Видно, что аналитические и численные предсказания заметно отличаются друг от друга. Согласие между ними, видимо, достигается лишь в пределе малой вязкости, а с ростом вязкости становится тем хуже, чем больше вязкость. Из рис. 4а видно также, что при фиксированном значении у? - /ц //ж
(равном 20 в данном случае) согласие тем лучше, чем «длиннее» цилиндр - чем больше параметр у = Н/И. При у —^ оо аналитические и численные результаты практически повторяют друг друга, а при конечном у располагаются по разные стороны от этой «кривой согласия» (отвечающей у = со) и отстоят от нее тем дальше, чем меньше у. Таким
Рис. 4. Зависимость коэффициента затухания колебаний от вязкости жидкости; а) р = 20 , б) /3 = 2
у. 1 аким образом,
предсказания зависимости коэффициента затухания от выбора параметра ^ при фиксированном значении вязкости жидкости в аналитических и численных моделях прямо противоположны друг другу: первая из них предсказывает уменьшение д с уменьшением у, а вторая - его увеличение. Наконец, при переходе в численной модели от у -1 к у - Ъ максимум на зависимости 5'(у) пропадает: либо перемещается за пределы исследованного интервала вязкостей, либо исчезает вовсе. Последняя возможность противоречит, однако, естественному предположению [1] о том, что функция 5{у) должна уменьшаться при больших к. При переходе от /? = 20 к /3 = 2 (рис. 46) ситуация не становится лучше. Здесь дело осложняется тем, что предположение анали-
Коренченко А.Е., Головня О.А., Бескачко В.П.
Численная модель крутильного вискозиметра, _заполненного ньютоновской жидкостью
тической теории о характере движения внутренних слоев жидкости (затухающие колебания) выполняется тем хуже, чем меньше отношение ¡5. При малых ¡3 внутренние области жидкости могут практически не участвовать в колебаниях (см. рис. 3), что может быть одной из причин расхождения аналитических и численных данных.
На данный момент трудно сказать, что является причиной несоответствия в предсказаниях аналитических и численных моделей. С одной стороны, аналитические модели основываются на приближениях, которые, как мы видели, оправдываются не безусловно. С другой стороны, при применении численных моделей встречаются свои подводные камни, препятствующие получению результатов необходимой точности. Среди них, помимо обычно обсуждаемых (таких как схема дискретизации, параметры сетки и т.п.), отметим специфические для данной задачи:
1) необходимость изучения асимптотических решений, справедливых при временах, достаточно удаленных от момента возбуждения колебаний,
2) необходимость получения максимально точной картины движения жидкости в окрестности твердых границ - в связи с расчетом вязких сил, действующих на твердое тело.
Преодоление каждой из отмеченных трудностей требует применения специальных вычислительных приемов и средств, которые не были использованы в настоящей работе в полной мере. Сейчас не ясно, насколько эти приемы и средства способны изменить описанные выше результаты. Для решения этого вопроса необходимы дальнейшие исследования.
Работа выполнена при поддержке РФФИ-Урал (№ 01-01-96424).
Литература
1. Швидковский Е.Г. Некоторые вопросы вязкости расплавленных металлов. - М.: ГИТТЛ, 1955.-206 с.
2. Коренченко А.Е., Бескачко В.П. Особенности установления колебаний в крутильном вискозиметре // Вестник ЮУрГУ. Серия «Математика, физика, химия». - 2002. - № 3. - Вып. 2. -С. 62-67.
3. Химмельблау Д. Прикладное нелинейное программирование. - М.: Мир, 1975- 534 с.
Поступила в редакцию 25 апреля 2003 г.