Синтез оптимального управления воздушным судном на основе уравнений нелинейной динамики Текст научной статьи по специальности «Физика»

Ovil Aviation High Technologies

Vol. 20, No. 03, 2017

ИНФОРМАТИКА, ВЫЧИСЛИТЕЛЬНАЯ ТЕХНИКА И УПРАВЛЕНИЕ

УДК 629.735:681.5

СИНТЕЗ ОПТИМАЛЬНОГО УПРАВЛЕНИЯ ВОЗДУШНЫМ СУДНОМ НА ОСНОВЕ УРАВНЕНИЙ НЕЛИНЕЙНОЙ ДИНАМИКИ

В.Ф. ДИЛЬ1, В.Н. СИЗЫХ2

1 Московский государственный технический университет гражданской авиации, Иркутский филиал, г. Иркутск, Россия

2Иркутский государственный университет путей сообщения, г. Иркутск, Россия

В статье рассматривается методика синтеза нелинейных систем управления воздушным судном (ВС) на основе оптимизации пилотажного контура по методу обратных задач динамики. Для синтеза алгоритмов управления применяется нелинейная математическая модель движения ВС на пилотажном и траекторном уровнях. Определены этапы методики синтеза на пилотажном уровне, включающие: назначение эталонных движений ВС по трем степеням свободы, выбор структуры алгоритмов управления и их параметров, определение степени близости текущего и эталонного движений с помощью квадратичного функционала с последующей организацией движения к экстремуму-минимуму градиентным методом. Через оптимизируемые параметры пилотажного уровня решается прямая задача динамики для траекторного уровня управления пространственным движением ВС. Основой для выполнения расчета траекторных параметров ВС является нелинейная математическая модель траекторного движения, входом для которого служат выходные параметры пилотажного уровня. При этом выходные параметры траекторного уровня определяются путем численного интегрирования входных сигналов с учетом динамических обдувочных коэффициентов ВС. Разработана структурная схема организации контуров управления пространственным движением ВС. Выполнено исследование функционирования пилотажного контура путем численного моделирования с использованием программ MathCad и Паскаль. Параметры эталонной модели определялись путем имитационного моделирования на языке Паскаль по реакции нелинейной модели ВС на ступенчатые отклонения аэродинамических рулей. Показано, что задача управления пространственным движением оптимальна в смысле реализации управляющего воздействия. Кроме того, в сравнении с [9] можно утверждать, что из-за обратимости преобразований энергий вращательного и поступательного движений изменяется только содержание прямых и обратных задач динамики.

Ключевые слова: автоматическое управление, пространственное движение, уравнения нелинейной динамики.

ВВЕДЕНИЕ

Задача синтеза систем автоматического управления (САУ) воздушного судна как нелинейных объектов остается актуальной. Одним из методов, позволяющим выполнить синтез системы управления для нелинейных объектов, является метод, основанный на концепции обратных задач динамики (ОЗД), глубоко проработанный рядом видных ученых нашей страны [1, 2, 11].

Сущность обратных задач динамики сводится к нахождению управляющих сил и моментов, под действием которых объект будет двигаться по назначенной траектории [1]. Применение методов ОЗД позволяет осуществить синтез САУ как для линейных, так и для нелинейных объектов. При этом синтезированный закон управления является оптимальным в смысле реализации назначенной траектории движения объекта [2]. Ввиду того, что синтез законов управления САУ методом ОЗД опирается на продувочные характеристики модели ВС, САУ становится адаптивной к изменяемым свойствам ВС. Эти качества метода ОЗД позволяют утверждать о перспективности синтеза САУ подобным образом.

АЛГОРИТМЫ УПРАВЛЕНИЯ ВРАЩАТЕЛЬНЫМ И ТРАЕКТОРНЫМ ДВИЖЕНИЯМИ ВОЗДУШНОГО СУДНА

Рассмотрим двухуровневый подход к процессу управления [5, 6, 10]. На верхнем уровне осуществляется траекторное движение, на нижнем (пилотажном) - угловое движение ВС. Решение задачи синтеза пилотажного уровня управления ВС организуется методом обратных задач динамики, а траекторный контур непосредственно реализуется путем решения прямой задачи динамики.

Vol. 20, No. 03, 2017

Ovil Aviation High Technologies

Вращательное движение ЛА описывается в связанной системе координат (СК) динамическими уравнениями Эйлера [3, 4, 6]

~МХ

Мсв = My = I <Ьу +

Mz d)z

(h — 1у)шуш2 Ох - Iz)шхшг

-0у - 1х)шхшу.

(1)

где I = й1ад(1х1у12) - главные моменты инерции ЛА.

В общем случае моменты в связанной системе координат представляют собой векторную сумму аэродинамических моментов, диссипативных и гироскопических моментов тяги авиационных двигателей. Управляющими являются аэродинамические моменты.

Аэродинамические моменты зависят от аэродинамической компоновки ВС, режимов полета, положений органов управления. Они выражаются через параметры ВС и аэродинамические коэффициенты Шг (г = х, у, £) [4, 6]:

Мг

= qSL[mxmymz] ,

(2)

где Ь = й1ад(1 I Ьа) - диагональная матрица характеристических линейных размеров ВС (I - размах крыла, Ьа - длина среднеквадратической хорды), 5 - эффективная площадь крыла, ц = рУ2/2 - скоростной напор, V - воздушная скорость, т£ - безразмерные коэффициенты аэродинамических моментов, действующих по осямсвязанной СК. Оставив в выражениях коэффициентов наиболее существенные составляющие, получим следующие их выражения:

тх = mxß + т%а + тхп8н + тхэ8э +

тхх1шх 2V

+

тху1шу 2V '

ßv^ а ^ <W ^ 8ЭХ ™уУ1иу ту = m,yß + т,уа + mvHoH + туэ5э +--—--+

mz = mZo + mza + mzBSB +

y " 2V 2V

(m"z + m* )Ьаш2

(3)

2V

Примечание. Далее выполним для удобства замену подстрочных индексов х, у, z на арабские цифры 1, 2, 3.

Динамические уравнения Эйлера с учетом примечания запишутся в виде

М1 d>1 Vs - ¡2)^2^3

М2 = 1 й>2 + (Ii - ¡3)^1^3

ш3\ .й)3. (h - Ii)toito2_

(4)

- управляющие моменты; ^ = 1, 2, 3.

Алгоритмы управления вращательным движением неманевренного ЛА наиболее просто синтезируются методом обратных задач динамики в сочетании с оптимизацией системы (4) по локальным функционалам, характеризующим энергию движения [1, 2].

Основные этапы построения методики синтеза алгоритмов сводятся к следующему.

1. Динамические характеристики проектируемой системы назначаются с помощью эталонных моделей для каждой степени свободы по реакции на «скачок» рулей [6, 8]. Назначенные траектории определяют требуемый пространственный маневр ЛА и главным образом характеризуются свойствами эталонных моделей вращательного движения.

Ovil Aviation High Technologies

Vol. 20, No. 03, 2017

2. Структура алгоритмов и их параметры выбираются из условия отслеживания с необходимой точностью фазовых траекторий эталонных моделей вращательного движения.

3. Степень близости угловых скоростей объекта и эталона оценивается с помощью квадратичных функционалов энергии ускорения.

Процесс движения к экстремуму - минимуму организуется по схеме градиентного метода или по дифференциальной схеме. Такие алгоритмы управления обладают свойствами естественной адаптивности, позволяют осуществить декомпозицию нелинейной системы (4) по степеням свободы за счет организации обратных связей с большими, но ограниченными коэффициентами усиления [2, 10].

Сформулируем задачу синтеза управления при стабилизации угловых скоростей на интервале времени [ г0, ^ ], полагая to = 0, момент времени нефиксированной величиной [2]. В текущий момент времени t = t0 = 0 состояние управляемого объекта характеризуется значениями угловых скоростей ф(0) = ф 0, и требуется найти такие управляющие воздействия ¿г (г = э, н, в), при которых

система (1.4) переходит из состояния ф(0) в новое состояние ф (в общем случае = ) определяет необходимую фигуру вращения) и находится в этом состоянии бесконечно долго. При этом необходимо, чтобы траектория управляемого движения по каждой степени свободы с требуемой точностью следовала за траекторией, определяемой эталонной моделью

который характеризует энергию ускорения вращательного движения в окрестности фазовых траекторий эталонной модели [2].

Задача синтеза решается в два этапа. На первом этапе синтезируются законы управления моментами М 3(ш), а на втором этапе через решение системы алгебраических уравнений (3) и (10) определяются искомые управляющие воздействия 81.

Законы управления моментами М3 (со) по каждой степени свободы получаем, применяя схему простого градиентного метода [1, 2]

<ü*s + YSlü*s + Ys0v*s = Ys0Vs, s = x,у,z,

(5)

где Уб0> Ys0 > 0, установившееся движение устойчиво: ^ со0, при t ^ ю.

Степень приближения управляемого процесса со3(0 ^ ш0 оценивается функционалом

3

(6)

3

(7)

Согласно (4), (6) компоненты градиента в (7) равны

dG(M)

1Щ)

(cb*j — Ö0j), j = 1,2,3 = x,y,z.

(8)

Отсюда следует, что

3

J=1

Vol. 20, No. 03, 2017

или после интегрирования по времени обеих частей

3

Ovil Aviation High Technologies

M,

(ш) = - ù)j),

(10)

j=i

где требуемые значения угловой скорости ш'2 вычисляются из модели эталонного движения (5). Структура матрицы коэффициентов ||?7у|| определяется характером перекрестных связей по угловым скоростям Ш] в (3) и с учетом примечания имеет вид

Гц Г12 0

Г21 г22 0 0 0 г33

Управляющие воздействия ¿1 определяются из уравнений (3, 10) с учетом примечания:

о 1 ( 3 ( ~ ~а \ Г ^у туН^х

5В = — (г33м3 - (r33 + mz г)ш3 - mz0 - т%а), 5Э =

m,

туэ - тхэтун

1 s

Su =— (Ах - тх5э)' т"

~ ~ V'

где mz0 = qSbamz0, = qSbam%' mz z = qSba

)ba

2V

' mzB = qSbamzB,

m,y = qSlm.y, m,y = qSlmy, m.yX =

u>r qSl2 ш

m m =-m

2V У ' У 2V У '

mSy = qSlmzэ' m,yu = qSlmzu' mx = qSlmx' = qSl

qSl2 й-, qSl2 côv ^5э 5э s,

I л ТП = 2У ■ X '"Lx

m -' = ' тихэ = qSlmuxэ' mxu = qSlmxu'

Ах = ГцШ*1 + Г12Ш*2 - (Гц +

m

Ых)ш1 - (r12 + тху^ ш2 - rn-xß -

(11)

Ау = Г21^1 + Г22^2 - (Г21 + ГПуХ)ы1 - (Г22 + ГПуУ) ^ - т^р - т$а.

Синтезированные законы управления рулями & (11) позволяют определить управляющие перегрузки (или выражаемые через них аэродинамические силы ) через аэродинамические коэффициенты [3, 4, 7]:

Cla = с10(а' Ю + АиС^а - 4в0в + с[э8э)

с2а = с20(а,' М) + с2?в (a)ôB + (о)<р + с2з (а)5з + с2,тщ (а)5тщ + с2пр (а)5пр; С3а = 4Н(а' ß)Su + cl(а' ß)ß + 4э0э;

'Fla c1a n1 1 ГР] 1 -Fla

F2a = qS c2a ' V-2 г 0 + DCB F2a ' G = mg'

-P3a- -c3a. П3. и 0 и - F3a -

DCB =

cosacosß -sinacosß sinß

sina cosa 0

cosasinß sinasinß cosß.

Fla F2a ■P3a

= mg

щ

V-2 V-3

Научный Вестник МГТУ ГА_Том 20, № 03, 2017

Civil Aviation High Technologies Vol. 20, No. 03, 2017

Траекторное движение ЛА описывается уравнениями [5, 8]

ет = ft£T, Ук = ftVk + д(п-Е2) (13)

где е2 - второй столбец матрицы направляющих косинусов £ c начальными условиями, определяемыми балансировочными значениями углов Эйлера: = dm, Y0> и положениями рулей Si0;

0 Шз -ш2

ft = -Шз 0 Ш1 - ортогональная, кососимметрическая матрица

- U>2 0

скоростей ш3 = ш5(Ь), вычисляется в пилотажном контуре; Ук = (Ук1Ук2Ук3) - вектор

земной скорости.

Полагая, что угловые скорости ш3 (б = 1...3) поступают из пилотажного контура, получим разрешимую систему дифференциальных уравнений, решение которой сводится к численному или аналитическому интегрированию этих уравнений.

Таким образом, реализуется следующая прямая задача динамики.

1. Известна математическая модель поступательного движения ВС (12, 13) и ее состояние в начальный момент времени:

Ук5(0) = ^о, ^(0) = ^0. (14)

2. Заданы управляющие перегрузки п5(а>(1)^, t Е [*;0, Ьк], - конечное время.

3. Решение сформулированной прямой задачи сводится к простому интегрированию системы (12)-(13) с начальными условиями (1.14) и матрицей направляющих косинусов е^].

4. Для определения текущих углов Эйлера достаточно вместо интегрирования кинематических уравнений Эйлера использовать 3 уравнения связи из матрицы направляющих косинусов

е21 = sind, е22 = cosycosd, £11 = cos^cosd. (15)

5. Далее по формулам

X \Vk1] гч \Vk1] AW1

Н = £ Vk2 , У2 = Vk2 - £ aw2

-Z ivj iVs\ \yks\ AWs.

определяются координаты пространственного движения ВС Ь = [X Н 2] в нормальной СК и составляющие вектора воздушной скорости в скоростной СК.

При отсутствии ветровых возмущений = 0 вектор земной скорости совпадает с вектором воздушной скорости Ук = V, и второе слагаемое уравнения (16) не учитывается. 6. Текущие значения углов атаки и скольжения а, @ определяются согласно (17)

а = агсгд ф2), 0 = агзт (у), V = ^У?2 + +У22 + У32, М = ^ (17)

На рис. 1 представлена структурная схема нелинейной САУ, в которой реализованы предложенные алгоритмы. Для контура вращательного движения входными параметрами являются эталонные и текущие значения угловых скоростей и параметры углового положения ВС. Выходными величинами для этого контура являются угловые отклонения рулевых поверх-

Vol. 20, No. 03, 2017

Civil Aviation High Technologies

ностей. Для контура траекторного управления входными величинами являются выходные параметры контура вращения, начальные значения траекторных параметров. Выходными величинами являются координаты пространственного движения ВС. Для реализации контуров управления требуются измерители угловых скоростей ДУС, измерители координат (акселерометры), измерители углового поля (ГВ, КС), измерители скорости, измерители сдвига ветра.

х Л

L

S

Рис. 1. Структура контуров управления нелинейной САУ Fig. 1. Structure nonlinear ACS control loops

ИССЛЕДОВАНИЕ ФУНКЦИОНИРОВАНИЯ КОНТУРОВ УПРАВЛЕНИЯ

Исследование осуществлялось путем численного моделирования контура вращательного движения с использованием программы MathCad. При этом была выбрана эталонная модель движения ВС, соответствующая динамическим свойствам маломаневренного ВС.

Параметры эталонной модели (5) рассчитывались по формулам

-1 = Р р - и

Yso — j,2 ' Ysi — y ' ps — 2 '

Постоянные времени в эталонной модели определялись путем имитационного моделирования на языке Паскаль по реакции нелинейной модели самолета ИЛ-76М на ступенчатые отклонения рулей высоты, направления и элеронов в 10 градусов.

Моделирование динамики движения ВС проводилось при различных значениях коэффициентов передачи контура управления угловой скорости, что соответствовало различным режимам полета. При всех выбранных значениях коэффициентов система оставалась устойчивой. При больших значениях коэффициента время сходимости было наименьшим.. Исследования показали справедливость разработанных ранее теоретических положений: чем больше значения

Civil Aviation High Technologies

Vol. 20, No. 03, 2017

коэффициентов усиления в законах управления, тем выше близость угловых скоростей вращательного движения к угловым скоростям эталонной модели (о>5(0 при г^^- сю). Результаты моделирования для случая средних значений коэффициентов усиления Г.. представлены соответственно на рис. 2-3.

®t,i ! ät,2

тЗЖ

-■г

Рис. 2. Отклонения рулей Fig. 2. Deviations rudder

Рис. 3. Изменение управляющих моментов Fig. 3. Change of control points

al

o3

Рис. 4. Подстройка угловых скоростей под эталон Fig. 4. Adjusting the angular velocity under standard

ЗАКЛЮЧЕНИЕ

Разработано алгоритмическое обеспечение пилотажного уровня управления ВС на основе метода обратных задач динамики и траекторного уровня на основе прямой задачи динамики. Выполнено моделирование вращательного движения самолета ИЛ-76 на основе алгоритмов управления с естественными свойствами адаптации под эталонную модель.

Сравнительный анализ результатов моделирования показывает следующее.

1. Контур управления вращательным движением остается устойчивым при всех исследуемых режимах полета ВС.

Vol. 20, No. 03, 2017

Civil Aviation High Technologies

2. Чем больше значения коэффициентов усиления в законах управления, тем выше близость угловых скоростей вращательного движения к угловым скоростям эталонной модели.

СПИСОК ЛИТЕРАТУРЫ

1. Крутько П.Д. Обратные задачи динамики управляемых систем. Нелинейные модели. М.: Наука. Гл. ред. физ.-мат. лит., 1988. 328 с.

2. Крутько П.Д. Обратные задачи динамики в теории автоматического управления. М.: Машиностроение. 2004. 576 с.

3. Ефремов А.В. Динамика полета / А.В. Ефремов [и др.]; под ред. Г.С. Бюшгенса. М.: Машиностроение, 2011. 776 с.

4. Бюшгенс Г.С., Студнев Р.В. Динамика продольного и бокового движения. М.: Машиностроение, 1979. 350 с.

5. Красовский А.А. Системы автоматического управления полетом и их аналитическое конструирование. М.: Наука, 1973. 240 с.

6. Красовский А.А. Динамика непрерывных самонастраивающихся систем. М.: Физма-тгиз, 1963. 467 с.

7. Воробьев В.Г., Кузнецов С.В. Автоматическое управление полетом самолетов. М.: Транспорт, 1995. 448 с.

8. Буков В.Н., Сизых В.Н. Метод и алгоритмы решения сингулярно-вырожденных задач аналитического конструирования регуляторов // Изв. АН. Теория и системы управления. 2001. № 5. С. 43-61.

9. Диль В.Ф., Сизых В.Н. Методика синтеза нелинейной системы управления воздушным судном на траекторном и пилотажном уровнях // Международный научно-исследовательский журнал. 2016. № 5 (47). Часть 3. С. 87-91. DOI: 10.18454/IRJ.2016.47.267

10. Пятницкий Е.С. Синтез иерархических систем управления механическими и электромеханическими объектами на принципе декомпозиции // Автомат. и телемех. 1989. № 1. C. 87-99.

11. Пупков К.А. Методы классической и современной теории автоматического управления. Том.3 Синтез регуляторов систем автоматического управления / К.А. Пупков [и др.] / под ред. К.А. Пупкова и Н.Д. Егупова. М.: Изд-во МГТУ им. М.Э. Баумана, 2004. 616 с.

СВЕДЕНИЯ ОБ АВТОРАХ

Диль Виктор Фридрихович, доцент кафедры Авиационных электросистем и пилотаж-но-навигационных комплексов Иркутского филиала МГТУ ГА, [email protected].

Сизых Виктор Николаевич, доктор технических наук, доцент, профессор кафедры Автоматизации производственных процессов Иркутского государственного университета путей сообщения, [email protected].

OPTIMAL AIRCRAFT CONTROL SYNTHESIS BASED ON THE EQUATIONS OF NON-LINEAR DYNAMICS

Viktor F. Dil1, Viktor N. Sizykh2

Moscow State Technical University of Civil Aviation, Irkutsk branch, Irkutsk, Russia

2Irkutsk State Transport University, Irkutsk, Russia

ABSTRACT

The article considers the technique of the synthesis of non-linear aircraft control systems by flight optimization using inverse dynamics problems. To synthesize control algorithms a non-linear model of aircraft flight and trajectory move-

Civil Aviation High Technologies

Vol. 20, No. 03, 2017

ment is used. The authors define method stages of flying level synthesis which include: selection of aircraft reference movements in accordance with three degrees of freedom, structuring the control algorithms and their parameters, defining the proximity of current and reference movements by means of a quadratic functional and further extremum-minimum movement organization by the gradient method. Through the optimized parameters of flying level the direct dynamics problem of trajectory level control of the aircraft spatial movement is solved. The basis for calculating the aircraft trajectory parameters is a non-linear model of the trajectory movement for which flying level output parameters serve as input data. The trajectory level output parameters are defined by numerical integration of input signals considering aircraft dynamic blow coefficients. The structure diagram of aircraft spatial movement control organization is developed. The flight contour functioning is examined using numerical modeling in MathCad and Paskal programs. Reference parameters were determined by Paskal simulation modeling according to the reaction of a non-linear aircraft model to the "bounces" of aerodynamical flight controls. It is shown that the spatial control problem is optimal in terms of input control realization. Besides, in comparison with [9] it is possible to state that due to energy reversibility of rotational and progressive movements only the content of direct and inversed problems of dynamics changes.

Key words: automatic control, spatial movement, equations of non-linear dynamics.

REFERENCES

1. Krutko P.D. Obratnye zadachi dinamiki upravljaemyh sistem. Nelinejnye modeli [The inverse problems of dynamics of the operated systems. Nonlinear models]. M., Science. Fizmatlit, 1988, 328 p. (in Russian)

2. Krutko P.D. Obratnye zadachi dinamiki v teorii avtomaticheskogo upravlenija [The in-versed problems of dynamics in the theory of automatic control.]. M., Mechanical engineering, 2004. 576 p. (in Russian)

3. Efremov A.V. et al. Dinamikapoleta [Dynamics of flight]. Ed. by Byushgens G.S. M., Mechanical engineering, 2011. 776 p. (in Russian)

4. Byushgens G.S., Studnev R.V. Dinamika prodol'nogo i bokovogo dvizhenija [Dynamics of longitudinal and lateral motion]. M., Mechanical engineering, 1979, 350 p. (in Russian)

5. Krasovsky A.A. Sistemy avtomaticheskogo upravlenija poletom i ih analiticheskoe kon-struirovanie [Systems of an automatic flight control and their analytical designing]. M., Science, 1973, 240 p. (in Russian)

6. Krasovsky A.A. Dinamika nepreryvnyh samonastraivajushhihsja sistem [The dynamics of continuous adaptive systems]. M., Fizmatgiz, 1963, 467 p. (in Russian)

7. Vorobyov V.G., Kuznetsov S.V. Avtomaticheskoe upravlenie poletom samoletov [Automatic control of aircraft flight]. M., Transport, 1995, 448 p. (in Russian)

8. Bukov V.N., Sizykh V.N. Metod i algoritmyi resheniya singulyarno-vyirozhdennyih zadach analiticheskogo konstruirovaniya regulyatorov [Method and algorithms of the solution of singular and degenerate problems of analytical designing of regulators]. Journal of Computer and Systems Sciences International, 2001, no. 5, pp. 43-61. (in Russian)

9. Dil V.F., Sizykh V.N. Metodika sinteza nelineynoy sistemyi upravleniya vozdushnyim sud-nom na traektornom i pilotazhnom urovnyah [Method of Synthesis of a nonlinear system of Aircraft Control at a trajectory and flight levels]. International Research Journal, May, 2016, no. 5 (47), part 3, pp. 87-91. DOI: 10.18454/IRJ.2016.47.267 (in Russian)

10. Pyatnitsky E.S. Sintez ierarhicheskih sistem upravleniya mehanicheskimi i elektrome-hanicheskimi ob'ektami na printsipe dekompozitsii [Synthesis of hierarchical control systems of mechanical and electromechanical objects on decomposition principle]. Automation and Remote Control, 1989, no. 1, pp. 87-99. (in Russian)

11. Pupkov K.A. Metody klassicheskoj i sovremennoj teorii avtomaticheskogo upravlenija. Tom.3 Sintez reguljatorov sistem avtomaticheskogo upravlenija [Methods of classical and modern automatic control theory. Vol. 3 Synthesis of regulators of automatic control systems]. Ed. by K.A. Pupkov and N.D. Egupov. M., Bauman MSTU Publishing House, 2004, 616 p. (in Russian)

Vol. 20, No. 03, 2017 Civil Aviation High Technologies

INFORMATION ABOUT THE AUTHORS

Viktor F. Dil, Associate Professor, Moscow State Technical University of Civil Aviation, Irkutsk Branch, [email protected].

Viktor N. Sizykh, Doctor of Technical Sciences, Associate Professor, Full Professor in Irkutsk State Transport University, [email protected].

Поступила в редакцию 18.04.2016 Received 18.04.2016

Принята в печать 27.04.2017 Accepted for publication 27.04.2017