# Error Analysis of Method for Calculation of Non-Contact Impact on Space Debris from Ion Thruster

A.P. Alpatov, A.A. Fokov, S.V. Khoroshylov, A.P. Savchuk

► To cite this version:

A.P. Alpatov, A.A. Fokov, S.V. Khoroshylov, A.P. Savchuk. Error Analysis of Method for Calculation of Non-Contact Impact on Space Debris from Ion Thruster. Mechanics, Materials Science & Engineering Journal, 2016, 10.13140/RG.2.1.3986.1361. hal-01325253

HAL Id: hal-01325253

https://hal.science/hal-01325253

Submitted on 2 Jun 2016

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Error Analysis of Method for Calculation of Non-Contact Impact on Space Debris from Ion Thruster

Alpatov A.P.¹, Fokov A.A.¹, Khoroshylov S.V.¹, Savchuk A.P.¹

1 – Institute of Technical Mechanics of National Academy of Sciences of Ukraine and State Space Agency of Ukraine, Dnipropetrovs'k, Ukraine

<figure><img src="image_2.png" /></figure>

DOI 10.13140/RG.2.1.3986.1361

**Keywords:** space debris removal, ion beam shepherd technology, spacecraft - space debris object system, contour of the central projection, simplified calculation of the impact, error analysis, simulation of the relative motion.

**ABSTRACT.** A simplified approach to determine the impact on a space-debris object (a target) from the ion thruster of a spacecraft (a shepherd), which was proposed before in the context the ion beam shepherd technology for space debris removal, was considered. This simplified approach is based on the assumption of the validity of the self-similar model of the plasma distribution in the thruster plume. A method for the calculation of the force impact using the information about the contour of the central projection of the object on a plane, which is perpendicular to the ion beam axis, was proposed within the framework of this model. The errors of this method, including the errors caused by an inaccuracy of its realization, are analyzed. The results of the analysis justify the admissibility of the application of the specified approach within the self-similar model of the plasma distribution. The preliminary conclusion has been made that this simplified approach can be used to control the relative motion of the shepherd - target system as well. This conclusion is based on the results of the simulation of the system motion, when the “real” value of the thruster impact is calculated by the direct integration of the elementary impacts over the target surface and the value of the same impact used in the control algorithms is determined using the information about the contour of the target. A number of factors such as the orbital motion of the system, external perturbations, and the attitude motion of the shepherd were neglected in the simplified model which was used for the simulation. These factors and errors in the interaction model are necessary to consider during a more detailed analysis of this approach. The analysis of the calculation errors presented in this paper can be used during implementation of the ion beam shepherd technology for active space debris removal.

**Introduction.** The technology for the removal of large debris from the low-Earth-orbit called ion beam shepherd is presented in the papers [1], [2]. This technology provides space debris de-orbiting due to the impact of the ion plume of the electric thruster (ET) of a spacecraft (a shepherd) located in close proximity to space debris object (a target). A certain distance between ion beam shepherd (IBS) and the target must be maintained to provide effectiveness of the impact of the ion beam. The ion beam impact can be determined by integration over the surface of a target when the mechanism of the ion interaction [3] with an elementary area of the surface and its relative position are known. This approach can be used for modeling of the system motion but not recommended to apply in a control loop due to its computational costs and the incomplete information about the position and shape of the target. A simplified approach to calculate the impact from the ion beam on the target, based on the information about the contour of its central projection on the image plane of the IBS camera was proposed in [4, 5]. The central projection of the target on the reference plane, which is perpendicular to the beam axis, is considered within the framework of this approach instead of its surface. This paper justifies the application of this approach to control the relative motion of the shepherd - target system by analyzing errors of the determination of the plasma beam impact on the target

**Calculation of the force impact by the direct integration over the target surface.** The model of the interaction of the ion beam with a space debris object, as well as the model of the ion plume should be considered for the calculation of the impact from the ET force on the target.

Neglecting the sputtering of the target material, the escaping ions from the target surface, and the electron pressure, the elementary force $dF_s$ transmitted to the target can be calculated as follows [6]:
$$d\mathbf{F}_s = m\mathbf{u}(-\mathbf{v} \cdot \mathbf{u})ds, \qquad (1)$$

where $m$ – is the ion mass;

$\mathbf{u}$ – is the vector of the ion velocity;

$ds$ – is the elementary area of the surface, which location will be characterized by the radius vector $\rho_s$ of its midpoint;

$\mathbf{v}$ – is the unit normal vector to the element of the surface;

$n$ – is the density of the plasma.

The force and torque transmitted from the ion beam to the target can be calculated by the integration of the elementary forces from equation (1) over the exposed surface $S$

$$\mathbf{F}_{\text{srf}} = \int_S d\mathbf{F}_s, \quad \mathbf{M}_{\text{srf}} = \int_S \rho_s \times d\mathbf{F}_s. \qquad (2)$$

There are ion beam models for the near and far regions of the plasma plume [7]. The far region of the plume presents the main interest in the context of the non-contact space debris removal, because this is where the plasma interacts with the target. Models with different degrees of complexity and accuracy were proposed for the description of the far region of an ET beam [8]. The so-called self-similar model of the beam is chosen for this study. Taking into account the fact that the Mach number at the beginning of the far region of the beam is much greater than 1, the character of the plasma distribution approaches to a cone. The plasma density can be determined for this case using the self-similar model at an arbitrary point as follows [8]:

$$n = \frac{n_0 R_0^2}{z^2 \operatorname{tg}^2 \alpha_0} \exp \left( -3 \frac{r^2}{z^2 \operatorname{tg}^2 \alpha_0} \right), \qquad (3)$$

where $r, z$ – are radial (distance from the midpoint of the surface element to the axis of the beam cone) and axial (the distance from the vertex of the cone along the axis of the beam) coordinates of the point;

$R_0$ – is the radius of the beam at the beginning of the far region (at the exit of the ET nozzle, $z = R_0 / \operatorname{tg}^2 \alpha_0$);

$n_0$ – is the plasma density at the beginning of the far region;

$\alpha_0$ – is the divergence angle of the beam.

Axial $u_z$ and radial $u_r$ velocity components of the plasma ions can be represented as follows:

$$u_z = u_{z0} = \text{const}, \quad u_r = u_{z0}r/z, \qquad (4)$$

where $u_{z0}$ is the axial component of the ion velocity at the beginning of the far region.
The target surface is divided into elements for determination of the integral force $\mathbf{F}_{srf}$ transmitted to the target by the ion beam. After that, the coordinates of the surface elements and the unit normal towards them are set with respect to the target reference frame (RF). Then the coordinates of the surface elements and their normals recalculated with respect to the IBS RF, or more precisely, to whose origin that is at the cone vertex of the ion beam. We can assume for simplicity reason that this RF coincides with the IBS RF. The elementary force is determined for each element of the surface using the equation (1) and taking into account its "visibility" from ET side and its position inside or outside of the beam cone. The integral force is determined using the expression (2).

**Calculation of the force impact using the contour of the target.** The simplified approach [5] for determination of the impact from the ion plume uses the information about the contour of the central projection of the target. According to this approach the force applied to an element of the surface is approximately equal to the force acting on the central projection of this area on a plane that is perpendicular to the axis of the beam cone, for example the plane of camera sensor placed next to the ET This assumption is explained by the fact that the cone cross-section increases with the square of the distance from the vertex of the cone and the plasma density decreases in inverse proportion to the square of the distance from the vertex of the cone.

The equation that determines the force $d\mathbf{F}_{\sigma}$ transmitted through the element of the surface can be written as follows:

$$d\mathbf{F}_{\sigma} = mn_c \mathbf{u}_c^2 \mathbf{e}_u d\sigma, \quad \mathbf{u}_c = u_{z0} \cdot \begin{bmatrix} x_c / f & y_c / f & 1 \end{bmatrix}^{\mathrm{T}}, \qquad (5)$$

$$n_c = \frac{n_0 R_0^2}{f^2 \mathrm{tg}^2 \alpha_0} \exp \left( -3 \frac{x_c^2 + y_c^2}{f^2 \mathrm{tg}^2 \alpha_0} \right), \qquad (6)$$

where T – is the transposition symbol;

$d\sigma$ – is the elementary area of the target projection on the camera image plane;

$\mathbf{e}_u$ – is the unit vector of the direction of the ion velocity $\mathbf{u}$;

$x_c, y_c$ – are coordinates of the point in the camera RF;

$f$ – is the focal length of the camera.

The full force $\mathbf{F}_{\mathrm{cnt}}$ transmitted from the ion plume to the target is calculated by this expression:

$$\mathbf{F}_{\mathrm{cnt}} = \int_{\Sigma} d\mathbf{F}_{\sigma}, \qquad (7)$$

where $\Sigma$ is the part of the image plane of the camera which is bounded by the contour.

**Analysis of the errors caused by the simplified calculation method.** Let us make a more detail analysis of the assumption about the equality of two elementary forces from the ion plume used in the simplified method. The first force acts on an element of the target surface and the second one is applied to the central projection of this area on a plane, that is perpendicular to the axis of the cone beam. First, we restrict ourselves to the case when the area $ds$ of the target surface is perpendicular to the direction of the ion beam. The equations (1), (3), (4) can be rewritten for this area as follows:
$$d\mathbf{F}_s = mn\mathbf{u}^2\mathbf{e}_u ds, \quad \mathbf{u}^2 = u_{z0}^2(1 + \operatorname{tg}^2\theta), \qquad (8)$$

$$n = m \frac{n_0 R_0^2}{z^2 \operatorname{tg}^2 \alpha_0} \exp\left(-3 \frac{\operatorname{tg}^2 \theta}{\operatorname{tg}^2 \alpha_0}\right), \quad ds = z^2 (\operatorname{tg}\theta / \cos\theta) d\theta d\varphi, \qquad (9)$$

where $\theta$ – is the angle between the direction of the beam distribution and the axis of the ion beam cone;

$\varphi$ – is the azimuth of the direction of the ion distribution, $\mathbf{e}_u = \mathbf{e}_u(\theta, \varphi)$.

The equations (8), (9) show that value $d\mathbf{F}_s = d\mathbf{F}_s(\theta, \varphi)$ is a function of variables $\theta$ and $\varphi$, and doesn't depend from the coordinates of the element $z$ of the surface. The "visible" area of the surface element varies in the case of its arbitrary orientation depending from the inclination angle of the element to the plane that is perpendicular to the direction of the ion beam. This variation of the area is taken into account in the equation (1) by the factor $(-\mathbf{v} \cdot \mathbf{u})$, that is equal to use the orthogonal projection of the area on the defined perpendicular plane. However, we do not deal with the orthogonal projection during the calculation of the force impact using the target contour but with the central one. In other words, a contour error is added during the calculations. Let's consider the following example on Fig. 1 to evaluate this error.

<figure><img src="image_3.png" /><figcaption>Fig. 1. Analysis of error caused by the replacement of the orthogonal projection by the central projection.</figcaption></figure>

We consider a part of the beam cone in the form of a small cone with the vertex angle $\alpha$ and name it elementary cone. Fig. 1a shows a cross-section of the cone by a vertical plane. The elementary cone "cuts" the element of the surface in the form of an ellipse (the segment AB in Fig. 1a). We assume
that the plane of the ellipse is inclined to the plane which is perpendicular to the axis $\zeta$ of the elementary cone at an angle $\beta$. The area $S_{cone}$ of the perpendicular cross-section of the elementary cone passing through the midpoint $M$ is the area of the central projection of the surface element on plane of this section. We also define the ratio $k_S$ of the area $S_{cone}$ to the area $S_{body}$ of the surface element. The difference between the value of this ratio and the value of $\cos \beta$ characterizes the error of the application of the equations (8) and (9). We introduce the so-called a coincidence coefficient defined by the function: $k_{coin} = 1 - (k_S - \cos \beta) / \cos \beta$.

The analysis of this function is not presented here due to its cumbersome, but it still shows the following: the independence of the coincidence coefficient of from the distance from a surface element to the cone vertex; this kind of error appears only when the values of $\beta$ are close to 90°; the acceptable results can be achieved by reducing the cone vertex angle $\alpha$ even when values of $\beta$ are close to 90°. Fig. 1 illustrates these conclusions by graphs of the variations of the coincidence coefficient from the following variables: the coordinate $\zeta_A$ (Fig. 1b); the cone vertex angle $\alpha$ (Fig. 1c); the inclination angle $\beta$ (Fig. 1d). One of these variables is varying while the other two variables are fixed in these graphs.

**Analysis of the errors caused by the implementation of the simplified approach.** Another kind of errors of the simplified approach can be caused by an inaccuracy of the camera placement. According to the approach the camera must be located so as its focal point coincides with the vertex of the imaginary cone of the ET beam. This requirement is difficult to meet from the engineering point of view and it is worth to consider an offset of the camera from the position which stipulated by the simplified method.

To analyze this offset we introduce the camera RF $O_c x_c y_c z_c$, whose origin $O_c$ is on the camera optical axis. The plane $x_c y_c$ coincides with the plane of the camera sensor. The axis $z_c$ is directed towards the target.

We consider the vector $d$, which connects the vertex of the imaginary cone of the beam and a point $P$ of the target. In the ideal case of the camera placement, the coordinates of the projection point $P$ on the plane of the camera sensor that included in the equation (5) are defined by following relations:

$$x_c = f \frac{d_1^{(cam)}}{d_3^{(cam)}}, \quad y_c = f \frac{d_2^{(cam)}}{d_3^{(cam)}}, \quad (10)$$

where $x_c, y_c$ – are the coordinates of the projection point in the camera RF, the subscript indicates the component number of the column $d^{(cam)}$ corresponding to the vector $d$.

Let us consider the case where the camera is mounted with a small offset relative to the ET beam. The target contour obtained using images from the camera with such offset is different from the one with the ideal placement of the camera. We will illustrate this considering the example of the contour determination.

The target is a circular cylinder. The height of the cylinder is 2,6 m. The diameter of the foundation of the cylinder is 2,2 m. The geometric center of the target is located on the axis of the beam cone on a distance of 7 m from its vertex The line 1 in Fig. 2 shows the contour corresponding to the case when the focal point of the camera coincidences with the vertex of the imaginary cone of the beam and the line 2 depicts the contour for the case when the camera offset is 0,2 m.
<figure><img src="image_4.png" /><figcaption>Fig. 2. Contours of a target projection.</figcaption></figure>

A method to correct coordinates of the contour points obtained using the equations (10) were proposed in [9] for the case when the cameras focal point is placed with offset from the imaginary vertex of the beam cone at the $\tilde{x}$ and $\tilde{y}$ along $O_c x_c$ and $O_c y_c$ axes, respectively. According to this method, the coordinates can be corrected as follows:

$$x_c^{\text{corr}} = x_c - f \frac{\tilde{x}}{z_{\text{nom}}}, \quad y_c^{\text{corr}} = y_c - f \frac{\tilde{y}}{z_{\text{nom}}}, \qquad (11)$$

where $x_c^{\text{corr}}$, $y_c^{\text{corr}}$ are the corrected coordinates of the point of the target projection;

$z_{\text{nom}}$ – is the nominal distance between the geometric center of the target and the vertex of the beam cone.

Line 3 in Fig. 2 shows the contour corrected according to the equation (11) for the case when $z_{\text{nom}} = 7$ m. This figure shows that the corrected contour is almost identical to the one for the no offset case (line 1).

### The possibility to use the simplified approach to control the shepherd-target relative motion.

The preliminary assessment of the possibility to use the simplified approach to determine the beam impact for the control of the shepherd-target relative motion is made on the basis of the simulation results. The "real" value of the beam impact is calculated by direct integration over the target surface and the values of the beam impact that used in the control algorithm are determined using the information about the target contour. The confirmation of such possibility is an indirect justification of the simplified approach for the calculation of the impact on a space debris object.

A simplified model of the motion was used for the preliminary evaluation where a number of factors were not taken into account, such as:

- the orbital motion of the IBS and target;

- external perturbations;

- the IBS attitude motion.
We start describing the reference frames that are required for a full model and then will apply the simplifying assumptions introduced above. The following reference frames are used:

$Ox_0y_0z_0$ is the inertial RF, the origin $O$ is located in the center of the Earth, axis $Oy_0$ is directed towards the north pole of the Earth, axis $Oz_0$ is the midpoint of the spring equinox;

$Sx_1y_1z_1$, $Tx_3y_3z_3$ are the orbital reference frames of the IBS and target, respectively, the origin of the RF $S$ and $T$ are located at their centers of mass, axis $z$ is directed along the radius vector connecting the center of the earth and the center of mass of the IBS or the target, axis $x$ is in the orbital plane and directed towards the orbital motion;

$Sx_2y_2z_2$, $Tx_4y_4z_4$ are the body RF of the IBS and the target respectively, the axes of which coincide with main central axes of inertia. For the ideal orientation of the IBS and target, these RF are parallel and coincide with the appropriate orbital reference frame.

We will also refer to this reference frames by zero, the first ..., the fourth according to the subscripts which are used in their notation. The transitions matrixes from $\nu$-th RF to $\mu$-th RF are denoted like $\Gamma_{\nu}^{\mu}$, $\nu, \mu = 0, 1, ..., 4$.

The IBS attitude position is determined by the following rotations of its body RF relative to its orbital RF on the angles: pitch $\vartheta_2$ (around the $y$ axis), roll $\varphi_2$ (around the $x$ axis) and yaw $\psi_2$ (around the $z$ axis). In order to avoid singularity of the kinematics relations for the case of the uncontrolled motion of the target, the attitude position of the target is more convenient to describe by four parameters of Rodrigues-Hamilton [10]. The parameters $\lambda_0, \lambda_1, \lambda_2, \lambda_3$ is used to specify the position of the body RF with respect to the orbital RF.

The dynamic equations of the attitude motion of the IBS and target is used in the Euler form [11]:

$$ \begin{cases} J_{vx} \dot{\omega}_{vx} + (J_{vz} - J_{vy}) \omega_{vy} \omega_{vz} = M_{vx}, \\ J_{vy} \dot{\omega}_{vy} + (J_{vx} - J_{vz}) \omega_{vx} \omega_{vz} = M_{vy}, \\ J_{vz} \dot{\omega}_{vz} + (J_{vy} - J_{vx}) \omega_{vx} \omega_{vy} = M_{vz}, \end{cases} \qquad (12) $$

where $J_{vx}$, $J_{vy}$, $J_{vz}$ are the moments of inertia of the IBS and target with respect to their principal central axes;

$M_{vx}$, $M_{vy}$, $M_{vz}$ – are projections of external torques acting on the IBS ($\nu = 2$) and target ($\nu = 4$).

These dynamic equations are supplemented with the kinematics relations:

$$ \begin{cases} \omega_{2x} = \phi_2 \cos \psi_2 + \vartheta_2 \cos \varphi_2 \sin \psi_2 + \omega_{20} \cos \varphi_2 \sin \psi_2, \\ \omega_{2y} = \vartheta_2 \cos \varphi_2 \cos \psi_2 - \phi_2 \sin \psi_2 + \omega_{20} \cos \varphi_2 \cos \psi_2, \\ \omega_{2z} = \psi_2 - \vartheta_2 \sin \varphi_2 - \omega_0 \sin \varphi_2, \end{cases} \qquad (13) $$

$$ \begin{cases} \omega_{4x} = 2(\lambda_0\lambda_1 - \lambda_1\lambda_0 + \lambda_2\lambda_3 - \lambda_3\lambda_2), \\ \omega_{4y} = 2(\lambda_0\lambda_2 - \lambda_2\lambda_0 + \lambda_3\lambda_1 - \lambda_1\lambda_3), \\ \omega_{4z} = 2(\lambda_0\lambda_3 - \lambda_3\lambda_0 + \lambda_1\lambda_2 - \lambda_2\lambda_1), \end{cases} \qquad (14) $$
$$
\left\{
\begin{aligned}
\phi_2 &= \omega_{2x} \cos \psi_2 - \omega_{2y} \sin \psi_2, \\
\theta_2 &= (\omega_{2x} \sin \psi_2 + \omega_{2y} \cos \psi_2) / \cos \phi_2 - \omega_0, \\
\psi_2 &= \omega_{2z} + (\omega_{2x} \sin \psi_2 + \omega_{2y} \cos \psi_2) \operatorname{tg}\phi_2,
\end{aligned}
\right.
\qquad (15)
$$

$$
\left\{
\begin{aligned}
2\lambda_0 &= -(\omega_{4x}\lambda_1 + \omega_{4y}\lambda_2 + \omega_{4z}\lambda_3), \\
2\lambda_1 &= \omega_{4x}\lambda_0 - \omega_{4y}\lambda_3 + \omega_{4z}\lambda_2, \\
2\lambda_2 &= \omega_{4y}\lambda_0 - \omega_{4z}\lambda_1 + \omega_{4x}\lambda_3, \\
2\lambda_3 &= \omega_{4z}\lambda_0 - \omega_{4x}\lambda_2 + \omega_{4y}\lambda_1,
\end{aligned}
\right.
\qquad (16)
$$

where $\omega_{vx}, \omega_{vy}, \omega_{vz}$ are the projections of the angular velocities of the IBS and target to their body RF, $v = 2, 4$;

$\omega_0$ – is the orbital angular velocity.

The components $\Gamma_{\nu}^{\nu-1}$, $\nu = 2,4$ of the transition matrixes from the body RF to the corresponding orbital ones are given by:

$$
\Gamma_2^1 = \begin{bmatrix}
c\vartheta_2 c\psi_2 + s\vartheta_2 s\varphi_2 c\psi_2 & -c\vartheta_2 s\psi_2 + s\vartheta_2 s\varphi_2 c\psi_2 & s\vartheta_2 c\varphi_2 \\
c\varphi_2 s\psi_2 & c\varphi_2 c\psi_2 & -s\varphi_2 \\
-s\vartheta_2 c\psi_2 + c\vartheta_2 s\varphi_2 c\psi_2 & s\vartheta_2 s\psi_2 + c\vartheta_2 s\varphi_2 c\psi_2 & c\vartheta_2 c\varphi_2
\end{bmatrix}, \qquad (17)
$$

$$
\Gamma_4^3 = \begin{bmatrix}
\lambda_0^2 + \lambda_1^2 - \lambda_2^2 - \lambda_3^2 & 2(\lambda_1\lambda_2 - \lambda_0\lambda_3) & 2(\lambda_1\lambda_3 + \lambda_0\lambda_2) \\
2(\lambda_1\lambda_2 + \lambda_0\lambda_3) & \lambda_0^2 + \lambda_2^2 - \lambda_3^2 - \lambda_1^2 & 2(\lambda_2\lambda_3 - \lambda_0\lambda_1) \\
2(\lambda_1\lambda_3 - \lambda_0\lambda_2) & 2(\lambda_2\lambda_3 + \lambda_0\lambda_1) & \lambda_0^2 + \lambda_3^2 - \lambda_1^2 - \lambda_2^2
\end{bmatrix}, \qquad (18)
$$

where $c\gamma$ and $s\gamma$ are cos $\gamma$ and sin $\gamma$, respectively.

As mentioned above, the orbital motion ($\omega_0 = 0$) is neglected in this study and it is assumed that $Ox_0y_0z_0$ is the inertial RF and the axes of reference frames $Sx_1y_1z_1$ and $Tx_3y_3z_3$ are parallel to the corresponding axes of $Ox_0y_0z_0$, and for the ideal orientation of the IBS and target their body RF coincide with the corresponding orbital one. The RF for this case are shown in Fig. 3, where $\rho$, $d$, $R$ are the radius vectors of a point $P$ in the corresponding reference frames; $r_{ST}$ is the radius vector of center of mass (CoM) of the target in the IBS body RF; $r_S$, $r_T$ are the radius vectors of CoM the IBS and target in $Ox_0y_0z_0$ RF.
<figure><img src="image_5.png" /><figcaption>Figure 3. Reference frames</figcaption></figure>

With these assumptions the following equations are correct:

$$d^{(2)} = (\Gamma_2^1)^T \cdot r_{ST}^{(0)} + (\Gamma_2^1)^T \Gamma_4^3 \cdot \rho^{(4)}, \quad r_{ST}^{(0)} = r_T^{(0)} - r_S^{(0)}, \qquad (19)$$

where the variables with superscript in parenthesis refer to the column of the vector projections to the axes of RF with the appropriate index.

These relations are necessary to determine the impact of the ion plume on the area of the target surface with the conditional center at a point P. The attitude angles and Rodrigues-Hamilton parameters included in the matrixes $\Gamma_2^1$, $\Gamma_4^3$ were determined integrating the equations (12) - (16) during the modeling of the system motion. Due to the accepted assumptions, the values of the components $\mathbf{r}_S$, $\mathbf{r}_T$, are calculated integrating the following equations:

$$m_S \cdot \dot{\mathbf{r}}_S^{(0)} = F_S^{(0)}, \quad m_T \cdot \dot{\mathbf{r}}_T^{(0)} = F_T^{(0)}, \qquad (20)$$

where $F_S^{(0)}, F_T^{(0)}$ are the vectors of the resultant force $\mathbf{F}_S, \mathbf{F}_T$, applied to the IBS and target respectively;

$m_S$ and $m_T$ – masses of the IBS and target. In this way, according to accepted simplified model, the system motion is described by the equations (12) - (20) under the condition $\omega_0 = 0$.

In the absence of external perturbations on the system, the force impact on the target is created only by the ion plume from the shepherd ET. This impact slows the orbital velocity of the target and allows to de-orbit it faster. The IBS concept assumes that two ETs are installed on the IBS. The second thruster is necessary to compensate the impact on the IBS motion from the main one, directed towards the target. However, in order to the IBS and target remain in the same orbit, the result of the force action on CoM of the IBS has to be the same as the result of the force impact on CoM of the target. This condition can be achieved by adjustment of the thrust of the second compensating thruster. The algorithm to adjust the thrust $\mathbf{F}_{E2}$ of compensating ET was chosen in the following simple form:

$$\mathbf{F}_{E2} = (m_S / m_T) \cdot \mathbf{F}_{\text{cnt}} - \mathbf{F}_{E1}, \qquad (21)$$
where $F_{E1}$ is the thrust of the main ET.

Assuming that the IBS attitude position is unchanged, this algorithm can be simplified as follows:

$$F_{E2\_x} = (m_s / m_T) \cdot F_{cnt\_x} - F_{E1\_x}, \qquad (22)$$

where the letter $x$ in the subscripts is used to denote the $x$-th component of the column, that corresponded to the thrust vector in the RF $Ox_0y_0z_0$ or the body RF of the IBS.

The algorithm to determine the target contour for solving modeling tasks is described in details in [5] and can be represented by the following steps:

- Approximate the target surface by basic elements;

- Project the central points of the surface elements on the image plane of the camera;

- Calculate the contour of the target projection by solving the problem of the polygon construction, which covers a set of points projected on a plane.

We assume that the origin $O_c$ of the camera RF coincides with the origin $S$ of the body RF of the IBS and the $z_c$ axis is directed opposite to the axis $x_2$ of the body RF of the IBS. Then, the coordinates of the point $P$ of the target is projected on the plane defined by the equations (10) or calculated using the relations (11) in the case of their correction.

The components of column $d^{(cam)}$ which are included in these equalities are given by

$$d^{(cam)} = \Gamma_2^{cam} \cdot d^{(2)}, \qquad (23)$$

where $\Gamma_2^{cam}$ is the transition matrix from the body RF of the IBS to the camera RF.

**Simulation results.** The assumption was made for the calculations that the only forces from the main and compensating ETs act on the IBS. The line of action of these forces passes through CoM of the IBS.

The following values of parameters were used during the simulation.

Parameters of the ETs are the following: the radius of the beam at the beginning of the far region $R_0 = 0,0805 \text{ m}$; the plasma density at the beginning of the far region $n_0 = 4,13 \cdot 10^{15} \text{ m}^{-3}$; the divergence angle of the beam $\alpha_0 = 7^\circ$; the axial component of the velocity of the plasma ions $u_{z0} = 71580 \text{ m/s}$; the ion mass $m = 2.18 \cdot 10^{-25} \text{ kg}$.

The IBS parameters were chosen as follows: the matrix of inertia is $diag(1283,4; 1379,5; 169,3) \text{ kg} \cdot \text{m}^2$; the mass is 500 kg.

The following parameters of the cylindrical target were used: the mass is 1000 kg; the height is 2,6 m; the diameter is 2,2 m. The inertia matrix of the target was calculated according to the formulas for the moments of inertia of a hollow cylinder.

The distance of 7 m from CoM of the IBS to the target along the axis $z$ of the camera RF was chosen as nominal. The focal length of the camera is 0,2 m.
The simulation results show that the difference between the nominal and the calculated distance from the target to the IBS is less than 1 sm after 800 s of the simulation when the initial location of CoM of the target is on the beam axis. The figures presented below show the results of the simulations for the case when the initial location of CoM of the target is 1 m away from the beam axis, the initial attitude position of the target is defined by the angle of 45° about the axis y of the body RF of the target, the camera offset from its nominal position is 0,2 m along the abscissa axis of the camera RF.

<figure><img src="image_6.png" /><figcaption>Fig. 4. The force acting on the IBS.</figcaption></figure>

Fig. 4 and Fig. 5 show the graphs of the forces acting on the IBS and target, respectively. Ideally these graphs should coincide with each other up to the scale. The variation of the target position (Fig. 6) from the nominal one for 800 s of the simulation is no more than 5 sm (this variation is denoted as $\delta^0 r_x$ in the figure) for the ideal camera position, and does not exceed 8 sm for the case when the camera has the position offset. The variations that obtained through the simulation (i.e. the errors of maintaining the nominal distance to the target) are not significant, taking into account to the real control algorithms will use the information about the distance between the target and IBS or its estimation unlike the applied simple algorithm (21), which does not use this information.

<figure><img src="image_7.png" /><figcaption>Fig. 5. The force acting on the target</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Fig. 6. The deviation of the target relative position</figcaption></figure>
Fig. 7 and Fig. 8 show the graphs of the torque acting on the target and the graph the angular velocity of the target, respectively.

<figure><img src="image_9.png" /><figcaption>Fig. 7. The torque acting on the target</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 8. The angular velocity of the target</figcaption></figure>

**Conclusions.** The errors of the simplified approach [5] to determine the impact from the ion thruster of a spacecraft (a shepherd) on a space debris object (a target) within the context the IBS technology [1, 2] has been analyzed. The approach is based on the information about the contour of the central projection of the target on a plane, which is perpendicular to the axis of the ion beam. The plane of the camera sensor, which is installed on the shepherd, has been considered as such a plane. The results of the error analysis justify the admissibility of the application of this simplified approach for the determination of the force impact of the ion thruster to the space debris within self-similar model of the plasma distribution [6]. The research results also allow to make a preliminary conclusion about the possibility to use this simplified approach for the control of the relative motion of the IBS- space debris object system. However, a more detailed analysis of this possibility is needed which require to consider the orbital motion, all range of acting disturbances, errors of the model of the interaction of the ion plume with the surface of the space debris object, and errors of implementation of the simplified approach.

**Acknowledgement.** The research leading to these results has received funding from the European Union Seventh Framework Program (FP7/2007-2013) under grant agreement n°607457.

**References**

[1] Bombardelli, C., Pelaez, J. *Ion Beam Shepherd for Contactless Space Debris Removal* // Journal of Guidance, Control and Dynamics. – 2011. – Vol. 34, #3. – Pp. 916–920.

[2] Bombardelli, C., Alpatov A.P., Pirozhenko, A.V., Baranov, E.Y., Osinovy, G.G., Zakrzhevskii, A.E. *Project "Space shepherd" with ion beam. Ideas and problems* // Space science and technology. – 2014. – V. 20, #2. – Pp. 55–60.

[3] Merino, M., Ahedo, E., Bombardelli, C., Urrutxua, H., Peláez, J. *Ion beam shepherd satellite for space debris removal* // Progress in Propulsion Physics. – 2013. – Vol. 4. – Pp. 789–802.

[4] Alpatov, A., Cichocki, F., Fokov, A., Khoroshylov, S., Merino, M., Zakrzhevskii, A. *Algorithm for determination of force transmitted by plume of ion thruster to orbital object using photo camera* // 66th International Astronautical Congress, Jerusalem, Israel. –2015. – Paper IAC-15-A6.5.5-x27732. – Pp. 1–9.

[5] Alpatov A., Cichocki F., Fokov A., Khoroshylov S., Merino M., Zakrzhevskii A. *Determination of the force transmitted by an ion thruster plasma plume to an orbital object* // Acta Astronautica.– 2016.–V. 119. – Pp. 241–251.
[6] Bombardelli C., Urrutxua H., Merino M., Ahedo E., Peláez J. *Relative dynamics and control of an ion beam shepherd satellite* // Spaceflight mechanics. – 2012. – Vol. 143. – Pp. 2145–2158.

[7] Merino M., Cichocki F., Ahedo E. *A collisionless plasma thruster plume expansion model* // Plasma Sources Science and Technology. – 2015. – Vol. 24(3), – Pp. 1–12.

[8] Bombardelli C., Merino M., Ahedo E., Peláez J., Urrutxua H., Iturri-Torreay A., Herrera-Montojoy J. *Ariadna call for ideas: Active removal of space debris ion beam shepherd for contactless debris removal* // Technical report. – 2011. – 90 p.

[9] Alpatov A.P., Zakrzhevskii A.E., Merino M., Fokov A.A., Khoroshylov S.V., Cichocki, F. *Determination of a force transmitted by a plume of an ion thruster to an orbital object* // Space science and technology. – 2016. – V. 22, #1. – Pp. 52–63.

[10] Kravets V.V., Kravets T.V., Kharchenko A.V. *Using quaternion matrices to describe the kinematics and nonlinear dynamics of an asymmetric rigid body* // Int. Appl. Mech. – 2009. – 44.#2. – Pp. 223-232.

[11] Lurie A.I. *Analytical Mechanics*. – M. Fizmatgiz. – 1961. – 824 p.