Allocation Optimization Strategy for High-Precision Control of Picosatellites and Nanosatellites

The high-precision control of picosatellites and nanosatellites has always plagued the astronautics field. Aiming to change the status quo of the actuators not being able to meet the high-precision attitude control of picosatellites and nanosatellites, this article formulates a control allocation strategy for picosatellites and nanosatellites using the solid propellant microthruster array (SPMA). To solve the problem of the diversity and complexity of ignition combinations brought about by the high integration of the SPMA, the energy consumption factor of the optimal allocation is established, and the relationships of the array’s energy consumption factor, the control accuracy, the number, and the ignition combinations of the thruster array are deduced. The optimization objective is introduced by Sherman-Morrison formula and singular value decomposition. Thus, the energy consumption problem is transformed into an integer programming problem, acquiring the control allocation strategy and the optimal thruster energy. Simulation results show that the proposed algorithm can effectively reduce the thrust energy consumption and improve the precision control, demonstrating the feasibility and efficiency of the proposed algorithm for picosatellites and nanosatellites.


Introduction
In the past years, a tremendous amount of work has been done on research and development of picosatellites and nanosatellites (PNS) because PNS can offer a wide variety of applications in space.PNS have significant scientific and cost advantages over using large, heavy satellites [1,2].Measurements of the lower thermosphere [3], earth observation [4], and technology demonstration are only a few examples of the potential of PNS.Till 2020, more than 400 PNS launches per year are estimated [5].All such missions need micropropulsion systems to maneuver the satellites.PNS have strict requirements for quality, size, power cost, and so on, and therefore, miniature and affordable propulsion systems with low power are urgently needed for PNS.Various systems have been researched to control the miniaturized satellites, such as cold gas microthrusters, electrospray microthrusters, field emission microthrusters, and mono/ bipropellant thrusters.The most commonly used actuators for PNS are magnetic coils and wheels [6][7][8][9].Wheels can provide continuous output without using fuel, whereas, when wheels are working at a low angular velocity, there are frictional torques that may have a strong influence on PNS.Wheels are excellent for using in conventional satellites because of its accurate output.For the PNS, miniaturized wheels would not only increase the cost but also reduce the accuracy.
This paper proposes solid propellant microthruster array (SPMA) as an interesting choice for pico-and nanosatellites.SPMA has several main advantages: (1) it has simple structure and high reliability, (2) it is free from the frictional forces inherent to moving parts, (3) there is no propellant leakage with high propellant stability, and (4) driving voltage is low [10][11][12][13][14].In SPMA, multiple independent addressing thrusters are highly integrated on the same chip (10 6  thrusters can be integrated on a 1 cm chip), which can deliver low thrust values and impulses in a short time.Simultaneously, our team has produced 10 × 10 scale microthruster arrays and 100 × 100 scale microthruster arrays with conventional fabrication methods [15].Nevertheless, the control allocation strategy for a large number of thrusters, misaligned and inappropriate positions for SPMA actuators, and irreproducible propellant and heat loss phenomena should be studied systematically before SPMA application.
Firstly, for a highly integrated actuator, a large number of thrusters have numerous combinations of ignition.Control torques produced by the thrusters in different locations are different and the thrusters in SPMA are not reproducible.Therefore, it is required that the actuators should not only satisfy the high-precision needs of satellite attitude control, but also have the smallest fuel consumption.For multiple and redundant thrusters in the SPMA, achieving the highest control accuracy with the least amount of fuel consumption is a dominant issue.
To solve the control allocation problems, there are several traditional methods: direct control allocation [16], daisy chaining [17], linear programming (LP) method [18], quadratic programming (QP) method [19], and pseudoinverse redistribution (PIR) method [20].Numerous works also have been carried out on the improvement of allocation algorithms [21], the application of the technique to engineering [22], and the extension to theoretical approaches in the recent years [23].Wang and Xie [24] designed an optimal thruster combinations table for the use of onboard realtime control allocation problems.In [25], a control allocation strategy is formulated as a min-max optimization problem, which ensures some robustness of system performances.Cristofaro and Johansen designed a fault-tolerant control allocation using unknown input observers [26].Kirchengast et al. [27] investigated the impact of input matrix factorization on two algorithms (direct allocation and redistributed pseudoinverse) for constrained control allocation.Liu et al. [28] proposed an adaptive control allocation algorithm by updating the number of operating actuators.Although in the work of Liu et al., the number of operating actuators is changed, the total number of actuators is small, and the control effects of actuators in array are the same, so the simplified model reduces the analysis.
To summarize, previous scientific works are generally based on few thrusters with a fixed number of operating actuators.There is no existing method that can solve control allocation problems for a large and changed number of actuators.In this paper, we propose a control allocation optimization approach for PNS by dynamically analyzing and updating the number, ignition combinations, and impulse of actuators.The contributions of this research are presented as follows: (1) For a large amount of thrusters like those in a highly integrated SPMA, the energy consumption of the thrusters is influenced by multiple factors such as the installation mode of SPMA, the ignition combinations, the unit impulse, and the precision control.
It is difficult to quantitatively analyze the energy consumption, ignition combinations, and control accuracy.Therefore, this paper proposes an array energy consumption factor (AECF) to evaluate the performance of thruster combinations, establishing the relationship between the number of thrusters and the total energy consumption.
(2) An allocation optimization algorithm based on the AECF is formulated that can reduce the thrust energy consumption and improve control resolution.

Control Allocation Model of SPMA
For the cubic satellite shown in Figure 1, O − xyz is the body coordinate system of satellite, the origin of the body coordinates locates at satellite's mass center, and principal axis coincides with inertial principal axes of satellite.The length of satellite is 2L R , and SPMA is installed in two symmetry planes perpendicular to the z-axis.As shown in Figure 1(a), the distance of each SPMA's edge from the satellite's symmetrical axis is L. The thruster arrays are installed on the diagonal line of satellite's surface, and the thruster serial number is shown in Figure 1(a).Supposing each SPMA has m × m, thrusters are represented as Γ m×m .The coordinates of the thruster in the SPMA can then be expressed as (i, j), i ∈ 1, m , j ∈ 1, m , i, j ∈ Z, and the position vector matrix of thruster array relative to the satellite centroid can be described as a matrix with the dimension of 1 × N(N = m 2 ): Δ ο indicates the unit distance of adjacent thrusters in the coordinate axis direction.The thruster control allocation model [29] is described as where η ∈ R n is the expected control command, and n is the degree of freedom (DOF) of the control requirement.F is the input control vector.A represents the thruster efficiency matrix.The control torque of the thruster with the coordinate (i, j) for the centroid is The impulse matrix composed of all the thrusters can be represented as Equation (3) shows an impulse matrix composed of N working thrusters.Currently, each impulse from the thrusters in the same array is equal, specifically F 11 = F ij = ⋯ = F N .The control torque at the satellite's centroid can be represented as The efficiency matrix A can be written as follows: For the satellite control system, the control dimension n = 3, A ∈ R 3×λ , and is the number of thrusters participating in control, λ ∈ 0, N .

Energy Consumption Model of SPMA
In this section, the energy consumption model of the SPMA according to the control allocation model is established, and energy consumption factor of SMPA is deduced.Thus, the relationship among the quantity, the ignition combinations, and the energy consumption of the thrusters can be studied.

Deduce Energy Consumption Factor. Rewriting (1) as
Supposing the control variance is σ 2 (representing the system error of actuator, namely, impulse error of thruster), I n×n is an N-order unit matrix and thus writing (7) in the form of substituting (8) into ( 9) and systemize as follows: The trace on both sides of (10) and then SPMA energy consumption model can be obtained: The left side of (11) represents the total required control impulse (thrust), namely, the energy consumption of the thruster in statistic.In this paper, tr is defined as the array energy consumption factor (AECF).Setting matrix K as The AECF can be expressed as Through the AECF's definition, it can be seen that the larger the AECF is, the larger the total energy consumption used for control will be.The smaller the AECF is, the smaller the total energy consumption used for control will be.Furthermore, according to (13) and the definition of matrix A, it can be seen that the AECF is influenced by three indexes: the number of working thrusters, the ignition combinations, and the unit impulse of thrusters in SPMA.

Energy Optimization Allocation Strategy.
When using the SPMA for high-precision control, because there are a large number of thrusters integrated in the array and under the condition of thruster redundancy configuration, allocation method that satisfies (1) is not unique.The total energy consumption will also be quite different when choosing a different number of thrusters to form ignition combinations.Therefore, to reduce the energy consumption of SPMA, the relationship between the number of working thrusters and value of AECF is studied; thus, an optimal energy allocation strategy based on AECF can be established.

International Journal of Aerospace Engineering
Setting A λ as the efficiency matrix of λ thrusters when being selected for ignition execution and randomly remove one thruster (supposing the coordinates of this thrust are (i, j)) from the λ thrusters, then the torque efficiency matrix of λ − 1 thrusters is A λ−1 .The following relationship ( 14) can be obtained using (7): variable K λ is introduced as Therefore, according to the Sherman-Morrison formula [30], ( 14) can be written as where According to the definition of the AECF in Section 3.1, we can obtain where From ( 17), it can be seen that AECF will add an item g ij when reducing a thruster.Therefore, whether g ij has a positive or negative value determines how the AECF is influenced by the number of working thrusters.
To determine the sign of g ij , A λ−1 is performed on singular value decomposition (SVD): where U and V T , respectively, represent a (λ − 1)-order orthogonal matrix and a 3-order orthogonal matrix and S is a (λ − 1) × 3-order diagonal matrix.Furthermore, K λ−1 can be written as S T S is a 3-order positive diagonal matrix as follows: where s 1 , s 2 , and s 3 are the eigenvalues of matrix A λ−1 .
Since the orthogonal transformation has no influence on the matrix trace, (21) can be obtained: Substituting ( 22) into ( 21), we can obtain International Journal of Aerospace Engineering We can also obtain Therefore, g ij can be written as Furthermore, we can determine that From ( 25) and ( 27), it can be seen that the AECF decreases with the increase of the number of working thrusters.Thus, when satisfying the same control variance σ 2 , increase of the number of working thrusters can reduce the demand of unit thruster energy consumption (thrust/ impulse), reduce the total energy consumption of SPMA, and improve the control resolution.

AECF Optimization Strategy
Based on the conclusion drawn from (27), the next step is optimizing the ignition combinations and achieving optimal energy consumption and precision control.
As deduced from Section 3, the smaller the AECF is, the lower the energy consumption is.From ( 25) and ( 27), we can see that the bigger the g ij value brought about by the addition of ignition thruster is, the smaller the AECF and the energy consumption are; thus, the key to optimal energy consumption is finding the optimal g ij value.To find this g ij value, a simplified optimization criterion is established for AECF optimization algorithm.Firstly, basement of chain control allocation is introduced to calculate the initial value of the AECF.

Initial Value of Chain
Basement.Chain basement algorithm considers each thruster in the thruster array as a basement and preferentially searches for the basement with a high control efficiency that meets the triaxial demand (force/torque).When one basement cannot complete the control force (torque), multiple sets of basements can be selected successively in a combined fashion to compensate.
The chain basement control allocation method for λ thrusters is expressed as

28
where C is the expected control torque given by (1); u i and a i are, respectively, the control basement and the direction vector of the corresponding basements.When the first group of basements cannot satisfy the control requirements, the next group of basements is added according to the basement efficiency until the expected control range is satisfied.
The basement efficiency can be defined from the following equation: where C ω is the control force (matrix) generated by the ω th thruster basement and C q is the expected control force (matrix).

AECF Optimization Principle.
Through the chain basement method, the initial χ thruster information A χ that satisfies the control requirements can be obtained.Equation (25) reflects the influence of interaction between the optimizing χ + 1 th thruster and the initial value A χ on AECF.s 1 , s 2 , and s 3 can be determined by SVD from A χ ; thus, finding the g ij value of the optimal χ + 1 th thruster (also, supposing the 5 International Journal of Aerospace Engineering coordinates for this thruster is (i, j)) can be changed to the following optimization problem.
where Φ is the array collection taking out the initial χ number of thrusters.This problem then turns into finding the appropriate v 1 , v 2 , and v 3 in order to reach the largest g ij value.
For microsatellites with small volumes and quality such as PNS, which have semidiameter of less than 1 m or far less than 1 m, it is easy to ascertain that each element a rirj ri ≤ 3, rj ≤ l in A χ is less than 1.In addition, according to the Frobenius theorem and the Gerschgorin estimation, the singular values s 1 , s 2 , and s 3 of matrix A χ are all less than 1.Therefore, in (25) Thus, the optimization objective can be further simplified as (32) and the calculation complexity can be markedly simplified.
For the SPMA, according to its layout collection Γ m×m in Figure 2, (31) can be further simplified: To optimize the overall energy consumption, adding one more thruster and making the incremental effect maximum, the energy impact factor of new ignition position is T .Then, ring the problem into the integer optimization concerning Finally, the AECF optimization can be easily solved using integer optimization algorithms such as branch and bound method.

Optimization of Three-Dimensional Configuration.
When SPMA is used to provide triaxial directional control output, SPMA should form a specific included angle with plane O − xyz as shown in Figure 3.
Therefore, (5) can be modified as Equation (34) is rewritten as The relevant descriptions are similar to the above equations.The horizontal axis in Figure 4 shows the number of added actuators requiring optimization, and the vertical axis  International Journal of Aerospace Engineering is the value of AECF; it can be seen from Figure 4 that as the number of optimized thrusters increases, the AECF constantly decreases, and the simulation results verify the validity of (27).According to the definition, as the AECF decreases, the thruster energy consumption required by control in Figure 4 decreases.When the optimizing thruster quantity is 26, the AECF is 478, reduced by more than half the amount compared to the initial AECF.The horizontal axis in Figure 5 is AECF, and the vertical axis is the desired impulse of thruster.From Figure 5, it is obvious that as the AECF gradually decreases, the desired impulse of thruster significantly decreases.Eventually when AECF reaches a value of 478, the desired impulse of the thruster falls sharply, by 47.6%.This means that the overall thruster array can save nearly half of the energy consumption while achieving the same control goal, and it can also provide a higher resolution of control as well as meeting a higher precision control.

Simulation Analysis
The horizontal axis in Figure 6 is the number of added optimized thrusters, and the vertical axis is the demanded thruster energy consumption.The red line is the thruster energy consumption demand using the AECF optimization algorithm proposed in this paper, and the blue line is the thruster energy consumption demand randomly optimized and obtained from the unignited available thrusters.We can see that as the number of optimized thrusters increases, the desired energy decreases, which agrees with the results of theoretical analysis.Significantly, adopting the proposed optimization algorithm can save more than 13% of the energy consumption compared to the stochastic optimization.From Figure 7, similar conclusions can be drawn that with the increase of the optimized number of array thrusters, the AECF decreases gradually; in the end of the AECF process, the energy consumption is only 44.72% of the original level, decreasing by 55.28%.The results prove the efficiency of proposed algorithm.By dynamically analyzing and processing the information of thrusters based on AECF, allocation problem of large-scale thrusters can also be effectively solved.
Using the chain basement algorithm from Section 4.1, the corresponding initial number of the execution thrusters range is [16,52].The simulation results are shown in Figures 8-11.
Figure 8 demonstrates the corresponding changes of AECF affected by the number of working thrusters, reflecting increases from 16 to 52, while the number of optimized thrusters increases from 0 to 12. From Figure 8, we can see that no matter how many thrusters are participating in the execution, the AECF will still decrease because of the optimization and addition of the number of the aiming thrusters in this layout.Figures 9 and 10 reflect the corresponding changes in desired energy.Taking Figure 9 as an example, with the increase of the number of working thrusters in the vertical axis, the corresponding energy demand gradually decreases.From Figure 10, we can see that when the number of the thrusters of initial participation is relatively small and then add the thrusters, the effects of energy demand will be more significantly improved.
The influence of the AECF variation on the unit impulse of a thruster can be seen more clearly in Figure 11.
The horizontal axis in Figure 11 is the change rate of AECF, and the vertical axis is the change in desired energy.It can be seen no matter how many initial thrusters there are, with the decrease of AECF, the unit thruster energy demand significantly decreases.Taking the red curve and  8 International Journal of Aerospace Engineering 16 thrusters with initial participation as an example, when the energy consumption factor AECF decreases by 3.5%, its energy demand decreases by 8.7%, and as the AECF decreases by 39%, its energy consumption decreases by 52%.The horizontal axis in Figure 13 is the number of addition of optimized thrusters, and the vertical axis is the change of AECF.With the gradual increase of optimized thrusters, the AECF gradually shrinks.Taking the number 2 cyan curve as an example, C 2 = [1.28,−0.67, −0.67] × 10 −4 N ⋅ s ⋅ m and the AECF is initially 2028.With the increase of the optimized of thrusters, the AECF decreases, and when the optimization increases by 19 thrusters, the AECF sharply reduces to 1144.
The horizontal axis in Figure 14 is the change rate of the AECF, and the vertical axis indicates the change in the energy demand.It can be seen that with the decrease of AECF, unit thruster energy demand decreases significantly.In addition, taking the number 2 cyan curve as an example, when the AECF optimization of this curve decreases by 13%, its energy demand decreases by 43%. Figure 15 also clearly reflects the relationship between the number of initial thrusters, the number of thrusters participating in optimization, and the thruster energy demand.

Verification of Attitude Control.
To estimate the performance of the proposed approach more completely, large angle maneuver with perturbations has been tested in this  9 International Journal of Aerospace Engineering section.A sliding mode controller is designed to perform the large angle maneuver efficiently.Comparison with wheels is also been described in the verification part.
Kinetic equation and the motion equation of the PNS are expressed as follows:

38
where I b is the rotational inertia of satellite, T c denotes the control torque, T d is the disturbance from space environment, q e is the error of attitude angular rate represented by quaternion, and H w is the angular momentum of wheel.
The torque (T R ) generated from the wheel can be obtained by the difference between the command torque (T c ) and the friction torque (T f ) as follows [31,32]: Friction torques generated by wheels are particularly critical to PNS.For numerical analysis, the friction torques are modeled by means of a viscous torque proportional to the wheel speed, coulomb torque, and a continuous Stribeck friction as follows:  The derivative of V can be obtained from both SPMA and wheel as follows: When appropriate parameters are chosen such as α 1 and α 2 , the stability of the closed loop system can be guaranteed.(5) The errors of sensors are zero for both SPMA and wheel.
The optimization process of the proposed ACEF algorithm is presented in Figure 16.The variations of the attitude angle and the attitude angular velocity using AECF are shown in Figure 17.
Figure 16(a) shows that the AECF before the optimization is 1686.1435,with the increase of the optimal number of array thrusters, AECF decreases gradually; eventually, AECF significantly decreases to 572.3267, reflecting the marked decrease of energy consumption demand.From the Figure 16(b), we can see that at the end of the AECF process, the energy consumption is only 39.78% of the original level, decreased by 60.22%, showing the ability to save energy consumption effectively.It can be determined from Figure 17  12 International Journal of Aerospace Engineering that the attitude precision controls (3σ) are 0.0362 °, 0.0325 °, and 0.0204 °, and the precisions of the attitude angular velocities (3σ) are 0.0013 °/s, 0.0025 °/s, and 0.0052 °/s.The variations of the attitude angle and the attitude angular velocity using wheel are shown in Figure 18.
Figure 18(a) shows the attitude angles controlled by wheel; Figure 18(b) presents the attitude angular velocities.Figures 18(c) and 18(d) show the angular rates of the wheel and the corresponding friction torque produced by wheel.As can be seen from Figures 18(c) and 18(d), during the stable state, the angular rates of wheel always change around 0. Thus, friction torques are produced by wheel as demonstrated in Figure 18(d), reducing the precision controls of PNS.
The control results of SPMA, SPMA with AECF, and the wheels are shown in Table 1.
Table 1 shows that SPMA using AECF achieves highest precision the best performance.The errors of the attitudes obtained by SPMA using the AECF are 65.94%, 46.49%, and 80.95% of the errors obtained by SPMA.The errors of attitude angular velocity obtained by SPMA using AECF are 61.90%,62.50%, and 50.98% of the errors obtained by SPMA.Meanwhile, the errors of attitudes obtained by SPMA using AECF are 51.49%,43.39%, and 18.46% of the errors obtained by wheels; the errors of attitude angular velocity obtained by SPMA using AECF are 0.68%, 1.68%, and 2.34% of the errors obtained by wheels, demonstrating the feasibility and effectiveness of the proposed algorithm.The proposed AECF algorithm can improve control resolutions and reach higher precision controls.

Conclusion
Aiming at complicated large-scale control allocation problems, this paper solves the problem of previous optimal allocation methods not being able to make the overall design of the number, energy consumption, and configuration of large-scale thrusters.Through the establishment of SPMA application model, the energy consumption factor is put forward.With the increase of AECF, the energy demand of thruster increases.The decrease of AECF can reduce the thruster energy consumption and alleviate the satellite loads.With the increase of the number of optimal thrusters, the AECF will be reduced.Moreover, the relationship between the ignition thrusters and AECF is provided.The simulation results indicate that through AECF optimizing allocation, the thruster energy demand can be efficiently reduced, control resolution and precision can be improved, and energy consumption can be reduced as well.To the best of our knowledge, this is the first research to apply highly integrated actuators to control the picosatellites and nanosatellites.This work only studies the control allocation problems.As SPMA is not reproducible, to achieve higher precision and reliability, optimal compensation strategy for misalignment and inappropriate positions, heat loss phenomena will be our future work for SPMA.
AECF change of desired energy Unoptimized AECF change of desired energy

Figure 9 :
Figure 9: Changes in energy demand.
Friction torque of wheel (N.m)(d) Frictional torques produced by wheel

Figure 18 :
Figure 18: Control results of wheel.
11 , d 12 , … , d 1m , d 21 , d 22 , … , d 2m , d i1 , d ij , … , d im , … , d mm , where d ij indicates the thruster position vector relative to the satellite centroid and can be expressed as d ij = x ij , y ij , z ij .The corresponding thruster-generated triaxial component of the unit impulse is expressed as e ij = e ijx , e ijy , e ijz T .