НАУЧНО-ТЕХНИЧЕСКИИ ВЕСТНИК ИНФОРМАЦИОННЫХ ТЕХНОЛОГИИ, МЕХАНИКИ И ОПТИКИ январь-февраль 2018 Том 18 №1 ISSN 2226-1494 http://ntv.i1mo.ru/
SCIENTIFIC AND TECHNICAL JOURNAL OF INFORMATION TECHNOLOGIES, MECHANICS AND OPTICS January-February 2018 Vol. 18 No 1 ISSN 2226-1494 http://ntv.ifmo.ru/en
УДК 532.529
ПРИМЕНЕНИЕ СХЕМЫ С НАСТРАИВАЕМЫМИ ДИССИПАТИВНЫМИ СВОЙСТВАМИ К РАСЧЕТУ ТЕЧЕНИЙ ГАЗА С РАЗВИТИЕМ НЕУСТОЙЧИВОСТИ НА КОНТАКТНОЙ ГРАНИЦЕ
Д.В. Садин"
a Военно-космическая академия имени А.Ф. Можайского, Санкт-Петербург, 197198, Российская Федерация Адрес для переписки: [email protected] Информация о статье
Поступила в редакцию 26.11.17, принята к печати 27.12.17 doi: 10.17586/2226-1494-2018-18-1-153-157 Язык статьи - русский
Ссылка для цитирования: Садин Д. В. Применение схемы с настраиваемыми диссипативными свойствами к расчету течений газа с развитием неустойчивости на контактной границе // Научно-технический вестник информационных технологий, механики и оптики. 2018. Т. 18. № 1. С. 153-157. doi: 10.17586/2226-1494-2018-18-1-153-157
Аннотация
Выполнено тестирование схемы с настраиваемыми диссипативными свойствами применительно к структурно-сложным задачам течений газа с развитием неустойчивости на контактной границе. Схема реализована с расщеплением на градиентные и деформационные члены уравнений Эйлера, аппроксимированные центральными разностями с TVD-ограничителями искусственной вязкости, а также на конвективные члены с TVD-ограничителями потоков. Схема дополнена TVD-методом Рунге-Кутта второго порядка по времени. Анализируется разрешающая способность и эффективность предлагаемой схемы в сравнении с некоторыми высокоточными схемами на примере решения задач импульсного сжатия и расширения газа, а также двойного маховского отражения. Схема с настраиваемыми диссипативными свойствами по соотношению цена (затраты машинного времени) - качество (разрешающая способность) находится на уровне, а для некоторых задач превосходит современные высокоточные схемы. Схема может быть рекомендована для численных исследований сложных ударно-волновых и вихревых течений с развитием неустойчивости на контактной границе. Ключевые слова
схема с настраиваемыми диссипативными свойствами, расчет течений газа, развитие неустойчивости, контактная граница
APPLICATION OF SCHEME WITH CUSTOMIZABLE DISSIPATIVE PROPERTIES FOR GAS FLOW CALCULATION WITH INTERFACE INSTABILITY EVOLUTION
D.V. Sadin"
a Mozhaisky Military Space Academy, Saint Petersburg, 197198, Russian Federation Corresponding author: [email protected] Article info
Received 26.11.17, accepted 27.12.17 doi: 10.17586/2226-1494-2018-18-1-153-157 Article in Russian
For citation: Sadin D.V. Application of scheme with customizable dissipative properties for gas flow calculation with interface instability evolution. Scientific and Technical Journal of Information Technologies, Mechanics and Optics, 2018, vol. 18, no. 1, pp. 153-157 (in Russian). doi: 10.17586/2226-1494-2018-18-1-153-157
Abstract
We have performed testing of the scheme with customizable dissipative properties as applied to problems for complicated gas flow structure with the evolution of interface instability. The scheme is implemented by the splitting into two phases. The first one uses the central differences of both deformation and gradient terms of the Euler equations with artificial viscosity TVD limiters. TVD type reconstructions for convective terms with flux limiters are used in the second phase. A two-order TVD Runge-Kutta algorithm is applied to march the solution in time. We address the issue of numerical resolution and efficiency of the proposed scheme in comparison with some high resolution schemes through numerical examples: the implosion problem, the explosion problem and the double Mach reflection problem. The scheme with customizable dissipative properties in terms of cost (cost of machine time) - quality (resolution) is at the level of modern high order schemes. For some numerical examples, the proposed scheme is superior. The scheme can be recommended for numerical studies of both complex shock wave and vortex flows with the evolution of interface instability. Keywords
scheme with customizable dissipative properties, gas flow calculation, instability evolution, interface
ПРИМЕНЕНИЕ СХЕМЫ С НАСТРАИВАЕМЫМИ ДИССИПАТИВНЫМИ СВОЙСТВАМИ.
Введение
Пространственные потоки газа, содержащие в начальные моменты времени поверхности раздела газов с различными термодинамическими свойствами (контактные разрывы), с течением времени приводят к сложным ударно-волновым и вихревым структурным течениям с развитием различных видов неус-тойчивостей (Рихтмайера-Мешкова, Кельвина-Гельмольца, Релея-Тейлора и др.). Численное моделирование указанных газодинамических явлений в определенном смысле является серьезным испытанием для применяемой дискретной модели. С одной стороны, дискретная модель должна обладать достаточной схемной или искусственной вязкостью для подавления осцилляций на разрывах, с другой стороны, целесообразным свойством является малая порожденная дискретизацией диффузия в областях гладкости и вихреобразования. В настоящее время имеется обширная литература по численному моделированию структурно сложных газодинамических течений. Интенсивно развиваются подходы, основанные на приближенном решении задачи Римана [1, 2], WENO (Weighted Essentially Non-Oscillatory)-реконструкции [3-6], компактных схемах [7-9], адаптивной искусственной диффузии [10-16], их комбинациях и др. Следует отметить, что стремление достичь высокой разрешающей способности схемы, как правило, сопряжено с заметным усложнением алгоритма и увеличением затрат вычислительных ресурсов. В связи с этим разработка алгоритмически простых и экономичных схем с точки зрения критерия (соотношения) цены (затрат машинного времени) и качества (разрешающей способности) не теряет своей актуальности.
Газодинамические течения в рамках сплошной среды в общем случае описываются уравнениями Навье-Стокса. Численное моделирование вязкого теплопроводного газа является корректным в случае, когда численная диффузия существенно меньше физической. Исходя из этого, в настоящей работе, как и в цитируемых работах [4, 13, 17], изучается уровень диссипативных свойств схем, аппроксимирующих градиентные, деформационные и конвективные члены уравнений Эйлера в консервативном виде. Такой подход используется для оценки разрешения сеток для корректного моделирования в рамках Навье-Стокса при заданном числе Рейнольдса [4].
Целью работы является исследование разрешающей способности схемы с настраиваемыми дисси-пативными свойствами [14, 15], дополненной дискретизацией по времени методом Рунге-Кутта второго порядка аппроксимации (customizable dissipative properties - CDP2) и ее чувствительности к развитию физической неустойчивости на контактной границе в сравнении с некоторыми высокоточными схемами.
Схема с настраиваемыми диссипативными свойствами строится путем расщепления по физическим процессам на два этапа. На первом из них отбрасываются конвективные члены, а градиентные и деформационные члены уравнений Эйлера аппроксимируются центральными разностями. Для монотонизации численного решения в схему вводится скалярная (как аддитивная добавка к давлению .Ри+1/2 + вП+ш) искусственная вязкость с ТУО-ограничителем вязкости уу типа Христенсена [10, 15]:
где B > 0 - коэффициент искусственной вязкости; у - показатель адиабаты газа; p, р, v - давление, плотность и проекция вектора скорости газа на соответствующую ось; k, n - индексация сетки по времени и пространству соответственно, полуцелый индекс относится к границе ячейки сетки. На втором этапе рассчитываются окончательные значения искомых функций путем взвешенной линейной комбинации проти-вопоточной и центральной аппроксимаций конвективных членов. TVD-ограничитель потоков уf является
весовым множителем для центральной разности, а величина 1 - у f - для аппроксимации против потока.
Обозначим описанные два этапа как разностный оператор q(1) = qk + tL (qk), где L - оператор пространственной аппроксимации. Используя TVD-метод Рунге-Кутта второго порядка [18], получим схему второго порядка по времени и пространству на гладких решениях O (h2 + т2) (h, т - шаги по пространству и времени; k - индекс временного слоя): qk+1 = 0,5(qk + q(1)) + 0,5tL(q(1)). Здесь q - вектор
консервативных искомых функций.
Настройка диссипативных свойств схемы для определенных классов задач выполняется заданием коэффициента искусственной вязкости (для большинства случаев и всюду ниже B = 1), а также выбором ограничителя искусственной вязкости уv (для всех задач ниже - ограничитель VAN LEER
Уv (r) = (r + И) / (1 + r)) и ограничителя потоков у f. Расчеты выполнялись с числом Куранта 0,4.
Метод расчета
Примеры решения задач
Для проверки разрешающей способности схемы CDP2 применительно к явлениям ударно-волновых течений газа с развитием неустойчивости на контактной границе приведем примеры трех широко цитируемых тестовых задач: импульсного сжатия (implosion problem), расширения (explosion problem) и двойного маховского отражения (double Mach reflection problem). Вычислительные свойства схемы CDP2 будем сравнивать со следующими высокоточными схемами: центральная схема с ограничителем (centered scheme with limiter - JT) [19], схема с кусочно-параболической реконструкцией (piecewise parabolic - PPM) [20], взвешенная существенно неосциллирующая схема (weighed essentially nonoscillatory - WENO) [3] и схема Году -нова-Колгана-Родионова - ГКР [13]. Газ полагается идеальным и калорически совершенным с показателем адиабаты у = 1,4. Начальные условия и решения приводятся в безразмерном виде.
Задача импульсного сжатия газа (implosion problem). Внутри квадрата (х, y)е(0;0,3)х(0;0,3) находится треугольная область с угловыми точками (0;0), (0;0,15) ,
(0,15;0) . Начальные условия для плотности р1 и давленияp газа соответственно: снаружи треугольника - p1o = po = 1 ; внутри треугольника - p1i = 0,125 и pt = 0,14 . Газ находится в состоянии покоя, начальные скорости равны нулю. На всех четырех границах заданы краевые условия отражения.
Рис. 1. Результаты решения задачи импульсного сжатия в момент времени 2,5: CDP2 (а); JT из [17] (б); WENO5 из [17] (в); PPM из [17] (г).
Плотность нанесена в виде 31 контурной линии от 0,35 до 1,1; скорости показаны стрелками
В результате распада начального разрыва формируется сходящаяся к началу координат ударная волна, а в обратном направлении распространяется волна разрежения. После отражения ударная волна начинает двигаться в противоположном направлении, испытывая с течением времени многократные взаимодействия с контактным разрывом, стенками и волнами разрежения. Через некоторое время на контактной границе развивается неустойчивость.
На рис. 1 приведены результаты расчетов структуры течения в момент времени tf = 2,5 на сетке
400x400. При расчете с использованием схемы CDP2 (ограничитель потоков MUSCL (Monotonie Upstream-Centred Scheme), у f = max [ min (2,2 r ,(1 + r )/2) ,0 ]) и WENO5 формируется крупномасштабная вихревая структура, которая затем отрывается от контактной поверхности и движется в диагональном направлении (рис. 1, а, в, соответственно). Центральная схема с ограничителем JT обладает большей схемной вязкостью (рис. 1, б). Течение на контактной границе, рассчитанное по схеме PPM, не обладает симметрией и упорядоченной структурой (рис. 1, г).
Задача импульсного расширения газа (explosion problem). Рассматривается цилиндрически симметричная двумерная задача. Центр круга радиусом 0,4 находится в начале координат. В начальный момент времени газ находится в состоянии покоя, плотность и давление р. = 1, p¡ = 1 заданы внутри круга и ро = 0,125, pt = 0,1 - снаружи. Ввиду симметрии решение ищется в квадранте (х, y) е (0;1,5)x (0;1,5) на сетке 400x400 в момент времени tf = 3,2 . На левой и нижней границе заданы условия отражения, на верхней и правой - экстраполяция искомых функций.
Эта задача является проверкой чувствительности схемы к развитию неустойчивости на границе раздела двух газов [17]. На рис. 2 представлены результаты расчетов, выполненные по различным схемам. Схемы CDP2 (ограничитель потоков MUSCL) и PPM более чувствительны к развитию неустойчиво -сти на контактной поверхности, чем WENO5 и JT.
На рис. 3 показаны фрагменты результатов расчетов [2;3]x[0;0,5] в виде 30 изолиний плотности
через равные интервалы от 1,5 до 22,9705. Разрешение сеток сравниваемых схем одинаковое - 1/480. Схемы CDP2 (ограничитель потоков VAN LEER) и WENO5 (рис. 3, а, б, соответственно) демонстрируют примерно одинаковые диссипативные свойства при существенно меньших затратах машинного времени
ПРИМЕНЕНИЕ СХЕМЫ С НАСТРАИВАЕМЫМИ ДИССИПАТИВНЫМИ СВОЙСТВАМИ...
предлагаемой дискретной модели. Схема ГКР не позволяет выявить развитие неустойчивости Кельвина-Гельмгольца на сетке указанного разрешения (рис. 3, в).
1,5
0,5 0,4 0,3 0,2 0,1 0
а б в г
Рис. 2. Результаты решения задачи импульсного расширения в момент времени 3,2: CDP2 (а); JT из [17] (б); WENO5 из [17] (в); PPM из [17] (г). Плотность нанесена в виде 27 контурных линий от 0,08 до 0,21; скорости показаны стрелками
0,5 0,4 0,3 0,2 0,1 0
2
2,4
а б в
Рис. 3. Результаты решения задачи двойного маховского отражения в момент времени 0,2: ОЭР2 (а); WENO5 из [4] (б); ГКР из [13] (в). Плотность нанесена в виде 30 контурных линий
от 1,5 до 22,9705
Заключение
Схема с настраиваемыми диссипативными свойствами продемонстрировала удовлетворительные вычислительные свойства для задач ударно-волновых течений газа с развитием неустойчивости на контактной границе. По соотношению цены (затрат машинного времени) и качества (разрешающей способности) схема СЭР2 находится на уровне, а для некоторых задач превосходит современные высокоточные схемы. Дальнейшие исследования направлены на применение схемы СЭР2 для численного моделирования многокомпонентных многофазных потоков.
Литература
Cocchi J.P., Saurel R., Loraud J.C. Treatment of interface
problems with Godunov-type schemes // Shock Waves. 1996.
V. 5. N 6. P. 347-357. doi: 10.1007/pl00003878
Toro E.F. Riemann Solvers and Numerical Methods for Fluid
Dynamics. Berlin, Springer-Verlag, 2009, 724 p.
Jiang G.-S., Shu C.-W. Efficient implementation of weighted
ENO schemes // Journal of Computational Physics. 1996.
V. 126. P. 202-228. doi: 10.1006/jcph.1996.0130
Shi J., Zhang Y.-T., Shu C.-W. Resolution of high order WENO
schemes for complicated flow structures // Journal of
Computational Physics. 2003. V. 186. N 2. P. 690-696. doi:
10.1016/S0021-9991(03)00094-9
Coralic V., Colonius T. Finite-volume WENO scheme for viscous compressible multicomponent flows // Journal of Computational Physics. 2014. V. 274. P. 95-121. doi: 10.1016/j.jcp.2014.06.003
Булат М.П., Волобуев И.А., Волков К.Н., Пронин В.А. Численное моделирование регулярного и маховского отражения ударной волны от стенки // Научно-технический вестник информационных технологий, механики и оптики. 2017. Т. 17. № 5. С. 920-928. doi: 10.17586/2226-1494-201717-5-920-928
Толстых А.И. О семействах компактных аппроксимаций 4-го и 5-го порядков с обращением двухточечных операторов для уравнений с конвективными членами // ЖВМ и МФ. 2010. Т. 50. № 5. С. 894-907.
References
Cocchi J.P., Saurel R., Loraud J.C. Treatment of interface
problems with Godunov-type schemes. Shock Waves, 1996,
vol. 5, no. 6, pp. 347-357. doi: 10.1007/pl00003878
Toro E.F. Riemann Solvers and Numerical Methods for Fluid
Dynamics. Berlin, Springer-Verlag, 2009, 724 p.
Jiang G.-S., Shu C.-W. Efficient implementation of weighted
ENO schemes. Journal of Computational Physics, 1996,
vol. 126, pp. 202-228. doi: 10.1006/jcph.1996.0130
Shi J., Zhang Y.-T., Shu C.-W. Resolution of high order
WENO schemes for complicated flow structures Journal of
Computational Physics, 2003, vol. 186, no. 2, pp. 690-696.
doi: 10.1016/S0021-9991(03)00094-9
Coralic V., Colonius T. Finite-volume WENO scheme for viscous compressible multicomponent flows. Journal of Computational Physics, 2014, vol. 274, pp. 95-121. doi: 10.1016/j.jcp.2014.06.003
Bulat M.P., Volobuev I.A., Volkov K.N., Pronin V.A. Numerical simulation of regular and Mach reflection of shock wave from the wall. Scientific and Technical Journal of Information Technologies, Mechanics and Optics, 2017, vol. 17, no. 5, pp. 920-928. (In Russian) doi: 10.17586/22261494-2017-17-5-920-928
Tolstykh A.I. On families of compact fourth- and fifth-order approximations involving the inversion of two-point operators for equations with convective terms. Computational
0
8. Shen Y.-Q., Zha G.-C. Generalized finite compact difference scheme for shock/complex flow field interaction // Journal of Computational Physics. 2011. V. 230. N 12. P. 4419-4436. doi: 10.1016/j.jcp.2011.01.039
9. Михайловская М.Н., Рогов Б.В. Монотонные компактные схемы бегущего счета для систем уравнений гиперболического типа // ЖВМ и МФ. 2012. Т. 52. № 4. С. 672-695.
10. Christensen R.B. Godunov Methods on a Staggered Mesh - An Improved Artificial Viscosity. Technical Report UCRL-JC-105269. 1990. 11 p.
11. Shankar S.K., Kawai S., Lele S. Numerical simulation of multicomponent shock accelerated flows and mixing using localized artificial diffusivity method // Proc. 48th AIAA Aerospace Sciences Meeting Including the New Horizons Forum and Aerospace Exposition. Orlando, USA, 2010. Art. 2010-0352.
12. Kurganov A., Liu Y. New adaptive artificial viscosity method for hyperbolic systems of conservation laws // Journal of Computational Physics. 2012. V. 231. N 24. P. 8114-8132. doi: 10.1016/j.jcp.2012.07.040
13. Тагирова И.Ю., Родионов А.В. Применение искусственной вязкости для борьбы с «карбункул»-неустойчивостью в схемах типа Годунова // Математическое моделирование. 2015. Т. 27. № 10. С. 47-64.
14. Садин Д.В. TVD-схема для жестких задач волновой динамики гетерогенных сред негиперболического неконсервативного типа // ЖВМ и МФ. 2016. Т. 56. № 12. С. 2098-2109. doi: 10.7868/S0044466916120152
15. Садин Д.В. Схемы с настраиваемыми диссипативными свойствами для численного моделирования течений газа и газовзвесей // Математическое моделирование. 2017. Т. 29. № 12. С. 89-104.
16. Wong M.L., Lele S.K. High-order localized dissipation weighted compact nonlinear scheme for shock- and interface-capturing in compressible flows // Journal of Computational Physics. 2017. V. 339. N 15. P. 179-209. doi: 10.1016/j.jcp.2017.03.008
17. Liska R., Wendroff B. Comparison of several difference schemes on 1D and 2D test problems for the Euler equations // SIAM Journal on Scientific Computing. V. 25. N 3. P. 995-1017. doi: 10.1137/S1064827502402120
18. Gottlieb S., Shu C.-W. Total variation diminishing Runge-Kutta schemes // Mathematics of Computation. 1998. V. 67. N 221. P. 73-85. doi: 10.1090/S0025-5718-98-00913-2
19. Jiang G.-S., Tadmor E. Nonoscillatory central schemes for multidimensional hyperbolic conservation laws // SIAM Journal on Scientific Computing. 1998. V. 19. N 6. P. 1892-1917. doi: 10.1137/S106482759631041X
20. Colella P., Woodward P. The piecewise parabolic method (PPM) for gas-dynamical simulations // Journal of Computational Physics. 1984. V. 54. N 1. P. 174-201. doi: 10.1016/0021-9991(84)90143-8
21. Woodward P., Colella P. The numerical simulation of two-dimensional fluid flow with strong shocks // Journal of Computational Physics. 1984. V. 54. N 1. P. 115-173. doi: 10.1016/0021-9991(84)90142-6
Авторы
Садин Дмитрий Викторович - доктор технических наук, профессор, профессор, Военно-космическая академия имени А.Ф. Можайского, Санкт-Петербург, 197198, Российская Федерация, Scopus ID: 6602924618, ORCID ID: 0000-0001-53354847, [email protected]
Mathematics and Mathematical Physics, 2010, vol. 50, no. 5, pp. 848-861. doi: 10.1134/S096554251005009X
8. Shen Y.-Q., Zha G.-C. Generalized finite compact difference scheme for shock/complex flow field interaction. Journal of Computational Physics, 2011, vol. 230, no. 12, pp. 44194436. doi: 10.1016/j.jcp.2011.01.039
9. Mikhailovskaya M.N., Rogov B.V. Monotone compact running schemes for systems of hyperbolic equations. Computational Mathematics and Mathematical Physics, 2012, vol. 52, no. 4, pp. 578-600. doi: 10.1134/S0965542512040124
10. Christensen R.B. Godunov Methods on a Staggered Mesh -An Improved Artificial Viscosity. Technical Report UCRL-JC-105269, 1990, 11 p.
11. Shankar S.K., Kawai S., Lele S. Numerical simulation of multicomponent shock accelerated flows and mixing using localized artificial diffusivity method. Proc. 48h AIAA Aerospace Sciences Meeting Including the New Horizons Forum and Aerospace Exposition. Orlando, USA, 2010, art. 2010-0352.
12. Kurganov A., Liu Y. New adaptive artificial viscosity method for hyperbolic systems of conservation laws. Journal of Computational Physic, 2012, vol. 231, no. 24, pp. 8114-8132. doi: 10.1016/j.jcp.2012.07.040
13. Tagirova I. Yu., Rodionov A. V. Application of the artificial viscosity for suppressing the carbuncle phenomenon in Godunov-type schemes. Mathematical Models and Computer Simulations, 2016, vol. 8, no. 3, pp. 249-262. doi: 10.1134/S2070048216030091
14. Sadin D.V. TVD scheme for stiff problems of wave dynamics of heterogeneous media of nonhyperbolic nonconservative type. Computational Mathematics and Mathematical Physics, 2016, vol. 56, no. 12, pp. 2068-2078. doi: 10.1134/S0965542516120137
15. Sadin D.V. Schemes with customizable dissipative properties as applied to gas-suspensions flow simulation. Mathematical Models and Computer Simulations, 2017.
16. Wong M.L., Lele S.K. High-order localized dissipation weighted compact nonlinear scheme for shock- and interface-capturing in compressible flows. Journal of Computational Physics, 2017, vol. 339, no. 15, pp. 179-209. doi: 10.1016/j.jcp.2017.03.008.
17. Liska R., Wendroff B. Comparison of several difference schemes on 1D and 2D test problems for the Euler equations. SIAM Journal on Scientific Computing, vol. 25, no. 3, pp. 995-1017. doi: 10.1137/S1064827502402120
18. Gottlieb S., Shu C.-W. Total variation diminishing Runge-Kutta schemes. Mathematics of Computation, 1998, vol. 67, no. 221, pp. 73-85. doi: 10.1090/S0025-5718-98-00913-2
19. Jiang G.-S., Tadmor E. Nonoscillatory central schemes for multidimensional hyperbolic conservation laws. SIAM Journal on Scientific Computing, 1998, vol. 19, no. 6, pp. 1892-1917. doi: 10.1137/S106482759631041X
20. Colella P., Woodward P. The piecewise parabolic method (PPM) for gas-dynamical simulations. Journal of Computational Physics, 1984, vol. 54, no. 1, pp. 174-201. doi: 10.1016/0021-9991(84)90143-8
21. Woodward P., Colella P. The numerical simulation of two-dimensional fluid flow with strong shocks. Journal of Computational Physics, 1984, vol. 54, no. 1, pp. 115-173. doi: 10.1016/0021-9991(84)90142-6
Authors
Dmitry V. Sadin - D.Sc., Full Professor, Mozhaisky Military Space
Academy, Saint Petersburg, 197198, Russian Federation, Scopus
ID: 6602924618, ORCID ID: 0000-0001-5335-4847,