CC BY
224
МОДЕЛЬ ДВИЖЕНИЯ СУДНА В ГОРИЗОНТАЛЬНОЙ ПЛОСКОСТИ
Дерябин В.В., к.т.н., доцент, Арктический морской институт им. В.И. Воронина - филиал Государственного университета морского и речного флота им. адм. С.О. Макарова в г. Архангельске, gmavitder@mail.ru, +79215506140
В статье предлагается модель движения судна в горизонтальной плоскости, которая учитывает влияние внешних факторов: ветра и волнения. Режим движения задаётся числом оборотов гребного винта и значением курса, для удержания на котором используется ПИД-регулятор. Интегрированием системы дифференциальных уравнений получаются координаты судна и его курс в заданный момент времени. На основе натурных наблюдений проводится проверка работоспособности модели, результаты которой позволяют говорить о возможности её использования при моделировании навигационных процессов и проектировании алгоритмов управления судном.
Ключевые слова: система координат, модель движение судна, вектор движения, горизонтальная плоскость.
Deryabin V., Ph.D., assistant professor, Arctic Marine Institute under the name of V Voronin - State University of Maritime and Inland
The paper proposes a model of the vessel in a horizontal plane , which takes into account the influence of external factors : wind and waves . Mode of movement given by the number of revolutions of the propeller and the value of the course, to retain that uses a specific regulator. Integration of the system of differential equations obtained coordinates of the vessel and its rate at a given time . Based on field observations conducted performance testing of the model, the results of which suggest the possibility of its use in modeling processes and navigation design algorithms to control the vessel.
Keywords: coordinate system model propulsion motion vector horizontal plane.
Определим сначала системы координат, в которых будет рассматриваться движение судна.
Неподвижная система координат Oxyz имеет своим началом точку, лежащую на уровне невозмущённой поверхности моря. Ось Ox совпадает с направлением меридиана на Север, а ось Oy - с направлением параллели на Восток. Ось Oz направлена вертикально вниз.
Подвижная система координат Gxj yj Zj связана с центром тяжести судна G . Ось Gxj лежит в диаметральной плоскости и параллельна основной плоскости. Направлена она в сторону носа. Ось Gyj ориентирована в сторону правого борта. Ось Gzj направлена вниз перпендикулярно основной плоскости.
Система координат теоретического чертежа O2 x2y2z2 имеет началом точку пересечения линии киля и плоскости мидель-шпангоута. Ось O2 x2 направлена по линии киля в сторону носа, O2 y2 лежит в основной плоскости и направлена в сторону правого борта. Ось
O2 z 2 направлена вертикально вверх.
При решении задач управления можно ограничиться моделью, имеющей три степени свободы. При этом получаются координаты центра тяжести судна в горизонтальной плоскости xg , yG и курс судна K . Вертикальная координата центра тяжести предполагается неизменной во времени, а углы крена и дифферента равными нулю. Связь между координатами произвольной точки в неподвижной и подвижной системах устанавливается следующим образом:
кц, к22, к66 - коэффициенты присоединённых масс, 12 - момент инерции судна относительно вертикальной оси, Т - упор винта; индексы у сил и моментов относятся: V - к силам и моментам вязкостного трения, р - к поперечной силе и моменту гребного винта, р - к
MODEL SHIP TRAFFIC ABOVE THE HORIZONTAL PLANE
Shipping n. a. Admiral Makarov, Arkhangelsk branch
x = xG + Xj cos K - yj sin K,
У = yG + xj sin K + yj cos K.
Связь координат в неподвижной системе с координатами в системе O 2 X2 y 2 z 2 записывается в виде:
(1)
x1 =-xg + x2> Уі =-yg + y 2 ’ z1 = zg - z 2-
(2)
где {xg > Уg > 2g ) - координаты центра тяжести в системе координат, в которой задаётся теоретический чертёж. Уравнения движения судна в подвижной системе координат имеют форму [5]:
т ■(1 + кп) ■ -V-1 = т ■(1 + к22) ■ Vyi ю + Fvxi + Т + FRA + FAA + FWA, dt
-Vyi
влиянию руля, A,W - характеризуют воздействие ветра и волнения.
Воздействие винта характеризуется упором, который рассчитывается следующим образом [3]:
т0 = Кт pD 2(п 2 D 2 + Va2)(1 -1) (4)
где Kт - коэффициент упора винта (в системе параметров, предложенных В.М. Лавреньтевым), Р - плотность воды, D - диаметр винта, п - число оборотов винта, Va - скорость осевого натекания воды на винт, ^ - коэффициент продольного засасывания. Индекс «0»
обозначает упор при осевом натекании. Коэффициент упора гребного винта Kт зависит от величины поступи JL = Va /д^2+~П202 ,
которая, в свою очередь, связана с относительной поступью / соотношением: /^ = / / дД + / 2 . Как правило, результаты модельных испытаний гребных винтов позволяют получить зависимости вида: кт = кт (/) - кривые действия. Здесь символом кт обозначен «обыч-
к +
гребного винта серии В4-55 с единичным шаговым отношением, построенная по данным [4], изображена на рис.1.
ный» коэффициент упора гребного винта, который связан с коэффициентом Кт следующим образом: Кт 1 + /2 . Кривая действия
Рис.1. Кривая действия гребного винта В4-55 с шаговым отношением 1
Воспользовавшись приведёнными выше соотношениями, можно перейти от кривой действия к универсальной диаграмме гребного винта (рис.2).
Рис.2. Универсальная диаграмма для гребного винта В4-55 с шаговым отношением 1 Для удобства вычислений на компьютере кривая KT = KT (Jl ) была аппроксимирована полиномом четвёртой степени:
KT = 2.4453(JL)4 + (-3.5788)(Jl)3 + 1.1793(Jl)2 + (-0.4991) JL + 0.4263 .
Скорость осевого натекания воды на винт может быть определена по следующей формуле:
К =)/ Voh +(Voyi - Ik ffl)2 cos(x8,) • (1 -у 0 (5)
где X - коэффициент скоса потока, в к - местный угол дрейфа, - коэффициент попутного потока, - расстояние от центра
тяжести до места установки винта. Вычисленная по формуле (5) скорость может принимать и отрицательные значения, при этом J^ также оказывается отрицательной. Несмотря на то, что кривая рис.2 построена только для неотрицательных значений поступи, мы будем использовать ту же аппроксимацию и для отрицательных значений в качестве первого приближения.
Коэффициент скоса потока определяется следующим образом [3]:
х =
Хо + (1 -Хо)
'I вк I 45
вк
,1'Эе Рк <Рк п
1 Л'эе Рк >Рк тах.
(6)
где коэффициент Хо можно принять равным 0.8 для судов с U-образными обводами, а max = 45° Коэффициент попутного потока определяется выражением [3]:
/ Рк max - 1 Рк ^
Рк max
о
,1'Эе Рк <Рк max
,1'Эе Рк >Ркmax.
(7)
где вкmax = 45° , а Щ - коэффициент попутного потока при нулевом угле дрейфа. Последний определяется по формуле [1]:
ф = 0.165 • ---Аф , где Сь - коэффициент общей полноты, V - объёмное водоизмещение судна, а поправка
, Рг = К_
Аф = 0.1 (Бг — 0.2) учитывается только при числах Фруда ’ больших 0.2.
Коэффициент продольного засасывания рассчитывается по рекомендациям [1]:
I = 0.6(1 + 0.67ф)ф
Местный угол дрейфа рассчитывается так:
вк = -asin
Voyi - 1k®
2
ox1
(8)
Зная осевую скорость Va , рассчитав поступь / , можно определить по универсальной диаграмме коэффициент упора Кт и рассчитать упор при осевом натекании по формуле (4).
В случае косого натекания изменяется коэффициент упора винта и появляется поперечная сила ^Ру1 и её момент М р21 [3]. Результаты испытаний гребных винтов серии В4-55 при круговом натекании набегающего потока [3] позволили составить диаграммы зависимости коэффициента увеличения результирующей силы и угла между её вектором и осью винта 0 от направления косого потока Хвк .
Величины Цр , 0 находятся по графикам в зависимости от направления скоса потока и коэффициента нагрузки гребного винта при осевом натекании От = 270 / pVa2Д . Здесь Д обозначает площадь диска винта. Будем предполагать, что данный коэффициент может изменяться в пределах 1 < От ^ 12 . Значения величин в промежуточных точках находятся линейной интерполяцией.
Находим величину результирующей силыК : К = и её проекции на оси подвижной системы координат - упор и поперечную
силу, а также момент вокруг вертикальной оси:
T = R COS0,
Fpyi = -R sin0,
MPz1 =~FPy1lK.
(9)
Настоящая модель движителя позволяет получить упор серийного гребного винта В4-55 в условиях произвольного угла натекания, боковую силу винта и её момент, появляющиеся в условиях косого натекания. Модель работает для режимов работы винта, которые соответствуют переднему ходу ( п > 0 ).
Силы и момент, действующие на судно со стороны руля, вычисляются по методу разделения нагрузок, предложенному в [3]:
FRx1 =-R cos(a -S) + Y sin(a -S) - R cos(a -S) + Y sin(a -8),
FRy1 =-Y cos(a -S) - R sin(a -S) - Y cos(a -S) - R sin(a -S),
MRz1 = - FRy1lR .
где р - сила лобового сопротивления, у - подъёмная сила, а - угол атаки руля, § - угол перекладки руля, ¡р - расстояние от центра
тяжести до руля. Символ ‘ относится к части руля, попадающей в струю винта, а ‘’ - находящейся вне влияния струи. Подъёмная сила и сила лобового сопротивления для каждой части руля рассчитываются следующим образом:
X = С XR (a ) AR ; X = C XR (a ) AR ,
Y = C YR (a ) AR ’ Y = C YR )
.„PVa
2
-Ar
(11)
где Cxr - коэффициент лобового сопротивления, Cyr - коэффициент подъёмной силы, Va - скорость натекания, Ar - площадь руля.
Коэффициенты лобового сопротивления и подъёмной силы зависят от относительного удлинения руля Я = Hr / bR , где Hr - высота руля, bR - средняя хорда (длина) руля. Коэффициент лобового сопротивления определяется соотношением [3]:
Cxr = Cxr0 + KD sin2 a + 2sin3 a, (12)
где Cxro = 0.01, коэффициент Kd зависит от относительного удлинения руля и выбирается по специальному графику.
с = 2пЯ a
Величину коэффициента подъёмной силы при докритических углах атаки можно определить по формуле Cyr = Т Ta .
2 + Я
Значение критического угла атаки оценивается по формуле a¿ ~ 33° /-Д . в случаях, когда угол атаки руля превышает значение ak , используются коэффициенты, характеризующие линейное изменение подъёмной силы при закритических углах атаки: Эти коэффициенты будем обозначать символами Cyr , Cyr . Они могут быть определены следующим способом:
Cyr = Cyr (ak) exp[-k *(a-ak)]
(13)
где к принимает значение некоторое отрицательное значение.
Несмотря на то, что струя винта в швартовном режиме «увеличивает» относительное эффективное удлинение руля, будем использовать эту характеристику, исходя только из геометрии руля, так как в режимах, характеризующихся малым коэффициентом нагрузки винта, подобного увеличения подъёмной силы не наблюдается.
Скорость натекания для части руля, попадающей в струю от винта, определяется соотношением:
(14)
где в - местный угол дрейфа в месте установки руля, а величины К0 и С^ определяются следующим образом [5]:
K о = 0.5
1 +
х ¡кр /0.5D , где ¡кр - расстояние от плоскости диска винта до оси баллера руля. Указанный коэффици-
ент вводят лишь в том случае, если х < 1.5 , в противном случае принимают К0 = 1. Величина СтУ = 0.00255т /(DVa )2 . Если
осевая скорость натекания принимает малые значения, например, 0.001 м/с или меньше, то СтУ рассчитывается для этого значения. Угол атаки для части руля, находящейся в струе винта, можно вычислить так:
a = S + asin
VoXl + (1 - lR®f fsin(-Хвr )
Va
(15)
Выражение для скорости натекания воды на часть руля, которая находится вне струи винта, имеет вид:
К =^0x1 + (Тау1 - ¡Р®ї ) <°2(ХРг )(1 -УУ + (^1 + (Тау1 - ¡Я®)2 ) )іп2(ХРг )
Угол атаки для части руля, находящейся вне винтовой струи, определяется следующим образом:
(16)
a = S + asin
VoXl +(oy1 - lRmJ fsin(-Xеr )
Va
(17)
2
2
2
2
В выражении (16) величина попутного потока для руляфгвычисляется по формулам (7), в которых величину ф необходимо понимать как коэффициент попутного потока руля, вычисляемый для судов с и-образными кормовыми шпангоутами по формуле [2]:
ф = 0.35 + 2(Сь — 0.8).
Предложенная модель рулевого комплекса применима при движении судна передним ходом, а также для швартовного режима переднего хода, то есть при неотрицательных значениях продольной составляющей относительной скорости. Форма руля в плане предполагается близкой к прямоугольной, влияние рудерпоста не учитывается.
Силы сопротивления неинерционной природы, действующие на корпус судна, определяются следующим образом [5]:
Fxi = Cxp (p/2)Vo2 ALa, ip
PALg l
Fyi = Cie (p/2)Vo2Alg ,
Mvz1
2
CMPV0 - CM0L 1 ® 1 ®-CMa~ (V0 + L® )sin
(18)
Формулы, определяющие коэффициенты Cxe, Cíe, Cmp , имеют вид:
CXp = -0.075sin{n - arcsin(CX0 /0.075)] • [1 - (в/tyx)]}, CYp = 0.5Cie sin 2в cos в + c2 sin2 в + c3 sin4 2Д, CMa = m1 sin 2в + m2 sin в + m3 sin3 2в + m4 sin4 2в,
(19)
где Сх0 - коэффициент сопротивления воды движению судна при нулевом угле дрейфа, фх - угол, определяемый по графику [5],
Сур , ^2» С3, т-1, Ш-2, Ш3, Ш4 - коэффициенты, определяемые по номограммам или рассчитываемые по формулам справочника.
Коэффициент сопротивления Сх0 = XV0 /(Р / 2)VoXlАьо , где XV0 - сопротивление трения при нулевом угле дрейфа, которое может быть определено как сопротивление эквивалентной пластины:
XV0 = Ср0 * (Р /2^0x1$, где £ - площадь смоченной поверхности судна, Ср0 - коэффициент трения эквивалентной пластины,
определяемый по формуле: Ср0 = 0.455/(1§Ке)2'58 . Число Рейнольдса Ке = V ^1 ' Ь¡V , где V - коэффициент кинематической вязкости воды.
Приведённая площадь погруженной части ДП судна вычисляется по формуле:
= LтО , где т'1 - осадка на миделе, О - приведённый коэффициент полноты погруженной части ДП судна, определяемый по материалам [5].
Величина См0 = 0.059^2 и представляет собой коэффициент момента при вращении на месте.
Коэффициент демпфирующего момента определяется следующим образом:
СМш = СМш0 + а1 I вш в I +а2 {1 — сов[(2я — 4 I в |) ео8 в + 0.11 вш2в I]} (20)
где СмШ0 - производная коэффициента демпфирующего момента при нулевом угле дрейфа, вычисляемая по формуле:
СМш0 = (0.739 + 8.7тг /¿)(1.611о2 — 2.873о + 1.33).
Коэффициенты а1, а2 определяются выражениями:
üi = 0.09 — Cm^0 — 0.0033
a2 = 0.008 — + 0.9
2 B
T-
-2— 0.05 L
L - 7 B
\ / \
- 20
2
-- 0.003
+ 0.4(ст- 0.9) + 0.05(в; - 0.9),
+ 0.45(ст- 0.955),
где в г - коэффициент полноты мидель-шпангоута.
Величина О носит название обобщённой кривизны траектории и определяется таким образом:
О = ш / + ш 2 , где ш = Ь / Р , и Р - радиус кривизны траектории центра тяжести судна. О принимает значения от -1 до 1.
Случаю О = 0 соответствует прямолинейное движение судна, а О = ±1 - вращение на месте.
Аэродинамические силы и моменты определяются следующими соотношениями [5]:
FAx1 = CAX (Pa / 2) VR AVB ,
FAy1 = CAí (pa / 2)VRAVL , MAz1 = CAM (pa / 2)VRAVLL
L
Табл.1. Параметры модели для т/х «Инженер Плавинский» с осадкой по летнюю грузовую марку
Гидродинамический коэффициент 0.12 с помощью номограммы [5]
Диаметр винта, м о 3.8 информация о винте
Момент инерции относительно вертикальной оси, КГ • М • М 1.5450 • Ю10 по формуле [3]
Коэффициент лобового сопротивления руля (часть, попадающая в струю винта) Кт 0.75 с помощью графика [3]
Коэффициент лобового сопротивления руля (часть, не попадающая в струю винта) К 02 0.65 с помощью графика [3]
Длина по ватерлинию, м 1 147 общие данные о судне
2 Площадь смоченной поверхности, м 5 3628 по данным теоретического чертежа
Площадь проекции надводной части корпуса на плоскость мидель-шпангоута, м2 А'в 340 по чертежу общего вида
Площадь проекции надводной части корпуса на диаметральную плоскость, м2 Ат 1050 по чертежу общего вида
Осадка на миделе, м Т т 7 общие данные о судне
Отстояние центра парусности от плоскости мидель-шпангоута, м ъ -40 по чертежу общего вида
Коэффициент полноты мидель-шпангоута Д., 0.9 при помощи теоретического чертежа
Средняя хорда руля, м 3 информация о рулевом устройстве
Г идродинамический коэффициент с2 0.6 при помощи номограммы [5]
Г идродинамический коэффициент с3 0.08 при помощи номограммы [5]
Коэффициент скоса потока при нулевом угле дрейфа Яо 0.8 по материалам [5]
Высота части руля, попадающей в поток винта Ьт 3.8 информация о рулевом устройстве
Высота части руля, расположенной вне струи от винта 2.2 информация о рулевом устройстве
Коэффициент присоединённой массы | 0.032 по материалам [4]
Коэффициент присоединённой массы к22 0.949 по материалам [4]
Коэффициент присоединённого момента инерции ^66 0.824 по материалам [4]
Расстояние от винта до плоскости мидель-шпангоута, м ^ К 75.5 по чертежу общего вида
Расстояние от оси баллера до плоскости мидель-шпангоута, м 1ц 76 по чертежу общего вида
Массовое водоизмещение, м т 14300000 исполнительный грузовой план
Г идродинамический коэффициент щ 0.045 при помощи номограммы [5]
Г идродинамический коэффициент пи 0.0086 при помощи формулы [5]
Г идродинамический коэффициент тъ 0.020 при помощи номограммы [5]
Г идродинамический коэффициент Ш4 -0.0025 при помощи номограммы [5]
Вспомогательный угол, рад Чх 1.6930 при помощи графика [5]
Приведённый коэффициент полноты погруженной части ДП с 0.962 по формуле [5]
Координаты ЦТ судна в системе координат 02х2У2г2 (0,0,4) исполнительный грузовой план
где ра - плотность воздуха, - скорость относительного ветра, Лу% - площадь проекции надводной части корпуса на плоскость
мидель-шпангоута, Лу^ - площадь проекции надводной части корпуса на диаметральную плоскость. Коэффициенты в выражении (21) определяются следующим образом:
Табл.2. Параметры регулятора и частота вращения гребного винта в экспериментах
Эксперимент №1 Эксперимент №2 Эксперимент №3
Обороты винта, об./мин. 130 120 130
Коэффициент к\ 1 1 1
Коэффициент кг) , с 180 100 100
Коэффициент к3 4 4 4
Интервал интегрирования, с 300 300 300
■1400 -1200 -1000 -800 -600 -400 -200 0 200
У,т
Рис.3 Траектории движения судна в первом эксперименте
4
■12000 -10000 -8000 -6000 -4000 -2000 0 2000
У,т
Рис.4 Траектории движения судна во втором эксперименте
Слх —-(0.03 + 0.08ео8 ак), СЛ¥ — -1.058т ак,
Слм — С
ЛУ
0.25 + — -^ - х Ь 2п 8
(22)
где СХК - курсовой угол относительного ветра, — - расстояние от центра боковой парусности до плоскости мидель-шпангоута.
Силы и момент, действующие на судно со стороны морского волнения, вычисляются исходя из линейной теории морских волн [6] в предположении, что корпус судна не искажает окружающую волновую поверхность. Потенциал и давление в точках корпуса определяются в таком случае следующими соотношениями:
Рис.5 Траектории движения судна в третьем эксперименте
1 —2п
h т ~az 2п , 2п .
у = -- — ge А cos(— х----------------1),
2 2п Ат
h -A~z . 2п , 2п (23)
p = pg— • e А sin(— х------1),
где h - высота волны, Т - истинный период, Я - длина волны. Распространение волны происходит вдоль оси O'x' . Значения сил Fwxl , FWy1 и момента MWzi находятся путём численного интегрирования давлений по смоченной поверхности корпуса.
Данная модель была реализована с параметрами для т/х «Инженер Плавинский», совершающего плавание с осадкой по летнюю грузовую ватерлинию. Числовые значения данных параметров приведены в таблице 1.
Для удержания судна на заданном курсе используется ПИД-регулятор. Режим работы рулевого устройства имеет ограничения: наибольшее значение угла перекладки руля составляет 35°, максимальное значение скорости перекладки руля - 5°/с.
Проверка адекватности имитационной модели выполнялась на основе данных, полученных в результате проведения трёх натурных наблюдений. На четырёхчасовом промежутке времени плавания фиксировались значения параметров внешних факторов (ветер и волнение), курса судна по гирокомпасу, продольной скорости по относительному лагу и числа оборотов винта. Также были известны «точные» координаты судна, полученные при помощи спутниковой навигационной системы GPS NAVSTAR. Критерием расхождения истинной и теоретической траекторий служит наибольшее значение модуля невязки в течение времени плавания.
Угол перекладки руля не фиксировался в ходе экспериментов, поэтому при проверке адекватности модели использовался ПИД-регу-лятор, стабилизирующий судно на «среднем» курсе, значения которого известны. Параметры регулятора и число оборотов гребного винта приведены в табл.2.
Ввиду ненадёжности имеющихся данных о течении, параметры последнего определялись путём сравнения векторов абсолютной и относительной скорости судна, получаемых на основе опытных данных.
В результате были получены следующие значения максимума модуля невязки на четырёхчасовом отрезке времени: 0.9, 1.0 и 1.7 мили. Теоретические и действительные траектории изображены на рисунках 3,4,5. На них красным цветом изображаются истинные координаты, а синим - полученные с использованием модели движения.
Предложенная модель движения судна является моделью, имеющей три степени свободы, то есть позволяет получить горизонтальные координаты судна и его курс в условиях воздействия внешних факторов. Она учитывает влияние ветра, волнения, движительно-рулевого комплекса. Зависимости, используемые при построении данной модели, получены на основе натурных испытаний моделей судов и гребных винтов в опытовых бассейнах. Подобная модель может быть использована для изучения процессов управления движением судна в условиях влияния внешних факторов.
Литература:
1. Басин А.М, Анфимов В.Н. Гидродинамика судна. - Л.: Речной транспорт, 1961. - 684 с.
2. Васильев А.В. Управляемость судов. Учебное пособие. - Л.: Судостроение, 1989. - 328 с.
3. Гофман А.Д. Движительно-рулевой комплекс и маневрирование судна. Справочник. - Л.: Судостроение, 1988. - 360 с.
4. Справочник по теории корабля: В трех томах. Том1. Гидромеханика. Сопротивление движению судов. Судовые движители/Под ред. Я. И. Войткунского. - Л.: Судостроение, 1985. - 768с.
5. Справочник по теории корабля: В трех томах. Том3. Управляемость водоизмещающих судов. Гидродинамика судов с динамическими принципами поддержания/Под ред. Я. И. Войткунского. - Л.: Судостроение, 1985. - 544с.
6. Чижиумов С.Д. Основы динамики судов на волнении. Учебное пособие. - Комсомольск на Амуре.: ГОУВПО «КнАГТУ», 2010. -110 с.
