2. 쿼드로터 시스템 모델링 및 상태방정식
그림 1의 쿼드로터 좌표계에서 각각의 모터에 대한 입력을 $u_{j}$$(j=1,\:2,\:3,\:4)$로 정의한다.
그림 1 쿼드로터 좌표계
Fig. 1 Coordinate system of quadrotor
쿼드로터의 고도 및 자세에 대한 동역학 방정식은 뉴턴-오일러 방정식을 통해 다음과 같이 구할 수 있다 (9),(10).
여기서 $z$는 $z$축 방향의 높이, $\phi$, $\theta$, $\psi$는 각각 roll, pitch, yaw이며 $I_{x}$, $I_{y}$,
$I_{z}$는 각 축에 해당하는 관성 모멘트, $g$는 중력가속도, $l$은 무게중심과 rotor 사이의 거리를 나타낸다. $m(t)$는 쿼드로터의
질량에 해당하며 공칭 질량 $m_{0}$와 물품 적재 등의 이유로 추가되는 질량에 대한 불확실성인 $\Delta m(t)$ 합으로, $m(t)=m_{0}+\Delta
m(t)$로 표현 된다. 질량에 대한 불확실성은 공칭 질량의 $20$% 이내의 값으로 $\Delta m(t)\le 0.2m_{0}$을 만족한다고 가정한다.
모델링 되지 않은 바람 등의 불확실한 외란은 $d_{i}(t)$$(i=z,\:\phi ,\:\theta ,\:\psi)$와 같이 나타내며 $\left
| d_{i}(t)\right |\le\overline{d}_{i}$$(i=z,\:\phi ,\:\theta ,\:\psi)$으로 bound 되어
있다 (11). $U_{i}$$(i=z,\:\phi ,\:\theta ,\:\psi)$는 고도, roll, pitch, yaw에 대한 입력이며, 이는 그림 1에서 각각의 모터 입력에 대해 다음과 같이 표현 가능하다.
여기서 $K_{t}$는 추력상수, $K_{y}$는 항력상수를 나타낸다. 식 (2)를 각 모터의 입력에 대해 정의하면 다음과 같다.
식 (1)의 동역학 방정식은 아래와 같은 상태방정식으로 표현 가능하며
여기서 상태변수 $x=[x_{1},\: x_{2},\: x_{3,\:}x_{4},\: x_{5},\: x_{6},\: x_{7},\: x_{8}]^{T}=[z,\:\dot
z ,\:\phi ,\:\dot\phi ,\:\theta ,\:$ $\dot\theta ,\:\psi ,\:\dot\psi]^{T}$이다. 쿼드로터
시스템에서의 파라미터는 표 1과 같다.
표 1 쿼드로터의 파라미터 (12)
Table 1 parameters of quadrotor (12)
Parameters
|
Value
|
$I_{z}$
|
$22.3\times 10^{-4}{N}\cdot{m}\cdot{s}^{2}$
|
$I_{y}$
|
$29.8\times 10^{-4}{N}\cdot{m}\cdot{s}^{2}$
|
$I_{z}$
|
$48\times 10^{-4}{N}\cdot{m}\cdot{s}^{2}$
|
$m_{0}$
|
$0.429{kg}$
|
$l$
|
$0.1785{m}$
|
$g$
|
$9.81$m/s²
|
$K_{t}$
|
$0.0161{N}$
|
$K_{y}$
|
$0.0159{N}\cdot{m}$
|
식 (4)에서 $z$축 높이에 대한 속도에 해당하는 상태 $x_{2}$의 경우 공칭 모델과 불확실성을 구분하면 다음과 같이 정리 가능하다.
여기서 $\Delta\widetilde m(t)=\Delta m(t)/(m_{0}+\Delta m(t))$, $\widetilde d_{z}(t)=-\Delta\widetilde
m(t)g+d_{z}(t)$ 이다. 쿼드로터의 오차 벡터를 아래와 같이 정의한다.
여기서 $x_{1d}$, $x_{3d}$, $x_{5d}$, $x_{7d}$는 각각 $z$축 높이, roll, pitch, yaw가 추종하려는 값이며
식 (6)을 바탕으로 오차에 대한 상태방정식은 다음과 같이 표현된다.
식 (7)은 비선형성을 포함하고 있다. 이는 다음과 같이 입력 $u_{z}$, $u_{\phi}$, $u_{\theta}$, $u_{\psi}$로 치환하여
궤환 선형화 시킬 수 있으며, 이 때 사용되는 제어 입력은 다음과 같다.
최종적으로 식 (8)을 식 (7)에 대입하여 정리되는 선형화 된 쿼드로터의 오차에 대한 상태방정식은 다음과 같다.
식 (9)에서 $z$축 높이, roll, pitch, yaw에 대해 각각의 상태 값을 받아와 독립적인 제어가 이루어진다. 이를 전체적인 개략도로 나타내면 그림 2와 같으며, $z$축 높이의 경우 초기 추력을 위해 초기에는 상수 값을 입력으로 한다.
그림 2 제어 개략도
Fig. 2 schematic diagram of control
3. $bold\gamma -$슬라이딩 모드 제어기 설계 및 분석
$bold z$ position : 이득 조절 요소 $\gamma_{1}$이 포함된 $\gamma -$슬라이딩 평면을 다음과 같이 정의한다.
여기서 $\gamma_{1}>0$은 이득 조절 요소, $k_{1}>0$은 제어 이득이다. $\gamma -$슬라이딩 평면 식 (10)을 미분하면 다음과 같다.
제어 입력 $u_{z}$에 대해 다음과 같이 정의한다,
여기서 $\beta_{1}>0$, $0<\epsilon_{1}<1$이며, ${sat}\left(s_{z}/\epsilon_{1}\right)$은 다음과
같다.
이때 $\epsilon_{1}$은 작을수록 ${sgn}(s_{z})$에 근사하게 된다. $\epsilon_{1}=\epsilon_{1}^{*}$으로
정의하면 $\left | s_{z}\right | >\epsilon_{1}^{*}$인 영역을 reaching phase, $\left | s_{z}\right
|\le\epsilon_{1}^{*}$인 영역을 boundary layer라고 한다 (14).
먼저 reaching phase에서 Lyapunov함수 $V\left(s_{z}\right)=(1/2)s_{z}^{2}$ 으로 정의한 후, 이를 식 (11)의 궤적에 따라 미분하면 다음과 같다.
여기서 질량에 대한 $\Delta\widetilde m(t)\le 0.2m_{0}$의 가정으로부터 다음과 같은 관계가 성립한다.
먼저 $\beta_{1}=\beta_{1}^{*}>0$인 유한한 값으로 정의하면 식(14)-(15)에서 얻을 수 있는 $\gamma_{1}$의 범위는 다음과 같다.
또한 식 (14)에서 $\xi_{2}$는 다음과 같은 조건을 만족해야 한다.
이때 $\xi_{2}$에 대한 조건은 $k_{1}$을 통해 조절 가능하며, 여기서 $k_{1}=k_{1}^{*}$로 정의한다. 최종적으로 식(16)-(17)의 조건을 만족하면 $\left | s_{z}\right |$는 유한한 시간에 boundary layer에 도달하게 된다.
Boundary layer에서의 제어입력은 식 (13)를 통해 다음과 같이 정의된다.
식 (18)의 제어 입력을 상태방정식에 대입하여 정리하면 다음과 같은 폐루프 시스템을 얻을 수 있다.
식 (19)은 다음과 같은 행렬의 형태로 표현 가능하며
여기서 $\xi_{z}=[\xi_{1},\:\xi_{2}]^{T}$, $\delta_{z}(t,\:s)=\left[0,\:\Delta\widetilde
m\left(\beta_{1}^{*}/\gamma_{1}\epsilon_{1}^{*}\right)s_{z}\right]^{T}$, $D_{z}(t)=[0,\:$$\widetilde
d_{z}(t)]^{T}$이며 $A_{z}(\gamma_{1})$은 다음과 같다.
$E_{z}(\gamma_{1})={diag}[1,\:\gamma_{1}]$으로 정의한 후, $A_{z}$에 대하여 다음과 같이 정의한다.
이때 $A_{z}$의 특성방정식은 다음과 같다.
식 (23)이 Hurwitz polynomial이 되기 위한 $\beta_{1}^{*}$, $k_{1}^{*}$, $\epsilon_{1}^{*}$의 범위는
Routh-Hurwitz 판별법을 이용하여 구할 수 있다 (14). 따라서 $\beta_{1}^{*}>0$, $k_{1}^{*}>0$, $\epsilon_{1}^{*}>0$이면 Hurwitz polynomial이며,
각 제어이득은 이를 만족하므로 식 (23)은 Hurwitz polynomial이다.
$E_{z}(\gamma_{1})$과 $A_{z}$는 다음과 같은 관계를 만족한다.
$A_{z}$가 Hurwitz 이므로 Lyapunov방정식 $A_{z}^{T}P_{z}+P_{z}A_{z}=-I$를 만족하며, 식 (24)를 $A_{z}$에 대한 식으로 이항 후 Lyapunov방정식에 대입하여 정리하면 다음과 같다.
여기서 $P_{z}(\gamma_{1})=E_{z}(\gamma_{1})P_{z}E_{z}(\gamma_{1})$이다.
Lyapunov 함수 $V(\xi_{z})=\xi_{z}^{T}P_{z}(\gamma_{1})\xi_{z}$로 정의한 후, 이를 식 (20)의 궤적에 따라 미분하면 다음과 정리 가능하다.
식 (26)는 다음과 같은 부등식을 만족한다.
식 (27)에서 마지막 항은 다음과 같이 행렬의 형태로 표현할 수 있다.
따라서 식 (28)을 바탕으로 $\xi_{1}$, $\xi_{2}$의 ultimate bound(UB)를 구하면 다음과 같다.
따라서 $z$축 높이에 대한 상태 $\xi_{1}$에 대하여 $\gamma_{1}$을 조절하면 시스템의 허용 한도 내에서 ultimate bound가
임의로 조절 가능하다.
Angle : roll에 대한 $\gamma -$슬라이딩 평면을 다음과 같이 정의한다.
여기서 $\gamma_{2}>0$은 이득조절 요소, $k_{2}>0$는 제어 이득이다. $\gamma -$슬라이딩 평면인 식 (30)을 미분하면 다음과 같이 정리된다.
roll에 대한 제어입력 $u_{\phi}$를 아래와 같이 정의하며
여기서 $\beta_{2}>0$, $0<\epsilon_{2}<1$을 만족하며, $\epsilon_{2}=\epsilon_{2}^{*}$으로 정의한다.
먼저 reaching phase에서 $s_{\phi}$에 대한 Lyapunov 함수 $V(s_{\phi})=(1/2)s_{\phi}^{2}$로 정의한
후, 식 (31)의 궤적에 따라 미분하면 다음과 같다.
식 (33)에서 $\beta_{2}=\beta_{2}^{*}$으로 정의하면 $\gamma_{2}$는 다음과 같은 범위를 가진다.
또한 $\xi_{4}$에 대하여 다음과 같은 조건을 만족해야 하며,
이는 $k_{2}$를 통해 조절 가능하며, 여기서 $k_{2}=k_{2}^{*}$으로 정의한다. 따라서 식(34)-(35)을 만족하면 유한한 시간에 boundary layer에 도달하게 된다.
Boundary layer에서의 제어 입력은 다음과 같이 정의된다.
식 (36)를 식 (9)에서 roll에 대한 상태방정식에 대입하여 정리하면 다음과 같이 표현 가능하다.
식 (37)은 다음과 같은 행렬의 형태로 표현 가능하며,
여기서 $\xi_{\phi}=[\xi_{3},\:\xi_{4}]^{T}$, $D_{\phi}(t)=[0,\: d_{\phi}(t)]^{T}$이며 $A_{\phi}(\gamma_{2})$는
다음과 같다.
또한 $E_{\phi}(\gamma_{2})={diag}[1,\:\gamma_{2}]$, $A_{\phi}$에 대하여 다음과 같이 정의한다.
이때 $A_{\phi}$의 특성방정식은 다음과 같다.
식 (41)에서 Hurwitz polynomial이 되기 위한 $\beta_{2}^{*}$, $k_{2}^{*}$, $\epsilon_{2}^{*}$의 범위는
Routh-Hurwitz 판별법을 이용하여 구할 수 있다. 따라서 $\beta_{2}^{*}>0$, $k_{2}^{*}>0$, $\epsilon_{2}^{*}>0$이면
Hurwitz polynomial이며, 각 제어이득은 이를 만족하므로 식 (41)은 Hurwitz polynomial이다.
$E_{\phi}(\gamma_{2})$과 $A_{\phi}$는 다음과 같은 관계를 만족한다.
$A_{\phi}$가 Hurwitz 이므로 $A_{\phi}^{T}P_{\phi}+P_{\phi}A_{\phi}=-I$ 를 만족하며, 식 (41)을 $A_{\phi}$에 대한 식으로 이항 후 Lyapunov방정식에 대입하여 정리하면 다음과 같다.
여기서 $P_{\phi}(\gamma_{2})=E_{\phi}(\gamma_{2})P_{\phi}E_{\phi}(\gamma_{2})$이다.
Lyapunov 함수 $V(\xi_{\phi})=\xi_{\phi}^{T}P_{\phi}(\gamma_{2})\xi_{\phi}$로 정의한 후, 이를
식 (38)의 궤적에 따라 미분하면 다음과 정리 가능하다.
식 (44)는 다음과 같은 부등식을 만족한다.
식 (45)의 마지막 식은 다음과 같이 행렬의 형태로 풀어 쓸 수 있다.
따라서 식 (46)을 바탕으로 $\xi_{3}$, $\xi_{4}$의 ultimate bound(UB)를 구하면 다음과 같다.
따라서 roll에 관련 된 상태 $\xi_{3}$, $\xi_{4}$에 대하여 $\gamma_{2}$을 조절하면 시스템의 허용 한도 내에서 ultimate
bound가 임의로 조절 가능하다.
동일한 방법으로 pitch, yaw의 제어 입력을 정의하면 다음과 같다.
4. 실험결과
$\gamma_{j}$$(j=1,\: 2,\: 3,\: 4)$의 변화에 따른 제어기의 성능을 검증하기 위해 고도 및 자세에 대하여 실험을 각각 진행하였다.
실험에 사용 된 드론은 Parrot사의 AR.Drone 이며, 그림 2는 쿼드로터 실험 환경을 나타낸다.
그림 3. 쿼드로터 실험 환경
Fig. 3. Quadrotor experiment environment
$bold z$ position : 불확실한 무게 변화 및 외란에 대해 $\gamma_{1}$의 변화에 따른 고도 변화를 비교한다. 이때 질량 변화에
대한 성능 검증을 위해 호버링한 상태에서 쿼드로터 무게의 $18%$에 해당하는 $0.078$kg의 질량을 추가하였다. 그림 4는 $\gamma_{1}$의 변화에 따른 고도를 나타내며, 초기 이득은 $\beta_{1}^{*}=4$, $k_{1}^{*}=1$, $\epsilon_{1}^{*}=0.8$로
설정하였다. 추종하려는 고도는 $1$m이다.
그림 4. $\gamma_{1}$에 따른 고도 변화
Fig. 4. Behavior of altitude with respect to $\gamma_{1}$
그림 4의 그래프를 보면 $\gamma_{1}$이 작을수록 무게변화가 없는 4초에서 8초 사이에 추종하려는 높이와의 오차가 작은 모습을 보인다. 또한 무게
변화가 생긴 8초 이후를 보면 $\gamma_{1}$이 작아짐에 따라 무게 변화에 대해 추종 하려는 고도와의 오차가 작아지는 것을 볼 수 있다. 보다
정확한 비교를 위해 $\gamma_{1}$의 변화에 따른 고도의 오차에 대한 그래프는 그림 4, 제어 입력에 대한 그래프는 그림 5와 같다. 및 입력에 대한 그래프를 나타내면 그림 5와 같다. 이때 오차에 대해서는 호버링한 상태인 4초부터 12초 까지를 나타내었다.
그림 5 $\gamma_{1}$에 따른 추종 오차
Fig. 5 Tracking error of altitude with respect to $\gamma_{1}$
그림 6 $\gamma_{1}$에 따른 제어 입력
Fig. 6 Control input with respect to $\gamma_{1}$
Angle : $\gamma_{j}$$(j=2,\: 3,\: 4)$의 변화에 따른 자세 각의 변화를 보기 위해 초기 이득 $\beta_{j}=4$,
$k_{j}=1$, $\epsilon_{j}=0.8$ $(j=2,\: 3,\: 4)$으로 두어 실험하였다. 그림 7은 각각 roll, pitch, yaw에 대한 실험 그래프를 나타낸다.
그림 7-9를 보면 $\gamma_{j}$$(j=2,\: 3,\: 4)$가 감소함에 따라 자세 각에 대한 오차가 줄어드는 것이 보인다. 보다 자세한 비교를 위해
그림 10-15에서 roll, pitch, yaw의 추종 오차에 대한 그래프와 제어 입력을 비교 하였다.
또한 $\gamma_{j}$$(j=1,\: 2,\: 3,\: 4)$에 따른 변화를 정량적으로 비교하기 위해 정상상태에서의 평균 제곱 오차를 구한다.
평균 제곱 오차를 구하는 식은 다음과 같으며, $i=(z,\:\phi ,\:\theta ,\:\psi)$, $j=(1,\:3,\: 5,\: 7)$이다.
그림 7 $\gamma_{2}$에 따른 roll의 각도 변화
Fig. 7 Behavior of roll angle with respect to $\gamma_{2}$
그림 8 $\gamma_{3}$에 따른 pitch의 각도 변화
Fig. 8 Behavior of pitch angle with respect to $\gamma_{3}$
그림 9 $\gamma_{4}$에 따른 yaw의 각도 변화
Fig. 9 Behavior of yaw angle with respect to $\gamma_{4}$
정상상태에서의 평균 제곱 오차를 구하기 위해 4초부터 12초까지의 구간에 대해 구하였으며, $n$은 해당 구간의 데이터 개수를 의미한다. 위의 식을
통해 얻은 평균 제곱 오차의 결과는 아래의 표 2와 같다.
그림 10 $\gamma_{2}$에 따른 roll의 추종 오차
Fig. 10 Tracking error of roll with respect to $\gamma_{2}$
그림 11 $\gamma_{3}$에 따른 pitch의 추종 오차
Fig. 11 Tracking error of pitch with respect to $\gamma_{3}$
그림 12 $\gamma_{4}$에 따른 yaw의 추종 오차
Fig. 12 Tracking error of yaw with respect to $\gamma_{4}$
따라서 $\gamma_{i}$$(i=1,\: 2,\:3,\:4)$의 조절에 따른 성능변화의 추이는 $z$축 높이, roll, pitch, yaw의
오차가 각각 $58.7%$, $78.4%$, $85.4%$, $81.8%$ 감소하였다는 것을 알 수 있다.
그림 13 $\gamma_{2}$에 따른 제어 입력
Fig. 13 Control input with respect to $\gamma_{2}$
그림 14 $\gamma_{3}$에 따른 제어 입력
Fig. 14 Control input with respect to $\gamma_{3}$
그림 15 $\gamma_{4}$에 따른 제어 입력
Fig. 15 Control input with respect to $\gamma_{4}$