CC BY
12
УДК 532.5.01; 532.5.032; 533.6.011
Г. Б. Сизых
Московский физико-технический институт (государственный университет)
О верификации численных расчетов вихревых течений методом проверки сохранения циркуляции
В общем пространственном случае для течения вязкого газа или жидкости получены удобные для верификации численных расчетов выражения скорости движения контура, по которому сохраняется циркуляция. Показано, что необходимо рассматривать эволюцию как минимум двух по-разному ориентированных контуров.
Ключевые слова: верификация, интегральный инвариант, циркуляция по контуру.
G. B. Sizykh
Moscow Institute of Physics and Technology (State University)
Verification of numerical calculation of eddy currents by checking the velocity circulation
We investigate viscous gas and fluid flows in space (three-dimensional) nonstationary cases. We consider the motion of a contour with constant circulation of the fluid rate. We find a new expression for the velocity of this contour. This expression is much simpler than the known one previously and can be used to check numerical calculations. It is shown that it is necessary to study the evolution of at least two differently oriented contours.
Key words: verification, integral invariant, circulation along the contour.
1. Введение
Обычно численный метод проверяется на известных точных решениях и путем сравнения с численными решениями других надежных, хорошо зарекомендовавших себя методов. Но желательна верификация и каждого конкретного расчета, для которого не известно ни точное решение, ни результат расчета надежного метода. Это можно сделать путем проверки законов сохранения, выполняющихся для выбранной математической модели исследуемого процесса. При этом необходимо использовать новые законы сохранения, не заложенные в основу численных схем.
В соответствии с теоремой Томсона [1], если массовые силы допускают потенциал, а идеальная жидкость баротропна, то циркуляция скорости по любому замкнутому «жидкому» контуру (который в любой момент времени состоит из одних и тех же частиц жидкости) остается постоянной. В случае вязкой жидкости указанная теорема неприменима, но можно сформулировать другое содержательное предложение. В частности, можно говорить о воображаемой среде, «частицы» которой движутся с некоторой скоростью и и составляют контуры с сохраняющейся во времени циркуляцией скорости жидкости V. Для плоских и незакрученных осесимметричных течений вязкой несжимаемой жидкости такие скорости найдены в [2, 3]. В общем пространственном случае течения вязкого газа или жидкости способ вычисления скорости И предложен в [4]. Но этот способ не удобен для использования, поскольку он не является локальным, и поэтому для расчета поля скорости и требуется интегрирование вдоль вихревых линий. В работе [5] впервые получено локальное выражение для скорости и и предложен способ вычисления, не требующий знания гидродинамических функций на конечном удалении от рассматриваемой точки (формула
(2) ниже). В [5] было предложено использовать это выражение для проверки расчетов, полученных конечно-разностными методами, путем проверки сохранения циркуляции по движущимся со скоростью U контурам. Однако попытки практически реализовать такой подход выявили две проблемы. Первая - действительно ли необходимо вычислять все компоненты у всех векторов, входящих в громоздкое выражение, предложенное в [5]? Может быть, в результате преобразований можно «избавиться» от необходимости вычислять некоторые компоненты этих векторов? Вторая - ошибки в каких именно параметрах можно обнаружить, отслеживая движение того или иного контура?
В данной работе удалось получить ответы на эти вопросы на основе новой формулы для скорости U. Основная идея работы состоит в использовании следующего свойства движущегося гладкого контура, т.е. контура, состоящего из воображаемых частиц, каждая из которых движется с некоторой скоростью. Если к скорости каждой частицы контура в точке расположения этой частицы прибавить произвольную скорость, направленную по касательной к контуру, то картина движения контура как геометрического места точек не изменится. Это свойство не было использовано в [5]. Поэтому полученные в данной работе формулы не являются модификацией формул работы [5] и представляют собой отдельный содержательный результат.
2. Скорость переноса простого вихревого контура
Сформулируем основной результат работы [5]. Он относится к стационарным и нестационарным течениям сжимаемых и несжимаемых, вязких (различных реологий) и невязких жидкостей. Здесь и далее, следуя [1], под жидкостью будем понимать как жидкость, так и газ, имея в виду, что жидкость может быть сжимаемой. Уравнение импульсов, описывающее течение, может быть представлено в виде
д V
— + (V ■ V)V = F, (1)
где V - скорость жидкости, F - равнодействующая всех сил, отнесенная к плотности. В общем случае здесь учитываются силы гидродинамического давления, вязкости и массовые силы.
В завихренном О = rot V = 0 течении гладкую замкнутую кривую без самопересечений Г назовем простым вихревым контуром (ПВК), если существует двухмерная односвязная область (поверхность) от, обладающая следующими свойствами:
а) от лежит в трехмерной области вихревого О = 0 течения, в которой компоненты скорости трижды непрерывно дифференцируемы;
б) вихревые линии пересекают от под острым углом к нормали;
в) Г С стг.
Будем обозначать e - единичный касательный вектор контура (вектор e задан в точках контура).
В работе [5] показано, что для произвольного ПВК и любого скалярного поля Л вектор-функция
U = V + + (О F)(° e) [О х e] + ЛО (2)
Q Q2 (Ъ2 - (О ■ e)2J
будет скоростью движения точек этого контура с сохранением циркуляции скорости жидкости V по нему.
3. Вариант формулы для скорости переноса простого вихревого контура
Во всех точках простого вихревого контура векторы е и О неколлинеарны. Поэтому любой вектор можно представить линейной комбинацией векторов е, О и [О х е]. В част-
ности,
Е = ае + 7 О + @ [О х е] , (3)
где коэффициенты а, 7 и @ в общем случае различны в разных точках контура. Умножая скалярно это векторное равенство на е и О, получим два скалярных уравнения, из которых выражаются коэффициенты а и 7:
а = (Е - е) - (О - е)2(Е О), (4)
О2 - (О • е)2
^ =(Е • °) - (О • е) (Е • е) (5)
7 = О2 - (о . е)2 • (5)
Запишем выражение (2) для скорости и с использованием формул (3) - (5). После несложных преобразований имеем
и = V + , (Е е) , [О х е] + *(О - е) П~* (О - О) е + А О.
/ - 2^ 1 J О2
(о2 - (О • е)2)
Поскольку функция Л может быть выбрана произвольно и поскольку скорость переноса контура определяется с точностью до векторного поля, коллинеарного касательному вектору е, каждая из скоростей
(Е е)
Ио = V + т-(-)-^ [О х е] (6)
(О2 - (О • е)2)
(Е е)
И1 = V - (V -е) е + т-^-^ [О х е] (7)
(О2 - (О • е)2)
будет скоростью движения точек контура Г с сохранением циркуляции скорости жидкости V по нему.
Обе формулы (6) и (7) менее громоздки по сравнению с формулой (2). Однако для проверки численных расчетов более полезной представляется формула (7).
Во-первых, из (7) следует необходимость «отслеживания» как минимум двух по-разному ориентированных контуров. Действительно, если обозначить через Ее касательную к контуру составляющую вектора Е, а через О^ и Vx - нормальные к вектору е (в точках контура) составляющие векторов О и V соответственно, то формулу (7) можно представить в виде
И1 = V± + х е] •
Таким образом, проверка закона сохранения циркуляции скорости V по выбранному контуру позволяет «отфильтровать» лишь ошибки вычисления Vx, Ее и О^. Для проверки других компонент V, Е и О требуется рассматривать движение другого контура, имеющего иную ориентацию в пространстве.
Во-вторых, формула (7) по сравнению с (6) не содержит составляющей, коллинеарной касательному вектору е. Это важно по следующей причине. В теории, наличие или отсутствие такой составляющей не меняет картину движения контура как геометрического места точек. Но в практических вычислениях составляющая, коллинеарная вектору е, приводит к дополнительным ошибкам, для уменьшения которых требуется уменьшать шаги пространственно-временной расчетной сетки. Поэтому формула (7) предпочтительнее, чем более простая на вид формула (6).
и
4. Скорость переноса простого вихревого контура по подвижной поверхности в пространстве
При проверке расчетов, полученных конечно-разностными схемами, естественно воспользоваться теми же сетками, которые используются в этих схемах. Как правило, трехмерную пространственную сетку можно представить набором двухмерных сеток, расположенных на, вообще говоря, подвижных поверхностях в пространстве. Поэтому было бы удобно при проверке расчета рассматривать движение контура, который все время остается на одной и той же поверхности. Однако скорость (7) не гарантирует, что контур, расположенный на одной из таких поверхностей, при движении со скоростью и останется лежать на этой поверхности. Это, в частности, следует из того, что вектор О^ не всегда нормален к поверхности, и поэтому вектор [О^ х е] может иметь ненулевую нормальную к плоскости составляющую. В этой связи представляется полезным получить формулу для скорости переноса контура, обеспечивающую расположение контура на поверхности во все время движения. Пусть п - единичный нормальный вектор к одной из таких поверхностей и. Рассмотрим контур, лежащий на неподвижной поверхности и. Если вновь воспользоваться тем, что скорость переноса контура определена с точностью до линейной комбинации векторов е и О, то после соответствующих преобразований формулы (7) получим
и2 = V. - £Ост + ^ [п х е] , (8)
1 'п 1гп
где векторы V0- и Оа - проекции векторов V и О на поверхность и; Уп = (V-п) и 1п = (О -п) - величины нормальных к поверхности составляющих векторов V и О;
= (Е-е) - величина проекции правой части уравнения движения (1) на касательный к контуру вектор е.
Векторы V0-, Ост и [п х е], входящие в (8), параллельны поверхности и. Поэтому, двигаясь со скоростью (8), контур будет оставаться на поверхности и.
Рассмотрим случай подвижной поверхности и. Движение поверхности можно описать движением частиц, составляющих эту поверхность с некоторой скоростью W( и). Вообще говоря, эта скорость различна для разных частиц поверхности. Положение поверхности в различные моменты времени определяется нормальной к поверхности и составляющей этой скорости Wn(и) = (W(и)•n)n = ^гап.
Скорость переноса контура определена с точностью до вектора, коллинеарного завихренности О. Поэтому после соответствующих преобразований формулы (8) получим
из = + Va + Уп)Ост + [п х е] , (9)
где обозначения имеют тот же смысл, что и для формулы (8). Двигаясь со скоростью Из, контур будет оставаться на подвижной поверхности и, если эта поверхность движется со скоростью ^гап.
5. Вязкая несжимаемая жидкость
Представляется естественным, что применение предложенных выше формул в первую очередь будет проводиться для наиболее простого случая вязких течений - для течений вязкой несжимаемой жидкости в потенциальном поле массовых сил. Представим полученные выше скорости переноса контура для такого течения. Динамическое уравнение системы уравнений Навье-Стокса запишем в форме [1, 6]:
д V (р V 2\
— + [О х V] = -г/ ^ О - V ( + П + — ) , (10)
где р - давление, р - плотность, П - потенциал массовых сил. Возможность исключения потенциальной составляющей ^^ ^ + П + ^г^ показана в работе [5]. Это значит, что в
формулу (2), а следовательно, и в формулу (9) в качестве F можно подставлять не всю правую часть (10), а только — v rot О. В результате приходим к следующему выводу. Пусть и = a(t) - подвижная поверхность с нормальным вектором n. Поверхность может быть не плоской, и в любой фиксированный момент времени нормаль n может быть различна в разных точках поверхности. Пусть поверхность а движется со скоростью ^ran так, что во всех точках этой поверхности во все рассматриваемые моменты времени нормальная составляющая завихренности Qra = (О -n) не обращается в ноль. Тогда, если в начальный момент времени контур с единичным касательным вектором e лежал на поверхности а, то при движении со скоростью
U3 = + V. + Ml-^ Ост — ,[n x e] (11)
он останется лежать на поверхности а, а циркуляция скорости жидкости V по этому контуру будет сохраняться. Смысл обозначений Va, Оа и Vn приведен после формулы (8). В случае неподвижной поверхности нужно в формуле (11) положить Wn = 0. Наконец, приведем выражение для скорости переноса простого вихревого контура, не привязанного ни к какой заданной подвижной или неподвижной поверхности. Такая скорость может быть вычислена по формулам (6), (7). Для вязкой несжимаемой жидкости, с учетом возможности исключения потенциальной составляющей, получаем, что каждая из скоростей
U = V — , . (e r0t0)2. [О x e]
m2 — (О -e)2J
и
Ui = V — (V-e) e — * . (e -r0t О)2. [О x e]
m2 — (О -e)2J
будет скоростью движения точек контура Г с сохранением циркуляции скорости жидкости V по нему.
6. Использование точных решений
Остановимся на одной проблеме применения полученных выше формул для проверки численных расчетов. Численные решения, как правило, получаются с некоторой ошибкой - с некоторым отклонением от точного решения. Поэтому точное сохранение численного значения циркуляции на численном решении невозможно. И для выбраковки неприемлемых решений необходимо решить вопрос о величине отклонения циркуляции от постоянного значения, превышение которой свидетельствует о «существенном» отличии между численным и точным решениями. Для этого представляется необходимым провести численные эксперименты для выяснения точности, с которой сохраняется циркуляция, когда численное решение имеет некоторую «приемлемую» точность. В рамках этого эксперимента необходимо иметь возможность вычислять отклонение численного решения от точного решения. Это возможно, если известно точное решение. Причем в рассматриваемом течении не должно быть точек с нулевой завихренностью. Последнее требование связано с тем, что полученные в данной работе формулы верны только для вихревых (О = 0) зон течения. В качестве таких точных решений уравнений Навье-Стокса можно использовать классические решения Куэтта [7], Пуазейля и Кармана [6]. А также решения, полученные в [8-10].
7. Заключение
Получены новые формулы для скорости движения точек контура с сохранением циркуляции скорости вязкого газа или жидкости V по нему. В качестве примера эти формулы
представлены для случая вязкой несжимаемой жидкости. Скорость переноса простого вихревого контура, вычисленная по формуле (7), по сравнению с ранее известными формулами [5], не имеет составляющей, коллинеарной касательному вектору е. Для вычисления по формуле (7) нет необходимости находить касательную к е составляющую завихренности и нормальную к е составляющую правой части уравнения движения (1). Эти обстоятельства и простота полученной формулы делают ее более удобной для верификации численных расчетов методом проверки сохранения циркуляции по контуру.
Показано, что при верификации необходимо рассматривать эволюцию как минимум двух по-разному ориентированных контуров.
Предложена формула (11), которая может быть использована для верификации полей течения, полученных конечно-разностными схемами, на тех же расчетных сетках.
Результаты работы предлагается использовать для проверки расчетов, полученных конечно-разностными схемами.
Литература
1. Кочин Н.Б., Кибель И.А., Розе Н.В. Теоретическая гидромеханика. Т. 1. М.: Физматгиз, 1963.
2. Голубкин В.Н., Сизых Г.Б. О некоторых общих свойствах плоскопараллельных течений вязкой жидкости // Изв. АН СССР. МЖГ. 1987. № 3. С. 176-178.
3. Брутян М.А., Голубкин В.Н., Крапивский П.Л. Об уравнении Бернулли для осесим-метричных течений вязкой жидкости // Ученые записки ЦАГИ. 1988. Т. 19, № 2. С. 98-100.
4. Марков В.В., Сизых Г.Б. Эволюция завихренности в жидкости и газе // Изв. РАН. МЖГ. 2015. № 2. С. 8-15.
5. Голубкин В.Н., Марков В.В., Сизых Г.Б. Интегральный инвариант уравнений движения вязкого газа // ПММ. 2015. Т. 79, вып. 6. С. 808-816.
6. Лойцянский Л. Г. Механика жидкости и газа. М.: Дрофа, 2003.
7. Couette M. Études sur le frottement des liquids // Ann. Chim. Phys. 1890. V. 21. P. 433510.
8. Аристов С.Н., Просвиряков Е.Ю. Неоднородные течения Куэтта // Нелинейная динамика. 2014. Т. 10, № 2. С. 177-182.
9. Аристов С.Н., Просвиряков Е.Ю. Волны Стокса в завихренной жидкости // Нелинейная динамика. 2014. Т. 10, № 3. С. 309-318.
10. Aristov S.N., Prosviryakov E.Yu. Unsteady Layered Vortical Fluid Flows // Fluid Dynamics. 2016. V. 51, N 2. P. 148-154.
References
1. Kochin N.E., Kibel' I.A., Roze N.V. Theoretical Hydromechanics. V. 1. New York: Interscience Publishers, 1964.
2. Golubkin V.N., Sizykh G.B. On some of the general properties of plane-parallel flows of a viscous fluid. Math. USSR Academy of Sciences. Fluid Dynamics. 1987. N 3. P. 176-178.
3. Brutyan M.A., Golubkin V.N., Krapivskiy P.L. About the Bernoulli's equation for axisymmetric flows of a viscous fluid. TsAGI Science Journal. 1988. V. 19, N 2. P. 98-100.
4. Markov V.V., Sizykh G.B. Évolution of vorticity in a liquid and gas. Bulletin of the Russian Academy of Sciences. Fluid Dynamics. 2015. N 2. P. 8-15.
5. Golubkin V.N., Markov V.V., Sizykh G.B. Integral invariant of the equations of motion of viscous gas. Journal of Applied Mathematics and Mechanics. 2015. V. 79, I. 6. P. 808-816.
6. Loytsyansky L.G. Mechanics of Liquids and Gases. Oxford: Pergamon Press, 1966.
7. Couette M. Studies on the friction of liquids. Annals of Chemistry and of Physics. 1890. V. 21. P. 433-510.
8. Aristov S.N., Prosviryakov E.Yu. Heterogeneous of the Couette flow. Nonlinear dynamics. 2014. V. 10, N 2. P. 177-182.
9. Aristov S.N., Prosviryakov E.Yu. Stokes waves in the swirling liquid. Nonlinear dynamics. 2014. V. 10, N 3. P. 309-318.
10. Aristov S.N., Prosviryakov E.Yu. Unsteady Layered Vortical Fluid Flows. Fluid Dynamics. 2016. V. 51, N 2. P. 148-154.
Поступила в редакцию 10.05.2016
