2.1. Definition of coordinate system

The study of kinematic analysis holds profound significance, serving as the foundational cornerstone of robotics. Given the intricate nature of the thrust-vectoring mechanism, the process of resolving the position solution for the vectored thruster represents a complex and crucial focal point. In order to comprehensively elucidate the dynamics governing the pose, encompassing both position and orientation, of the upper moving platform, it is imperative to establish coordinate systems within the thrust-vectoring mechanism to articulate relative motion. As illustrated in Figure 3, we have introduced two coordinate systems within the thrust-vectoring mechanism.

Figure 3.

Schematic of coordinate systems in the thrust-vectoring mechanism.

As depicted in Figure 3, the reference frame $O_{xyz}$ is affixed to the central point $O$ of the fixed base. In this frame, the $x$ -axis aligns itself parallel to the direction of $A_{3}A_{2}$ , while the $z$ -axis extends normal to the fixed base. Conversely, the moving reference frame $P_{ijk}$ is anchored to point $P$ , located on the lower surface of the moving platform. Within this frame, the $i$ -axis runs parallel to $C_{3}C_{2}$ , while the $k$ -axis is oriented perpendicular to the lower surface of the moving platform.

The pose state ${\boldsymbol{q}}$ of the moving platform within the parallel manipulator can be succinctly represented in relation to the fixed coordinate system. Given that the thrust-vectoring mechanism encompasses one translational degree-of-freedom and two rotational degrees-of-freedom, the generalized vector ${\boldsymbol{q}}$ emerges as a pivotal construct. Within this vector, each element corresponds to one of the six variables judiciously selected to describe both the position and orientation of the upper moving platform with precision. The essential motion parameters encompassed in defining the pose state are as follows:

(1) \begin{equation} \begin{array}{c} {\boldsymbol{q}}=\left[x_{P},y_{P},z_{P},\alpha, \beta, \gamma \right]^{\mathrm{T}}\end{array} \end{equation}

where $\alpha$ , $\beta$ , and $\gamma$ signify the $i$ -axis, $j$ -axis, and $k$ -axis rotation angles of the moving reference frame $P_{ijk}$ in relation to the fixed reference frame $O_{xyz}$ , respectively. $x_{P}$ , $y_{P}$ , and $z_{P}$ represent the displacements of point $P$ within the fixed reference frame. Hence, the pose state of the moving platform can be unequivocally ascertained through the six parameters delineated in Equation 1. When the moving coordinate system on the moving platform of the parallel mechanism undergoes rotation in relation to the fixed coordinate system, the rotation matrix of the moving coordinate system within the mechanism in the fixed coordinate system can be aptly expressed as [Reference Liu, Hu, Xu, Wang and Du2]:

(2) \begin{equation} {\boldsymbol{R}}=\left[\begin{array}{c@{\quad}c@{\quad}c} c\beta c\gamma & {s}{\alpha }s\beta c\gamma -c\alpha s\gamma & {c}{\alpha }s\beta c\gamma +s\alpha s\gamma \\[5pt] c\beta s\gamma & {s}{\alpha }s\beta s\gamma +c\alpha c\gamma & {c}{\alpha }s\beta s\gamma -s\alpha c\gamma \\[5pt] -s\beta & s\alpha c\beta & c\alpha c\beta \end{array}\right] \end{equation}

where $s(\!\cdot\!)$ represents the sine function, and $c(\!\cdot\!)$ represents the cosine function.

In cases where the origin points of the moving coordinate system and the fixed coordinate system within the parallel mechanism do not align, and their respective poses exhibit disparities, the position of point $S$ within the fixed coordinate system ${ }^{{\boldsymbol{O}}}{{\boldsymbol{S}}}{}$ can be aptly articulated by the moving coordinate ${ }^{{\boldsymbol{P}}}{{\boldsymbol{S}}}{}$ :

(3) \begin{equation} { }^{{\boldsymbol{O}}}{{\boldsymbol{S}}}{}={\boldsymbol{R}}\cdot { }^{{\boldsymbol{P}}}{{\boldsymbol{S}}}{}+{ }^{{\boldsymbol{O}}}{{\boldsymbol{P}}}{} \end{equation}

where ${ }^{{\boldsymbol{O}}}{{\boldsymbol{P}}}{}=[x_{P},y_{P},z_{P}]^{\mathrm{T}}$ denotes the position vector representing the origin point $P$ of the moving coordinate system within the fixed coordinate system.

2.2. Inverse position analysis of the designed thrust-vectoring mechanism

When the positional conditions and structural parameters of the moving platform are provided, determining the variation values in each kinematic pair under these positional constraints constitutes the process referred to as the inverse kinematics solution. As depicted in Figure 3, three kinematic branches interconnected between the fixed base and the moving platform form an equilateral triangle. The distance from the origin $O$ of the fixed coordinate system to the center of the rotating hub is designated as $r_{1}$ , representing the circumradius of the equilateral triangle formed by the three points $A_{1}$ , $A_{2}$ , $A_{3}$ . Similarly, the circumradius of the equilateral triangle formed by the vertices $C_{1}$ , $C_{2}$ , $C_{3}$ on the moving platform is denoted as $r_{2}$ .

It’s worth noting that the plane formed by the points $C_{1}$ , $C_{2}$ , $C_{3}$ shares the same orientation as the plane formed by the points $D_{1}$ , $D_{2}$ , $D_{3}$ . The key distinction lies in the origin of the dynamic coordinate system, which is established by the points $C_{1}$ , $C_{2}$ , $C_{3}$ and differs from the center point formed by the points $D_{1}$ , $D_{2}$ , $D_{3}$ solely in the $k$ -axis direction. Consequently, to streamline the kinematic analysis process, this paper exclusively scrutinizes the moving coordinate system formed by the points $C_{i}$ , from which the orientation of the entire moving platform can be deduced.

As illustrated in Figure 3, the angle between the rotation arm $A_{i}B_{i}$ and the plane formed by $A_{1}$ , $A_{2}$ , $A_{3}$ is defined as $\theta _{Ai}$ , signifying the rotation angle of the drive motor. Given the positional and orientational specifications of the moving platform within the thrust-vectoring mechanism, the kinematic inverse analysis furnishes the rotation angle of the servo motor. We further define the angle formed by the plane $O_{xy}$ and the connecting rod $B_{i}C_{i}$ as $\theta _{Bi}$ .

The structural design of the thrust-vectoring mechanism showcases an arrangement of three equidistant kinematic chains connecting the two platforms, with each chain plane orthogonal to the rotational axis of its respective revolute pair. Thanks to these structural and kinematic attributes, each revolute pair possesses a single degree-of-freedom, culminating in the emergence of three constraint equations at the central point of the thrust-vectoring mechanism. These constraint equations, pertinent to each joint within the kinematic branch chain, can be succinctly expressed in mathematical terms as follows:

(4) \begin{equation} \left\{\begin{array}{l} x_{1}=0\\[5pt] x_{2}=\sqrt{3}y_{2}\\[5pt] x_{3}=-\sqrt{3}y_{3} \end{array} \right.\end{equation}

Through the utilization of Equation 3, the position transformation relationship for point $C_{1}$ within the fixed coordinate system can be established as follows:

(5) \begin{equation} \begin{array}{c}{\,}^{O}\left[x_{C1},y_{C1},z_{C1}\right]^{\mathrm{T}}={\boldsymbol{R}}\cdot{\,}^{P}\left[0,-r_{2},0\right]^{T}+^{O}\left[x_{P},y_{P},z_{P}\right]^{\mathrm{T}}\end{array}\end{equation}

Consequently, the position of point $C_{1}$ within the fixed coordinate system can be derived as follows:

(6) \begin{equation}\begin{array}{c@{\quad}c} {}^{O}\left[x_{C1},y_{C1},z_{C1}\right]^{\mathrm{T}}= & { }^{O}\left[-r_{2}\left({s}{\alpha }s\beta c\gamma -c\alpha s\gamma \right),-r_{2}\left({s}{\alpha }s\beta s\gamma +c\alpha c\gamma \right),-r_{2}s\alpha c\beta \right]^{\mathrm{T}}\\[5pt] & +^{O}\left[x_{P},y_{P},z_{P}\right]^{\mathrm{T}}\end{array} \end{equation}

Given the structural characteristics of the thrust-vectoring mechanism and as per the constrained Equation 4, the $x$ -axis component of point $C_{1}$ in the fixed coordinate system consistently remains zero. Therefore, the determination of the $x$ -axis direction of point $C_{1}$ can be succinctly expressed as:

(7) \begin{equation} \begin{array}{c} { }^{O}{x_{C1}}{}=-r_{2}\left({s}{\alpha }s\beta c\gamma -c\alpha s\gamma \right)+^{O}x_{P}=0 \end{array} \end{equation}

Consequently, the coordinates of the origin point $P$ within the moving coordinate system ${ }^{O}{x_{p}}{}$ can be expressed as:

(8) \begin{equation} \begin{array}{c} { }^{O}{x_{p}}{}=r_{2}\left({s}{\alpha }s\beta c\gamma -c\alpha s\gamma \right) \end{array} \end{equation}

In accordance with the analysis detailed in Equation 8, the alteration in the position of the central point $P$ along the $x$ -axis direction is associated with the orientation angles $\alpha$ , $\beta$ , $\gamma$ and the radius $r_{2}$ .

Similarly, the position coordinates of points $C_{2}$ and $C_{3}$ within the fixed coordinate system can be formulated as follows:

(9) \begin{equation} \left\{\begin{array}{l} { }^{O}{{\boldsymbol{C}}}{}_{2}={\boldsymbol{R}}\cdot { }^{P}{{\boldsymbol{C}}}{}_{2}+{ }^{{\boldsymbol{O}}}{{\boldsymbol{P}}}{}\\[5pt] { }^{O}{{\boldsymbol{C}}}{}_{3}={\boldsymbol{R}}\cdot { }^{P}{{\boldsymbol{C}}}{}_{2}+{ }^{{\boldsymbol{O}}}{{\boldsymbol{P}}}{} \end{array} \right.\end{equation}

By employing the coordinate transformation equation as detailed in Equation 9, and taking into account the constraint conditions for points $C_{2}$ and $C_{3}$ as expressed in Equation 4, we can deduce the connection between the alteration in the position of the origin point $P$ along the $y$ -axis direction and the orientation angles $\alpha$ , $\beta$ , $\gamma$ as well as the radius $r_{2}$ as follows:

(10) \begin{equation} \begin{array}{c} { }^{O}{y}{}_{p}=\frac{r_{2}}{2}\left(c\beta c\gamma -{s}{\alpha }s\beta s\gamma -c\alpha c\gamma \right) \end{array} \end{equation}

The orientation angle $\gamma$ arising from the rotation of the moving platform around the $z$ -axis of the fixed coordinate system can be expressed as:

(11) \begin{equation} \begin{array}{c} \gamma =\arctan \left[\dfrac{s\alpha s\beta }{c\beta +c\alpha }\right] \end{array} \end{equation}

With the known parameters ( $z_{P}$ , $\alpha$ , $\beta$ ) of the moving platform, the pose parameters ( $x_{P}$ , $y_{P}$ , $\gamma$ ) can be determined by utilizing Equation 8 through $11$ . The inverse kinematics of the thrust-vectoring mechanism is examined based on these pose parameters of the moving platform. As shown in Figure 3, set $| \overrightarrow{A_{{i}}B_{{i}}}| =l_{ab}$ , $| \overrightarrow{B_{{i}}C_{{i}}}| =l_{bc}$ , and the vector $OC_{i}$ can be expressed as:

(12) \begin{equation} \begin{array}{c} \overrightarrow {OC_{{i}}}=\overrightarrow {OP}+\overrightarrow {PC_{{i}}}=\overrightarrow {OA_{{i}}}+\overrightarrow {AB_{{i}}}+\overrightarrow {BC_{{i}}} \end{array} \end{equation}

On the branch of $A_{1}-B_{1}-C_{1}$ , Equation 12 can be expressed as:

(13) \begin{equation} \begin{array}{c} \left[\begin{array}{c} x_{P}\\[5pt] y_{P}\\[5pt] z_{P} \end{array}\right]+{\boldsymbol{R}}\left[\begin{array}{c} 0\\[5pt] -r_{2}\\[5pt] 0 \end{array}\right]=\left[\begin{array}{c} x_{c1}\\[5pt] y_{c1}\\[5pt] z_{c1} \end{array}\right]=\left[\begin{array}{c} 0\\[5pt] -r_{1}\\[5pt] 0 \end{array}\right]+\left[\begin{array}{c} 0\\[5pt] l_{ab}{\cos } \theta _{A1}\\[5pt] l_{ab}{\sin } \theta _{A1} \end{array}\right]+\left[\begin{array}{c} 0\\[5pt] l_{bc}{\cos } \theta _{B1}\\[5pt] l_{bc}{\sin } \theta _{B1} \end{array}\right] \end{array} \end{equation}

Consider the constraint $x_{1}=0$ , the first equation can be eliminated, and the system of equations can be simplified as:

(14a) \begin{equation} \begin{array}{c} \begin{cases} l_{ab}{\cos } \theta _{A1}+l_{bc}{\cos } \theta _{B1}=b_{11}\\[5pt] l_{ab}{\sin } \theta _{A1}+l_{bc}{\sin } \theta _{B1}=b_{21} \end{cases} \end{array} \end{equation}

with:

(14b) \begin{equation} \begin{array}{c} \begin{cases} b_{11}=y_{P}-r_{2}\left({s}{\alpha }s\beta s\gamma +c\alpha c\gamma \right)\\[5pt] b_{21}=z_{P}-r_{2}\left(s\alpha c\beta \right) \end{cases} \end{array} \end{equation}

Similarly, on the branches of $A_{2}-B_{2}-C_{2}$ and $A_{3}-B_{3}-C_{3}$ , the following equations are found:

(15a) \begin{equation} \begin{array}{c} \begin{cases} l_{ab}{\cos } \theta _{A2}+l_{bc}{\cos } \theta _{B2}=b_{12}\\[5pt] l_{ab}{\sin } \theta _{A2}+l_{bc}{\sin } \theta _{B1}=b_{22} \end{cases} \end{array} \end{equation}

(16a) \begin{equation} \begin{array}{c} \begin{cases} l_{ab}{\cos } \theta _{A3}+l_{bc}{\cos } \theta _{B3}=b_{13}\\[5pt] l_{ab}{\sin } \theta _{A3}+l_{bc}{\sin } \theta _{B3}=b_{23} \end{cases} \end{array} \end{equation}

with:

(15b) \begin{equation} \begin{array}{c} \begin{cases} b_{12}=r_{1}-2\left[y_{p}+\dfrac{\sqrt{3}}{2}r_{2}\left(c\beta s\gamma \right)+\dfrac{1}{2}r_{2}\left({s}{\alpha }s\beta s\gamma +c\alpha c\gamma \right)\right]\\[10pt] b_{22}=z_{P}+\dfrac{\sqrt{3}}{2}r_{2}\left(-s\beta \right)+\dfrac{1}{2}r_{2}\left(s\alpha c\beta \right) \end{cases} \end{array} \end{equation}

(16b) \begin{equation} \begin{array}{c} \begin{cases} b_{13}=r_{1}-2\left[y_{p}-\dfrac{\sqrt{3}}{2}r_{2}\left(c\beta s\gamma \right)+\dfrac{1}{2}r_{2}\left({s}{\alpha }s\beta s\gamma +c\alpha c\gamma \right)\right]\\[10pt] b_{23}=z_{P}-\dfrac{\sqrt{3}}{2}r_{2}\left(-s\beta \right)+\dfrac{1}{2}r_{2}\left(s\alpha c\beta \right) \end{cases} \end{array} \end{equation}

Figure 4.

Rotation angle of moving platform $\alpha = {15}^{\circ}$ (0.2617 rad).

Based on the kinematic analysis of the thrust-vectoring mechanism discussed above, the rotation angles $\theta _{Ai}$ of the three driving servo motors can be computed. However, owing to the properties of trigonometric functions, the calculated rotation angles $\theta _{Ai}$ and $\theta _{Bi}$ are not unique. For example, given the pose of the moving platform with $\alpha =15^{^{\circ}}$ , $\beta =0^{^{\circ}}$ , $z_{P} =105\text{ mm}$ , there are two different methods for actuating the servo motors, both of which can satisfy the pose requirements of the moving platform. In Figure 4(a), the rotation angle $\theta _{A1}$ at servo motor $A_{1}$ is $143.7678^{^{\circ}}$ , and the rotation angle $\theta _{B1}$ is $66.5273^{^{\circ}}$ ; as shown in Figure 4(b), the corresponding angles are $\theta _{A1}=17.3618^{^{\circ}}$ , $\theta _{B1}=94.6023^{^{\circ}}$ . To account for the practical constraints and spatial limitations of the thrust-vectoring mechanism, in this design, the rotation angle $\theta _{Ai}$ has been selected to be greater than $90^{^{\circ}}$ , that is $\theta _{Ai}\gt 90^{^{\circ}} (i=1,2,3)$ , thereby circumventing issues such as non-unique rotation angles and abrupt changes. At this point, the calculation results for angles $\theta _{Ai}$ and $\theta _{Ai}$ are given by Equation 17.

(17a) \begin{equation} \begin{array}{c} \begin{cases} \theta _{A1}=\pi -{\arcsin }\dfrac{{b_{11}}^{2}+{b_{21}}^{2}+{l_{ab}}^{2}-{l_{bc}}^{2}}{2l_{ab}\sqrt{{b_{11}}^{2}+{b_{21}}^{2}}}-{\arctan }\dfrac{b_{11}}{b_{21}}\\[10pt] \theta _{B1}={\arcsin }\dfrac{{b_{11}}^{2}+{b_{21}}^{2}+{l_{bc}}^{2}-{l_{ab}}^{2}}{2l_{ab}\sqrt{{b_{11}}^{2}+{b_{21}}^{2}}}-{\arctan }\dfrac{b_{11}}{b_{21}} \end{cases} \end{array} \end{equation}

(17b) \begin{equation} \begin{array}{c} \begin{cases} \theta _{A2}=\pi -{\arcsin }\dfrac{{b_{12}}^{2}+{b_{22}}^{2}+{l_{ab}}^{2}-{l_{bc}}^{2}}{2l_{ab}\sqrt{{b_{12}}^{2}+{b_{22}}^{2}}}-{\arctan }\dfrac{b_{12}}{b_{22}}\\[10pt] \theta _{B2}={\arcsin }\dfrac{{b_{12}}^{2}+{b_{22}}^{2}+{l_{bc}}^{2}-{l_{ab}}^{2}}{2l_{ab}\sqrt{{b_{12}}^{2}+{b_{22}}^{2}}}-{\arctan }\dfrac{b_{12}}{b_{22}} \end{cases} \end{array} \end{equation}

(17c) \begin{equation} \begin{array}{c} \begin{cases} \theta _{A3}=\pi -{\arcsin }\dfrac{{b_{13}}^{2}+{b_{23}}^{2}+{l_{ab}}^{2}-{l_{bc}}^{2}}{2l_{ab}\sqrt{{b_{13}}^{2}+{b_{23}}^{2}}}-{\arctan }\dfrac{b_{13}}{b_{23}}\\[10pt] \theta _{B3}={\arcsin }\dfrac{{b_{13}}^{2}+{b_{23}}^{2}+{l_{bc}}^{2}-{l_{ab}}^{2}}{2l_{ab}\sqrt{{b_{13}}^{2}+{b_{23}}^{2}}}-{\arctan }\dfrac{b_{13}}{b_{23}} \end{cases} \end{array} \end{equation}

2.3. Forward position analysis of the designed thrust-vectoring mechanism

When the rotation angles of each joint are specified, the process of determining the pose of the moving platform is referred to as the forward kinematics of the thrust-vectoring mechanism. Research on forward kinematics analysis is valuable for subsequent investigations into control and structural optimization. In this paper, the forward kinematics process of the thrust-vectoring mechanism is employed to compute the relative pose parameters of the moving platform using the known rotation angles $\theta _{A1}$ , $\theta _{A2}$ , and $\theta _{A3}$ of the servo motors.

From Eq. (15), the position of point $C_{i}$ can be described as:

(18) \begin{equation} \left\{\begin{array}{l} {\boldsymbol{C}}_{1}={ }^{O}\left[0,-\left(r_{1}-l_{ab}c\left(\theta _{A1}\right)+l_{bc}c\left(\theta _{B1}\right)\right),l_{ab}s\left(\theta _{A1}\right)+l_{bc}s\left(\theta _{B1}\right)\right]^{T}\\[5pt] {\boldsymbol{C}}_{2}={ }^{O}\left[\sqrt{3}\left(r_{1}-l_{ab}c\left(\theta _{A2}\right)+l_{bc}c\left(\theta _{B2}\right)\right)/2,\left(r_{ab}c\left(\theta _{A2}\right)+l_{bc}c\left(\theta _{B2}\right)\right)/2,l_{ab}s\left(\theta _{A2}\right)+l_{bc}s\left(\theta _{B2}\right)\right]^{T}\\[5pt] {\boldsymbol{C}}_{3}={ }^{O}\left[-\sqrt{3}\left(r_{1}-l_{ab}c\left(\theta _{A3}\right)+l_{bc}c\left(\theta _{B3}\right)\right)/2,\left(r_{1}-l_{ab}c\left(\theta _{A3}\right)+l_{bc}c\left(\theta _{B3}\right)\right)_{1}/2,l_{ab}s\left(\theta _{A3}\right)+l_{bc}s\left(\theta _{B3}\right)\right]^{T} \end{array} \right.\end{equation}

As the moving platform takes the form of an equilateral triangle with a circumradius of $r_{2}$ , the connection between the three points $C_{1}$ , $C_{2}$ , and $C_{3}$ on the moving platform can be expressed as $|C_{1}C_{2}|=|C_{1}C_{3}|=|C_{2}C_{3}|=\sqrt{3}r_{2}$ . Meanwhile, the position coordinates of point $P$ in the fixed coordinate system can be articulated as:

(19) \begin{equation} \begin{array}{c} \left.\begin{cases} { }^{O}{x_{P}}{}={ }^{O}{x_{C1}}{}+{ }^{O}{x_{C2}}{}+{ }^{O}{x_{C3}}{}\\[5pt] { }^{O}{y_{P}}{}={ }^{O}{y_{C1}}{}+{ }^{O}{y_{C2}}{}+{ }^{O}{y_{C3}}{}\\[5pt] { }^{O}{z_{P}}{}={ }^{O}{z_{C1}}{}+{ }^{O}{z_{C2}}{}+{ }^{O}{z_{C3}}{} \end{cases} \right. \end{array} \end{equation}

Simultaneously, the normal vector ${\boldsymbol{n}}$ of the moving platform can be described as:

(20) \begin{align} {\boldsymbol{n}}=\overrightarrow {C_{1}C_{2}}\times \overrightarrow {C_{1}C_{3}} & =\left| \begin{array}{c@{\quad}c@{\quad}c} i & j & k\\[5pt] x_{2}-x_{1} & y_{2}-y_{1} & z_{2}-z_{1}\\[5pt] x_{3}-x_{1} & y_{3}-y_{1} & z_{3}-z_{1} \end{array}\right| \nonumber\\[5pt] &= ai+bj+ck=\left(a,b,c\right) \end{align}

Hence, the alteration in the attitude angle $(\alpha, \beta, \gamma )$ of the moving platform can be defined as:

(21) \begin{equation} \left\{\begin{array}{l} \alpha =\mathrm{c}^{-1}\left(\dfrac{{\boldsymbol{n}}_{1}\cdot {\boldsymbol{m}}_{1}}{\left| {\boldsymbol{n}}_{1}\right| \cdot \left| {\boldsymbol{m}}_{1}\right| }\right)\\[10pt] \beta =\mathrm{c}^{-1}\left(\dfrac{{\boldsymbol{n}}_{2}\cdot {\boldsymbol{m}}_{2}}{\left| {\boldsymbol{n}}_{1}\right| \cdot \left| {\boldsymbol{m}}_{1}\right| }\right)\\[10pt] \gamma =\mathrm{c}^{-1}\left(\dfrac{{\boldsymbol{n}}_{3}\cdot {\boldsymbol{m}}_{3}}{\left| {\boldsymbol{n}}_{1}\right| \cdot \left| {\boldsymbol{m}}_{1}\right| }\right) \end{array} \right.\end{equation}

where ${\boldsymbol{n}}_{1}$ , ${\boldsymbol{n}}_{2}$ , and ${\boldsymbol{n}}_{3}$ denote the direction vectors of the $i$ , $j$ , and $k$ axes of the moving coordinate system, respectively. Meanwhile, ${\boldsymbol{m}}_{1}$ , ${\boldsymbol{m}}_{2}$ , and ${\boldsymbol{m}}_{3}$ represent the direction vectors of the $x$ , $y$ , and $z$ axes of the fixed coordinate system. Based on the Equations 18 through 21, the position vector $[x_{P},y_{P},z_{P}]^{\mathrm{T}}$ and the attitude vector $[\alpha, \beta, \gamma ]^{\mathrm{T}}$ of the moving platform in the fixed coordinate system can be determined through mathematical calculations.

4.1. PINN-MPC for the vectored thruster

According to Equation 17, the rotation angle of revolute joint $\theta _{Ai}$ can be concisely represented as $\theta _{Ai}=f_{i}(\alpha, \beta, \gamma (\alpha, \beta ))$ . Subsequently, the derivative of $\theta _{Ai}$ can be derived as:

(22) \begin{equation} \dot{\theta }_{Ai}=\left(\frac{\partial f_{i}}{\partial \alpha }+\frac{\partial f_{i}}{\partial \gamma }\frac{\partial \gamma }{\partial \alpha }\right)\dot{\alpha }+\left(\frac{\partial f_{i}}{\partial \beta }+\frac{\partial f_{i}}{\partial \gamma }\frac{\partial \gamma }{\partial \beta }\right)\dot{\beta } \end{equation}

To facilitate analysis and control design, let ${\boldsymbol{x}}=[\alpha, \beta ]^{\mathrm{T}}$ represent the state vector encompassing the attitude angles, and ${\boldsymbol{u}}=[\dot{\theta }_{A1},\dot{\theta }_{A2},\dot{\theta }_{A3}]^{\mathrm{T}}$ represent the input vector consisting of the rates of change of the revolute joint angles. The state equation can then be expressed in a compact form as:

(23) \begin{equation} \begin{array}{c} {\boldsymbol{u}}=\left[\begin{array}{c} \dot{\theta }_{A1}\\[10pt] \dot{\theta }_{A2}\\[10pt] \dot{\theta }_{A3} \end{array}\right]=\left[\begin{array}{c@{\quad}c} \left(\dfrac{\partial f_{1}}{\partial \alpha }+\dfrac{\partial f_{1}}{\partial \gamma }\dfrac{\partial \gamma }{\partial \alpha }\right) & \left(\dfrac{\partial f_{1}}{\partial \beta }+\dfrac{\partial f_{1}}{\partial \gamma }\dfrac{\partial \gamma }{\partial \beta }\right)\\[10pt] \left(\dfrac{\partial f_{2}}{\partial \alpha }+\dfrac{\partial f_{2}}{\partial \gamma }\dfrac{\partial \gamma }{\partial \alpha }\right) & \left(\dfrac{\partial f_{2}}{\partial \beta }+\dfrac{\partial f_{2}}{\partial \gamma }\dfrac{\partial \gamma }{\partial \beta }\right)\\[10pt] \left(\dfrac{\partial f_{3}}{\partial \alpha }+\dfrac{\partial f_{3}}{\partial \gamma }\dfrac{\partial \gamma }{\partial \alpha }\right) & \left(\dfrac{\partial f_{3}}{\partial \beta }+\dfrac{\partial f_{3}}{\partial \gamma }\dfrac{\partial \gamma }{\partial \beta }\right) \end{array}\right]\left[\begin{array}{c} \dot{\alpha }\\[10pt] \dot{\beta } \end{array}\right]=\boldsymbol{\varphi }\left({\boldsymbol{x}}\right){\boldsymbol{x}} \end{array} \end{equation}

which can also be written as:

(24) \begin{equation} \begin{array}{c} \dot{{\boldsymbol{x}}}=\left(\boldsymbol{\varphi }^{\mathrm{T}}\left({\boldsymbol{x}}\right)\boldsymbol{\varphi }\left({\boldsymbol{x}}\right)\right)^{-1}\boldsymbol{\varphi }^{\mathrm{T}}\left({\boldsymbol{x}}\right){\boldsymbol{u}}={\boldsymbol{B}}\left({\boldsymbol{x}}\right){\boldsymbol{u}} \end{array} \end{equation}

PINNs uniquely integrate empirical data with the governing physical laws, expressed as Partial Differential Equations (PDEs), into the training process. In order for the PINN to accurately fit the state Equation 24, it is imperative to devise a well-suited loss function for the optimized network. This loss function, denoted as $\mathcal{L}$ , serves as a composite metric that combines both supervised and unsupervised losses. The supervised loss $\mathcal{L}_{data}$ , captures the discrepancy between the network’s predictions and actual data measurements, while the unsupervised loss $\mathcal{L}_{phys}$ , enforces the satisfaction of the governing PDEs. And the loss function $\mathcal{L}$ can be expressed as:

(25) \begin{equation} \begin{array}{c} \mathcal{L}=w_{data}\mathcal{L}_{data}+w_{phys}\mathcal{L}_{phys} \end{array} \end{equation}

where $w_{data}$ and $w_{phys}$ are weighting factors that judiciously balance the contributions of the data-driven loss $\mathcal{L}_{data}$ and the physics-informed loss $\mathcal{L}_{phys}$ , respectively. These weights are tuned to ensure that both empirical data and physical laws are adequately represented in the training process.

The supervised loss $\mathcal{L}_{data}$ typically measures the difference between the network’s output and the corresponding target values from the available data. This can be achieved using various metrics such as mean squared error or mean absolute error. On the other hand, the unsupervised loss $\mathcal{L}_{phys}$ ensures that the network’s predictions satisfy the governing PDEs. This is achieved by incorporating the PDEs into the loss function using techniques like automatic differentiation. The PDEs define relationships between different variables and their derivatives, and the unsupervised loss penalizes any violations of these relationships by the network’s predictions. The specific form of the PDE violation measure depends on the nature of the PDEs and how they are encoded into the loss function. It typically involves computing the residuals of the PDEs using the network’s predictions and then summing or integrating these residuals over the relevant domain. Therefore, the losses of $\mathcal{L}_{data}$ and $\mathcal{L}_{phys}$ can be defined as:

(26) \begin{equation} \mathcal{L}_{data}=\frac{1}{N_{data}}\sum _{i=1}^{N_{data}}\left\| {\boldsymbol{x}}-\hat{{\boldsymbol{x}}}\right\| ^{2} \end{equation}

(27) \begin{equation} \mathcal{L}_{phys}=\frac{1}{N_{phys}}\sum _{i=1}^{N_{phys}}\left\| {\boldsymbol{B}}\left({\boldsymbol{x}}\right){\boldsymbol{u}}-\dot{\hat{{\boldsymbol{x}}}}\right\| ^{2} \end{equation}

where $N_{data}$ and $N_{phys}$ respectively represent the number of data points applied to $\mathcal{L}_{data}$ and $\mathcal{L}_{phys}$ , and $\hat{{\boldsymbol{x}}}$ and ${\boldsymbol{x}}$ represent network predicted values and true values, respectively. By minimizing the composite loss function $\mathcal{L}$ , the PINN is trained to not only fit the available data but also adhere to the governing physical laws expressed by the PDEs. This ensures that the network’s predictions are both accurate and physically meaningful, making it suitable for use in systems like vector thruster controllers.

A precisely defined cost function is pivotal for evaluating the performance of the control system in MPC. The cost function, commonly denoted by $J$ , aims to simultaneously achieve two primary objectives: minimizing deviations from the reference trajectory and optimizing the utilization of control inputs within the prediction horizon. The mathematical expression representing this cost function can be refined and polished as follows:

(28) \begin{equation} \begin{array}{c} \min J=\int _{0}^{Nt}\left(\left\| \hat{{\boldsymbol{x}}}\left(s\right)-{\boldsymbol{x}}_{ref}\right\| _{Q}^{2}+\left\| \hat{{\boldsymbol{u}}}\right\| _{R}^{2}\right)\text{ds } \end{array} \end{equation}

where $N$ is the length of the predicted time domain, ${\boldsymbol{x}}_{ref}$ represents the target trajectory, $Q$ and $R$ are weight matrices judiciously employed to balance state bias and the relative importance of control inputs. Meanwhile, $\hat{{\boldsymbol{x}}}$ in the optimization equation is obtained through the trained PINN model. The design of this cost function considers both the deviation between the system state and the reference trajectory, as well as the utilization of control inputs. By minimizing the cost function $J$ , the MPC controller can find an optimal sequence of control inputs that enable the system to operate along the desired trajectory while minimizing the consumption of control inputs. This optimization approach helps improve the performance and efficiency of the control system.

The control constraints of the vector thruster principally consist of two parts, the state ${\boldsymbol{x}}$ and input ${\boldsymbol{u}}$ , which can be mathematically expressed as follows:

(29a) \begin{equation} \begin{array}{c} {\boldsymbol{u}}_{max}=-{\boldsymbol{u}}_{min}=\left[5^{^{\circ}}/s,5^{^{\circ}}/s,5^{^{\circ}}/s\right]^{\mathrm{T}} \end{array} \end{equation}

(29b) \begin{equation} \begin{array}{c} {\boldsymbol{x}}_{max}=-{\boldsymbol{x}}_{min}=\left[15^{^{\circ}},15^{^{\circ}}\right]^{\mathrm{T}} \end{array} \end{equation}

Equation 29a specifies the maximum and minimum permissible values for each component of the input vector u, ensuring that the control inputs do not exceed a certain angular velocity limit. Equation 29b defines the maximum and minimum allowable values for each component of the state vector x, limiting the range of possible attitude angles for the moving platform. By adhering to these constraints, the control system ensures that the vector thruster operates within safe and effective parameters, preventing excessive or unsafe maneuvers.

4.2. Simulation analysis

For the training of this network, the loss function $\mathcal{L}$ has been calibrated using $N_{data}=1000$ data points for the data-driven loss $\mathcal{L}_{data}$ and $N_{phys}=2000$ . for the physics-informed loss $\mathcal{L}_{phys}$ . The weights assigned to these components are $w_{data}=1$ and $w_{phys}=0.1$ , respectively, ensuring a judicious balance between data and physics in the loss computation. This careful balancing of the loss function terms is crucial, as it allows the neural network to effectively learn from the available data while adhering to the underlying physical principles governing the system dynamics. By properly tuning these weights, the model can achieve an optimal tradeoff between fitting the empirical observations and satisfying the theoretical equations, resulting in enhanced accuracy and robustness.

The controller operates with a sampling period of $t=0.1$ second, and a prediction horizon of $N=10$ steps is employed to optimize the control decisions. The weighting matrices in the MPC algorithm have been calibrated as $Q=\text{diag}(200,200)$ and $R=\text{ diag}(0.1,0.1,0.1)$ , reflecting the relative importance of minimizing deviations from the target trajectory and controlling the magnitude of control inputs. The judicious selection of these weight matrices is essential for achieving the desired control performance, as it determines the relative emphasis placed on tracking accuracy and control effort minimization. By carefully balancing these competing objectives, the controller can navigate the tradeoff between precise trajectory following and efficient resource utilization, thereby maximizing overall system effectiveness. For the convenience of calculation, the Pytorch toolkit is utilized for derivative calculation.

To assess the performance of this control method, a simulation experiment has been designed. We opt for the following trajectory:

(30) \begin{equation} {\boldsymbol{x}}_{ref}= \left\{\begin{array}{l}\alpha =15^{^{\circ}}\cos \left(\dfrac{\pi }{50t}\right)\\[10pt] \beta =15^{^{\circ}}\sin \left(\dfrac{\pi }{50t}\right) \end{array}\right. \end{equation}

As depicted in Figure 12, subgraph (a) presents the trajectory tracking outcomes in the $\alpha$ and $\beta$ planes, offering a clear visualization of the platform’s movement. Subgraph (b) exhibits the integral value of the input ${\boldsymbol{u}}$ , which corresponds to the rotation angle of the rotary joint, providing insights into the dynamic behavior of the system. Additionally, subgraphs (c) and (d) showcase the tracking results and tracking errors in the $\alpha$ and $\beta$ directions, respectively, enabling a comprehensive evaluation of the tracking performance. During the initial tracking phase, which can be seen in (a) and (c), due to the mismatch between the initial attitude and the desired target value’s, the attitude angle must promptly converge toward the trajectory of the target value’s evolution and subsequently achieve a stable tracking outcome. This initial transient response is a natural consequence of the system’s dynamics and highlights the controller’s ability to rapidly compensate for deviations and synchronize with the desired trajectory.

Figure 12.

The movement state change of the thrust-vectoring mechanism.

As depicted in (d), after 20 s, the tracking task is essentially in a stable stage. During this stabilized phase, the maximum tracking error for $\alpha$ is $0.103^{^{\circ}}$ , with an average error of $0.064^{^{\circ}}$ . Analogously, the maximum error for $\beta$ is $0.107^{^{\circ}}$ , averaging at $0.058^{^{\circ}}$ . These outcomes demonstrate the precision of the tracking performance and a satisfactory convergence toward the desired trajectory, indicating effective control and accuracy in achieving the tracking objectives. The low magnitude of the tracking errors, both in terms of maximum deviations and average values, underscores the controller’s capability to maintain precise and consistent trajectory following once the initial transient has subsided.

To enhance the understanding and analysis of the control effectiveness of this control method on vector thrusters, a motion trajectory is delineated for the moving platform. This trajectory, as detailed in Equation 31, represents an 8-shaped path, allowing for a comprehensive assessment of the control method’s capabilities. By studying this trajectory, we can gain insights into how the control method performs in maneuvering the platform through complex and challenging patterns. The inclusion of such intricate trajectories in the simulation is crucial, as it enables a thorough evaluation of the controller’s robustness and adaptability under various operational scenarios, which is essential for real-world applications where unpredictable and dynamic conditions may be encountered.

(31) \begin{equation} {\boldsymbol{x}}_{ref}=\left\{\begin{array}{l} \alpha =-15^{^{\circ}}\sin \left(\frac{2\pi }{50t}\right)\\[5pt] \beta =15^{^{\circ}}\cos \left(\frac{\pi }{50t}\right) \end{array} \right.\end{equation}

When examining the figure-8 trajectory tracking results in Figure 13, the controller demonstrates outstanding performance. Despite the increased complexity of this motion pattern, the controller maintains exceptional tracking precision during the stable phase (after 20 s). For the $\alpha$ angle component, the maximum tracking error is $0.252^{\circ}$ , with an average error of $0.128^{\circ}$ ; for the $\beta$ angle component, the maximum error is $0.161^{\circ}$ , with an average error of $0.065^{\circ}$ . These data reflect the control method’s ability to achieve excellent tracking performance even when handling intricate trajectories, ultimately allowing the moving platform’s attitude to converge accurately to the desired values.

Figure 13.

The movement state change of the thrust-vectoring mechanism.

Notably, although the tracking errors are slightly elevated compared to the simpler sinusoidal trajectory, their magnitudes remain within an acceptable and reasonable range. This situation confirms the controller’s remarkable adaptability, as its performance only degrades moderately when faced with more challenging motion patterns, rather than suffering severe deterioration. The controller’s capability to flexibly handle various complex trajectories fully demonstrates its robust and reliable nature.

Overall, through the simulated figure-8 trajectory trial, we clearly observe the control strategy’s outstanding performance. Even under stringent conditions, the controller can rapidly and accurately track the target trajectory, achieving the desired attitude control effect. This lays a solid foundation for the widespread application of this method in practical engineering scenarios.

Figure 14.

The movement state change of the thrust-vectoring mechanism.

When the moving platform of the parallel mechanism tracks the infinite motion trajectory defined by Equation 32, the control results obtained are presented in Figure 14. This scenario introduces an additional layer of complexity, as the trajectory’s non-periodic nature poses a more demanding challenge for the controller to maintain accurate tracking over an extended duration.

(32) \begin{equation} {\boldsymbol{x}}_{ref}=\left\{\begin{array}{l} \alpha =15^{^{\circ}}\cos \left(\dfrac{\pi }{50t}\right)\\[10pt] \beta =15^{^{\circ}}\sin \left(\dfrac{2\pi }{50t}\right) \end{array}\right. \end{equation}

In Figure 14, the tracking precision is rigorously evaluated after the initial 20-second period. Specifically, the maximum tracking error for the α component is measured to be a modest $0.117^{^{\circ}}$ , with an average error of only $0.063^{^{\circ}}$ . In contrast, the $\beta$ component exhibits a slightly higher maximum error of $0.202^{^{\circ}}$ , averaging at $0.127^{^{\circ}}$ . Comparing the results presented in Figure 14 with those in Figure 13, it becomes evident that the tracking error tends to be more pronounced for angle components undergoing more rapid changes, as opposed to those experiencing slower variations. This behavior is expected and aligns with theoretical predictions, as abrupt changes in the desired trajectory pose a greater challenge for the controller to respond and compensate promptly.

Despite the added complexity of the non-periodic trajectory, it is noteworthy that the tracking errors remain well within an acceptable range, effectively demonstrating the proposed method’s capability in accurately tracking the target trajectory of the vector propulsion device, even for complex, aperiodic motions. This underscores the practical feasibility and reliability of the control strategy, rendering it a viable option for real-world applications where the system may encounter diverse and unpredictable motion profiles. The consistent and robust performance across various trajectory complexities highlights the controller’s versatility and resilience, essential qualities for practical deployment in dynamic and demanding environments.

The controller’s ability to maintain precise tracking, even when faced with the compounded challenges of non-periodicity and rapid trajectory changes, is a testament to the efficacy and sophistication of the proposed control approach. This remarkable performance not only validates the theoretical foundations but also instills confidence in the method’s potential for practical implementation, where the ability to handle diverse and unpredictable conditions is paramount.