Ovil Aviation High TECHNOLOGIES
Vol. 19, No. 06, 2016
УДК 532.517.4, 621.4
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ПАРАМЕТРОВ ТЕЧЕНИЯ ОДИНОЧНОЙ ВЕТРОЭЛЕКТРИЧЕСКОЙ УСТАНОВКИ
СВ. СТРИЖАК1
1 Институт Системного Программирования Российской Академии Наук, г. Москва, Россия
Известно, что на территории РФ планируется строительство нескольких крупных ветропарков. Задачи, связанные с проектированием и с оценкой эффективности работы ветропарков, являются актуальными на сегодняшний день. Одно из возможных направлений в проектировании связано с математическим моделированием. Метод крупных вихрей (вихреразрешающее моделирование), разработанный в рамках направления вычислительной гидродинамики, позволяет в деталях воспроизводить нестационарную структуру течения и определить различные интегральные характеристики.
В данной статье проведен расчет работы одиночной ветроэлектрической установки с помощью метода крупных вихрей и метода плоских сечений вдоль лопасти турбины. Для постановки задачи рассматривалась расчетная область в форме параллелепипеда и использовалась адаптированная неструктурированная сетка. Математическая модель включала в себя основные уравнения неразрывности и количества движения для несжимаемой жидкости. Крупномасштабные вихревые структуры рассчитывались при помощи интегрирования фильтрованных уравнений. Расчет был проведен с использованием модели Смагоринского для определения значения турбулентной подсеточной вязкости. Геометрические параметры ветроэлектрической установки задавались исходя из открытых источников в интернете.
Все физические величины в расчетной области определялись в центре расчетной ячейки. Аппроксимация слагаемых в исходных уравнениях была выполнена со вторым порядком точности по времени и пространству. Уравнения для связи скорости и давления решались с помощью итерационного алгоритма PIMPLE.
Общее количество рассчитываемых физических величин на каждом временном шаге равнялось 18. В связи с этим требовались ресурсы вычислительного кластера.
В результате расчета течения в следе для трехлопастной турбины получены осредненные и мгновенные значения скорости, давления, подсеточной кинетической энергии и турбулентной вязкости, компоненты тензора подсеточных напряжений. Полученные результаты, качественно совпадающие с известными результатами экспериментов и численных расчетов, свидетельствуют о возможности адекватно рассчитать параметры течения для одиночной ветроэлектрической установки.
Ключевые слова: ветропарки, ветроэлектрические установки, решатель, метод крупных вихрей, модель Смагоринского, расчетная область, вихревой след, профиль сечения, угол атаки, профиль скорости.
ВВЕДЕНИЕ
В 2016 г. в средствах массовой информации появились сообщения о строительстве новых ветропарков на территории РФ [1]. Дочернее предприятие ОАО «Атомэнергомаш» ГК «Ро-сатом» АО «Ветрогенерирующая Компания» («ВетроОГК») в 2016 г. завершило первый этап работ по ветропарку «Шовгеновский» в Республике Адыгея. Были выполнены работы по вет-ромониторингу. К 2018 г. планируется построить 75 установок мощностью по 2 МВт каждая. Всего в Адыгее для реализации проектов по развитию ветроэнергетики подготовлены три площадки. Планируется, что ветровые электростанции суммарной мощностью 441 МВт будут построены в Шовгеновском, Гиагинском и Кошехабльском районах региона. На Кубани первые ветропарки планируют построить в Темрюкском районе и на Тамани. ЗАО «ВетроОГК» собирается возвести два объекта суммарной мощностью 364 млн кВтч в год. Станцию в районе поселка Бегового, расположенного на полуострове Тамань, могут ввести в эксплуатацию в 2018 г. Для объекта на Таманском полуострове найдена площадка для строительства.
В связи с развитием ветроэнергетики в РФ, проектированием новых ветрогенераторов, ветроэлектрических установок (ВЭУ) и их эксплуатацией в различных климатических условиях на обширной территории РФ (Северо-Запад, Север, Юг, Дальний Восток), строительством крупных ветропарков возникает потребность в решении различных научных прикладных задач. К актуальным задачам в ветроэнергетике можно отнести:
Vol. 19, No. 06, 2016
Ovil Aviation High TECHNOLOGIES
- определение оптимальной формы лопасти турбины с целью повышения эффективности ВЭУ и снижения уровня акустического шума при работе ВЭУ для определения безопасной зоны проживания людей;
- оценку влияния вибраций лопастей турбины и их разрушения в процессе эксплуатации ВЭУ;
- оценку влияния неровностей, шероховатости на поверхности лопастей турбины в процессе ледообразования на КПД лопастей;
- оценку образования льда различной формы, срыва льда с поверхности лопастей и определение безопасной зоны для нахождения людей;
- изучение влияния формы лопастей, расположения вихрегенераторов в пограничном слое с целью управления отрывом потока и повышения КПД ВЭУ;
- оценку эффективности работы ВЭУ при их групповом расположении с учетом влияния рельефа местности и топографии;
- учет влияния атмосферного пограничного слоя на работу ветропарка и оценку влияния вихревых следов от ВЭУ на локальный климат региона.
Известно, что ранее проводились исследования по изучению структуры течения за одиночной ВЭУ в гидроканале [2, 3] и в аэродинамической трубе [4]. Также выполнялись аналитические исследования [5] и численные расчеты для одиночной ВЭУ с использованием различных подходов на базе URANS [6] и LES [7-10].
МЕТОДЫ И МЕТОДОЛОГИЯ ИССЛЕДОВАНИЯ
В рамках web-лаборатории UniHUB ИСП РАН на вычислительном кластере установлен программный комплекс в составе открытой библиотеки OpenFOAM. Данный программный комплекс позволяет рассчитывать параметры для ВЭУ, оценивать мощность работы ветропарка при изменении параметров атмосферного пограничного слоя с учетом влияния стратификации среды, учета вращения Земли, учета реального рельефа местности (шероховатость поверхности) с использованием метода крупных вихрей (вихреразрешающее моделирование) и различных моделей для подсеточной турбулентной вязкости, моделей напряжения сдвига на поверхности. Программный комплекс включает в себя решатели и служебные утилиты:
- решатель для моделирования атмосферного пограничного слоя с учетом различных моделей турбулентности для подсеточной кинетической энергии в постановке метода крупных вихрей - ABLSolver;
- решатель для моделирования работы одиночной лопастной турбины - pisoturbineFoam;
- решатель для моделирования течения в ветропарке при произвольном количестве лопастных турбин - windPlantSolver;
- утилита для инициализации начальных данных для физических полей - setFieldsABL;
- набор специализированных граничных условий для задания параметров атмосферы.
Дополнительно существует возможность подключения отдельного модуля FAST в решатель windPlantSolverFAST для решения задач аэроупругости лопастей турбины c учетом задания геометрических параметров турбины.
ПОСТАНОВКА ЗАДАЧИ
Рассматривалась расчетная область в форме параллелепипеда с координатами в крайних точках: (-318; -318; -318) (2646; 318; 318). В точке с координатами (5,0; 0,0; -90,0) была расположена одиночная ВЭУ. Геометрические параметры ВЭУ задавались исходя из открытых источников в интернете. Для ВЭУ задавалось: количество лопастей равным 3, радиус лопасти 63,0 м, радиус основания цилиндрической втулки, к которой крепятся лопасти, равен 1,5 м, высота установки 87,6 м, частота вращения ротора турбины равна 9,15 оборотов в минуту.
Ovil Aviation High TECHNOLOGIES
Vol. 19, No. 06, 2016
Для построения базовой сетки использовались утилита blockMesh, утилита refineMesh использовалась для построения адаптированной неструктурированной сетки вокруг трехлопастной турбины. Общее количество расчетных ячеек составило 2 979 840. Как минимум 40 вычислительных ячеек приходилось на диаметр ротора для разрешения концевых и корневых вихрей.
На входе расчетной области задавалось значение скорости ветра Ио = (8,0; 0,0; 0,0) м/а На верхней и нижней границах для значения скорости задавалось условие проскальзывания потока, на боковых границах - условие цикличности, на выходной границе задавалось условие продолжения решения. Для давления на боковых границах задавалось условие цикличности, на остальных границах задавалось условие равенства нулю градиента. На рис. 1 показана расчетная область и адаптированная неструктурированная сетка.
а б
Рис. 1. а) расчетная область, б) расчетная сетка Fig. 1. a) numerical domain, b) numerical grid
Требовалось определить параметры потока в следе за трехлопастной турбиной, с целью оценки величины действующих сил, распределения физических параметров по длине оси симметрии расчетной области.
МАТЕМАТИЧЕСКАЯ МОДЕЛЬ
Математическая модель включала в себя основные уравнения неразрывности и количества движения для несжимаемой жидкости. Крупномасштабные вихревые структуры рассчитывались при помощи интегрирования фильтрованных уравнений [11, 12]. Для получения фильтрованных уравнений использовался коробочный фильтр. Мелкие вихри, размер которых не превышал шага расчетной сетки, моделировались с помощью модели Смагоринского для под-сеточной вязкости [12].
SU:
IT = 0, (1)
j
где Uj = Uj - Uj - вектор значения разрешимой скорости после процедуры фильтрации уравнений,
dUL + А(U.U.) = -„... оu -др.- AipDkfPb-,1 _/Ф'
-1 gi-(§} + f,. (2)
dt dxj j j j oX, dXj ~ "J' j
где в ук - альтернативный символ, О] - вектор вращения Земли, р - значение разрешимого давления после процедуры фильтрации, И? - часть подсеточного тензора напряжений. Слагаемое
Том 19, № 06, 2016_Научный Вестник МГТУ ГА
Vol. 19, No. 06, 2016 Ovil Aviation High TECHNOLOGIES
- Pb
для расчета параметра для плавучести содержит отношение плотностей —, g. - составляю-
Po
щая ускорения свободного падения, — градиент давления в горизонтальном направлении.
Учитывая характерные размеры лопастей турбины, при которых число Рейнольдса Re может достигать 107^108, бывает сложно разрешать с помощью метода крупных вихрей все особенности течения. Для этого потребуются слишком большие расчетные сетки. Известно, что подход на базе Actuator Line Model (метод плоских сечений) не требует слишком подробных сеток вокруг лопастей турбины. Данный подход позволяет воспроизводить различные типы вихрей: в следе, концевые, корневые и граничные вихри. В рамках подхода Actuator Line Model лопасть турбины разбивается (дискретизируется) на отдельные плоские сечения с заданной геометрией в форме известных профилей с заданной толщиной и круткой профиля (рис.2). Для каждого профиля рассчитываются значения подъемной силы и силы сопротивления. Сила, действующая на поток, равна по модулю силе действующей на лопасти турбины. Далее проводится процедура нормирования силы, и величина силы проецируется на поток путем добавления отдельного слагаемого в уравнения для количества движения. Результирующая сила F определялась по методике [6]
j actuator
fturbine (r ) = Lexp S П
ч2
r
-
(3)
Данные для аэродинамических коэффициентов Сх(а) и Су(а) при изменении угла атаки в диапазоне от -180° до 180° для тел, описывающих ВЭУ: «СуНпёег1», «СуНпёег2», аэродинамических профилей «Би40_Л17», «БШ5_Л17», «ОШ0_Л17», «БШ5_Л17», «БШ1_Л17», «КЛСЛ64_Л17» заданы в табличном виде с шагом в 5°. Значения для толщины профиля и его крутка также задавались.
Рис. 2. Фотография ВЭУ и схема для метода плоских сечений Fig. 2. Photo of wind turbine and scheme for Actuator Line Model
Количество сечений для лопасти турбины задавалось равным 64. Распределение положения сечений по длине лопасти выбиралось равномерным. Величина параметра epsilon (ширина Гауссовского распределения) в формуле (3) задавалось равным 10,0.
Все физические величины в расчетной области определялись в центре расчетной ячейки [13]. Аппроксимация слагаемых в исходных уравнениях была выполнена со вторым
Ovil Aviation High TECHNOLOGIES
Vol. 19, No. 06, 2016
порядком точности по времени и пространству. Уравнения для связи скорости, давления, решались с помощью итерационного алгоритма PIMPLE. Полученная система алгебраических уравнений решалась итерационным методом сопряженных градиентов с предобуславлива-телем для скорости, давления, тензора напряжений и параметров для модели подсеточной вязкости. Общее количество рассчитываемых физических величин (скалярные, векторные, тензорные) для выбранной модели для подсеточной турбулентной вязкости равно 18. В связи с этим требовались ресурсы вычислительного кластера. Величина напряжения сдвига на поверхности может определяться по модели Schumann. На поверхности компоненты тензора равны нулю, кроме значений Т13, Т23.
РЕЗУЛЬТАТЫ РАСЧЕТОВ
В ходе расчета были получены поля осредненных и мгновенных величин (скорость, давление, подсеточная вязкость, тензор напряжения и другие). Общее время счета задачи составило 1000 с, шаг по времени задавался dt = 0,033 с. На рис. 3 представлены результаты расчета значения модуля скорости U от времени.
| |
в г
Рис. 3. Значения модуля скорости в различные моменты времени: а - 150 с, б - 300 с, в - 500 с, г - 1000 с Fig. 3. U_0 magnitude value at different time: a) 150 s, b) 300 s, c) 500 s, d) 1000 s
Видно, что вихревой след за ВЭУ меняется во времени и распространяется до выходной границы расчетной области. Результаты расчетов показывают, что конечное время в 1000 с оказывается достаточным для первой оценки картины течения. Большое количество временных шагов требует больших вычислительных ресурсов.
На рис. 4 представлены зависимости распределения полного и осредненного давления и скорости от количества итераций шагов по времени. Шаг записи по времени dt = 0,033 с. Видно, что в месте расположения ВЭУ значения величины давления резко меняются.
На рис. 5 представлены зависимости кинетической энергии турбулентности и турбулентной вязкости от количества итераций по времени. Их распределение имеет близкий характер. На рис. 6 представлены зависимости для средних значений нелинейных величин для пульсации компонент скорости от количества итераций по времени. Данные зависимости характеризуют анизотропию течения. Наибольшие значения имеют величины с компонентами по осям OZ и OZ, что характерно для течения с вращением потока.
Vol. 19, No. 06, 2016
Ovil Aviation High TECHNOLOGIES
ю-,
-U_Magnltude
U Mean_Mag nltude pMean
500 1000 li>00 2000 2i)U0 SJ00
Рис. 4. Зависимость полного и осредненного давления и скорости от количества итераций по времени Fig. 4. Total and Mean values of pressure, velocity from number of iterations
asoo
0.ÛÛ5
и S00 1000 I SOU 2000 2SUU 3000 S5UU
Рис. 5. Зависимость кинетической энергии турбулентности и турбулентной вязкости от количества итераций по времени Fig. 5. Subgrid scale turbulent kinetic energy and viscosity from number of iterations
-o.i
SOU 1000 I £.00 2UUU 2Й00 3UUU S5UU
Рис. 6. Зависимость средних значений нелинейных величин пульсации компонент скорости
от количества итераций по времени Fig. 6. Subgrid scale mean values of nonlinear fluctuation components of velocity
from number of iterations
Научный Вестник МГТУ ГА_Том 19, № 06, 2016
Ovil Aviation High TECHNOLOGIES Vol. 19, No. 06, 2016
На рис. 7 представлены зависимости осредненных величин компонент тензора напряжения от количества итераций по времени. Распределение компонент по осям OХ и ОХ, ОY и OY имеет близкий характер и дает максимальный вклад в общее значение модуля тензора напряжения.
Рис. 7. Зависимость для осредненных величин компонент тензора напряжения от количества итераций по времени Fig. 7. Subgrid scale mean values of stress tensor components from number of iterations
ОБСУЖДЕНИЕ ПОЛУЧЕННЫХ РЕЗУЛЬТАТОВ И ЗАКЛЮЧЕНИЕ
Полученные результаты, качественно совпадающие с известными данными экспериментов и численных расчетов [2-10], свидетельствуют о возможности адекватно рассчитать параметры течения для одиночной ВЭУ. Расчеты проведены с использованием ресурсов вычислительного кластера, с использованием 12-96 ядер, в рамках web-лаборатории UniHUB. В дальнейшем планируется провести исследования динамики вихревых структур с учетом использования лагранжевой динамической модели Смагоринского [14], влияния изменений величины теплового потока у поверхности Земли.
СПИСОК ЛИТЕРАТУРЫ
1. Росатом до 2020 г. построит ветропарки мощностью 612 МВТ на Юге России: 22 июля 2016 г. [Электронный ресурс]. URL: http://kuban.rbc.ru/krasnodar/freenews/5792012f9a 79473230837559?from=newsfeed (дата обращения: 27.10.2016).
2. Naumov I.V., Rahmanov V.V., Okulov V.L., Velte K.M, Meyer K.E, Mikkelsen R.F. Diagnostic of flow behind the model of windmill turbine // Thermophysics and aeromechanics. 2012. Vol. 19. No. 3. Pp. 267-278.
3. Okulov V.L., Naumov I.V., Mikkelsen R.F., Kabardin I.K., Sorensen J.N. A regular Strouhal number for large-scale instability in the far wake of a rotor. J. Fluid Mech. 2014. Vol. 747. Pp. 369-380.
4. Krogstad P.A., Lund J.A. An experimental and numerical study of the performance of a model turbine. Wind Energ. 2012. 15. Pp. 443-457.
5. Окулов В.Л., Соренсен Ж.Н., ван Куик Г.А.М. Развитие теории оптимального ротора. / К 100-летию вихревой теории гребного винта профессора Н.Е. Жуковского. М. - Ижевск: НИЦ «Регулярная и хаотическая динамика», 2013. 120 с.
6. Sorensen J.N., Shen W.Z. Numerical Modeling of Wind Turbine Wakes, Journal of Fluids Engineering 124, 2002, pp. 393-399.
7. Churchfield M.J., Moriarty P.J., Vijayakumar G., Brasseur J.G. Wind Energy-Related Atmospheric Boundary Layer Large-Eddy Simulation Using OpenFOAM. Conference Paper NREL/CP-500-48905 August 2010. Pp. 1-26.
Vol. 19, No. 06, 2016
Ovil Aviation High TECHNOLOGIES
8. Calaf M., Meneveau C., Meyers J. Large eddy simulations of fully developed wind-turbine array boundary layers, Phys. Fluids 22. 015110 (2010).
9. Churchfield M.J., Lee S., Michalakes J., Moriarty P.J. A numerical study of the effects of atmospheric and wake turbulence on wind turbine dynamics. Journal of Turbulence. Vol. 13. No. 14. 2012. Pp. 1-32.
10. Munters M., Meneveau C., Meyers J. Turbulent Inflow Precursor Method with Time-Varying Direction for Large-Eddy Simulations and Applications to Wind Farms. Boundary-Layer Me-teorol. 2016. DOI 10.1007/s10546-016-0127-z
11. Sagaut P. Large eddy simulation for incompressible flows: an introduction. Springer. Berlin. 2002. 426 p.
12. Ferziger J.H., Peric M. Computational Methods for Fluid Dynamics. Springer-Verlag, Berlin et al.: Springer. 2002. 423 p.
13. Weller H.G., Tabor G., Jasak H., Fureby C. A tensorial approach to computational continuum mechanics using object oriented techniques. Computers in Physics. 1998. Vol. 12. No. 6. Pp. 620-631.
14. Meneveau C., Lund T.S., Cabot W.H. A Lagrangian dynamic subgrid-scale model of turbulence. J. Fluid. Mech 1996. 319. Pp. 353-385.
СВЕДЕНИЯ ОБ АВТОРАХ
Стрижак Сергей Владимирович, кандидат технических наук, инженер ИСП РАН, e-mail: strijhak@yandex.ru.
MATHEMATICAL MODELING OF FLOW PARAMETERS FOR SINGLE WIND TURBINE
Sergei V. Strijhak1
1The Institute for System Programming of the Russian Academy of Sciences, Moscow, Russia
ABSTRACT
It is known that on the territory of the Russian Federation the construction of several large wind farms is planned. The tasks connected with design and efficiency evaluation of wind farm work are in demand today. One of the possible directions in design is connected with mathematical modeling. The method of large eddy simulation developed within the direction of computational hydrodynamics allows to reproduce unsteady structure of the flow in details and to determine various integrated values. The calculation of work for single wind turbine installation by means of large eddy simulation and Actuator Line Method along the turbine blade is given in this work. For problem definition the numerical method in the form of a box was considered and the adapted unstructured grid was used.
The mathematical model included the main equations of continuity and momentum equations for incompressible fluid. The large-scale vortex structures were calculated by means of integration of the filtered equations. The calculation was carried out with Smagorinsky model for determination of subgrid scale turbulent viscosity. The geometrical parameters of wind turbine were set proceeding from open sources in the Internet.
All physical values were defined at center of computational cell. The approximation of items in equations was executed with the second order of accuracy for time and space. The equations for coupling velocity and pressure were solved by means of iterative algorithm PIMPLE. The total quantity of the calculated physical values on each time step was equal to 18. So, the resources of a high performance cluster were required.
As a result of flow calculation in wake for the three-bladed turbine average and instantaneous values of velocity, pressure, subgrid kinetic energy and turbulent viscosity, components of subgrid stress tensor were worked out. The received results matched the known results of experiments and numerical simulation, testify the opportunity to adequately calculate the flow parameters for a single wind turbine.
Key words: wind farms, wind turbine, solver, large eddy simulation, Smagorinsky model, numerical domain, vortex wake, airfoil profile, angle of attack, velocity profile.
Civil Aviation High TECHNOLOGIES
Vol. 19, No. 06, 2016
REFERENCES
1. Rosatom do 2020 goda postroit vetroparki moshnostiyu 612 Mvt yf yuge Rossii 22 iulya, 2016 [Rosatom will have built to 2020 windparks with capacity of 612 MW in the South of Russia: 22 July 2016]. Available at: http://kuban.rbc.ru/krasnodar/freenews/5792012f9a79473230837559? from=newsfeed. (Accessed: 27.10.2016). (in Russian)
2. Naumov I.V., Rahmanov V.V., Okulov V.L., Velte K.M, Meyer K.E, Mikkelsen R.F. Diagnostic of flow behind the model of windmill turbine. Thermophysics and aeromechanics. 2012, vol. 19, no. 3, pp. 267-278.
3. Okulov V.L., Naumov I.V., Mikkelsen R.F., Kabardin I.K., Sorensen J.N. A regular Strouhal number for large-scale instability in the far wake of a rotor. J. Fluid Mech. (2014). Vol. 747, pp. 369-380.
4. Krogstad P.A., Lund J.A. An experimental and numerical study of the performance of a model turbine. Wind Energ, 2012, 15, pp. 443-457.
5. Okulov V.L., Sarensen, J.N, van Kuik G.A.M. Razvitie teorii optimalnogo rotora [The development of the optimal rotor theory]. K 100-letiu vichrevoi teorii grebnogo vintaprofessora N.E. Zhu-kovskogo [On the 100th anniversary of the vortex theory of screw propeller professor N.E. Zhukovsky.] M. - Izhevsk: NITS "Regulyarnaya i chaoticheskya dinamika". 2013. 120 p. (in Russian)
6. Sorensen, J.N., Shen, W.Z. Numerical Modeling of Wind Turbine Wakes, Journal of Fluids Engineering 124, 2002, pp. 393-399.
7. Churchfield M.J., Moriarty P.J., Vijayakumar G., Brasseur J.G. Wind Energy-Related Atmospheric Boundary Layer Large-Eddy Simulation Using OpenFOAM. Conference Paper NREL/CP-500-48905 August 2010. pp.1 - 26.
8. Calaf M., Meneveau C., Meyers J. Large eddy simulations of fully developed wind-turbine array boundary layers, Phys. Fluids 22. 015110 (2010).
9. Churchfield M.J., Lee S., Michalakes J., Moriarty P.J. A numerical study of the effects of atmospheric and wake turbulence on wind turbine dynamics. Journal of Turbulence, vol. 13, no. 14. 2012, pp. 1-32.
10. Munters M., Meneveau C., Meyers J. Turbulent Inflow Precursor Method with Time-Varying Direction for Large-Eddy Simulations and Applications to Wind Farms. Boundary-Layer Me-teorol. 2016. DOI 10.1007/s10546-016-0127-z
11. Sagaut P. Large eddy simulation for incompressible flows: an introduction. Springer. Berlin. 2002. 426 p.
12. Ferziger J.H., Peric M. Computational Methods for Fluid Dynamics. Springer-Verlag, Berlin et al.: Springer. 2002. 423 p.
13. Weller H.G., Tabor G., Jasak H., Fureby C. A tensorial approach to computational continuum mechanics using object oriented techniques. Computers in Physics. 1998, vol. 12, no. 6, pp. 620-631.
14. Meneveau C., Lund T.S., Cabot W.H. A Lagrangian dynamic subgrid-scale model of turbulence. J. Fluid. Mech 1996. 319. pp. 353-385.
INFORMATION ABOUT THE AUTHOR
Strijhak Sergei Vladimirovich, PhD, Engineer of Institute for System Programming of the Russian Academy of Sciences, strijhak@yandex.ru.