Планирование оптимального эксперимента по одновременному определению вязкости и плотности ньютоновской среды Текст научной статьи по специальности «Физика»

УДК 532.137.3

ПЛАНИРОВАНИЕ ОПТИМАЛЬНОГО ЭКСПЕРИМЕНТА ПО ОДНОВРЕМЕННОМУ ОПРЕДЕЛЕНИЮ ВЯЗКОСТИ И ПЛОТНОСТИ НЬЮТОНОВСКОЙ СРЕДЫ

И.В. Елюхина

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

Введение. Экспериментальные данные, получаемые разными авторами с помощью основного метода исследования вязкостных свойств агрессивных и высокотемпературных сред, а именно: метода крутильных колебаний Швидковского Е.Г. Г1], являются достаточно противоречивыми (см., например, данные по вязкости жидкого железа (рис. 1), полученные в подавляющем большинстве

именно данным методом при декларируемой

77,П 0,07-

0,06-

0,05-

0,04

0,03

1500

1600

1700

Рис. 1. Зависимость вязкости от температуры для жидкого железа 1, 2 ...16 - см. номера литературных источников в [2]

авторами погрешности менее 5 %).

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

В одном из таких методов заключение о получении надежной оценки вязкости жидкости выводится из сравнения известного из независимых источников значений плотности и оцененного при одновременном определении вязкости и плотности среды из вискозиметрического уравнения. Идея метода появилась уже достаточно давно [3, 4], но до сих пор, насколько это известно, не реализована на практике.

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

Р'

2 2 \ Ро + Чо

Ц =-2 уМ/З

Р +4'

2 А

К + 1-

ч-

2 2 Л Ро +4о 2 2 р +ч

К-(Ь+12) = 0

МЮ

Ь2=а

М к

__у МвпН)

н V Ьх М

0)

(2)

где / = л/—Т ; р-З/т - коэффициент затухания колебаний; q - 2к !т - циклическая частота; р0=$0/т0; - 2лг/г0; 8 и т - логарифмический декремент и период колебаний заполненного жидкостью цилиндра соответственно; т0 и 30 - период и логарифмический декремент собственных установившихся затухающих колебаний пустой системы; К - момент инерции всей подвесной системы без жидкости относительно оси цилиндра; М - масса среды; у - кинематическая

вязкость; р = К\[к!у ; к - р + / -д ; Я - внутренний радиус цилиндра; Я - полувысота цилиндра;

вп =//п2-£/у; ¡ип - характеристические числа, определяемые из уравнения -

функция Бесселя первого рода первого порядка; , - функции трения, отражающие роль сил трения на боковой поверхности цилиндра и его торцах соответственно; а = 4 - при учете контакта жидкости с крышкой; а = 2 - в случае свободной поверхности.

Неизвестные параметры - вязкость у и плотность р жидкости - определим методами параметрической идентификации [5] из условия минимума функции качества вида

/(V,р) - 4с, • Ле2+ с2 ■, (3)

где с2 - весовые коэффициенты; Ые, 1т - действительная и мнимая части от функции F.

Детальное обсуждение возможности одновременной оценки вязкости и плотности ньютоновской жидкости из измеряемых в эксперименте параметров крутильных колебаний цилиндра проведено в работе Ньюводта Дж., Сенгерса Дж. и Кестина Дж. [4]. Но полученные результаты являются лишь хорошим начальным приближением в силу многих факторов, в частности:

1) исследования проведены в приближении двух частных случаев большого и малого цилиндров и в этом смысле не являются универсальными;

2) принятая для исследования точность измерения периода и декремента не всегда может быть реализована практически;

3) при определении ошибок измерения плотности и вязкости были приняты во внимание только ошибки прямых измерений периода и декремента затухания подвесной системы;

4) обычно в вискозиметрической практике реализуются значения радиуса цилиндра Я ~ 0,5-1,5 см с целью повышения качества эксперимента, например, для обеспечения изотер-мичности при исследовании металлических расплавов. В [4] данный диапазон значений радиуса является нерабочим при точности измерений Дг~10~5 и Д£~1(Г6, где знак означает, что значения даны в относительных единицах.

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

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

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

локальное направление дна оврага, т.е. кривой, где- (у(/?),/?) =0, в пределах точности численных методов можно считать совпадающим с мгновенным значением производной

* до/ др V , * *\ *

ру = /\А =—---- в точке (V ,/? ), где индекс соответствует значению величины в некоторой

/оу Ъу р исследуемой точке).

Существование оврага связано с видом вискозиметрического уравнения (1), где плотность входит только в торцевые слагаемые в функции трения и не содержится в Ц (2), отражающем роль сил трения на боковых поверхностях вискозиметра. В случае только одного торцевого слагаемого р[, увеличивается приблизительно в два раза, и несколько увеличивается чувствительность оцениваемых характеристик у'х и рх к измеряемым величинам х. Здесь у'х =ду/ ~ — х-г, т0, 3, М, И или К и

/ох дх у

р'х = к р: ■ (4)

Проанализируем изменение чувствительности ру в зависимости от изменения параметров установки и колебаний. Очевидно, что значение р'у зависит в значительной степени от отношения Н/К, а при варьировании параметрами т0 и К наблюдается слабый поворот дна оврага.

Елюхина И.В.

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

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

у о ооб

О 0045

0 003

0 0015

Ч j

I I

I I !

6 002

2 5

7.5

10 Р

а)

б)

Рис. 2. Зависимость положения оси оврага от параметров установки: а) от массы: Я = 1; 1 - М = 50; 2 - /И = 100; 3- М= 160; б) от радиуса: М= 40; 1 - Я = 0,7; 2 - Я = 1; 3- Я = 1,3; /7=6; у= 0,003; «ь = 5; К- 120 (значения параметров здесь и ниже приведены в СГС)

Пределы изменения параметров Мв настоящей работе были приняты следующими:

1) М выбирается из условия: 2Н ~(1...10)Л; 5) 30 ~ 0.001...0.004;

2) Л - 0.5...1.5 см; 6) р~ 1...20 г/см3;

3) К ~ 20...200 гсм2; 7) V/ ~ 0.001...0.05 ст.

4)г0~1...10 с;

Ошибки оцениваемых параметров свяжем с ошибками измеряемых в эксперименте величин

с помощью соотношения, имеющего, например, для вязкости, вид

______

dv — dv delv- — • At" +

dr

Эт*

dv тт dv — 9v —-• ATn + — ■ Д£ + -—АЛТ + ---AM +--A/?,

ЭАГ

ЭМ

ЭЛ

(5)

1/

0.004

где Ax-относительная точность измерения какого-либо параметра. Для расчетов принимаем точность измеряемых коэффициентов, обычно реализуемую в экспериментах: А г= 0.0001, ДТо = 0.0001, Д£ = 0.001, A£ = 0.01, ДМ =0.001, АЛ = 0.01.

При существенном изменении значения производной v'x Ах (например, на 1 %) происходит заметный сдвиг оврага и, следовательно, минимума функции качества, соответствующего оцениваемым вязкости и плотности. На рис. 3 проиллюстрировано данное смещение на примере

ошибки в измерении К . Здесь овраг 1 соответствует К -150 (истинному значению), 2 - £ = 152, 3-К = 154;/? = 6, v = 0.003, М=50, Л = 1, г0=5.

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

чувствительность плотности к вязкости p'v изменялась в среднем в пределах от 1 до 100, и в связи с этим изменялись значения рх (4).

Результаты расчетов v'x для двух торцевых поверхностей приведены в таблице. Значения чувствительности по периоду колебаний т в пределах 1 % совпадают со значениями чувствительности по т0, и поэтому этот параметр отдельно не

Рис. 3. Смещение оврага в результате ошибки рассматривается измеряемых параметров

mm +

1

0 00375

0 0035

0.00325

0,003

0 00275

0 0025

Таблица

Чувствительность вязкости к измеряемым в эксперименте параметрам

ду/ /драг Э у/ Уд 5 ду/ /дк ду/ /дМ ду/ /дя ду/ /дт0

Интервал изменения менее 1 2-2 5 1 5-2 5 2-3 5 до 1000

Обратим внимание на ге параметры, которые имею! самую низкую точность измерений Это, прежде всего, момент инерции К и радиус тигля Л, точность измерения которых не превышает 1 % Из таблицы видно, что максимальная ошибка в определении вязкости, вносимая данными параметрами может достигать 4-6 % Отклонение остальных параметров (исключая пока из рассмотрения г и т0) не оказывает существенного влияния на точность оценивания вязкости и плотности

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

¡х{у,р,г,д) = \Кс(Р)\ и ¡2{у,р,т,6) = \тЩ, (6)

установлено, что чувствительности вязкости к измеряемым в эксперименте параметрам, рассчитанные через функции (6) различны на несущественную величину, кроме чувствительности к периодам т и т0, которая может быть в сотни раз выше для мнимой части по сравнению с действительной При определении из вискозиметрической функции (1) только вязкости среды требуется одно уравнение (из двух (6)), наиболее подходящее для этой цели в конкретной задаче, а при оценке двух параметров его необходимо дополнить вторым соотношением (например, в [4] проводится совместное решение для и /2 (6))

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

1 Пусть мы решаем задачу, аналогичную рассмотренной в работе [4], те, учитывая отклонение от истинных значений только периода и декремента затухания колебаний, определяем области одновременной оценки вязкости и плотности при ошибке менее 1 %

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

- «по периоду колебаний» (иначе говоря, у'Т принимает высокие значения только при расчете через /2 (6)) Если выбрать две экспериментальные точки, соответствующие различным наборам параметров М ,Я,К,т0,80,р,у и провести расчет только по действительным частям уравнений,

то можно получить соответствующие решения с ошибкой менее 0,1 % и по вязкости, и по плотности Достаточно простым для практической реализации представляется выбор точек с различными значениями масс или периодов г0 Лучший поворот оси оврага, исходя из анализа чувствительности, обеспечивается вариацией массы, поэтому остановимся именно на этом варианте

Рассмотрим планирование данного вискозиметрического эксперимента в терминах теории чувствительности На множестве всех допустимых выборок измерений определим точки, минимальное число которых (а именно, 2) обеспечивает однозначное, устойчивое и надежное решение рассматриваемой задачи Воспользуемся условием (см , например, [6])

—> шах , (7)

т е

Ъ8Х/ Ъ8г/ (дд{/ _дд2/ ч /3у /дук /Эр /Эр'

—»шах , где - номера измерении.

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

Елюхина И.В.

Планирование оптимального эксперимента по одновременному _определению вязкости и плотности ньютоновской среды

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

1) использования при расчетах действительной и мнимой частей функции качества для одной экспериментальной точки;

2) непрерывного изменения массы, т.е. построения функций вида г = г(М), 8 = 8(М). Здесь, к тому же, существует большая возможность накопления ошибок измерения массы, периодов г, т0 и декремента.

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

2. Теперь учтем зависимость ошибок оцениваемых параметров от ошибок всех измеряемых в эксперименте величин, т.е. воспользуемся формулой (4). Напомним, что в большей степени ошибки в одновременном измерении р и У связаны с ошибками в опытных значениях параметров установки: момента инерции К и радиуса тигля Я, а V1 = у\. В связи с этим кривые, определяющие 1%-ную границу в [4], описывают границу в среднем до 5 % при измерении V и до 5-Ру (т.е. до нескольких сотен) % при оценке р.

Учтем, что в экспериментальных значениях характеристик колебаний заключены некие эффективные значения радиуса и момента К . Достаточно очевидным в данном направлении представляется метод, основанный на определении 4-х параметров (у,р, Я, К) из нескольких действительных частей вискозиметрических уравнений, построенных для различных экспериментальных точек(т,т0,8,80,М ) с одинаковым набором (Я,К).

Выбор точек осуществим как в терминах матрицы Якоби, так и путем численного анализа поведения функций чувствительности для различных выборок измеряемых параметров. Математическое моделирование вискозиметрического эксперимента для ньютоновской жидкости показало, что использование комплекса параметров с максимально возможно различными массами и периодами т0 (т.е. четыре точки с Мшп, Мш при т0тш и М^, М^ при т0гшх) позволяет/>е-гиить данную задачу оптимально в виду поворота оврага в различных плоскостях пространства (V,

Теперь обратимся к численной реализации указанных выше алгоритмов. Функция качества

= (индекс ] соответствует номеру экспериментальной точки) имеет че-

з

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

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

1) периода (приводящие к ошибкам в значениях вязкости в несколько десятков и в значениях плотности в несколько сотен процентов);

2) момента инерции подвесной системы без жидкости и радиуса тигля (приводящие к ошибкам при определении плотности в отдельных областях до нескольких сотен процентов).

Повторим еще раз, что в общем случае определение неизвестных параметров должно проводиться по действительным частям вискозиметрических уравнений, построенных для четырех точек (K9R9T0mm9SQ9Mmm —(K,R9T0msJl9SQ9Mmn ~~>t29S2)9 (К,/?,г0^пД,Мтах -*ТЪ90Ъ) и ( К, R, T0mdX, S0, Mmax -» т4, ), совместно с определением характеристик вискозиметра R иК .

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

Работа выполнена при поддержке РФФИ-Урал (№ 01-01-96424) и программы поддержки научного творчества молодежи в вузах Челябинской области.

Литература

1. Швидковский Е.Г. Некоторые вопросы вязкости расплавленных металлов. - М.: ГИТТЛ, 1955. -206 с.

2. Островский О.И., Григорян В.А. О структурных превращениях в металлических расплавах // Известия вузов. Черная металлургия. - 1985. - № 5. - С. 1-12.

3. Бескачко В.П., Вяткин Г.П., Уткин Е.А., Щека А.И. Моделирование экспериментов по измерению вязкости методом Швидковского // Расплавы. - 1990. - № 2. - С. 57-64.

4. Nieuwoudt J.C., Sengers J.V., Kestin J. On the theory of oscillating-cup viscometers / Physica. -1988. - 149A. - P. 107-122.

5. Yelyukhina I.V., Toropov E.V. Estimation of heat and mass transfer coefficient by parametric identification methods // Mater, of Internat, seminar «Modeilling, advanced process technology, expert and control system of heat and mass transfer phenomena». - Ekaterinburg, 1996.

6. Елюхин В.А., Холпанов Л.П. Статистическое оценивание параметров в задачах идентификации // Теорет. основы хим. технологии. - 1990. - Т. 24, № 6. - С.784-793.