CC BY
40
УДК 004.67 + 532.57 + 533.6.08
М. Н. Карчевский \ И. Е. Полетаев 2, Г. С. Сухоруков 2
1 Институт теплофизики им. С. С. Кутателадзе СО РАН пр. Акад. Лаврентьева, 1, Новосибирск, 630090, Россия
2 Новосибирский государственный университет ул. Пирогова, 2, Новосибирск, 630090, Россия
karchevskymi@gmail.com,poletaev.igor.phys@gmail.com, grigoriy.sukhorukov@gmail.com
АЛГОРИТМЫ РАСПОЗНАВАНИЯ И СЛЕЖЕНИЯ ЗА ПУЗЫРЯМИ ДЛЯ ИЗМЕРЕНИЯ ПАРАМЕТРОВ КАВИТАЦИИ НА ГИДРОКРЫЛЕ
Работа посвящена модификации алгоритмов расчета поля скорости в объеме потока по изображениям на примере решения задач, связанных с кавитацией. Кавитация порождает целый комплекс физических процессов, протекающих в жидкости и твердом теле, результатом которых является кавитационная эрозия - образование, перенос потоком и коллапс облаков кавитационных пузырей, пластическая деформация твердого тела, образование трещин и т. п. В связи с этим разработка эффективного метода определения гидродинамических свойств кавитирующих течений является крайне актуальной задачей с точки зрения совершенствования и оптимизации современного гидрооборудования. Компьютерный анализ пар мгновенных последовательных изображений пузырей, зарегистрированных с малым промежутком времени между кадрами, позволяет рассчитать их смещение за время между кадрами, распределение пузырей на поверхности пластины и динамику роста. По полученным данным проводится статистический анализ, делаются выводы о структуре течения, а также приводится сравнение результатов с теоретическими исследованиями. В ходе работы усовершенствованы существующие и созданы новые алгоритмы обнаружения пузырей на изображениях и расчета поля скоростей. Повышены точность и пространственное разрешение измеренного поля скоростей, а также апробированы теоретические выводы.
Ключевые слова: кавитация, трассеры, распознавание Particlelmage Velocimetry (PIV), Particle Tracking Velocimetry (PTV), поле скорости, алгоритмы, цифровая обработка изображений.
Введение
Важной задачей для проектирования эффективных гидроагрегатов и обтекаемых тел является изучение процессов, происходящих во взаимодействующих с ними потоках жидкости и газа. Для проведения фундаментальных исследований потоков, содержащих сложные вихревые структуры, требуется развивать методы вычислительных экспериментов. Данная работа посвящена модификациии применению алгоритмов обработки экспериментальных данных для решения проблемы кавитации [1-3] - процесса парообразования в областях с резкими перепадами давлений. В большинстве случаев кавитация несет негативные последствия, так как газы в пузырях химически агрессивны, а также могут иметь высокую температуру, что негативно сказывается на поверхности материалов, взаимодействующих с потоками. В связи с этим изучение гидродинамики кавитирующих течений является крайне актуальной задачей с точки зрения совершенствования и оптимизации современного гидрооборудования.
В литературе описано большое количество экспериментов и теорий, посвященных обтеканию одиночных двумерных гидрокрыльев. При этом основные результаты экспериментальных исследований представляют собой визуальные наблюдения кавитационных каверн,
Карчевский М. Н, Полетаев И. Е., Сухоруков Г. С. Алгоритмы распознавания и слежения за пузырями для измерения параметров кавитации на гидрокрыле // Вестн. Новосиб. гос. ун-та. Серия: Информационные технологии. 2016. Т. 14, № 1. С. 23-38.
ISSN 1818-7900. Вестник НГУ. Серия: Информационные технологии. 2016. Том 14, № 1 © М. Н. Карчевский, И. Е. Полетаев, Г. С. Сухоруков, 2016
построенные карты режимов течений с качественным описанием их характеристик и точечные измерения давления и скорости [4; 5]. Тем не менее, несмотря на большое количество обширных работ по исследованию кавитирующих течений, до сих пор детальная количественная информация, необходимая для оптимизации существующих математических моделей, крайне ограничена, даже для упрощенных условий, а результаты работ различных авторов зачастую расходятся. Это связано с проблематичностью, а иногда и невозможностью проведения экспериментов в натурных условиях и высокой стоимостью испытаний на дорогостоящих модельных стендах. Для изучения пространственной структуры и динамики кавита-ционных каверн, а также описания и анализа реализующихся режимов течений применялась высокоскоростная фотосъемка с разрешением несколько мегапикселей.
Основными методами обработки вышеописанных экспериментов измерения распределения скорости течений и турбулентных характеристик вблизи гидрокрыльев являются так называемые алгоритм кросскорреляционного анализа трассерных картин PIV (Particlelmage Velocimetry) и алгоритм слежения за каждой отдельно взятой частицей PTV (Particle Tracking Velocimetry). PIV обрабатывает изображения целиком без их предварительного анализа, PTV оперирует объектами и их координатами. В нашем случае такими объектами являются кави-тационные пузыри.
Постановка эксперимента
Бесконтактный способ анализа различного рода потоков [6] является наиболее эффективным, так как он позволяет получить информацию о структуре и распределении скорости потока без воздействия на него. Для получения данных в экспериментах мы использовали так называемый метод цифровой трассерной визуализации (ЦТВ), который на данный момент является наиболее современным, эффективным и динамично развивающимся методом количественного исследования потоков жидкостей и газов.
Для определения поля скорости в поток добавляются полимерные светоотражающие частицы (d ~ 10 мкм). При освещении лазерным ножом объема потока ПЗС-камера фиксирует мгновенную трассерную картину. По серии таких картин определяется смещение частиц между кадрами, а зная время между экспозициями, можно определить поле скорости в структуре потока. В подавляющем большинстве реализаций методов цифровой трассерной визуализации необходимо только два последующих кадра, отснятых с небольшим временным промежутком. Таким образом, «мгновенное» поле скорости как результат обработки пары трассерных картин - это среднее значение скорости за время между кадрами. Для успешной обработки изображений с трассерами необходимо, чтобы частицы не покидали область, засвеченную лазерным ножом. К тому же требуется, чтобы смещение частиц было достаточно малым для корректной работы алгоритмов. Эти условия накладывают ограничение сверху на временной межкадровый интервал (до 100 нс для гиперзвуковых потоков), что позволяет называть результирующее поле скорости мгновенным (рис. 1).
Рис. 1. Принципиальная схема эксперимента с применением методов цифровой трассерной визуализации (слева)
и фотография экспериментальной установки (справа)
Рис. 2. Эксперимент: кавитационные пузыри, отснятые высокоскоростной ПЗС-камерой
Ввиду того, что требовалось проанализировать явление кавитации, в схему метода были привнесены изменения: в качестве трассеров выступали кавитационные пузыри (рис. 2), алгоритм распознавания которых описан ниже.
Алгоритм распознавания пузырей на изображениях
Идентификация пузырей на изображениях -весьма нетривиальная задача. В данном эксперименте пузыри освещены не полностью из-за сложного оптического доступа, с течением времени они могут резко менять форму и скорость, а также зачастую наплывать друг на
друга. Для решения этой проблемы был создан следующий эффективный алгоритм нахождения пузырей. В данной серии экспериментов использовались камеры с высокой частотой съемки (порядка 20 кГц). Из-за ограничений внутренней памяти камеры не удавалось получить более 10 000 последовательных изображений. Такой большой набор данных позволил применить различные временные фильтры, а также вероятностную фильтрацию.
На первом шаге данной фильтрации вычислялось среднее поля интенсивности I (х, у) по всем ансамблям данной серии экспериментов:
I (х, у )=N Ьп (х, у).
^ п=1
Далее для каждого отдельного пикселя на изображении смотрелась разница интенсивности со средним изображением:
I(х уЛ0,7^у\>1 ^у)
^У): [7(х,у),I(х,у)<I(х,у).
Если интенсивности сильно отличаются, то значение пикселя оставлялось исходя из предположения о том, что с высокой долей вероятности это движущийся объект является пузырем. После применялся медианный фильтр и бинаризация изображения по порогу Р:
( ) I0,7 (х, у)> Р
7(уН1,I(х,у)<Р •
Итого, все оставшиеся области анализировались на предмет их схожести с эллипсом и аппроксимации параметров эллипса по точкам области методом Левенберга - Маркова [7] (метод оптимизации, направленный на решение задач о наименьших квадратах, являющийся альтернативой методу Ньютона).
Для повышения точности определения параметров проводилась повторная аппроксимация параметров эллипса, взятая с большим числом точек. Дополнительные точки брались из предположения, что центр эллипса, найденный на предыдущем шаге, является истинным центром, а значит, бинаризацию данной области можно производить с меньшим порогом. Результатом работы алгоритма являлся массив координат центров пузырей и их радиусы.
Алгоритм Particle Image Velocimetry
Алгоритм PIV [8-10] обрабатывает два ансамбля изображений, которые строятся по следующему принципу: вся серия изображений разделяется на двухкадровые подсерии, аналогично проведению эксперимента с двумя кадрами. Далее, первый ансамбль заполняется изображениями в первый момент времени, второй - следующий момент времени.
Первым шагом этого алгоритма является расчет коэффициента корреляции Пирсона [8] между фиксированной областью изображения из первого ансамбля кадров и различными областями изображения из второго ансамбля по формуле
R.
j ) =
INi (( j) - i ( j)) • (in Q + Ar, j + As) -12 (/ + Ar, j + As)) _
ai (, j)a2 ( + Ar, j + As) '
£ ( j )=N 11 in (j);
(i)
( j )=
Здесь Щ (i, j) - интенсивность области (i, j) на n-м кадре из k-го ансамбля изображений.
Таким образом, получены корреляционные плоскости.
Как уже говорилось, одной из проблем PIV-анализа является низкое качество полученных корреляционных плоскостей. Для устранения не коррелируемого шума соседние корреляционные функции, имеющие схожую структуру, перемножаются. В качестве итоговой корреляционной плоскости берется их среднее геометрическое. Определение положения перемножаемых плоскостей зависит от качества исходной области. В обрабатываемых нами изображениях наилучшее отношение качества к времени обработки достигалось при перемножении четырех соседних плоскостей (верхняя, нижняя, левая, правая).
Далее для определения наиболее вероятного смещения выполняется поиск координаты максимума корреляционной плоскости. Как известно, форма корреляционного пика описывается функцией Гаусса [6; 11]. Аппроксимация пика функцией Гаусса проводится при помощи метода fit с параметром gaussl - этот метод является стандартным в среде разработки MatLab.
Рис. 3. PIV-обработка: принципиальная схема расчета корреляционной плоскости
Затем, после нахождения наиболее вероятного смещения оно делится на указанный временной период между кадрами в парах, тем самым получая значения скорости в исследуемой области. Далее, из других восстановленных параметров, а также поворота и растяжения корреляционного пика на плоскости можно получить дополнительную информацию о потоке в данной точке [6], такую как пульсации и завихренности.
Пульсации скорости можно определить (3) из параметров дисперсии функций Гаусса (2), описывающих кривую пересечения корреляционного пика с плоскостями X—z и Y—z, где (X, 7) - дискретные координаты максимума корреляционного пика.
Рис. 4. PIV-обработка: восстановление параметров пика корреляционной плоскости
1 (х-а.)
/ (к) =-е 2°2 , к = х, у. (2)
ок>12%
Пульсации скорости:
Рх =у1°~х ; Ру . (3)
О завихренности можно говорить в том случае, если максимум корреляционной плоскости повернут. Для этого находятся Х-координаты Х^- и Х^Х максимумов корреляционного пика по сечениям У -1 и У +1 соответственно. Далее, вычисляется угол между осью У и осью найденных максимумов (4), это значение нормируется.
Коэффициент завихренности, необходимый для определения локальных характеристик потока:
4 2
Ф = — агссоэ^ . (4)
% - УУах1 )2+4
Алгоритм Particle Tracking Velocimetry
Для измерения поля скорости, наряду с подходом PIV, используется методика слежения за каждой отдельно взятой частицей - Particle Tracking Velocimetry (PTV) [6; 11-13]. После обработки изображений алгоритмом можно получить так называемое нерегулярное поле скорости: каждой частице на каждом кадре в соответствие ставиться свой вектор, в отличие от регулярных полей, где последнее строится с векторами в узлах регулярной сетки (рис. 5).
Рис. 5. Пример нерегулярного (слева) и регулярного поля скорости (справа)
РТУ использует два кадра, отснятых с небольшим интервалом времени - А?. На первом этапе в каждом из двух кадров определяются положения частиц (рис. 6, 1). В случае эксперимента с кавитацией это происходит посредством алгоритма распознавания пузырей на изображениях, в других случаях применяются методы корреляции с маской частицы или метод простой бинаризации изображения по заданному порогу [6]. Затем для каждой частицы на первом кадре ведется поиск пары на втором кадре. Таким образом, с учетом межкадровой задержки для каждой частицы определяется ее вектор скорости. После проводится фильтрация данных, в результате чего отбрасываются ошибочные векторы (подробнее о фильтрах ниже). Итого, на выходе у РТУ после обработки для двух кадров - мгновенное нерегулярное поле скорости потока.
Для поиска пар частиц (далее - точек) используется так называемый релаксационный метод [6; 12], основывающийся на двух предположениях.
1 2
Рис. 6. Иллюстрация к работе метода РТУ. Синие и красные точки на изображениях - центры частиц (пузырей) на г-м и (г + 1)-м кадрах соответственно
1. Максимальная скорость Утах определена. Известна максимальная скорость точек в потоке и, соответственно, максимальное расстояние - Тт, на которое она может переместиться за время между кадрами А??.
2. Непрерывное распределение скоростей. Локальная группа точек движется примерно в одном направлении.
Реализация алгоритма: для каждой точки на первом кадре ведется поиск пары на втором кадре в пределах окружности с радиусом максимально возможного смещения Тт = Утах • А? (рис. 6, 2). Каждому из возможных исходов присваивается одна и та же начальная вероятность:
р(0) = р(о) (5)
у ' N +1 1 7
где рук) - вероятность перехода г-й точки в]-ю на к-м шаге итерации, Р^) - вероятность отсутствия пары у точки, N- количество возможных пар для данной (г-й) точки (см. рис. 6, 3). Сумма вероятностей всех исходов при такой нормировке, очевидно, равна единице:
X Р^к)= 1, Vk. (6)
N
Предположение 2 приводит к тому, что можно ввести порог Тп, в пределах которого соседние точки показывают схожее движение, т. е. векторы смещения соседних точек А1 (находящихся в пределах Тп) не должны отличаться больше, чем на величину Тд (рис. 6, 4) -
условие жесткости а. Посчитанная на предыдущем (к - 1)-м шаге итерации сумма _1) вероятностей переходов всех тех точек, которые находятся в пределах Тп и у которых в ожи-
даемой области Тт найдена возможная пара, - дает определенный вклад в ненормализованную вероятность текущего события * Р^) на (к)-й итерации (перехода 7-й точки из первого
кадра в у-ю точку во втором кадре), а именно: сумма 0(к -1) умножается на коэффициент сходимости вероятности В, а вероятность перехода, посчитанная на предыдущем шаге итерации Ррк -1), умножается на коэффициент сходимости вероятности А, или
* у=а . р-) + в. ек-1), ек-1)р1к-1) . (7)
Здесь ^щ! =
ЧА' - Т -
- условие жесткости; А1 - отклонение от возможного смещения со-
>1 Д,|| > Т
седа; коэффициенты - А < 1, В > 1. Окончательно, с учетом нормировки вероятностей на
единицу, расчетные формулы приобретают вид
* Р(к) р(к-1) р(к)__у_ р(к)__17__(8)
у р(к-1)+у *р(к)' Гу р(к-1)+у *р(к)• ^
7 у а 7 у а
Таким образом, после определенного числа итераций для каждой точки определяется ее максимально вероятная пара на следующем кадре, после строится нерегулярное поле скорости для текущего кадра.
Приведенный метод обладает достаточной точностью без применения предварительной обработки. Однако имеет некоторые ограничения ввиду использованных предположений. Оптимальная настройка параметров стандартного РТУметода, используемых в формулах (5)-(8), более подробно описана в [6].
Улучшение стандартных алгоритмов
Перед тем, как перейти к основной части статьи, обобщим информацию о РТУ и Р1У в виде следующей таблицы.
Характеристики стандартных РТУи Р1У
Particle Image Velocimetry (Single Pixel) Particle Trace Velocimetry
Достоинства • не требует выделенных объектов; • средние характеристики потока; • пульсации, завихренности • высокое пространственное разрешение; • слабая чувствительность к неравномерному засеву; • возможен анализ нестационарных потоков
Недостатки • высокая точность достигается только при анализе большого числа изображений; • низкая скорость обработки; • невозможен анализ нестационарных потоков • высокая ресурсоемкость; • зависимость производительности от параметров самого потока. Сложная настройка параметров алгоритма
Ограничения • ламинарность течения. На всем ансамбле изображений поток должен быть стационарным • максимальная скорость определена. Известна максимальная скорость частицы в потоке; • непрерывное распределение скоростей. Локально группа частиц движется примерно в одном направлении
В различного рода экспериментах с кавитацией могут получаться изображения разной степени однородности. Как уже говорилось, изображения подобных потоков хорошо обраба-
тываются алгоритмом Р1У. Обработка же подобных изображений алгоритмом РТУ затруднена ввиду сильной зависимости эффективности обработки от изначально заданных параметров. Очевидно, что в данном случае параметры завихренности потока различны для разных областей изображения.
Решением данной проблемы является деление изображений на области со схожим коэффициентом завихренности по данным предварительного РГУ-анализа. Первоначально изображения делятся на одинаковые области, размеры которых зависят от размеров исходного изображения. Далее, для каждой области рассчитывается средний параметр завихренности, а также модуль отклонения от медианного значения. Затем проводится объединение и перерасчет областей со схожими параметрами. В итоге, алгоритм РТУ теперь будет обрабатывать каждый из получившихся массивов субкадров, выбирая настройки индивидуально для каждого из массивов.
Улучшение PTV
Как следует из вышеизложенного, для правильной работы PTV необходима его предварительная настройка (а именно 6 параметров, описанных в алгоритме Particle Tracking Velocimetry), что является достаточно сложной задачей ввиду того, что от кадра к кадру, очевидно, картина потока может сильно меняться. Для решения этой проблемы и упрощения работы экспериментатора сделали следующее.
Первое нововведение: коррекция радиуса поиска соседей Tn в зависимости от средней плотности частиц р на субкадрах:
T =-
а • Tm
Р
Здесь а - экспериментально определяемый параметр Тт - максимальное смещение, которое в свою очередь теперь определяется предварительной РГУ-обработкой. Теперь благодаря тому, что не приходится обрабатывать «лишних соседей», алгоритм стал работать быстрее при
практически нулевой потере точности в . Как видно из рис. 7, в среднем скорость работы возросла на 18 %, погрешность осталась практически неизменной.
Второе нововведение: автоматическое определение оптимального критерия жесткости Ту и числа итераций
в зависимости от коэффициента завихренности ф:
Ту =р(1 + ф)Тт;
N = N (ф)е[3,5].
Здесь р - экспериментально определяемый параметр. Как видно из рис. 8, при добавлении каждой из итераций время работы увеличилось в среднем на 3-4 минуты соответственно, но качество результатов при этом, очевидно, улучшилось.
Все вместе это позволяет теперь наиболее эффективно использовать вычислительные мощности и понизить погрешность алгоритма.
Рис. 7. Тест производительности 1 (1 000 кадров; 3 итерации; Tm = 10; Tq = 3)
Очевидная потребность извлечь из имеющихся экспериментальных данных максимум информации подтолкнула нас к созданию принципиально нового алгоритма, позволяющего повысить качество его результатов и сократить время работы, благодаря уменьшению необходимого числа итераций. На рис. 9 приведена более подробная иллюстрация этой идеи.
Коротко реализацию алгоритма, описанного выше, можно представить в виде следующих шагов.
• В зависимости от коэффициента завихренности ф предварительно проводится определенное число итераций.
• Далее, на трех последовательных субкадрах выбирается некоторая частица на промежуточном (рис. 10, 1). Она, очевидно, является возможной парой и возможным предшественником для некоторых точек со следующего и предыдущего субкадров.
• Рассмотрим все полученные тройки частиц, у которых первая и последняя находятся по разные стороны от промежуточной. Получится некоторый набор комбинаций (рис. 10, 2).
• Отсеем заранее неверные направления возможных смещений, используя предварительную Р1У-обработку (рис. 10, 3).
• Соединим любого предшественника со всеми возможными парами для второй частицы, обозначив при этом наименьшее расстояние от 2-й частицы до полученных отрезков как Н соответственно (рис. 10, 4). После среди всех к^ ищется минимальное - Нтш (рис. 10, 5).
• Теперь, вспоминая о том, что чем ближе тройка частиц находится к одной прямой, тем вероятнее, что эта одна и та же частица, добавим к вероятности р^1' для 1-й и Р^1' для 2-й точки перейти в 2 и 3 соответственно специально рассчитанную индивидуальную добавку х,, [7]:
Рис. 8. Тест производительности 2 (1 000 кадров; р = 20; Тт = 10; Тп = 3)
* р(^+1) _ р( 1 О ■*■ 1 О
* р^1) _ р(*
хр _ тах ^
Л
Л
1
2 V 2
Т
V ч
• После все вероятности нормализуются по формулам (8).
Алгоритм проходит по всем частицам каждого кадра, тем самым существенно уменьшая погрешность результатов (рис. 10, 6).
Подводя промежуточные итоги, представим следующий график, иллюстрирующей время работы по построению поля скорости РТУ методом и погрешности результатов в зависимости от средней плотности частиц на субкадр до и после применения всех вышеописанных улучшений. Для слабозавихренных изображений (ф«0,2) удалось ускорить работу алгоритма в среднем на 150-170 % и сократить погрешность на 28 %.
Улучшение стандартных алгоритмов. Фильтры
Автоматическая обработка изображений методом РТУ неизбежно приводит к появлению ошибочных векторов смещения. Ввиду того, что негативный вклад данных векторов в итоговый результат может быть весьма существенным, их удаление в процессе обработки является важной задачей. С этой задачей справляются фильтры промежуточных результатов, запускаемые между различными этапами обработки. Первый фильтр является стандартным для РГУ-анализа, он основан на алгоритме скользящего среднего на нерегулярной сетке [14]
Рис. 9. Идея улучшения стандартного РТУ
Рис. 10. Иллюстрация к работе алгоритма улучшения РТУ. Черные, синие и красные объекты на изображениях - центры частиц (пузырей) на (г - 1)-м, г-м и (г + 1)-м кадрах соответственно
и предположении о локальной гладкости векторного поля. Значение каждого вектора сравнивается со значением соседних векторов. Вектор V(г0) рассматривается
как ошибочный, если выполняется соотношение
|V (Г0 )- V (Г0 )|> к, где к - пороговое значение:
к = а • | V, (Го)- V (Го )|,
пеЩго )1
0 < а < 1.
Здесь V(г0) - это вектор средней локальной скорости, посчитанный по соседним векторам:
'('0 ) =
I п
V №
1Ы(г0) п п
( '0 )
N ('0
№
('0) ■
Рис. 11. Тест производительности 3 (1 000 кадров; Тт = 10; Тп = 30; ф0 = 0,2)
Описанный метод может быть применен несколько раз итерационно. Каждая итерация удаляет векторы, наиболее отличные от среднего вектора, рассчитанного с применением весовой функции:
( ) ( II' - Г0II ^
(Г0 ) = ехР .
Сложность применения алгоритма скользящего среднего к нерегулярным данным состоит в том, что нельзя подобрать оптимальный размер области осреднения, которая, безусловно, будет содержать соседние векторы. В случае если рассматривать достаточно большую область, которая заведомо включает соседние векторы, метод будет работать неточно, так как среднее значение скорости может сильно отличаться от локального значения. Для решения этой проблемы было предложено рассматривать большую область с учетом весовой функции по Гауссу.
Второй и третий фильтры являются оригинальными. Один из них использует данные Р1У и отсекает ошибочные векторы, находящиеся вне плоскости возможного смещения. Другой же ищет частицы, у которых векторы скорости попали в одну и ту же точку, а после корректирует их возможное перемещение в зависимости от корреляции между средним вектором скорости, полученным РТУ-методом, и исследуемым ошибочным вектором.
Улучшение стандартных алгоритмов. Оптимизация Р1У
В некоторых случаях Р1У-обработка является очень ресурсоемкой, ввиду большого числа обрабатываемых изображений (в некоторых экспериментах их количество может достигать сотен тысяч). Обработка может длиться несколько дней даже на достаточно мощных вычислительных системах, не говоря уже о стандартных устройствах, используемых для обработки. В связи с этим актуальной задачей является оптимизация и ввод дополнительных настраиваемых параметров, регулирующих отношение времени выполнения к погрешности.
В связи с тем, что разрешение изображений различается от случая к случаю, в целях экономии расчетных мощностей целесообразно деление изображений на области. Стандартный алгоритм Р1У в качестве области, для которой рассчитывается корреляционная плоскость, берет один пиксел. В случае же больших изображений данный подход приводит к существенному возрастанию времени обработки. Предварительное деление изображений на малые области частично решает эту проблему, однако приводит к увеличению погрешности. Размеры этих областей в пикселах задаются вручную после первоначальной оценки качества изображений.
Также для ускорения PIV-обработки вводится понятие пространственной частоты обработки (далее Spatial Frequency) - отношение пропускаемых областей к областям, для которых рассчитываются корреляционные плоскости. В процессе тестирования, после установки параметра Spatial Frequency равным 2 и 3, удалось уменьшить время обработки в 4 и 9 раз с увеличением погрешностей на 5,1 и 11,7 % соответственно.
Еще одним настраиваемым параметром является размер корреляционной плоскости. Он определяется из предполагаемых параметров потока - средняя скорость в пикселах (или областях, в случае больших изображений), и оборудования эксперимента - частота кадров. Рекомендованные для выбора размеры также зависят от направления и однородности потока: размер оси перпендикулярной вектору скорости для однородных потоков составляет 7 или 9 пикселей (или областей), размер оси параллельной вектору скорости должен быть больше, чем произведение периода между кадрами и средней скорости в пикселях (областях). Оптимальная настройка размеров корреляционной плоскости существенно снижает время, затрачиваемое на обработку.
Результаты
Для итогового анализа изображения был создан комбинированный алгоритм, включающий как РГУ-, так и РТУ-обработку. Обобщая все, что написано выше, можно сформулировать принцип работы итогового комбинированного алгоритма: на первом этапе обработкой РГУ получается поле скорости, далее находятся центры изображений, затем на основе данных РГУ изображения делятся на субкадры, эти данные, в свою очередь, передаются к обработке алгоритму РТУ, в итоге имеются два поля скорости, которые, после фильтрации ошибок и регуляризации поля скорости РТУ, объединяются в одно в зависимости от расчетных погрешностей алгоритмов РГУ и РТУ соответственно.
В процессе разработки тесты проводились не на реальных фотографиях, а на двух типах синтетических изображений потоков: однородно завихренных и смешанного вида (рис. 12). Для создания таких изображений была создана программа на языке С++ с использованием графической библиотеки ОрепСУ. Параметры завихренности, как сказано выше, определяются коэффициентом (1), локально для каждой области изображения (обычно изображение делилось на 15-20 частей), и после общий параметр завихренности ф0 для изображения вычислялся как среднее арифметическое ф для каждой из областей.
Результаты итоговой обработки каждого из массивов синтетических изображений представлены на рис. 13. Как видно, алгоритм существенно быстрее справляется с менее завихренными изображениями ввиду малого количества необходимых итераций для достижения таких значений погрешностей.
, I
i'I - 4
■ ■
* .
| г • ^
• > * • л
1 ■ • г
* 4
■ - < . : * - ■ *
! :
б
а
Рис. 12. Типы синтетических изображений: а - однородное и ф0 < 0,2; б - с неоднородностями и ф0 > 0,2
Средняя плотность частиц р. шт/субкадр Средняя плотность частиц р, шт/субкадр
а б
Рис. 13. Тест производительности 4 (1 000 кадров 64 х 64; Тт = 10: а - слабая завихренность ф0 = 0,15; Тд ~ 2,5; Ыи = 1; б - сильная завихренность ф0 = 0,5; Тд ~ 3; Ми = 3)
Рисунок 14 хорошо иллюстрирует то, насколько эффективнее стал работать итоговый алгоритм после всех улучшений. Результаты обработки исходных изображений представлены на рис. 15.
Рис. 14. Тест производительности 5 (1 000 кадров 64 х 64; Тт = 10; ф0 = 0,4; Тд ~ 3)
'ЫМ
а Ш 190 2И Я» МО
Рис. 15. Итоговое поле скорости эксперимента с кавитацией (2 000 кадров, 1024 х 256 р1х; Т= 26000 8; Тт = 3,4; р = 11; Тп = 18; ф0 = 0,24; Тд =0,8; = 2)
Заключение
В заключение отметим, что все изначально поставленные цели и задачи были выполнены: было разработано программное обеспечение, способное обрабатывать различного рода экспериментальные изображения, в том числе изображения кавитационных потоков.
Также, благодаря различного рода улучшениям, удалось существенно повысить качество результатов, ускорить время их получения и автоматизировать определение необходимых оптимальных настроек для работы алгоритмов. Благодаря этому в среднем синтетические изображения с равной завихренностью стали обрабатываться на 60 % быстрее, а средняя точность результатов увеличилась на 22 %. Итогом проделанной работы стало создание фреймворка, объединяющего улучшенные РГУ- и РТУ-алгоритмы, а также включающего алгоритм поиска центров кавитационных пузырей. В будущем данный фреймворк будет использоваться для анализа результатов эксперимента с кавитацией и без нее.
В перспективе возможно применение итогового комбинированного метода для анализа объемных потоков - подобные эксперименты уже ставятся, однако обработка их результатов - это более сложная задача. Еще одним направлением для разработок является улучшение качества поиска пузырей на изображениях/
Список литературы
1. Каупин Ф., Герберт Э. Кавитации в воде // Труды физики. 2006. С. 1-15.
2. Лабертиукс К. Р., Цессио С. Л. Кавитационные потоки. Часть 1. Кавитация, возникающая на моделях без изменения масштаба // Жидкая механика. В. 2001. Т. 431. С. 1-41.
3. ПромтовМ. Кавитация. Тамбов, 2012. С. 1-17.
4. Леру Ж.-Б., Коути-Дельгоша О., Астолфи Ж. А. Экспериментальное и численное исследование механизмов, связанных с нестабильностью кавитацией на двумерных подводных крыльях // Физика жидкости. 2005. Т. 17, № 5. С. 1-20.
5. Сан Т. И., Фэиф Г. М. Структура турбулентно пузырьковой струи-ГГ. Свойства фазовых профилей // Многофазный поток. 1986. Т. 12. С. 115-126.
6. Ахметбеков Ю. Алгоритмы и программное обеспечение для обработки информации измерительных систем оптической диагностики двухфазных потоков: Дис. ... канд. техн. наук. Новосибирск, 2011.
7. Гилл П. Э., Мюррей У., Райт М. Х. Практическая оптимизация. М.: Мир, 1985. 509 с.
8. Линдкей Р., Мерзкирч В. Новый метод РГУ для измерений в многофазных потоках и его применение в двухфазных потоках // Эксперименты в жидкостях. 2002. Т. 33. С. 814-825.
0 16 50 71 100 150 1ть 200 25С
9. Скарано Ф. Обзор РТУ в сверхзвуковых поток // Р1У. Берлин: Гейельберг: Спрингер, 2008.
10. Станислас М. Основные результаты третьего международного РТУ чемпионата // Эксперименты в жидкостях. 2008. Т. 45. С. 27-71.
11. Ахметбеков Ю., Алексеенко С., Дулин В., Маркович Д. Применение новой методики РЫБ в сочетании с РТУ и РТУ для исследования турбулентных струй с пузырями // Тр. VII Междунар. симпозиума по Р1У. Рим, 2007.
12. Эстрада-Перез C. E. Исследование пристенных турбулентных канальных потоков с улучшенным алгоритмом РТУ // У1 Междунар. симпозиум по РТУ. Пасадена, Калифорния, США, 2005.
13. Такехара К., Этох Т. Исследование по идентификации частиц в РТУ. Корреляционный метод маски частиц // Визуализация. 1999. Т. 1. С. 313-323.
14. Хост-Масден А., МакКласки Д. Р. Точность и надежность измерений Р1У // У11 Междунар. симпозиум. Лиссабон, 1994.
Материал поступил в редколлегию 05.10.2015
M. N. Karchevsky 1, I. E. Poletaev 2, G. S. Sukhorukov 3
1 Institute of Thermophysics of SB RAN 1 Lavrentiev Ave., Novosibirsk, 630090, Russian Federation
2 Novosibirsk State University 2 Pirogov Str., Novosibirsk, 630090, Russian Federation
karchevskymi@gmail.com,poletaev.igor.phys@gmail.com, grigoriy.sukhorukov@gmail.com
DETECTION AND MONITORING BUBBLES FOR MEASURING PARAMETERS
OF CAVITATION AT HYDROFOILS
This work is devoted to the modification of algorithms for calculating the velocity field in a volume of flow per images at example of solution of problems associated with cavitation. Cavitation generates a set of physical processes occurring in the liquid and solid body, the result of which is the cavitation erosion - formation, transferal by flow and the collapse of cavitation bubbles clouds, formation and propagation of waves of compression and tension in the liquid, the plastic deformation of solids, the accumulation of fatigue, formation of cracks, etc. In this regard, the development of an effective method for determining the hydrodynamic properties of cavitating flows is a very urgent problem from point of view of improvement and optimization of modern hydraulic equipment. Computer analysis of consecutive pairs of instant images of bubbles, registered with the small time interval between frames, makes it possible to calculate the displacement of the particles in the time between frames, the distribution of bubbles on the surface of the plate. By the data obtained is conducted statistical analysis to draw conclusions about the structure of the flow, also it's gave a comparison with the results of theoretical research. During the work, there were made the improvements existing and creating of new algorithms for detection of bubbles in the images and the calculation of the field of velocities. There's increased accuracy and spatial resolution of the measured velocity field and tested theoretical conclusions.
Keywords: cavitation, tracers, recognition Particle Image Velocimetry (PIV), Particle Tracking Velocimetry (PTV), the velocity field, the algorithms.
References
1. Caupin F., Herbert E. Cavitation in water. C. R. Phys, 2006, p. 1-15.
2. Laberteaux K. R. &Ceccio S. L. (2001) Partial cavity flows. Part 1. Cavities forming on models without spanwise variation // J. Fluid Mech. 2001, vol. 431, p. 1-41.
38
M. H. KapMeBCK^R, E. noneTaeB, f. C. CyxopyKOB
3. Promtov M. Cavitation. Tambov, 2012, p. 1-17.
4. Leroux J. B., Coutier-Delgosha O. and Astolfi J. A. (2005) A joint experimental and numerical study of mechanisms associated to instability of partial cavitation on two-dimensional hydrofoil // Phys. Fluids, Vol. 17, No. 5, p. 1-20.
5. Sun T.Y. and Faeth G.M. Structure of turbulent bubbly jets-II. Phase property profiles [Journal] // Int J Multiphase Flow. 1986. T. 12. p. 115-126.
6. Akhmetbekov Y. Algorithms and software information processing means for measuring systems of optical diagnostics of two-phase flows. Sci. Dis. Novosibirsk, 2011, p. 4-120.
7. Gill P. E., Murray W., H. Wright M. Practical optimization. M.: World, 1985. 509 c.
8. Lindken R. and Merzkirch W. A novel PIV technique for measurements in multiphase flows and its application to two-phase bubbly flows [Journal] // ExpFluids. 2002 . Vol. 33. p. 814-825.
9. Scarano F. Overview of PIV in supersonic Flows [Part of book] // Particle Image Velocimetry. Berlin: Heidelberg: Springer, 2008.
10. Stanislas M. M ain results of the third international PIV challenge. [Journal] // Exp. Fluids. 2008 r. 1: T. 45. p. 27-71.
11. Akhmetbekov Y., Alekseenko S., Dulin V., Markovich D. Application of a novel PLIF technique combined with PIV and PTV to bubble turbulent jets study. Proceedings of «7th International Symposium on PIV». Rome, 2007.
12. Estrada-Perez C. E. Near-wall study of turbulent channel flows with an improved PTV algorithm. 2005.
13. Takehara K. and Etoh T. A study on particle identification in PTV. Particle Mask Correlation Method. [Journal] // J Visualization. 1999. T. 1. p. 313-323.
14. Host-Masden A. and McCluskey D. R. On the Accuracy and Reliability of PIV measurements [Conference]. 1994.
