Relaxation oscillations in the Darie wind power plant model

Full Text

1. Предварительные сведения

В данной статье проводится математическое моделирование колебаний ветроэнергетической установки Дарье (ВЭУ). В статьях [1-6] рассматривалась математическая модель ветрогенератора, состоящая из двух уравнений. В [1; 2] были предложены новые типы конструкции ветроэнергетической установки (рис 1.1.). Исследуются движения электромеханической системы, моделирующей стационарные режимы генератора с ветротурбинным движителем с помощью качественных методов теоретической механики в статьях [3-5]. Рассмотрены вопросы существования и устойчивости установившихся режимов и описаны области их притяжения [3-5]. Анализ критических режимов, основанный на использовании техники траекторий уток, был проведен в [7]. В работе [6] исследовано уравнение стационарных режимов при различных значениях внешнего сопротивления динамической модели. Также описаны различные аспекты конструктивной теории малых ветроэлектростанций [8; 9]. Математическая модель малой ветроэнергетической установки Дарье представляет собой трехтемповую дифференциальную систему:

$\begin{array}{rrr}J\dot{\Omega }& =& -kI+M\left(\Omega \right),\\ L\dot{I}& =& k\Omega -(R+r)I,\\ \dot{R}& =& \epsilon F(\Omega ,I,R),\end{array}$ (1.1)

где $J$ — момент инерции турбины, $L$ — индуктивность якоря, $I$ — сила тока в обмотке якоря, $\Omega$ — коэффициент быстроходности, $R$ — внешнее сопротивление, $M$ — относительный момент аэродинамических сил. Параметр $k$ — коэффициент электромеханического взаимодействия, $r$ — малое внутреннее сопротивление якоря. Точка означает производную по времени $t$. График экспериментальных данных на рис. 2.1 задает функцию $M\left(\Omega \right)$ [3]. Величины $L$ и $\epsilon$ рассматриваются как малые параметры. Это означает, что в данной системе переменная $I$ является самой быстрой, $R$ — самой медленной из трех переменных.

Рис. 1.1.: Ветрогенератор Дарье

Fig. 1.1. Wind turbine Darrieus

2. Основные результаты

2.1. Аппроксимация функции $\mathrm{M\Omega }$

Рассмотрим график экспериментальных данных $M\Omega$, изображенный на рис. 2.1. Относительный момент аэродинамических сил и коэффициент быстроходности на рис. 2.1 являются безразмерными ($M=\frac{M_a}{\mathrm{0,5}\rho Sb{V}^2}$, $\Omega =\frac{b\omega }{V}$, где $\omega$ — угловая скорость, $b$ — расстояние от эффективного давления лопастей для оси вращения, $V$ — скорость воздуха, $M_a$ — момент аэродинамических сил, $\rho$ — плотность воздуха, $S$ — площадь лопастей) [3; 5; 7]. Для построения функции $M\left(\Omega \right)$ зададим ее как кусочно-линейную функцию по промежуткам:

$\Omega \in [0;0.86];[0.86;1.8];[1.8;4];[4;8.64]$

и аппроксимируем на каждом из участков отдельно.

Рис. 2.1: Экспериментальные данные

Fig. 2.1. Experimental data

В результате получается следующее задание функции $M\left(\Omega \right)$:

$M\left(\Omega \right)=\left\{\begin{array}{r}0.02906976745\Omega +\mathrm{0.12,}\\ где\Omega \in [0;0.86];\\ -0.0638297872\Omega +\mathrm{0.199893617,}\\ где\Omega \in [0.86;1.8];\\ 0.643181818\Omega -\mathrm{1.0727272725,}\\ где\Omega \in [1.8;4];\\ -0.28556034485\Omega +\mathrm{2.6422413795,}\\ где\Omega \in [4;8.64].\end{array}\right.$ (2.1)

Рисунок 2 показывает график функции $M\left(\Omega \right)$, заданной (2.1).

Рис. 2.2: График функции $\mathbf{M}\mathbf{\left(}\mathbf{\Omega }\mathbf{\right)}$, заданной (1)

Fig. 2.2. Graph of the function $\mathbf{M}\mathbf{\left(}\mathbf{\Omega }\mathbf{\right)}$, given by (1)

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

2.2 Медленная кривая

Система (1.1) является сингулярно возмущенной, так как малый параметр $L$ умножается на производную силы тока. Система (1.1) имеет двумерное инвариантное многообразие [7; 10; 11]

$I=h(\Omega ,I,R).$ (2.2)

Функция $h$ может быть найдена из уравнения инвариантности

$L\frac{\partial h}{\partial R}\epsilon F(\Omega ,I,R)+L\frac{\partial h}{\partial \Omega }\frac{1}{J}(-kh+M(\Omega \left)\right)=k\Omega -(R+r)h.$ (2.3)

Пренебрегая членами порядка $O\left(L\right)$, получается

$h=\frac{k\Omega }{R+r}.$ (2.4)

Движение по инвариантному многообразию описывается дифференциальной системой, где $k=\mathrm{0.5,}r=\mathrm{0.1,}J=1$:

$\begin{array}{rrr}J\dot{\Omega }& =& -\frac{k^2\Omega }{R+r}+M\left(\Omega \right),\\ \dot{R}& =& \epsilon F(\Omega ,I,R).\end{array}$ (2.5)

Полученная система является быстро-медленной, так как правая часть уравнения для медленной переменной $R$ умножается на малый параметр $\epsilon$. В связи с этим мы можем использовать аппарат теории релаксационных колебаний для анализа этой системы [12; 13]. Медленная кривая системы (2.5) задается уравнением:

$-\frac{k^2\Omega }{R+r}+M\left(\Omega \right)=0.$ (2.6)

Из (6) получим уравнение медленной кривой в явной форме, для этого выразим $R\left(\Omega \right)$:

$R=R\left(\Omega \right)=\frac{k^2\Omega }{M\left(\Omega \right)}-r.$ (2.7)

График медленной кривой представлен на рис. 2.3, при $k=\mathrm{0.5,}r=\mathrm{0.1,}J=1$.

Рис. 2.3: Графики медленной кривой

Fig. 2.3. Slow curve charts

Для исследования устойчивости медленной кривой (2.6) необходимо найти производную функции $R\left(\Omega \right)$ на каждом из промежутков отдельно:

$\frac{\partial R\left(\Omega \right)}{\partial \Omega }=\frac{k^{2}M\left(\Omega \right)-M^{\text{'}}\left(\Omega \right)k^2\Omega }{M^2\left(\Omega \right)}.$ (2.8)

На первом участке медленной кривой при $\Omega \in [0;0.86]$: пусть $\Omega =\mathrm{0.5,}R=0.8291270529$, тогда $\frac{\partial R\left(\Omega \right)}{\partial \Omega }(0.5;0.8291270529)=-0.2400000001<0$, следовательно, участок является устойчивым.

На втором участке медленной кривой при $\Omega \in [0.86;1.8]$: пусть $\Omega =\mathrm{1,}R=1.737373042$, тогда $\frac{\partial R\left(\Omega \right)}{\partial \Omega }(1;1.737373042)=-0.1998936100<0$, следовательно, участок является устойчивым.

На третьем участке медленной кривой при $\Omega \in [1.8;4]$: пусть $\Omega =\mathrm{3,}R=0.77533156489$, тогда $\frac{\partial R\left(\Omega \right)}{\partial \Omega }(3;0.7753315648)=0.3575757572>0$, следовательно, участок является неустойчивым.

На четвертом участке медленной кривой при $\Omega \in [4;8.64]$: пусть $\Omega =\mathrm{6,}R=1.514849187$, тогда $\frac{\partial R\left(\Omega \right)}{\partial \Omega }(6;1.514849187)=-0.4403735633<0$, следовательно, участок является устойчивым.

В результате устойчивыми участками медленной кривой являются промежутки:

$\Omega \in [0;0.86];[0.86;1.8];[4;8.64],$

а неустойчивым участком:

$\Omega \in [1.8;4].$

Точками срыва являются точки экстремумов медленной кривой, следовательно, исходя из рис. 2.3 имеются две точки срыва с координатами: $A(1.8;5.194117653)$ и $B(4;0.5666666669)$. Точки $A$ и $B$ являются точками смены устойчивости (рис. 2.3).

2.3 Релаксационные колебания

Пусть $F(\Omega ,I,R)=\epsilon (\alpha -\Omega )$. Тогда система уравнений (2.5) будет выглядеть следующим образом:

$\begin{array}{rrr}J\dot{\Omega }& =& -\frac{k^2\Omega }{R+r}+M\left(\Omega \right),\\ \dot{R}& =& \epsilon (\alpha -\Omega ).\end{array}$ (2.9)

Зададим точку $({\Omega }_{*},R_{*})$ как положение равновесия системы уравнений (2.9):

$\begin{array}{rrr}-\frac{k^2{\Omega }_{*}}{R_{*}+r}+M\left({\Omega }_{*}\right)& =& \mathrm{0,}\\ \epsilon (\alpha -{\Omega }_{*})& =& 0.\end{array}$ (2.10)

Выполнив преобразования, получим координаты положения равновесия:

$\begin{array}{rrr}{\Omega }_{*}& =& \alpha ,\\ R_{*}& =& \frac{k^2{\Omega }_{*}}{M\left({\Omega }_{*}\right)}-r.\end{array}$ (2.11)

Исследуем на устойчивость положения равновесия при различных значениях параметра $\alpha \mathrm{.}$

Матрица Якоби системы (2.9) в общем виде будет выглядеть следующим образом:

$\begin{array}{r}A=\left(\begin{array}{ll}\frac{\partial f}{\partial \Omega }& \frac{\partial f}{\partial R}\\ & \\ \frac{\partial g}{\partial \Omega }& \frac{\partial g}{\partial R}\end{array}\right).\end{array}$ (2.12)

Рассмотрим 3 случая, когда положение равновесия находится на:

1. на первом устойчивом участке медленной кривой при $\Omega \in [0;1.8]$;
2. на втором устойчивом участке медленной кривой при $\Omega \in [4;8.64]$;
3. на неустойчивом участке медленной кривой при $\Omega \in [1.8;4]$.

Рассмотрим первый случай. Так как момент аэродинамических сил задается кусочно-линейным образом и $\Omega \in [0;1.8]$, то возьмем уравнение момента аэродинамических сил для второго участка.

Матрица Якоби в точке $(\alpha ,\frac{k^2\alpha }{\alpha }-r)$, при $\Omega \in [0.86;1.8],k=\mathrm{0.5,}r=\mathrm{0.1,}J=1$:

$A{|}_{(\alpha ,\frac{k^2\alpha }{\alpha }-r)}==\left(\begin{array}{ll}-\frac{0.199893617-0.0638297872\alpha }{\alpha }-0.0638297872& \frac{{(0.199893617-0.0638297872\alpha )}^2}{k^2\alpha }\\ & \\ -1& 0\end{array}\right).$ (2.13)

Характеристическое уравнение будет равно:

${\lambda }^2+\left(\frac{0.19989361700}{\alpha }\right)\lambda +0.016296966936\alpha -0.10207333629+\frac{0.15982983247}{\alpha }=0.$ (2.14)

Найдем дискриминант и корни уравнения (2.14) при $\alpha =1$:

$\begin{array}{r}{\lambda }_1=-0.099946808500-0.25310886707i;\\ {\lambda }_2=-0.099946808500+0.25310886707i.\end{array}$ (2.15)

При полученных собственных значениях положение равновесия является асимптотически устойчивым фокусом (рис. 2.4).

Рис. 2.4: Медленная кривая (тонкая линия) и траектория системы (жирная линия) при $\mathbf{\epsilon }\mathbf{=}\mathbf{0.01}\mathbf{;}\mathbf{k}\mathbf{=}\mathbf{0.5}\mathbf{;}\mathbf{r}\mathbf{=}\mathbf{0.1}\mathbf{;}\mathbf{J}\mathbf{=}\mathbf{1}\mathbf{;}\mathbf{\alpha }\mathbf{=}\mathbf{1}\mathbf{;}$ начальная точка: $\mathbf{\Omega }\mathbf{\left(}\mathbf{0}\mathbf{\right)}\mathbf{=}\mathbf{7.5,}\mathbf{R}\mathbf{\left(}\mathbf{0}\mathbf{\right)}\mathbf{=}\mathbf{2.3}$

Fig. 2.4. Slow curve (thin line) and system trajectory (thick line) at $\mathbf{\epsilon }\mathbf{=}\mathbf{0.01}\mathbf{;}\mathbf{k}\mathbf{=}\mathbf{0.5}\mathbf{;}\mathbf{r}\mathbf{=}\mathbf{0.1}\mathbf{;}\mathbf{J}\mathbf{=}\mathbf{1}\mathbf{;}\mathbf{\alpha }\mathbf{=}\mathbf{1}\mathbf{;}$; starting point: $\mathbf{\Omega }\mathbf{\left(}\mathbf{0}\mathbf{\right)}\mathbf{=}\mathbf{7.5,}\mathbf{R}\mathbf{\left(}\mathbf{0}\mathbf{\right)}\mathbf{=}\mathbf{2.3}$

Рассмотрим второй случай. Так как момент аэродинамических сил задается кусочно-линейным образом и $\Omega \in [4;8.64]$, то возьмем уравнение момента аэродинамических сил для четвертого участка. Действуя аналогичным образом, как и в предыдущем случае, находим корни характеристического уравнения при $\alpha =4.5$:

$\begin{array}{r}{\lambda }_1=-0.29358237550-1.2454650976i;\\ {\lambda }_2=-0.29358237550+1.2454650976i.\end{array}$ (2.16)

При полученных собственных значениях положение равновесия является асимптотически устойчивым фокусом (рис. 2.5).

Рис. 2.5: Медленная кривая (тонкая линия) и траектория системы (жирная линия) при $\mathbf{\epsilon }\mathbf{=}\mathbf{0.01}\mathbf{;}\mathbf{k}\mathbf{=}\mathbf{0.5}\mathbf{;}\mathbf{r}\mathbf{=}\mathbf{0.1}\mathbf{;}\mathbf{J}\mathbf{=}\mathbf{1}\mathbf{;}\mathbf{\alpha }\mathbf{=}\mathbf{4.5}\mathbf{;}$ начальная точка: $\mathbf{\Omega }\mathbf{\left(}\mathbf{0}\mathbf{\right)}\mathbf{=}\mathbf{7.5,}\mathbf{R}\mathbf{\left(}\mathbf{0}\mathbf{\right)}\mathbf{=}\mathbf{2.3}$

Fig. 2.5. Slow curve (thin line) and system trajectory (thick line) at $\mathbf{\epsilon }\mathbf{=}\mathbf{0.01}\mathbf{;}\mathbf{k}\mathbf{=}\mathbf{0.5}\mathbf{;}\mathbf{r}\mathbf{=}\mathbf{0.1}\mathbf{;}\mathbf{J}\mathbf{=}\mathbf{1}\mathbf{;}\mathbf{\alpha }\mathbf{=}\mathbf{4.5}\mathbf{;}$; starting point: $\mathbf{\Omega }\mathbf{\left(}\mathbf{0}\mathbf{\right)}\mathbf{=}\mathbf{7.5,}\mathbf{R}\mathbf{\left(}\mathbf{0}\mathbf{\right)}\mathbf{=}\mathbf{2.3}$

Наконец, рассмотрим третий случай. Так как момент аэродинамических сил задается кусочно-линейным образом и $\Omega \in [4;8.64]$, то возьмем уравнение момента аэродинамических сил для третьего участка. Действуя аналогичным образом, находим корни характеристического уравнения при $\alpha =3$:

$\begin{array}{r}{\lambda }_1=0.17878787875-0.97308003575i;\\ {\lambda }_2=0.17878787875+0.97308003575i.\end{array}$ (2.17)

Рис. 2.6: Медленная кривая (тонкая линия) и предельный цикл (жирная линия) при $\mathbf{\epsilon }\mathbf{=}\mathbf{0.01}\mathbf{;}\mathbf{k}\mathbf{=}\mathbf{0.5}\mathbf{;}\mathbf{r}\mathbf{=}\mathbf{0.1}\mathbf{;}\mathbf{J}\mathbf{=}\mathbf{1}\mathbf{;}\mathbf{\alpha }\mathbf{=}\mathbf{3}\mathbf{;}$ начальная точка: $\mathbf{\Omega }\mathbf{\left(}\mathbf{0}\mathbf{\right)}\mathbf{=}\mathbf{7.5,}\mathbf{R}\mathbf{\left(}\mathbf{0}\mathbf{\right)}\mathbf{=}\mathbf{2.3}$

Fig. 2.6. Slow curve (thin line) and limit cycle (thick line) at $\mathbf{\epsilon }\mathbf{=}\mathbf{0.01}\mathbf{;}\mathbf{k}\mathbf{=}\mathbf{0.5}\mathbf{;}\mathbf{r}\mathbf{=}\mathbf{0.1}\mathbf{;}\mathbf{J}\mathbf{=}\mathbf{1}\mathbf{;}\mathbf{\alpha }\mathbf{=}\mathbf{3}\mathbf{;}$ starting point: $\mathbf{\Omega }\mathbf{\left(}\mathbf{0}\mathbf{\right)}\mathbf{=}\mathbf{7.5,}\mathbf{R}\mathbf{\left(}\mathbf{0}\mathbf{\right)}\mathbf{=}\mathbf{2.3}$

При полученных собственных значениях положение равновесия находится на неустойчивом участке медленной кривой и является неустойчивым фокусом. В этом случае в системе наблюдается релаксационный цикл (рис. 2.6). При этом траектория сначала движется вдоль устойчивого участка медленной кривой, далее особая точка проходит точку экстремума (точку срыва) и при этом теряет устойчивость, что соответствует бифуркации Андронова — Хопфа, т. е. при прохождении параметра $\alpha$ через бифуркационное значение устойчивое положение равновесия дифференциальной системы становится неустойчивым. Далее срывается и переходит в $\epsilon$-окрестность другого устойчивого участка, следует вдоль него до второй точки срыва, после которого возвращается на первый устойчивый участок (рис. 2.6) [12].

Суммируем вышесказанное. В зависимости от изменении значения параметра $\alpha$ особая точка меняет свое положение на медленной кривой. Если положение равновесия лежит на устойчивой части медленной кривой, то оно является асимптотически устойчивым. С физической точки зрения такая ситуация представляется предпочтительной, поскольку соответствует стационарному режиму работы системы. Случай, когда положение равновесия находится на неустойчивом участке и возникает релаксационный цикл, представляется неприемлемым с физической точки зрения, так как в системе возникают колебания с достаточно большой амплитудой.

Выводы

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

About the authors

Alena S. Kirsanova

Author for correspondence.
Email: askirsanova99@gmail.com
ORCID iD: 0009-0005-1799-1423

postgraduate student of the Department of Differential Equations and Control Theory

References

Supplementary files

Supplementary Files

Action

1. JATS XML

Download

2. Fig. 1.1. Wind turbine Darrieus

Download (245KB)

Indexing metadata

3. Fig. 2.1. Experimental data

Download (60KB)

Indexing metadata

4. Fig. 2.2. Graph of the function , given by (1)

Download (69KB)

Indexing metadata

5. Fig. 2.3. Slow curve charts

Download (53KB)

Indexing metadata

6. Fig. 2.4. Slow curve (thin line) and system trajectory (thick line) at ; starting point:

Download (40KB)

Indexing metadata

7. Fig. 2.5. Slow curve (thin line) and system trajectory (thick line) at ; starting point:

Download (36KB)

Indexing metadata

8. Fig. 2.6. Slow curve (thin line) and limit cycle (thick line) at starting point:

Download (39KB)

Indexing metadata