CC BY
149
Особенности построения расчетной схемы для моделирования динамики стабилизатора расхода в пакете OpenFOAM
В.Г. Мельникова <vg-melnikova@yandex.ru> О.С. Коцур <oskotsur@gmail.com> Г.А. Щеглов <shcheglov_ga@jbmstu.ru> Федеральное государственное бюджетное образовательное учреждение высшего образования «Московский государственный технический университет имени Н.Э. Баумана (национальный исследовательский университет)» (МГТУим. Н.Э. Баумана), 105005, Россия, Москва, 2-я Бауманская ул., д. 5, стр. 1
Аннотация. Представлены результаты методического исследования, направленного на тестирование возможностей свободно распространяемого пакета OpenFOAM по построению расчетной схемы для моделирования динамики гидравлических агрегатов с помощью метода контрольного объема и подвижных сеток. Рассмотрены выбор наиболее подходящих технологий динамических скользящих сеток; построение расчетной схемы для сопряженной задачи FSI о взаимодействии подвижного регулирующего элемента (золотника) с течением рабочей жидкости; расчет переходного режима движения золотника для прототипа стабилизатора расхода новой конструкции; анализ нестационарных течений во внутренней полости агрегата; анализ устойчивости и производительности счета. На примере упрощенной осесимметричной модели стабилизатора расхода несжимаемой среды описаны основные этапы подготовки расчетной схемы для моделирования задачи об установлении золотника в положение равновесия под действием гидродинамической силы и силы реакции пружины. Приведены результаты предварительного стационарного расчета при фиксированном золотнике, которые использовались в качестве начального условия. Подробно описаны методы реализации движения сетки, отслеживающей движение золотника и технология скользящих сеток GGI, которая использовалась для моделирования перекрытия отверстия. Приведены результаты расчета переходного режима: поля скоростей и давления, графики обобщенных сил, действующих на золотник, график перемещений золотника от времени. Расчеты показали механизм возникновения нестационарных струйных и вихревых течений в проточной части стабилизатора. Проведенное методическое исследование позволяет сделать вывод о том, что пакет с открытым исходным кодом OpenFOAM в версии extend может быть успешно использован в качестве альтернативы коммерческим пакетам программ вычислительной гидродинамики (CFD), поскольку в нем имеются все необходимые средства для
построения расчетных схем с деформируемыми сетками и расчета переходных режимов в агрегатах гидравлической автоматики, имеющих в конструкции подвижные элементы.
Ключевые слова: численное моделирование; вычислительная гидродинамика; FSI; OpenFOAM; открытое программное обеспечение; стабилизатор расхода
DOI: 10.15514/ISPRAS-2017-29( 1 )-4
Для цитирования: Мельникова В.Г., Коцур О.С., Щеглов Г.А. Особенности построения расчетной схемы для моделирования динамики стабилизатора расхода в пакете OpenFOAM. ТрудыИСПРАН,том29,вып. 1,2017г., стр. 53-70. DOI: 10.15514ASPRAS-2017-29(1 >4
1. Введение
Для наиболее полного удовлетворения специальным эксплуатационным требованиям при создании новых агрегатов гидравлической автоматики аэрокосмической техники важно проведение математического моделирования режимов их работы. В частности, для обеспечения синхронности выполнения операций исполнительными органами важной задачей является расчет параметров процессов, протекающих в стабилизаторах расхода (CP) золотникового типа. Данные агрегаты, за счет малых перемещений золотника, должны поддерживать с заданной точностью расход рабочей жидкости (РЖ) при значительных (более 10 МПа) перепадах давления между входом и выходом. При этом в проточной части регулируемых сечений возникают гидродинамические силы, вызывающие заклинивание золотника, что искажает расходно-перепадную характеристику СР. Движение РЖ здесь является существенно нестационарным из-за чего применение приближенных аналитических моделей гидродинамического нагружения золотника для выбора проектных параметров конструкции CP не позволяет получить удовлетворительного согласования с данными экспериментов [1]. При численном моделировании течения РЖ в регуляторах методом конечного объема в современных программных комплексах вычислительной гидродинамики (CFD) подвижностью золотника часто пренебрегают для уменьшения вычислительных затрат [2, 3]. Однако, чтобы в условиях больших перепадов давления в полной мере проанализировать нелинейное нестационарное нагружение золотника, необходимо рассматривать сопряженную задачу его движения при взаимодействии с потоком среды во внутренней полости агрегата - т.н. задачу FSI (Fluid-Structure Interaction). Полноценное FSI-моделирование работы CP в программах CFD требует построения расчетных схем с использованием подвижной (перестраиваемой) сетки, описывающей такие особенности работы гидравлического агрегата как перекрытие отверстий и зазоров, резкое изменение объемов внутренних полостей при перемещениях золотника. Наиболее приспособленными для работы с подобными сложными расчетными схемами являются универсальные
коммерческие программы CFD, в частности. Fluent и CFX, входящие в комплекс ANSYS [4, 5]. Однако высокая стоимость владения пакетом, закрытость исходного кода и ограничения коммерческой лицензии, которые не защищают промышленные предприятия от возможных санкций со стороны иностранных правообладателей, делают задачу поиска альтернативных программных средств с открытым исходным кодом безусловно актуальной. В настоящее время среди программ с открытым исходным кодом наибольшими возможностями по расчету полей различной природы методом контрольного объема обладает пакет OpenFOAM, который предоставляет пользователю большую свободу, позволяя модифицировать существующие и реализовывать новые численные модели с использованием обширной библиотеки алгоритмов [6]. Но универсальность данного пакета и полный доступ ко всем настройкам имеет в качестве побочного эффекта сложность построения и настройки расчетных схем. Для отработки и отладки технологий работы с подвижными деформируемыми и скользящими сетками в данной работе рассматривается тестовая модель CP, основанная на новой конструкции золотника, обеспечивающей высокую точность управления расходом РЖ в широком диапазоне изменений перепадов давления на регуляторе (от 5 до 35 МПа) [7]. Сечение геометрической модели прототипа CP показано на рис. 1. Агрегат состоит из разъемного корпуса со штуцерами входа и выхода, гильзы с отверстиями золотникового дросселя, установленной в корпусе, а также из размещенного в гильзе подпружиненного золотника, который выполнен в форме стакана с двумя скользящими опорами на концах. Одна из опор имеет острую кромку, которая перекрывает проходные сечения отверстий в гильзе как показано на рис. 2.
1 2 з * 5
Рис. 1. Исследуемый стабилизатор расхода золотникового типа: 1 - штуцер входа, 2 - золотник, 3 - направляюгций штифт, 4 - пружина, 5 - опора пружины, б - штуцер
выхода, 7-гильза
Fig. 1. The studied flow rate regulator valve prototype: 1 - inlet fitting, 2 - plunger, 3-pin, 4 - spring, 5 - adjustment washer, б - outlet fitting, 7 - cylinder
Puc. 2. Увеличенная зона с отверстием золотникового дросселя
Fig. 2. Scaled-up zone with the piston throttle hole
Кольцевой канал, образованный внутренней поверхностью гильзы и поверхностью золотника между опорами, сообщен с внутренней полостью золотника системой отверстий. Часть отверстий расположена напротив отверстий золотникового дросселя и способствует взаимной компенсации энергии истекающих струй РЖ. Направляющий штифт служит для точного позиционирования отверстий в гильзе напротив отверстий в золотнике. Другие отверстия предназначены для уменьшения воздействия на золотник радиальной составляющей гидродинамической силы путем снижения перепада давления между его внутренней полостью и кольцевым каналом. Таким образом, золотник CP разгружается от неравномерного распределения давления на цилиндрических поверхностях золотника, возникающих из-за нестационарных эффектов течения РЖ. Вынесение пружины из зоны регулирования дополнительно снижает нежелательное возмущающее действие гидродинамических сил [7].
В процессе работы при изменении перепада давления на агрегате золотник перемещается и изменяет размер проходного сечения отверстий золотникового дросселя (рис. 2), тем самым сохраняя постоянной величину расхода. Например, при уменьшении давления на входе расход уменьшается относительно требуемого значения. Соответственно, гидродинамическая сила, действующая на золотник, уменьшается, а сила упругости поджатой пружины превосходит гидродинамическую силу в данном положении, и золотник смещается влево, приоткрывая отверстие дросселя (разжимая пружину) до того момента, пока не будет достигнуто новое положение равновесия. Результаты экспериментальных исследований CP показали, что описанную конструкцию возможно применять в широком диапазоне перепадов давления, однако отклонение величины массового расхода от номинала превышает заданное поле допуска [1]. Снижение статизма характеристики массового
расхода возможно путем выбора наилучшей геометрии золотника, для чего требуется проведение численного моделирования режима работы СР. Целью работы является тестирование возможностей свободно распространяемого пакета ОрепРОАМ по построению расчетной схемы для моделирования динамики стабилизатора расхода с помощью метода контрольного объема и подвижных сеток. В качестве основных задач были определены: выбор наиболее подходящих технологий динамических скользящих сеток [8]; построение расчетной схемы для сопряженной тестовой задачи Р81; расчет переходного режима движения золотника для прототипа СР; анализ нестационарных течений во внутренней полости агрегата; анализ устойчивости и производительности счета.
2. Постановка задачи и метод решения
Для ускорения методических расчетов с подвижной сеткой задача численного моделирования переходного режима работы СР рассматривается с учетом следующих допущений: РЖ считается несжимаемой вязкой средой, а ее течение считается ламинарным. Тепловые эффекты не рассматриваются. Учет сжимаемости и турбулентности течения жидкости, а также тепловых эффектов может быть выполнен средствами ОрепРОАМ в дальнейшем. С целью отработки методики построения расчетной схемы со скользящими сетками модель СР была упрощена до осесимметричной. Четыре отверстия золотникового дросселя при этом фактически заменены кольцевым каналом. Также кольцевыми каналами заменены отверстия в золотнике. Таким образом, золотник является в данной схеме набором из четырех тел, движущихся как единое целое. В силу сделанного допущения в области течения штифт был заменен неподвижной сферой и из рассмотрения исключены геометрическая модель пружины, канавки лабиринтных уплотнений опор золотника, а также зазоры между золотником и гильзой. Считается, что золотник имеет одну поступательную степень свободы вдоль продольной оси агрегата. Действие пружины на золотник заменено линейной вязкоупругой восстанавливающей силой. Такие упрощения искажают реальную картину течения РЖ в СР, однако сохраняют все основные особенности расчетной схемы, которые должны учитываться при построении подвижной сетки, давая возможность проводить методические расчеты с малыми затратами ресурсов.
Математическая задача Р81 включает две подзадачи: динамики РЖ и динамики золотника. Задача динамики РЖ включает уравнение неразрывности и уравнение сохранения импульса несжимаемой среды, записанные относительно неизвестных полей скоростей и давления с граничными и начальными условиями, рассмотренными ниже. Динамика золотника описывается уравнением его малых колебаний с учетом линейного демпфирования и действия нестационарной позиционной силы, вычисляемой путем интегрирования распределения давления РЖ по поверхности золотника и проецирования полученной силы на продольную ось агрегата. Нулевое
начальное положение золотника выбрано таким образом, что отверстие золотникового дросселя частично перекрыто. Начальная скорость золотника принята нулевой.
Численное моделирование производится методом контрольного объема в пакете OpenFOAM-extend-3.2 [6]. На этапе препроцессинга проводится генерация сетки, задание граничных условий и параметров расчета, получение начальных условий для динамического расчета.
Расчет переходного режима производится в течение заданного промежутка времени при известном постоянном перепаде давления. При этом золотник в процессе движения переходит от начального положения в положение равновесия, обусловленное равенством гидродинамической силы и силы упругости пружины. На этапе постпроцессинга исследуются параметры движения золотника и течения РЖ.
2.1 Подготовка расчетной сетки
В пакете OpenFOAM для работы с динамической сеткой имеется несколько классов. Объект класса dynamicTopoFvMesh позволяет изменять топологию сетки при больших деформациях методом mesquiteMotionSolver (аналог метода remeshing в ANSYS). Родительский класс dynamicFvMesh, позволяет производить добавление и удаление слоев сетки методом layerAdditionRemoval (аналог метода layering в ANSYS), а также деформацию сетки методом displacementLaplacian (аналог метода mesh smoothing в ANSYS). Для описания перемещений золотника использован последний из указанных методов. Однако из-за сильного вырождения ячеек прямое использование деформируемой сетки в зоне перекрытия отверстия золотникового дросселя невозможно. Решением этой проблемы стало разделение расчетной области на две подобласти -неподвижную (область течения между штуцером входа и гильзой) и подвижную (области течения во внутренней полости золотника, внутри гильзы и в выходном штуцере), - взаимодействующих в зоне отверстия с помощью скользящего интерфейса. Хотя объект класса dynamicFvMesh имеет встроенные интерфейсы для скользящих сеток slidinglnterface и attachDetach of boundaries, в данной расчетной схеме был использован обобщенный интерфейс GGI, описанный ниже.
Подготовка сетки в OpenFOAM производится в несколько этапов, как показано на рис. 3. На первом этапе для каждого из доменов по отдельности при помощи утилит blockMesh и snappyHexMesh строится пространственная расчетная сетка с использованием фоновой блочной сетки расчетной области и электронной геометрической модели обтекаемого тела в формате стереолитографии (*.stl). Для повышения качества сетки в области со сложной формой границы используется утилита surfaceFeatureExtract, позволяющая выделить особенности границы в отдельный файл. На втором этапе производится объединение расчетных доменов в один с помощью утилиты mergeMeshes. Для получения сетки осесимметричной задачи на третьем этапе проводится 58
вращение передней грани трехмерной стеки вокруг продольной оси агрегата на Г, в результате чего получается сетка в виде сектора цилиндра толщиной в одну ячейку, как показано на рис. 4. На четвертом шаге при помощи утилит гниоРгисИ и сгсгПсРгПсИ сетка разбивается на участки для задания граничных условий.
Рис. 3. Этапы создания расчетной сетки в OpenFOAKi
Fig. 3. OpenFOAKi mesh creation pipeline
Полученная сетка, показанная на рис. 4, имеет четыре уровня сгущения относительно базовой. Лучше всего разрешена область вокруг отверстия золотникового дросселя и вблизи поверхности обтекаемого тела. Более крупные ячейки находятся во внутренней полости золотника. Самая грубая сетка - в зонах, где отсутствуют подвижные элементы: в полостях штуцера входа и выхода и в зоне расположения пружины.
Рис. 4. Двухмерная расчетная сетка осесимметричной задачи Fig. 4. 2-dimensional axisymmetric mesh
2.2 Подготовка начальных условий
Для получения начальных полей давления и скорости РЖ перед расчетом переходного режима выполняется стационарный расчет с неподвижным золотником в заданном начальном положении. При этом в качестве граничных условий на всех поверхностях кроме входа, выхода и боковых граней задано условие жесткой непроницаемой стенки (для давления - zeroGradient, для скорости - нулевое значение). Для боковых граней задано условие циклической симметрии (wedge). На входе в CP задано полное давление 22,6 МПа, а на выходе - атмосферное давление (0,1 МПа). Таким образом перепад давления на агрегате составил 22,5 МПа.
Вычисления проводятся в ламинарной постановке с помощью решателя simpleFoam для несжимаемого течения, использующего алгоритм SIMPLE [9]. Параметры РЖ соответствуют авиационному маслу с плотностью 860 кг/м3 и кинематической вязкостью 21,5 ТО6 м2/с. Пример распределений модуля скорости и давления представлен на рис. 5. Видно, что в течении вблизи отверстия золотникового дросселя имеют место струйные эффекты.
Рис. 5. Полученные поля давлений (снизу) и скоростей (сверху) при стационарном расчете. Зеленым цветом выделены поверхности золотника, по которым интегрировалась гидродинамическая сипа
Fig. 5. Pressure (bottom) and velocity (top) fields at steady-state solution. Green lines -plunger surfaces where hydrodynamic force was calculated
График проекции на продольную ось суммарной гидродинамической силы, действующей на золотник, приведенный на рис. 6, показывает, что наличие струй нарушает монотонную сходимость стационарного решения: гидродинамическая сила пульсирует с периодом около 55 итераций, что затрудняет получение начальных условий. В качестве стационарного решения принимались осредненные по нескольким периодам значения рассчитываемых
полей. На основе среднего значения гидродинамической силы был проведен анализ сеточной сходимости, результаты которого показаны на рис. 7. В результате анализа была выбрана сетка из 138 ООО ячеек, на которой в дальнейшем проводились расчеты нестационарной задачи.
0 100 200 300 400 500 600 700 800
Итерация / Iteration
Рис. б. График гидродинамической сипы, действующей на золотник Fig. б. Plunger hydrodynamic force
aj is
СП 5 го 5 !- I 10
О) CT ^ -
< CL 5
и
...................±.
• 37012 ■ 38584
♦ 93005 á. 130000
* 138717
7Z
.............i..............:..
; ; I
О 50000 100000 150000 200000 250000 300000 350000
Количество ячеек / Number of cells
Рис. 7. График сходимости значения гидродинамической сипы по сетке Fig. 7. Mesh convergence
2.3 Подготовка граничных условий
Расчет переходного режима производится с помощью нестационарного решателя pimpleDyMFoam для несжимаемого течения, использующего гибридный алгоритм PISO-SIMPLE и возможность движения сетки [10]. На каждом шаге по времени делается около 20 итераций по давлению и 1-2 итерации по скорости. Основные трудности при построении расчетной схемы задачи связаны с заданием граничных условий на подвижном золотнике с
учетом действия суммарной гидродинамическои силы, а также моделирование перекрытия отверстия золотникового дросселя.
Положения подвижных и неподвижных границ течения показаны на рис. 8, а описание соответствующих им граничных условий представлено в табл. 1. Для течения несжимаемой жидкости значение давления задается отнесенным к плотности - кинематическое давление. На участке channelinlet, где входная скорость потока неизвестна, а известно только значение давления, задана комбинация граничных условий totalPressure и pressurelnletVelocity. В качестве параметра условия totalPressure задается значение полного давления. Граничные условия pressurelnletVelocity и pressurelnletOutletVelocity предназначены для входных/выходных отверстий, где определено значение давления. Условие movingWallVelocity задает условие прилипания на движущейся стенке.
с ha п л е l_wa 11, су I inderJ| ed_ws 11 channeLoutlet plunger_moving_wall cylindeioslip_wall
Рис. 8. Схема используемых названий участков граничных поверхностей при расчете.
Fig. 8. The scheme of using names of boundary surfaces patches. Табл. 1. Задание граничных условий OpenFOAKi Table 1. Boundary conditions in OpenFOAKi
Название участка Тип P U pointDisplacement
channel inlet patch totalPressure (26227) pressurelnletVelocity fixedValue (0 0 0)
channel outlet ggi ggi ggi fixedValue (0 0 0)
plunger inlet ggi ggi ggi; slip
plunger moving wall wall zeroGradient movingWallVelocity (0 0 0) sixDoFRigidBody
cylinder slip wall wall zeroGradient fixedValue (0 0 0) slip
channel wall, cylinder fixed wall, plunger side, cylinder finger wall zeroGradient fixedValue (0 0 0) fixedValue (0 0 0)
cylinder outlet patch fixedValue (1163) pressurelnletOutletVel ocity; fixedValue (0 0 0)
sidel, side 2 wedge wedge wedge wedge
В расчетной схеме отслеживание подвижных границ золотника (р1ш^ег_пкш^_\уа11 на рис. 8) реализовано с помощью объекта класса (ЗупаписРуМевЬ [11, 12]. Перестроение сетки производится на основании
перемещения граничных поверхностей подвижных тел. Модель расчета перемещения узлов сетки displacementLaplacian широко применяется для задач FSI в OpenFOAM, поскольку она не изменяет топологию сетки. Однако при больших перемещениях тела возможна потеря качества сетки и вырождение ее ячеек [13]. Используемая модель диффузии - directional. Величина смещения узлов сетки рассчитывается на основании решения уравнения Лапласа.
7-(yVdm) = 0
где dm - вектор перемещения узлов, м; у - коэффициент «диффузии» узлов сетки при движении. Смещение границ производится с использованием объекта класса sixDoFRigidBodyDisplacement (движение границ как абсолютно жесткого тела) с учетом различных механических связей, а также поверхностных гидродинамических сил [14].
Участок сетки, на котором происходит перекрытие отверстия при движении золотника, моделируется в версии OpenFOAM-extend с помощью технологии Generalized Grid Interface (GGI) [8, 15]. GGI позволяет установить связь между двумя смежными областями сетки, ячейки которых не совпадают на границе, в том числе и для сеток в зоне отверстия, где одна область сетки является неподвижной, а другая - подвижной (неподвижная и подвижная подобласти на рис 4), скользящей по границе первой. Сопряжение между поверхностями channel outlet и plunger inlet в зоне контакта (рис. 8) определяется набором перекрывающихся (где протекает поток) и неперекрывающихся (где задано условие стенки) ячеек как показано на рис. 9.
Рис. 9. Принцип работы GGI
Fig. 9. GGI principle of work
Для каждого временного шага зона перекрытия разбивается на участки (face Sets) между пересечениями ячеек двух смежных доменов. Через эти участки передаются рассчитываемые параметры потока. Для каждого участка оценивается его весовой коэффициент по отношению к содержащей его ячейке
в каждом из смежных доменов (зоны inlet и outlet) и вводятся условия совместности:
п m
к=1 к=1 п m
^ W^1 = 1, ^ = 1;
к=1 к=1
где ф - численный поток; Wk - веса; i - номер ячейки inlet; / - номер ячейки outlet; m - число ячеек inlet соприкасающихся с j-й ячейкой; п - число ячеек outlet соприкасающихся с /-й ячейкой; площадь к ячейки зоны outlet; S'k -площадь к ячейки зоны inlet; S0I - площадь зоны пересечения outlet и inlet.
3. Результаты расчета
Расчет переходного режима проводился для перепада давления 22,5 МПа. Рассматривался золотник массой (с учетом приведенной массы пружины) 0,046 кг. Свободная длина пружины выбрана равной 0,060 м. Параметры вязкоупругой восстанавливающей силы, моделирующей пружину: жесткость -400 Н/м, коэффициент демпфирования - 40 Н с/м. Собственная частота колебаний данной колебательной системы в вакууме (без демпфирования) составила 14,8 Гц (период колебаний 0,068 с). Демпфирование обуславливает апериодическое свободное движение золотника в вакууме. В начальный момент времени золотник расположен так, что отверстие дросселя перекрыто на 30 % (0,0006 м), длина пружины при этом равна - 0,035 м. В результате расчета были получены картины нестационарного течения, а также характеристики гидродинамического нагружения золотника. На рис. 10 представлены зависимости от времени перемещений золотника, а также обобщенных гидродинамической и восстанавливающей сил, действующих на золотник. Движение золотника включает переход в положение равновесия - затухающие колебания с периодом около 0,0055 с (частота около 180 Гц) и вибрацию вблизи положения равновесия с периодом около 0.0025 с (частота около 400 Гц), вызываемые нестационарностью течения в перекрытом дросселе. Максимальное смещение золотника составило 0,44 мм, при ширине всего отверстия 2 мм. В момент окончания расчета (0.016 с) отверстие дросселя перекрыто на 11 % (0,0002 м), а длина пружины составляет 0,0346 мм. В процессе перекрытия отверстия скорость в струе возрастает, и она притягивается к наклонной поверхности золотника как показано на рис. 11. На рис. 12 представлены характерные фазы течения, показывающие изменение модуля скорости в конце переходного режима за период времени порядка 0,001 с. Видно, что вблизи наклонной кромки золотника струйное течение
эволюционирует, порождая вихревые структуры в проточной части золотника. Обнаруженные струйные эффекты и вихревые течения могут возникать и при пространственном течении в проточной части реального агрегата, приводя к существенным изменениям распределения давления по поверхности золотника, возникновению дисбаланса радиальной составляющей гидродинамических сил и могут быть причиной различия расчетных и экспериментальных характеристик СР.
3° г-----—
5 0
£ 2
£ ÎD
a. 0.1
4-Î* ai
С <D E * 1 0.2
<u ч
u с
tv ai
Q. s 0 3
1Л 1
ТЭ a>
L_ 3"
£L> ai aj S 0.4
С ai
ZI Q_
Cl ai С 0
\
y /
y
20
оГ
и
о 10
LL.
X о
§ -10
и
-20
" V/'MM
0.005 о. оа 0.015 Время, c/Time, s
0.02
0.005 0 01 0 015 Время, с/Time, s
О 02
Рис. 10. Результаты расчета переходного режима: а - перемещение золотника; b -обобщенные силы, действующие на золотник: 1 - гидродинамическая сила, 2 -восстанавливающая сипа Fig. 10. Transient solution graphs: a - plunger displacement: b - forces applied on plunger: 1 - hydrodynamic forces, 2 - restoring force
U Magnilude, ffi/s
10J.9
Puc. 11. Поле скоростей вблизи отверстия золотникового дросселя в начальный (а) и конечный моменты времени (Ь)
Fig. 11. Velocity field in zone of the piston throttle hole at initial (a) and final (b) time steps
Puc. 12. Поле скоростей в золотнике в разные моменты времени: 1-0.012 с. 2-0.0122 с. 3-0.0123 с. 4 -0.0126 с. 5-0.0127 с, 4-0.0129 с.
Fig. 12. Velocity.field at time moments: 1 - 0.012 s, 2-0.0122 s, 3 - 0.0123 s, 4 - 0.0126 s,
5-0.0127 s. 4-0.0129 s.
Для решения системы уравнений для давления применяется метод GAMG (Geometric agglomerated algebraic multigrid solver). Для сглаживания используется метод Гаусса - Зейделя (предсглаживающих итераций не делалось, постсглаживающих итераций - 2). Система уравнений для прогнозов скоростей решается методом BiCGStab (Biconjugate gradient stabilized method) с предобуславливателем, основанным на DiLU-факторизации. Итерации выполнялись до достижения точности Ю-7. Шаг по времени выбирался автоматически из условия Сотах < 2. Характерная величина шага по времени -Ю-7 с. Время расчета переходного режима длительностью 0,016 с составило 128 ч на 1 ядре (Intel(R) Xeon(R) CPU Х5670, 2.93GHz). Для проведения расчетов был использован открытый облачный сервис UNIHUB разработки ИСПРАН [16].
4. Заключение
Проделанная работа позволяет сделать вывод о том, что пакет с открытым исходным кодом OpenFOAM в версии extend может быть использован в качестве альтернативы коммерческим пакетам программ вычислительной гидродинамики (CFD) для расчета переходных режимов в агрегатах гидравлической автоматики. В данном пакете имеются все необходимые
средства для построения расчетных схем с деформируемыми сетками, описывающими движение золотников и перекрытие отверстий золотниковых дросселей. Подобные расчетные схемы также можно использовать при моделировании в OpenFOAM двигателей внутреннего сгорания, гидравлических дросселей и других элементов пневмогидросистем, имеющих в конструкции подвижные элементы: поршни, клапаны и пр. Расчеты модельной задачи показали, что наличие струй вблизи наклонной кромки золотника нарушает монотонную сходимость стационарного решения при вычислении начальных условий и существенно влияет на динамику движения золотника. В целом можно сделать вывод, что технология построения расчетной схемы успешно отлажена в двумерной постановке и готова к переносу на трехмерный случай. Дальнейшая работа предполагает полноценное моделирование реальной конструкции CP в условиях, приближенных к экспериментальному исследованию [1]. Сравнение результатов расчета с экспериментом позволит произвести верификацию расчетной схемы и использовать ее в дальнейшем для оптимизации внутренней геометрии CP с целью уменьшения статизма его расходной характеристики. Работа поддержана грантом РФФИ (проект № 17-08-01468 А).
Список литературы
[1]. Копков Г. А., Кучин А. П., Новиков А. Е., Иванов М. Ю., Реш Г. Ф., Антонов Д. С. Стабилизаторы расхода для синхронизации перемещения исполнительных органов систем летательных аппаратов. Научно-технический юбилейный сборник. КБ химавтоматики: В 3 томах. Под ред. В. С. Рачука. Воронеж: Кварта, 2012. Т.1., стр. 219-223.
[2]. Салман М. И., Попов Д. Н. Компьютерное исследование и расчет гидродинамических нагрузок на золотник. Наука и образование. МГТУ им. Н. Э. Баумана. Электрон, журн. 2012. № 10. С. 79-92. DOI: 10.7463/1112.0491484. Режим доступа: http://technomag.bmstu.ru/doc/491484.html (дата обращения: 02.01.2017).
[3]. Широкова К. А., Целшцев В. А., Целшцев Д. В., Галлямов Ш. Р. Численное моделирование потоков в струйно-золотниковом гидроусилителе. Вестник УГАТУ. 2008. Т.11, №2 (29). С. 55-59.
[4]. ANSYS Fluent: CFD Simulation: сайт. Режим доступа: http://www.ansys.com/ products/fluids/ansys-fluent (дата обращения: 29.12.2016).
[5]. ANS YS CFX: Turbomachinery CFD Simulation: сайт. Режим доступа: http://www.ansys.com/products/fluids/ansys-cfx (дата обращения: 29.12.2016).
[6]. OpenFOAM. Free CFD Software: сайт. Режим доступа: http://openfoam.org/ (дата обращения: 29.12.2016).
[7]. Регулятор расхода: пат. 2548613 Рос. Федерация: МПК G05D 7/01 (2006.01). Дергачев А. А., Иванов М. Ю., Копков Г. А., КучинА. П., Новиков А. Е., Реш Г. Ф., Синявин В.Г.; заявитель и патентообладатель: Открытое акционерное общество «Военно-промышленная корпорация «Научно-производственное объединение машиностроения»-№2014102669/28; заявл. 29.01.14; опубл. 20.04.15, Бюл. № 11.
[8]. Beaudoin M., Jasak H. Development of a Generalized Grid Interface for Turbomachinery simulations with OpenFOAM. Open Source CFD International Conference 2008, Berlin, Germany.
[9]. Ferziger J.H., Peric M. Computational Methods for Fluid Dynamics, Springer, 3rd Ed., 2001. 431 p.
[10]. OpenFOAM guide. The PIMPLE algorithm in OpenFOAM - OpenFOAMWiki: сайт. Режим доступа: http://openfoamwiki.net/index.php/OpenFOAM_guide/The_PIMPLE_algorithm_ in_OpenFOAM (дата обращения: 29.12.2016).
[11]. OpenFOAM dynamicFvMesh. С++ Source Code Guide. Режим доступа: http://www.openfoam.com/documentation/cpp-guide/html/a05056.html (дата обращения:
29.12.2016).
[12]. Gonzlez А.О. Mesh motion alternatives in OpenFOAM. Tech. rep.. Institution of Applied Mechanics at Chalmers Technical University, Gteborg, Sweden, 2009. http://www.tfd.chalmers.se/~hani/kurser/OS_CFD_2009/Andreu01iverGonzalez/ProjectRep ort_Corrected.pdf (дата обращения: 29.12.2016).
[13]. Kotsur О., Scheglov G., Leyland P. Verification of Modelling of Fluid-structure Interaction (FSI) Problems Based on Experimental Research of Bluff Body Oscillations in Fluids. 29th Congress of the International Council of the Aeronautical Sciences, 7-12 September 2014, Saint-Petersburg, Russia. Paper ICAS2014-2.6.3-2014_0953.
[14]. OpenFOAM sixDoFRigidBodyMotion. С++ Source Code Guide. Режим доступа: http://cpp.openfoam.org/v3/a09731.html (дата обращения: 29.12.2016).
[15]. Jasak Н. NUMAP-FOAM Summer School, Zagreb 2-15 Sep 2009 General Grid Interface Theoretical Basis and Implementation. Режим доступа: http://powerlab.fsb.hr/ ped/kturbo/OpenFOAM/SummerSchool2009/lectures/TurboGGI.pdf (дата обращения:
03.01.2017).
[16]. UniHUB: сайт. Режим доступа: http://www.unicluster.ru/unihub.html (дата обращения: 29.12.2016).
Numerical simulation of the flow rate regulator valve using
OpenFOAM
V.G. Melnikova <vg-melnikova@yandex.ru>
O.S. Kotsur <oskotsur@gmail.com> G.A. Shcheglov <shcheglov_ga@bmstu.ru> Federal state budgetary institution of higher professional education «Bauman Moscow State Technical University (National research university of technology) », ul. Baumanskaya 2-ya, 5/1, Moscow, 105005, Russia
Abstract. The results of methodical investigation, aimed on testing of the performances of free opensource CFD toolbox OpenFOAM in the field of simulation of hydraulic units' dynamics using Finite Volume Method and dynamic meshes are presented. The following key features are reviewed: the choice of appropriate algorithms managing dynamic sliding meshes; design of the model case for the FSI-problem of interaction between a moving regulating element (plunger) and power fluid; plunger dynamics simulation for the prototype of the flow rate regulator of new design; transient analysis of regulator's internal flow; analysis of stability and computational efficiency. As the example simplified axisymmetric regulator model with incompressible power fluid is considered. The main steps of model case preparation are
described. The model case is designed to simulate the problem of plunger equilibration under hydrodynamic forces and spring reaction. Results are given also for a preliminary steady-state simulation with fixed plunger, they have been used as initial conditions. Detailed description is given for the methods of mesh motion simulation, which follows plunger, as well as technology of sliding meshes (GGI), used for the piston throttle hole overlapping. Results are given for transient simulation: velocity, pressure fields, forces graphs, acting on the plunger, its displacement over time. The mechanism of jet streams and vertical flows creation in regulator's flow channel is described. Methodical investigation, held in this paper, confirms that OpenFOAM in the "extend" version can be successfully used as an alternative for commercial CFD codes, as it contains all necessary tools to create and simulate cases incorporating dynamic meshes. It provides means for simulation of transient problems of hydraulic units which have moving parts.
Keywords: CFD; numerical simulation; OpenFOAM; open source CFD codes; flow rate regulator.
DOI: 10.15514/ISPRAS-2017-29(1 )-4
For citation: Melnikova V.G., Kotsur O.S., Shcheglov G.A. Numerical simulation of the flow rate regulator valve using OpenFOAM. Trudy ISP RAN/Proc. ISP RAS, vol. 29, issue 1, 2017, pp. 53-70 (in Russian). DOI: 10.15514/ISPRAS-2017-29(l)-4
References
[1]. Kopkov G.A., Kuchin A.P., Novikov A.E., Ivanov M.Ju., Resh G.F., Antonov D.S. Flow regulators for synchronization of executive devices of aircraft systems. Nauchno-tehnicheskij jubilejnyj sborrtik. KB himavtomatiki [Anniversary proc. DD Chemical automatics], 2012, vol.1, pp. 219-223 (in Russian).
[2]. Salman M.I., Popov D.N. Computer study and calculation of hydrodynamic loads on the valve. Nauka i obrazovanie. MGTU im. N.E. Baumana. Jelektron. zhurn [Science and education. BMSTU. Electronic journal], 2012, no. 10, pp. 79-92. (in Russian). DOI: 10.7463/1112.0491484. Available at: http://technomag.bmstu.ru/doc/491484.html (Accessed 2 January 2017)
[3]. Shirokova K.A., Celishhev V.A., Celishhev D.V., Galljamov Sh.R. Numerical simulation of the flows in hydraulic amplifier with a sliding valve. Vestnik UGATU [Proc. of USATU], 2008, vol. 11, no. 2 (29), pp. 55-59 (in Russian).
[4]. ANSYS Fluent: CFD Simulation. Available at: http://www.ansys.com/products/fluids/ ansys-fluent (Accessed 29 December 2016).
[5]. ANSYS CFX: Turbomachinery CFD Simulation. Available at: http://www.ansys.com/ products/fluids/ansys-cfx (Accessed 29 December 2016).
[6]. OpenFOAM. Free CFD Software. Available at: http://openfoam.org/ (Accessed 29 December 2016).
[7]. Dergachev A. A., Ivanov M. Ju., Kopkov G. A., Kuchin A. P., Novikov A. E., Resh G. F., Sinjavin V. G. Regulator rashoda [Flow rate regulator]. Patent RF, no. 2548613, 2006.
[8]. Beaudoin M., Jasak H. Development of a Generalized Grid Interface for Turbomachinery simulations with OpenFOAM. Open Source CFD International Conference 2008, Berlin, Germany.
[9]. Ferziger J. H, Peric M., Computational Methods for Fluid Dynamics, Springer, 3rd Ed., 2001.
69
[10]. OpenFOAM guide. The PIMPLE algorithm in OpenFOAM - OpenFOAMWiki. Available at: http://openfoamwM.neiyindex.php/OpenFOAM_guide/The_PIMPLE_algori1iim_in_OpenF OAM (Accessed 29 December 2016).
[11]. OpenFOAM dynamicFvMesh. C++ Source Code Guide. Available at: http://www.openfoam.com/documentation/cpp-guide/html/a05056.html (Accessed 29 December 2016).
[12]. Gonzlez A.O. Mesh motion alternatives in OpenFOAM, tech. rep., Institution of Applied Mechanics at Chalmers Technical University, Gteborg, Sweden, 2009. Available at: http://www.tfd.chalmers.se/~hani/kurser/OS_CFD_2009/Andreu01iverGonzalez/ProjectRep ort_Corrected.pdf (Accessed 29 December 2016).
[13]. Kotsur O., Scheglov G., Leyland P. Verification of Modelling of Fluid-structure Interaction (FSI) Problems Based on Experimental Research of Bluff Body Oscillations in Fluids. 29th Congress of the International Council of the Aeronautical Sciences, 7-12 September 2014, Saint-Petersburg, Russia. Paper ICAS2014-2.6.3-2014_0953.
[14]. OpenFOAM sixDoFRigidBodyMotion. C++ Source Code Guide. Available at: http://cpp.openfoam.org/v3/a09731.html (Accessed29 December 2016).
[15]. Jasak H. NUMAP-FOAM Summer School, Zagreb 2-15 Sep 2009 General Grid Interface Theoretical Basis and Implementation. Available at: http://powerlab.fsb.hr/ped/ kturbo/OpenFOAM/SummerSchool2009/lectures/TurboGGI.pdf (Accessed 03 January 2017).
[16]. UniHUB. Available at: http://www.unicluster.ru/unihub.html (Accessed 29 December 2016).
