Consider a multi-output nonlinear system in the form of
$$\begin{align*}
\dot{x} &= f(x, u)\\
y &= h(x),
\end{align*}$$
where $x \in \mathbb{R}^n$ is the state, $u \in \mathbb{R}^m$ is the control input, $y \in \mathbb{R}^p$ is the output, $f$ and $h$ are smooth vector fields. The control input $u: \mathbb{R} \rightarrow \mathbb{R}^m$ is assumed to be an analytic time function. In particular, we will restrict our interest to the class of systems of the following form:
$\eqalignno{
\dot{x} &= A_ix_i+g_i(x_1,\ldots,x_i;u;y_{i+1},\ldots,y_p)\\
y &= C_ix_i & 1 \leq i \leq p,
}$
where $x=[x_1^{\mathrm T}, x_2^{\mathrm T}, \ldots, x_p^{\mathrm T}]^{\mathrm T} \in \mathbb{R}^n$, $x_i=[x_{i1}, x_{i2}, \ldots, x_{in}]^T \in \mathbb{R}^{n_i}$, $y=[y_1, \ldots, y_p]^{\mathrm T} \in \mathbb{R}^p$
$A_i = \begin{bmatrix}
0 & 1 & \cdots & 0 \\
\vdots & \vdots & \ddots &\vdots\\
0 & 0 & \cdots & 1\\
0 & 0 & \cdots & 0
\end{bmatrix} \in \mathbb{R}^{n_i \times n_i}, C_i = [1, 0, \ldots, 0] \in R^{1 \times n_i}$
$g_i = \begin{bmatrix}
g_{i1} (x_{[1,i-1]}x_{i1},u,y_{[i+1, p]} )\\
g_{i1} (x_{[1,i-1]}x_{i1},x_{i2},u,y_{[i+1, p]} )\\
\vdots\\
g_{i1} (x_{[1,i-1]}x_{[i1,in_i]},u,y_{[i+1, p]} )
\end{bmatrix}$
with $x_{[1,i-1]}=[x_1^{\mathrm T}, \ldots, x_{i-1}^{\mathrm T}]^{\mathrm T}$ and $y_{[i+1,p]}=[y_{i+1},\ldots,y_p]^{\mathrm T}$, and $g$ is a smooth vector field.
The system has the form
$$\begin{align*}
\dot{x}_1 &= x_2 + 0.01x_1u \\
\dot{x}_2 &= -x-1 + (1-x_1^2)+x_3u \\
\dot{x}_3 &= x_4 + 0.01x_2x_3 \exp(u) \\
\dot{x}_4 &= -x_3 + (1 - x_3^2)x_4 + u \\
y_1 &= x_1 \\
y_2 &= x_3
\end{align*}$$
with $u= 2 \sin {3t}$.