CC BY
437
Сер. 10. 2010. Вып. 4
ВЕСТНИК САНКТ-ПЕТЕРБУРГСКОГО УНИВЕРСИТЕТА
УДК 629.12.035
Д. В. Никущенко, Е. Н. Надымов, Р. А. Шушков
РАСЧЕТ ГИДРОДИНАМИЧЕСКИХ ХАРАКТЕРИСТИК ПОДВОДНЫХ АППАРАТОВ С ВЫСТУПАЮЩИМИ ЧАСТЯМИ, РУЛЯМИ И СТАБИЛИЗАТОРАМИ
При создании подводных аппаратов (ПА) одной из важных задач является оценка их мореходных качеств, прежде всего ходкости и управляемости [1—3], точность которой тесно связана с качеством прогнозирования гидродинамических сил и моментов, действующих на ПА при движении его по криволинейной траектории. В настоящее время основным способом определения остается испытание масштабных моделей. Однако при этом возникает ряд сложностей как технического, так и фундаментального характера, в частности трудоемкости изготовления моделей и проведения измерений, учета масштабного эффекта и т. п. В связи с этим все большую популярность приобретает использование численных методов. На ранних стадиях проектирования компьютерное моделирование позволяет проработать целый ряд конструктивных решений, например выполнить оценку эффективности органов управления и маневренных качеств. Таким образом, задача разработки и совершенствования надежных численных методов прогнозирования гидродинамических характеристик (ГДХ) является актуальной проблемой.
Классификация сил и моментов, действующих на подводный объект (ПО) при его движении, показана на рис. 1. Все силы и моменты можно разделить на следующие группы: инерционные, гидродинамические и реактивные, развиваемые движительно-рулевым комплексом (ДРК) ПО. Гидродинамические силы, действующие на корпус ПО, в свою очередь, делятся на силы и моменты при поступательном движении, а также силы и моменты при вращении. Первые из них называются позиционными силами, вторые - вращательными [2, 3].
При определении составляющих гидродинамических сил и моментов часто используют метод разделения ПО на голый корпус и выступающие части (ВЧ): ограждение
Никущенко Дмитрий Владимирович — кандидат технических наук, доцент кафедры теории корабля факультета кораблестроения и океанотехники Санкт-Петербургского государственного морского технического университета. Количество опубликованных работ: 46. Научные направления: теория турбулентности, динамика судов. E-mail: ndmitry@list.ru.
Надымов Евгений Николаевич — студент кафедры вычислительных методов механики деформируемого тела факультета прикладной математики-процессов управления Санкт-Петербургского государственного университета. Научный руководитель: доктор физико-математических наук,
проф. В. А. Павловский. Научные направления: гидродинамика, теория турбулентности. E-mail: johnnypmpu@gmail.com.
Шушков Роман Александрович — студент кафедры вычислительных методов механики деформируемого тела факультета прикладной математики-процессов управления Санкт-Петербургского государственного университета. Научный руководитель: доктор физико-математических наук,
проф. В. А. Павловский. Научные направления: гидродинамика, теория турбулентности. E-mail: hunterfordreams@inbox.ru.
© Д. В. Никущенко, Е. Н. Надымов, Р. А. Шушков, 2010
Рис. 1. Классификация сил и моментов
рубки, носовые и рубочные горизонтальные рули, кормовое оперение и пр. [4]. При этом влияние оперения на обтекание корпуса учитывается косвенно, что не позволяет производить точный расчет ГДХ системы в целом. Кроме того, ГДХ как корпуса, так и ВЧ оказываются под сильным влиянием движителя, в связи с этим учет поля течения в области его расположения также является одной из важных проблем.
Методы определения позиционных и вращательных составляющих ГДХ ПО можно разделить на теоретические (приближенные), экспериментальные и численные (рис. 2). Основными экспериментальными способами являются испытания модели в опытовых бассейнах, аэродинамических трубах (для позиционных составляющих), на ротативных установках или планарных механизмах (для вращательных составляющих).
Рис. 2. Способы определения гидродинамических характеристик ПО
В настоящее время в связи с развитием и ростом производительности вычислительных машин численные методы получают все большее распространение. Они становятся реальной альтернативой экспериментальным методам и имеют перед ними ряд преимуществ, например:
• численный эксперимент позволяет получать более полную и широкую информацию о течении;
• стоимость численных экспериментов по сравнению с экспериментальными методами непрерывно понижается;
• численные методы позволяют моделировать условия течения, не воспроизводимые или трудно воспроизводимые при экспериментальных испытаниях на моделях;
• благодаря численным методам время проектирования существенно уменьшается, численное исследование можно произвести очень быстро, т. е. за короткие сроки просчитать сотни вариантов и выбрать оптимальную конструкцию.
Численные методы базируются на моделях как вязкой, так и невязкой жидкости. Использование первой приводит к необходимости моделирования турбулентных течений. Следует отметить, что до сих пор не разработана универсальная модель турбулентности, применимая для любых видов течений.
Существует три основных подхода к моделированию турбулентных течений: прямое численное моделирование (Direct Numerical Simulation, DNS), метод крупных вихрей (Large Eddy Simulation, LES), решение осредненных по времени уравнений Навье-Стокса (Reynolds-Averaged Navier-Stockes’ equations, RANS) [5, 6]. Однако их применение требует сравнительно больших вычислительных ресурсов и времени подготовки исходных данных.
При использовании подходов на основе модели невязкой жидкости применяют два основных метода: прямое решение уравнений Эйлера и методы гидродинамических особенностей (МГО). При использовании последних воздействие поверхности исследуемого объекта на поток заменяется воздействием непрерывного распределения источников/стоков, дипольных или вихревых слоев. При этом задача сводится к решению уравнения Лапласа в заданной области с граничными условиями Неймана (задано либо граничными условиями Дирихле (задано у>).
За прошедшее время методы теоретического определения сил и моментов, действующих на корпус ПО и его элементы, не получили такой степени развития, чтобы широко войти в практику инженерных расчетов при проектировании ПО. Основная причина этого - то, что теоретические методы расчета ГДХ корпуса, кормового оперения и других ВЧ являются полуэмпирическими и не учитывают ряд важнейших физических факторов [7]. Разработка надежного, качественного, эффективного и физически адекватного метода в современных условиях возможна только на базе применения компьютерных технологий [8]. Эти технологии позволяют разрабатывать расчетные методы определения ГДХ морских подвижных объектов в двух направлениях - инженерном и исследовательском. Так как экспериментальный и отчасти исследовательские методы (основанные на модели вязкой жидкости) требуют трудоемкой подготовки и, в конечном счете, большого времени исполнения, разумно для анализа различных вариантов использовать упрощенные подходы, к которым относится и метод дискретных вихрей [9].
При решении поставленной задачи методом дискретных вихрей жидкость считается идеальной и невесомой, а движение ее вне вихревых поверхностей - потенциальным: v = УФ, где Ф(г, t) - потенциал вызванных скоростей в данной точке пространства.
Поверхность каждого моделируемого объекта заменяется непрерывным распределением вихревых особенностей (присоединенных вихрей), тогда в поток с нее должны сходить свободные вихри, обусловленные изменением интенсивности присоединенных (т. е. местной скорости потока). В соответствии с уравнением неразрывности
V - V = У-УФ = 0, вне вихревых поверхностей (т. е. поверхности тела Б и вихревой пелены Бр (см. рис. 1)) потенциал скоростей Ф является гармонической функцией (т. е. удовлетворяет уравнению Лапласа).
Таким образом, задача об обтекании тела потоком невязкой жидкости сводится к отысканию потенциала вызванных скоростей Ф(г, Ь), удовлетворяющего уравнению Лапласа
ДФ = 0,
при выполнении следующих граничных условий:
1) на поверхности объекта Б - условие непротекания:
дФ (г)
= (Ь) ■ п(Ь),
Б
дп(Ь)
где ух(г) - скорость движения точек поверхности Б (или скорость набегающего потока при обращенном движении);
2) на бесконечности - условие убывания возмущений: |УФ (гм, г) | ^ 0, при гм ^ го;
3) на пелене свободных вихрей Бр - условия совместности течения:
тШт = д-Ш'
4) на линии схода вихревой пелены - условие Чаплыгина-Жуковского о конечности скоростей (т. е. градиент Ф ограничен в окрестностях точек Мо € Б в любой момент времени г).
В результате выполнения условий непротекания и убывания возмущений задача сводится к решению интегрально-дифференциального уравнения
7-7Г— X / 7Г~ (~~= _у°° ' пмА*)> (!)
4п дпм0 ^ ] дпм \гмм0 )
г=1 Бъ
где N - общее количество вихревых поверхностей; - интенсивности вихревых слоев; М0 - точка, принадлежащая вихревой поверхности; М - точка вне вихревых поверхностей.
Детальное рассмотрение вопросов существования и единственности решения приведенного уравнения содержится в монографии И. К. Лифанова [10].
Окончательно, задача сводится к решению уравнения (1) относительно неизвестных интенсивностей присоединенных вихревых особенностей ^г, по которым возможно определять интенсивности свободных вихрей, если геометрия вихревых пелен известна. Отсутствие нагрузки на пелене свободных вихрей (краевое условие 3) приводит к следующему условию: точки вихревой пелены должны перемещаться по направлению местной скорости. То есть для точки (^0), сошедшей с линии Ьр в момент го, в момент г должно выполняться равенство
—Мр^_ ’ °') = V (мр (г, г0), г) + \р (мр (г, г0), г) + v00 (мр (г, г0), г),
где V - вызванная скорость в точке Мр от системы присоединенных вихрей; \р - вызванная скорость в точке Мр от системы свободных вихрей на временном слое г.
Поставленная задача может решаться в связанной системе координат для случая как обращенного движения (т. е. когда неподвижное тело обтекается потоком жидкости), так и абсолютного движения (т. е. когда тело движется с поступательной скоростью V0(£) и угловой ш(Ь) в изначально неподвижной жидкости).
Для решения уравнения (1) непрерывное распределение вихрей заменяется дискретным, для чего поверхности крыльев системы продольными и поперечными линиями разбиваются на элементарные площадки - панели. Воздействие каждой панели на поток заменяется воздействием вихревой рамки постоянной интенсивности. Тогда уравнение (1) заменяется системой алгебраических уравнений относительно интенсивностей вихревых рамок:
N т и
Ги аЧ2 = -4п,шп1 — ЕЕЕ д6рг4-8+1> 1 = 1, 9 =1, 2>;Н> (2)
3=1 р=1 1=1 в = 1
здесь тпз = vnj /vто - безразмерная нормальная составляющая скорости невозмущенного потока в точке Мз; Г* = Г+/(^ТО|В) - безразмерная циркуляция г-й вихревой рамки; ХN - регуляризирующая переменная; I - номер рамки с нулевой циркуляцией; 3^ = Гт—1 — Г- - циркуляция рамки, сходящей с линии Ьр (равная разности циркуляций прилегающих рамок); а.з - координатная функция, определяемая взаимным расположением точек М^ и М.. на поверхности тела Б:
Е(гй + 1 — Гк) ' (гк — Г])
( г/ \ / м-
к=1 (гк - ) (гк+1 - Гк) - [(гк+1 - Гк){гк - Гз )]
х [ (г&+1 -Гк)- (гк -гд-) (гй+1 - Гк) ■ (гк -Гд-)
I |гй + 1-Г^'| \*к-Гз\
Вихревую рамку с нулевой циркуляцией необходимо вводить в силу того, что в противном случае система линейных алгебраических уравнений получается вырожденной. В итоге система (2) является переопределенной, поэтому в работе [11] было предложено вводить регуляризирующую переменную ЛN.
После вычисления значений циркуляций вихревых рамок можно рассчитать скорости потока во всей области течения. Затем, используя интеграл Коши-Лагранжа или теорему Н. Е. Жуковского «в малом» [12], определить гидродинамические силы во всех элементарных площадках поверхностей тел системы, после чего найти силы и моменты, действующие со стороны жидкости как на всю систему целиком, так и на каждое тело в составе системы.
На основе изложенного алгоритма разработан программный комплекс, предназначенный для расчета гидродинамических характеристик морских ПО вместе с системой крыльев (к числу последних относятся органы управления, стабилизации и ограждения выдвижных устройств).
Расчет крыльевых систем имеет ряд характерных особенностей [13], таких как малое удлинение, влияющее на взаимодействие между потоками, сходящими с торцевых кромок, и большая относительная толщина, что приводит к изменению эпюры распределения давлений по поверхности крыла [14]. Численная схема, использованная для создания данного комплекса, позволяет учитывать указанные особенности обтекания ПО.
Расчет характеристик объекта в целом выполняется при помощи двух подпрограмм: WingSim и Sub_3D. Первая, представленная в виде ехе-файла и использующая
графический интерфейс (систему окон), определяет геометрию крыльевой системы, расположение крыльев относительно корпуса. Во второй, представленной в виде Фортран-проекта, рассчитываются силы и моменты, действующие на изолированный корпус. Общую структуру расчетного комплекса иллюстрирует рис. 3.
Обработка Расчет крыльев
WingSim
данных
Обработка Расчет корпуса
Sub_3D
данных
.0 Настройка
5
to CD О СО .0 расчета > „ Получение SubObject
§ ГДХ ПА
d
Рис. 3. Структура программы SubObject
Алгоритм расчета объекта, состоящего из корпуса и системы крыльев, представляет собой итерационный процесс, где в качестве начального приближения берутся ГДД крыльев и корпуса по отдельности без взаимного влияния. Далее последовательно на каждой итерации вычисляются ГДД (силы и моменты) на корпусе во внешнем поле скоростей, вызванных наличием крыльев, затем ГДД крыльев с учетом внешнего поля скоростей от корпуса.
Итерационный процесс продолжается до тех пор, пока не будет достигнута некоторая, заданная наперед точность вычислений, или процесс не будет остановлен самим пользователем.
Таким образом, процедуру расчета задачи обтекания систем объектов, состоящих из корпуса и крыльев, можно представить в виде следующих шагов:
• создание расчетной геометрии корпуса и крыльев: их взаимное расположение, параметры разбиения объектов, задание системы координат;
• задание параметров расчета: количества итераций, угла атаки, характерных размеров, способов вычисления сил и моментов;
• последовательный итерационный расчет ГДД объектов системы;
• визуализация полученных результатов.
По изложенному алгоритму был создан программный расчетный комплекс SubObject, который позволяет выполнять моделирование как по отдельности систем крыльев и корпуса ПА, так и с учетом взаимного влияния. Он основан на концепции объектно-ориентированного программирования и реализован на языке C++ в среде Borland C++ Builder.
Разработанная математическая модель и предложенный программный комплекс были использованы для решения ряда практических задач. Ниже излагаются результаты исследования двух схематизированных моделей подводных объектов.
Процесс решения приводится на рис. 4. Сначала задаются геометрия корпуса, расположение крыльев, затем происходит выбор параметров системы (угол атаки, характерные размеры, положение центра масс). После этого выполняется расчет ГДД ПА, представленный на рис. 5.
6S
Рис. 4. Задание геометрии ПА
Коэффициенты гидродинамических сил и моментов определяются следующим образом:
М,
Сі =
иру2)У 3
Нру2)у'
і = 1,2, 3,
где ¥і - продольная сила, нормальная сила и поперечная сила, действующие на модель; Мі - поперечный, поворотный и продольный моменты, действующие на модель; ^(ру2) - скоростной напор потока; V - полное водоизмещение ПА в масштабе модели.
Первая модель представляла собой корпус жесткого дирижабля «Акрон» (рис. 6, а), являющийся телом вращения с относительным удлинением 5,9. Геометрические характеристики корпуса: объемное водоизмещение V = 7808, 47 м3, длина Ь = 80 м, ширина В = 13, 5 м, высота Н = 13, 5 м; геометрические характеристики крыльев: размах крыла I = 4, 3 м, хорда крыла Ь = 15, 3 м. Моделировалось обтекание модели напором воздуха со скоростью 50 м/с, число Рейнольдса Яв = 0, 98 • 107. Определение гидродинамических характеристик объекта производилось в диапазоне углов атаки а Є [0..14]0.
На рис. 7 и 8 приведены результаты расчета изолированного корпуса в сравнении с данными экспериментального исследования в аэродинамической трубе ЦАГИ [3], а также корпуса с оперением. Как видно из рис. 7, погрешность произведенных вычислений составляет не более 15%. На рис. 8 показано влияние оперения на изменение гидродинамических характеристик дирижабля. Можно видеть, что в результате нормальная сила увеличивается на 25%, а продольный момент уменьшается на 9%.
Рис. 5. Вычисление ГДХ ПА в программе SubObject
Второй моделью послужил схематизированный корпус ПА вместе с системой горизонтальных и вертикальных стабилизаторов, а также ограждением рубки (см. рис. 6, б). Геометрические характеристики корпуса: объемное водоизмещение V = 0,1106 м3, длина Ь = 2, 5 м, ширина В = 0, 282 м, высота Н = 0, 282 м; вертикальные стабилизаторы: размах крыла I = 0,151 м, хорда крыла Ь = 0,112 м; горизонтальные стабилизаторы: размах крыла I = 0, 099 м, хорда крыла Ь = 0,112 м; ограждение рубки: высота Н = 0,101, длина I = 0, 326 м. Моделировалось исследование в аэродинамической трубе со скоростью 40 м/с, число Рейнольдса Яв = 0, 98 • 107. Расчет производился в диапазоне углов атаки а Е [0.. 10]0. Силы и моменты определялись в связанной системе координат, начало которой находилось в центре тяжести модели.
На рис. 9 а, б показано изменение интегральных характеристик (голого корпуса, корпуса со стабилизаторами и ограждением рубки) в зависимости от угла атаки, а также влияние системы крыльев на изменение сил и моментов. Из него видно, что нормальная сила голого корпуса на 50% меньше, чем корпуса вместе с хвостовым оперением и ограждением рубки. Продольный момент по стабилизируется при добавлении крыльев, его график практически параллелен оси абсцисс для всего диапазона углов атаки.
Рис. 6. Корпус жесткого дирижабля «Акрон» вместе с хвостовым оперением (а), схематизированный корпус модели ПО вместе с хвостовым оперением и ограждением рубки (б)
Рис. 7. Коэффициент нормальной силы (а) и продольного момента (б) дирижабля «Акрон» - голый корпус 1 - эксперимент; 2 - расчет; 3 и 4 - йЯТ и «стандартная» к — є модели.
Рис. 8. Коэффициент нормальной силы (а) и продольного момента (б) дирижабля «Акрон» - корпус + крылья 1 - расчет; 2 - эксперимент.
В целом можно сказать, что разработанный расчетный комплекс представляет собой эффективное средство моделирования обтекания систем разнородных элементов, в частности оперенного корпуса ПО. Данный алгоритм дает возможность производить расчет при различных комбинациях корпуса и ВЧ, позволяя оперативно оценивать мореходные качества ПО в зависимости от геометрии корпуса и расположения ВЧ. При этом определение ГДХ возможно как для каждого объекта по отдельности, так и для всей системы в целом.
Су а Ст б
Рис. 9. Коэффициент нормальной силы (а) и продольного момента (б) модели подводного аппарата 1 и 3 - голый корпус (эксперимент и расчет); 2 и 4 - корпус + крылья + ограждение рубки
(эксперимент и расчет).
Литература
1. Белоцерковский С. М., Ништ М. И. Отрывное и безотрывное обтекание тонких крыльев идеальной жидкостью. М.: Наука, 1978. 352 с.
2. Гурьев Ю. В., Никущенко Д. В. Опыт и перспективы использования компьютерных технологий для решения задач управляемости подводных объектов // Материалы науч.-техн. конференции «XII Макеевские чтения». СПб., 2005. С. 21-26.
3. Федяевский К. К. Материалы к аэродинамическому расчету воздушных кораблей. Ч. II // Труды Центр. аэрогидродин. ин-та. 1933. Т. 178. 72 с.
4. Никущенко Д. В., Рогожина Е. А. Влияние рулей на гидродинамические характеристики оперения подводных аппаратов // XLII Крыловские чтения: тез. докл. СПб., 2006. С. 67.
5. Ferzieger J. L., Peric M. Computational methods for fluid dynamics. Berlin: Springer-Verlag, 2002. 423 p.
6. Mathieu J., Scott J. An introduction to turbulent flow. Cambridge: Cambridge Univ. Press, 2000. 374 p.
7. Sato T., Izumi K., Miyata H. Numerical Simulation of Maneuvering Motion // Twenty-Second Symposium on Naval Hydrodynamics. Fukuoka, Japan, 2000. P. 724-737.
8. Sung C., Jiang M., Rhee B. at al. Validation of the flow around a turning submarine // Twenty-Forth Symposium on naval Hydrodynamics. Fukuoka, Japan, 2002. P. 669-681.
9. Sarpkaya T. Computational methods with vortices. The 1988 Freeman Scholar Lecture // J. of Fluids Engineering. 1989. Vol. 111. P. 5-52.
10. Лифанов И. К. Метод сингулярных интегральных уравнений и численный эксперимент в математической физике, аэродинамике, теории упругости и дифракции волн. М.: Янус, 1995. 519 с.
11. Лифанов И. К., Михайлов А. А. К расчету безотрывного и отрывного обтекания тел // Труды ВВИА. 1986. № 1313. С. 137-145.
12. Жуковский Н. Е. О присоединенных вихрях // Собр. соч. Т. IV: Аэродинамика. М.; Л.: Гос. изд-во техн.-теор. лит-ры, 1949. C. 69-91.
13. Белоцерковский С. М., Гиневский А. С. Моделирование турбулентных струй и следов на основе метода дискретных вихрей. М.: Физматлит, 1995. 368 с.
14. Никущенко Д. В., Калинин О. С., Кра>сиков В. И., Рогожина Е. А. Определение гидродинамических характеристик изолированного оперения подводного аппарата на основе метода дискретных вихрей // 6-я Междунар. конференция по морским интеллектуальным технологиям «Моринтех-2005». СПб., 2005. С. 111-113.
Статья рекомендована к печати проф. Ю. М. Далем. Статья принята к печати 10 июня 2010 г.
