High-safety and high-efficiency trajectory planning method and system for active observation of non-cooperative target

Abstract

The invention provides a high-safety and high-efficiency trajectory planning method for active observation of a non-cooperative target for a non-cooperative target observation task with relatively large navigation information error and relatively low updating frequency, and relates to the technical field of spacecraft orbit optimization. A state gain reachable set is calculated by a navigation initial error through a matrix differential equation, spacecraft structure size compensation is carried out on the state gain reachable set, the boundary of a security strategy is used as a path inequality constraint, path fuel consumption and path length are used as an optimization model of a composite performance index, and optimization variables are observation path states and corresponding control. And comprehensively establishing a trajectory planning model, and solving the optimization problem by using model prediction control so as to obtain a corresponding path state. The method is suitable for track planning of active observation of the non-cooperative target in aerospace orbit optimization.

Description

Technical Field

[0001] The invention relates to the technical field of spacecraft orbit optimization, and in particular to a trajectory planning method for active observation of non-cooperative targets using a state gain reachable set and fast model predictive control. Background Art

[0002] In recent years, the rapid development of on-orbit service (OOS) technology has not only made it possible to maintain and repair failed or faulty spacecraft, but also significantly increased the sensitivity of non-cooperative targets and their potential threats. It is particularly noteworthy that the United States successfully completed the on-orbit docking and capture mission through the "Orbit Express" (OE) in 2007, and achieved on-orbit takeover and life extension using the "Mission Extension Vehicle" (MEV) in 2021. These technological achievements fully demonstrate the potential of on-orbit service. Therefore, it is particularly important to accurately identify the main functions of non-cooperative spacecraft and their threat level.

[0003] Traditional on-orbit identification methods rely on the orbiting trajectory generated by the periodic solution of the HCW equation. The service spacecraft uses the navigation error ellipsoid of the target as the safety boundary, orbits and observes the target payload and structure, and thus calculates its attitude and relative orbit state. However, this method has two major limitations:

[0004] First, affected by the velocity measurement error, the navigation error ellipsoid will expand over time. In non-cooperative observation scenarios where positioning and velocity measurement errors are large and the navigation information update frequency is low, the navigation error ellipsoid tends to significantly underestimate the possible state range of the target, thereby increasing the risk of spacecraft center-of-mass collision. In addition, the navigation error ellipsoid cannot reflect the size and structural characteristics of the two spacecraft, which further increases the probability of spacecraft structure collision.

[0005] Second, since the closed trajectory generated by the periodic solution of the HCW equation is used to achieve the orbiting observation of the target, the observation mission cycle is usually at least one target orbit cycle, which is much longer than the ground orbit determination cycle. This scheme has poor information concealment and is easy to expose the mission object and intention, putting us in a passive situation. An effective way to shorten the mission cycle is to actively advance the orbiting observation target. It is necessary to solve the corresponding mathematical programming model to solve the observation trajectory sequence and its corresponding control sequence. Classical methods for solving trajectory planning problems include optimization algorithms such as pseudo-generalization, sequential convex optimization and heuristic algorithms, but these methods are all less effective. Summary of the invention

[0006] The present invention aims at non-cooperative target observation tasks with large navigation information errors and low update frequency, and thus proposes a new spacecraft observation method, namely, a trajectory planning method for active observation of non-cooperative targets with high security and high efficiency.

[0007] To achieve the above object, the present invention provides the following technical solutions:

[0008] The present invention proposes a trajectory planning method for active observation of non-cooperative targets with high security and efficiency, and the method comprises the following steps:

[0009] Step S1: determining the target navigation estimation state and the target navigation error covariance matrix according to the navigation information output by the non-cooperative target spacecraft;

[0010] Step S2: determining the state gain reachable set of the non-cooperative target spacecraft according to the target navigation estimated state and the target navigation error covariance matrix;

[0011] Step S3: using the state gain reachable set of the non-cooperative target spacecraft to perform size compensation on the structures of the non-cooperative target spacecraft and the service spacecraft, and determining the safety observation ellipsoids of the non-cooperative target spacecraft and the service spacecraft;

[0012] Step S4: During the process of actively observing non-cooperative targets, the service spacecraft must avoid entering the safe observation ellipsoid;

[0013] Step S5: Select the diagonal plane of the non-cooperative target spacecraft as the planning plane of the observation trajectory, and use the four intersection points of the two diagonals in the planning plane and the safety observation ellipsoid as observation nodes, and form a reference closed trajectory that meets the observation task by sequentially connecting the four observation nodes;

[0014] Step S6: Select the discretized observation trajectory sequence and its corresponding control sequence in the reference closed trajectory as optimization variables, select the distance-energy mixed quadratic index as the optimization index, use the discretized inequality as the path constraint condition, and calculate the corresponding optimized trajectory sequence and the corresponding control sequence by using the trajectory planning solver.

[0015] Furthermore, the output navigation information is expressed in the form of a navigation error ellipsoid;

[0016] The navigation error ellipsoid is X 0 =ε(x 0 ,p -1 X 0 ), where p is the confidence level, and the center of the navigation error ellipsoid x 0 Represents the target navigation estimation state, the shape of the navigation error ellipsoid X 0 represents the target navigation error covariance matrix.

[0017] Furthermore, the target navigation estimated state x 0 and the target navigation error covariance matrix X 0 Substitute the following differential equation:

[0018]

[0019] And solve it to get the state gain reachable set;

[0020] Among them, A represents the coefficient matrix in the HCW dynamics equation, A ú Represents the transposed matrix of the A matrix.

[0021] Furthermore, the above state gain reachable set is expressed as:

[0022]

[0023] Among them, [t 0 ,t f ] is the exercise time range.

[0024] Furthermore, the above safety observation ellipsoid is expressed as:

[0025]

[0026] Among them, e At represents the state transition matrix for time interval t, represents the transpose of the state transfer matrix, and I represents the 6×6 identity matrix.

[0027] Furthermore, the expression that the service spacecraft needs to avoid entering the safe observation ellipsoid is:

[0028]

[0029] Furthermore, the trajectory planning method is applicable to the scenario where the non-cooperative target and the service spacecraft are located near the same near-circular orbit, and the navigation and guidance system measures that the non-cooperative target is located within the error ellipsoid ε(x 0 ,X 0 )Inside.

[0030] The trajectory planning method for active observation of non-cooperative targets with high security and efficiency described in the present invention can be fully implemented using computer software. Therefore, correspondingly, the present invention also provides a trajectory planning system for active observation of non-cooperative targets with high security and efficiency. The system includes a storage device, which is used to execute the above-mentioned trajectory planning method and steps for active observation of non-cooperative targets with high security and efficiency.

[0031] The present invention also provides a computer-readable storage medium, on which a computer program is stored. When the computer program is executed by a processor, any one of the above-mentioned methods for trajectory planning with strong security and high efficiency for active observation of non-cooperative targets is executed.

[0032] The present invention also provides a computer device, which includes a memory and a processor, wherein the memory stores a computer program, and when the processor runs the computer program stored in the memory, the processor executes any one of the above-mentioned methods for trajectory planning for active observation of non-cooperative targets with high security and high efficiency.

[0033] The beneficial effects of the present invention are:

[0034] 1. The present invention proposes a trajectory planning method for active observation of non-cooperative targets with high efficiency and strong security, which can enhance the security of observation tasks and significantly shorten the observation task cycle. The principle is: the state gain reachable set is calculated from the initial navigation error through a matrix differential equation, and the spacecraft structure size is compensated for it. The boundary of the safety strategy is used as a path inequality constraint, and the optimization model with path fuel consumption and path length as composite performance indicators is established. The optimization variables are the observed path state and its corresponding control, and the trajectory planning model is established comprehensively, and then the model predictive control is used to solve the optimization problem for the problem to obtain the corresponding path state.

[0035] Furthermore, based on the state gain reachable set theory, the present invention proposes an analytical safety ellipsoid model, which effectively improves the safety of the collision avoidance strategy while ensuring real-time performance, and solves the problem that the traditional navigation error ellipsoid has poor timeliness and cannot describe the structural size.

[0036] Furthermore, compared with the traditional passive observation method, the active observation trajectory designed by the present invention realizes the comprehensive acquisition of target information in a single circle, solving the problems of low information acquisition efficiency and obvious task intention caused by the long cycle of the current passive observation task.

[0037] Furthermore, the applicable scenario of the present invention is: the non-cooperative target and the service spacecraft are located near the same near-circular orbit, and the navigation and guidance system measures that the non-cooperative target is located in the error ellipsoid ε(x 0 ,X 0 ) within a certain period of time, and ensure the active observation of non-cooperative targets in a shorter mission cycle while ensuring safety.

[0038] The present invention is applicable to trajectory planning for active observation of non-cooperative targets in aerospace orbit optimization. BRIEF DESCRIPTION OF THE DRAWINGS

[0039] Figure 1It is a flow chart of a trajectory planning method for active observation of non-cooperative targets with high efficiency and strong security according to the present invention;

[0040] Figure 2 It is a schematic diagram of the safe collision avoidance strategy based on the navigation error ellipsoid and its reachable set according to the present invention;

[0041] Figure 3 is a schematic diagram of the active observation trajectory of the non-cooperative target spacecraft according to the present invention;

[0042] Figure 4 is the non-cooperative observation simulation scenario parameter of the present invention;

[0043] Figure 5 It is the safe collision avoidance ellipsoid and path point distribution of the present invention, wherein Figure (a) is the drift initial value path point distribution, and Figure (b) is the periodic initial value path point distribution;

[0044] Figure 6 is a comparison of the computational efficiency of the security strategy described in the present invention;

[0045] Figure 7 It is the drift initial value non-cooperative target active observation trajectory described in the present invention;

[0046] Figure 8 Actively observe relative trajectories for drift initialization non-cooperative targets;

[0047] Fig. 9 Actively observe trajectories for non-cooperative targets with periodic initialization values;

[0048] Fig.10 Actively observe relative trajectories for periodic initialization non-cooperative targets;

[0049] Fig.11 Observe control sequences for non-cooperative targets. DETAILED DESCRIPTION

[0050] In the following description, specific details such as specific system structures and technologies are provided for the purpose of illustration rather than limitation, so as to provide a thorough understanding of the embodiments of the present application. However, it should be clear to those skilled in the art that the present application may also be implemented in other embodiments without these specific details. In other cases, detailed descriptions of well-known systems, devices, circuits, and methods are omitted to prevent unnecessary details from obstructing the description of the present application.

[0051] The specific embodiments of the present invention are further described in detail below in conjunction with the accompanying drawings. The following embodiments will help those skilled in the art to further understand the present invention, but do not limit the present invention in any form. It should be noted that, for those of ordinary skill in the art, several changes and improvements can be made without departing from the concept of the present invention, which all belong to the protection scope of the present invention.

[0052] Implementation method 1, see Figure 1 This embodiment is described. This embodiment proposes a new spacecraft observation method for non-cooperative target observation tasks with large navigation information errors and low update frequency, that is, a trajectory planning method for active observation of non-cooperative targets with strong security and high efficiency, which can enhance the security of observation tasks and significantly shorten the observation task cycle. The principle is: the state gain reachable set is calculated from the initial navigation error through a matrix differential equation, and the spacecraft structure size is compensated for it. The boundary of the safety strategy is used as a path inequality constraint, and the optimization model with path fuel consumption and path length as composite performance indicators is used. The optimization variables are the observed path state and its corresponding control, and a trajectory planning model is comprehensively established, and then the model predictive control is used to solve the optimization problem for the problem to obtain the corresponding path state.

[0053] like Figure 1 As shown, the path planning method includes the following steps:

[0054] Step S1: determining the target navigation estimation state and the target navigation error covariance matrix according to the navigation information output by the non-cooperative target spacecraft;

[0055] Step S2: determining a state gain reachable set of the non-cooperative target spacecraft according to the target navigation estimated state and the target navigation error covariance matrix;

[0056] Step S3: using the state gain reachable set of the non-cooperative target spacecraft to perform size compensation on the structures of the non-cooperative target spacecraft and the service spacecraft, and determining the safety observation ellipsoids of the non-cooperative target spacecraft and the service spacecraft;

[0057] Step S4: During the process of actively observing non-cooperative targets, the service spacecraft must avoid entering the safe observation ellipsoid;

[0058] Step S5: Select the diagonal plane of the non-cooperative target spacecraft as the planning plane of the observation trajectory, and use the four intersection points of the two diagonals in the planning plane and the safety observation ellipsoid as observation nodes, and form a reference closed trajectory that meets the observation task by sequentially connecting the four observation nodes;

[0059] S6: Select the discretized observation trajectory sequence [x (k) ;x (k+1);...;x (k+N) ] and its corresponding The control sequence [u (k) ;u (k+1) ;...;u (k+N) ] is used as the optimization variable, the distance-energy mixed quadratic index is selected as the optimization index, the discretized inequality is used as the path constraint, and the corresponding optimized trajectory sequence and the corresponding control sequence are calculated by using the trajectory planning solver.

[0060] Implementation method 2: This implementation method specifically describes a trajectory planning method for active observation of non-cooperative targets with high security and high efficiency proposed in the above implementation method 1.

[0061] First, the relative orbital dynamics model is described:

[0062] Both the service spacecraft and the non-cooperative target spacecraft are operating in circular orbits, and limited by the effective distance of the observation payload, the relative distance between the two is much smaller than the orbit radius. In order to describe the relative motion between the two, this embodiment chooses to establish an LVLH coordinate system with virtual points near the two as the origin. In this coordinate system, the relative motion between the two can be approximately described by the linearized HCW equation:

[0063]

[0064] Among them, u x (t),u y (t),u z (t) represents the acceleration component caused by all external forces except gravity; represents the average angular velocity of the orbit where the origin of the reference system is located, r c is the average orbital radius of the orbit where the origin of the coordinate system is located, and μ is the gravitational constant of the earth.

[0065] The above linear time-invariant system (1) can be rewritten into the following state space form:

[0066]

[0067] Where x(t) is the state vector, u(t) is the control vector, u(t)=[u x (t),u y (t),u z (t)] u′ ∈R 3 , A represents the coefficient matrix of the dynamic equation, and B represents the control matrix of the dynamic equation.

[0068] In the linear time-invariant system (2), given the initial state x(t 0) and control input u(t), the corresponding analytical solution is:

[0069]

[0070] Among them, the state transfer matrix is the homogeneous solution matrix of the above dynamic equation (2).

[0071] Assume the sampling period is Δt and adopt the zero-order hold assumption: u (k) =u(t), kΔt≤t<(k+1)Δt, the above formula (3) can be discretized into:

[0072] x (k+1) =A d x (k) +B d u (k) , k=1,...,N (4)

[0073] Among them, A d With B d They correspond to the discrete state matrix and control matrix of the continuous thrust spacecraft respectively:

[0074] A d =e AΔt ,

[0075] This implementation method establishes a relative orbital dynamics model, which makes it possible to describe the motion laws of all spacecraft as the operation basis of this implementation method. At the same time, the state gain reachable set designed subsequently is the result of the motion of the navigation error ellipsoid in the relative orbital dynamics model.

[0076] Then a safe collision avoidance strategy based on the navigation error ellipsoid and its reachable set is constructed;

[0077] The minimum safe area for non-cooperative target observation is usually determined by the navigation error sphere generated by the Kalman filter. The conservatism of this method of statically estimating the target state range will significantly decrease as the frequency of navigation information updates decreases. In order to overcome the limitations of traditional methods, this implementation method studies the expansion characteristics of the shape of the navigation error sphere under the influence of celestial dynamics based on matrix differential equations, and accurately gives the law of change of the navigation error sphere over time (state gain reachable set). The safety strategy based on the state gain reachable set has the advantages of not decreasing conservatism over time and effectively avoiding spacecraft structure collisions after expanding the size of the spacecraft, thereby significantly improving the overall safety of the observation mission.

[0078] Specifically:

[0079] The navigation and guidance system relies on the navigation information output by the Extended Kalman Filter (EKF), which is presented in the form of a navigation error ellipsoid. The navigation error ellipsoid consists of the following elements: 0 represents the estimated state of the error Gaussian distribution; the shape is given by the symmetric positive definite covariance matrix X 0 Description, used to quantify the uncertainty of the estimated state; and the corresponding confidence p. In other words, the navigation information indicates that the target state x(t 0 ) is located within the navigation error ellipsoid with a confidence probability p 0 =E(x 0 ,p - 1 X 0 ):

[0080] X 0 = {x(t 0 )∈R 6 |(x(t 0 )-x 0 ) ú X 0 -1 (x(t 0 )-x 0 )≤p} (5)

[0081] Under the influence of the dynamic system, the initial value of the target is estimated ellipsoid X 0 =E(x 0 ,p -1 X 0 ) will shift position and expand rapidly over time. When the navigation information is updated less frequently, the original navigation error ellipsoid X 0 This will seriously deviate from the position of the target's possible state set and seriously underestimate the size of the set, significantly increasing the risk of collision.

[0082] Therefore, the safe collision avoidance strategy with the state gain reachable set as the boundary in this implementation can effectively ensure that the trajectory of the service spacecraft avoids the center of mass of the target spacecraft, thereby significantly reducing the risk of center of mass collision. However, in order to eliminate the risk of structural collision as much as possible, it is also necessary to compensate the structural size of the state gain reachable set. Traditional reachable set calculations rely on numerical traversal to solve the boundary, which has low computational efficiency and is difficult to adapt to the application requirements of space platforms with scarce computing resources. This implementation is based on the Kalman filter to represent the navigation error information as a shape matrix and a center function, and focuses on studying and analytically solving the covariance matrix X through matrix differential equations. 0 The time-varying characteristics of the state gain are integrated with the central function into a state gain reachable set. This method realizes the real-time solution of the reachable set, meets the needs of practical applications on low-computing platforms, and significantly improves its practicality.

[0083] The state gain reachable set is defined as: given the motion time range [t 0 ,t f ], the system equations satisfy Initial state set X 0 The terminal state x(t f ) is called the state gain reachable set, which is defined as:

[0084]

[0085] For a non-cooperative target spacecraft, the navigation system estimates its initial state set as X 0 =E(x 0 ,p -1 X 0 ), whose shape function X 0 is a symmetric positive definite matrix. Since a symmetric positive definite matrix is ​​still a symmetric positive definite matrix after linear transformation, the shape function X after the integral accumulation transformation is 0 (t) is still a symmetric positive definite matrix.

[0086] The above properties ensure that the state gain can reach the set R(t f ,X 0 ) is still an ellipsoid set, and its shape function X 0 (t) and the central function x 0 (t) Under the action of gravity, the following differential equations are satisfied:

[0087]

[0088] When the initial condition is X 0 (t 0 )=X 0 , x 0 (t 0 )=x 0 When , the above differential equation (7) can be directly solved analytically by integration:

[0089]

[0090] Therefore, the specific form of the state gain reachable set is as follows:

[0091]

[0092] When the sampling period is Δt, the discrete step length is k∈Z + When , the discrete state gain reachable set is:

[0093]

[0094] in, Therefore, the state gain reachable set can be simply written as:

[0095]

[0096] The state gain reachable set proposed in this implementation solves the problem that the navigation error ellipsoid estimates the target center of mass state and becomes increasingly distorted over time, and accurately gives the upper limit of the range of the target state at any time. However, due to the large physical structure size of the spacecraft, relying solely on the reachable set for collision judgment may lead to the risk of structural collision. Therefore, it is necessary to correct the structural size of the reachable set boundary. Define the size correction coefficient as ρ P +ρ I +ρ, where ρ P and ρ I Respectively represent the structural envelope of the service spacecraft and the target spacecraft, and ρ is the fault tolerance distance. By correcting the structural size of the reachable set boundary, the potential risk area is effectively covered, thereby significantly improving the conservatism and accuracy of the collision avoidance strategy. Figure 2 shown.

[0097] The result of extending the sphere outside the boundary of the ellipsoid is still an ellipsoid. Although this process is relatively simple in concept, there is currently a lack of mathematical tools that can analytically represent this process. In order to solve the extended ellipsoid in real time, this implementation method appropriately sacrifices some accuracy while ensuring conservatism, and approximates a slightly more conservative analytical solution. To this end, this embodiment introduces an auxiliary function δ * , which is used to describe any direction vector l∈R 6 The distance from the intersection point with the ellipsoid shell E(q,Q) to the center of the ellipsoid:

[0098] δ * (l,E(q,Q))=(l,q)+(l,Ql) 1 / 2 (12)

[0099] Where (·,·) represents the vector inner product. Auxiliary function δ * (l,E(q,Q)) describes the geometric properties of the ellipsoid in direction l. Figure 2 As shown, the distance difference Δδ between the outer extension ellipsoid and the outer shell of the state gain reachable set in any direction is * Both are ρ P +ρ I +ρ:

[0100]

[0101] in, are the center and shape matrices of the externalized ellipsoid. The externalization process of the ellipsoid does not change the center position, that is: Therefore, the distance difference Δδ *is the difference in shell size between the old and new ellipsoids caused by the shape matrix:

[0102]

[0103] Where I∈R 6×6 is the identity matrix. Since It is necessary to satisfy formula (14) in any direction l. Solving the matrix directly is very complicated and inefficient. 1 / 2 is a strictly concave function, so we can use the Jensen inequality Scaling the above inequality:

[0104]

[0105] Therefore, the shape function of the collision avoidance ellipsoid can be approximated into a more conservative analytical form that can be calculated in real time.

[0106]

[0107] The analytical approximate solution The error tolerance distance ρ in is set as the optimization variable, and the accuracy is improved by reducing the approximation error Δd:

[0108]

[0109] The optimal fault tolerance distance ρ is obtained * Finally, taking into account the timeliness of navigation information, the size of the spacecraft structure and the real-time nature of the calculation, this implementation method designs the following safe collision avoidance strategy:

[0110]

[0111] According to the above safety collision avoidance strategy, during the active observation of non-cooperative targets, the service spacecraft needs to avoid entering the safety observation ellipsoid; it can be expressed as:

[0112]

[0113] Finally, spacecraft active observation safety trajectory planning based on fast model predictive control (FMPC);

[0114] Although the above-mentioned safe collision avoidance strategy (17) effectively characterizes all possible states and structural dimensions of non-cooperative targets, clearly gives the trajectory forbidden zone in the observation mission, and directly applies it to the trajectory planning model as a path constraint, it effectively reduces the potential collision risk. However, in the observation mission of non-cooperative spacecraft, in addition to considering safety factors such as collision, it is also necessary to improve the mission concealment to avoid the disclosure of one's own intentions as much as possible. The current space situational awareness relies on ground orbit determination information. Once the non-cooperative party solves the six orbital numbers of our service spacecraft, it can infer our mission intentions and observation objects, causing us to fall into a very passive position. Therefore, the traditional free companion flight observation method that requires at least one orbital cycle obviously cannot meet the concealment requirements and is not suitable for non-cooperative target observation tasks with high concealment requirements. To solve this problem, this implementation method designs an observation trajectory that can obtain all the key information of the target in a single circle, and tracks the trajectory by active thrust, which significantly shortens the mission cycle. Finally, based on the fast model predictive control (FMPC) strategy (the computational complexity is approximately O(N 3 ), where N is the total discrete step length), realizes the rapid solution of high-complexity trajectory planning and provides an efficient and reliable solution for the covert observation of non-cooperative targets.

[0115] Specifically:

[0116] The active observation trajectory of non-cooperative target spacecraft needs to meet multiple complex constraints, including nonlinear inequality path constraints for dynamically avoiding the safe collision avoidance ellipsoid, equation path constraints for reaching a specific observation position, and model constraints for satisfying conditions such as dynamic models and control amplitudes in real time. In order to achieve accurate identification of non-cooperative targets, the sensor field of view of the service spacecraft is required to completely cover the target, and its focal length must meet the recognition accuracy requirements, which means that the observation trajectory must be stably maintained near the target trajectory. Generally, the semi-major axis of the navigation error ellipsoid of high-value orbit satellites such as GEO measured on the ground is greater than or equal to 50m, while the diameter of the structural dimensions of the traditional satellite body and key payloads does not exceed 6m. According to a simple estimate, if the field of view angle of the service spacecraft sensor is ≥3.5°, it can completely cover the target satellite body and its key payloads.

[0117] For traditional cube satellites, a large field of view sensor can cover the three body planes of the target at one time by aligning it with its body diagonal. Based on this, this implementation method selects the body diagonal plane of the non-cooperative target as the planning plane of the observation trajectory, and designs the four intersection points of the two diagonals in the plane and the safe collision avoidance ellipsoid as observation nodes. By sequentially connecting these four observation nodes, a reference closed trajectory that meets the observation task is formed. Along the mission trajectory generated by the iterative planning of the initial observation trajectory, the service spacecraft can fully obtain the target's key payload and platform sensor information during a single orbit. The specific observation schematic diagram is shown below. Figure 3 shown.

[0118] Assume that the attitude angle of the non-cooperative target spacecraft relative to the LVLH coordinate system is (The rotation order is zyx), then the direction vector of the diagonal of the spacecraft body is:

[0119]

[0120] in, is the direction vector of the diagonal of the cube body, and and Located on two diagonal surfaces:

[0121]

[0122] The above completely derives and constructs four key observation nodes. The service spacecraft needs to make the sensor optical axis coincide with the diagonal direction of the target body at these four nodes to fully obtain all the information of the non-cooperative target body and key payloads. By combining the safety ellipsoid equation (17) and the two diagonal direction vectors p in the diagonal plane of the target body i The specific coordinates of these four key observation nodes can be solved using the quadratic programming method.

[0123] according to Figure 4 The scene parameters listed in the plot are obtained Figure 5 The safe collision avoidance ellipsoid and path point distribution, Figure 5 (a) and (b) show different navigation error initial value matrices X 0 (t 0 ) The difference in initial values ​​leads to changes in the shape and size of the safety ellipsoid, which directly affects the position distribution of the four key observation nodes and further affects the expected path x in the secondary index. p In addition, the safe collision avoidance ellipsoid as a trajectory restricted area directly affects the coefficient matrix in the inequality path constraint The value of .

[0124] Therefore, the initial value trajectory is designed as a closed curve that sequentially connects the four observation nodes. In order to track the above desired state and minimize energy consumption at the same time, the optimization objective function is selected as follows:

[0125]

[0126] Among them, l:R 9 →R is a convex quadratic index function, the specific form is:

[0127]

[0128] In the indicator function (21), Q∈R 6×6 , R∈R 3×3 are all symmetric positive definite matrices, x p (t) is the expected trajectory, u(t)∈R 3 is the control vector.

[0129] In addition, avoiding structural collisions between spacecraft during the observation process is a constraint that must be considered in path planning. The above fully considers the constraints such as information update cycle and spacecraft structure size, and constructs a sufficiently conservative safety collision avoidance ellipsoid (17). Therefore, the observation path planning should satisfy the inequality constraint of avoiding the safety ellipsoid at each discrete time k:

[0130]

[0131] Among them, x (k) To service the current status of the spacecraft, and They are the center state and shape matrix of the safe collision avoidance ellipsoid:

[0132]

[0133] The nonlinear and non-convex inequality constraint (22) is difficult to handle in trajectory planning. This implementation method linearizes the constraint by a first-order Taylor expansion, which not only ensures the linearity of the simplified constraint but also realizes the convexity of the problem, thereby significantly reducing the computational complexity. The constraint function (22) is denoted as: By using the result of the current iteration Performing a first-order Taylor expansion at , the original constraint can be convexified into the following linear constraint:

[0134]

[0135] The direct linearization of non-convex constraints ignores high-order small quantities, which inevitably leads to some states that meet constraint (24) entering the collision avoidance ellipsoid. In order to reduce the impact of linearization error on safety, the corresponding safety margin ρ has been designed above to compensate for the error. Combined with formulas (22) and (23), the safety constraint after linear simplification can be expressed as:

[0136]

[0137] Constraint (25) can be simplified into matrix form:

[0138]

[0139] in:

[0140]

[0141] In addition, since the thrust of the service spacecraft is limited in amplitude, the control input needs to satisfy the following inequality constraints:

[0142] u min ≤u (k) ≤u max , k=0,1,...,N (27)

[0143] For ease of calculation, the above constraints can be converted into matrix form:

[0144] F u u (k) ≤f u (28)

[0145] in:

[0146]

[0147] Combining the above task analysis and constraint transformation process, the trajectory planning problem of actively observing non-cooperative targets can be summarized as the following optimization problem:

[0148]

[0149] The optimization variable finally solved for the above optimization problem is the trajectory sequence x (1) ,K,x (N) and the control sequence u (0) ,K,u (N -1) , N is the prediction time domain. The quadratic programming problem is defined by the following parameters:

[0150] x p (k) ,A d ,B d ,Q,R,F x ,F u ,f x ,f u (30)

[0151] Finally, the trajectory planning solver based on fast model predictive control (FMPC) is used to calculate the corresponding optimized trajectory sequence and the corresponding control sequence.

[0152] Implementation method 3: This implementation method specifically describes the corresponding optimized trajectory sequence and the corresponding control sequence calculated by the trajectory planning solver based on fast model predictive control (FMPC) described in the above implementation method 2;

[0153] The convex quadratic programming problem (29) proposed in the above-mentioned second embodiment is a classic optimization problem. The fast model predictive control (FMPC) method adopted in this embodiment is an efficient explicit MPC, and its problem solver is an interior point method - the original barrier method.

[0154] First, define the optimization variable sequence z and the expected sequence z p for:

[0155]

[0156] With the help of the above definition, the original problem can be rewritten into a compact matrix form:

[0157]

[0158] Among them, H is the quadratic term matrix of the indicator function; P and h (C and b) correspond to the parameter matrices of the inequality constraints (equality constraints) respectively. The specific matrix structures are as follows:

[0159]

[0160]

[0161] The original barrier method embeds the inequality constraints in the optimization problem (32) into the objective function by introducing a barrier function, thereby transforming the original problem into an optimization problem containing only equality constraints without changing the KKT optimality condition. The new optimization problem can be efficiently solved by constructing a residual sequence and using the Newton iteration method. Therefore, the optimization problem (32) can be equivalently transformed into:

[0162]

[0163] The barrier function φ has the following expression:

[0164]

[0165] in, is the i-th row vector of the inequality constraint coefficient matrix P.

[0166] In the process of iteratively solving the new optimization problem (33), the weight coefficient γ of the barrier function will be proportionally reduced to 1 / 10 of the original value in each iteration, gradually reducing the degree of violation of the inequality constraint. Finally, when γ→0, the solution of the new optimization problem (33) will gradually converge to the optimal solution of the original optimization problem (32).

[0167] The active observation trajectory planning problem for non-cooperative targets (32) is a complex constrained optimization problem with many constraints and high dimensions. It is almost impossible to directly construct the initial trajectory and control sequence in the feasible domain. Although the standard Newton iteration method has a very high solution efficiency, it requires that the initial value of the iteration strictly satisfies all constraints, which puts strict requirements on the construction of the initial value. In contrast, the initial point infeasible Newton method alleviates the difficulty of guessing the initial value of the iteration by relaxing the strict condition that the initial value must satisfy the equality constraint, while retaining the high computational efficiency of the Newton iteration method.

[0168] To use this solver to quickly solve the optimization problem (33), we first need to define the Lagrangian function L(z,ν) corresponding to the problem:

[0169] L(z,ν)=(zz p )H(zz p )+γφ(z)+ν(Cz-b) (35)

[0170] Where ν∈R 6N is the Lagrange multiplier corresponding to the equality constraint Cz-b = 0. From the minimum-maximum principle, it can be seen that the zeros of the two partial differential equations of Lagrange equation (35) about (z, ν) are necessary conditions for the optimal solution:

[0171]

[0172] Among them, γP ú d is the gradient of the implicit expression γφ(z) of the inequality constraint Pz<h; vector d∈R 9N+6 The specific expression is: P i represents the i-th row of the inequality coefficient matrix P; (z * ,ν * ) are the optimal independent variable sequence and the optimal dual variable sequence respectively.

[0173] Let equation group (36) be F(z * ,ν * )=0, the general solution to this nonlinear system of equations is the Newton iteration method, which selects the initial value in the solution neighborhood and gradually reduces the residual of the function F(z,ν) so that the approximate solution (z k ,ν k ) gradually converges to the true solution (z * ,ν * ).

[0174] To solve the equation group (36) using the Newton iteration method, we first need to k ,ν k ), the function F(z,ν) is expanded by the first order Taylor:

[0175]

[0176] Among them, ▽F(z k ,ν k ) means that the function F(z,ν) is k ,ν k ) is the gradient matrix at the position. In order to make the equation gradually approach the zero residual condition, take F(z k+1 ,ν k+1 )=0, then the residual and update direction satisfy:

[0177]

[0178] Among them, Δz and Δν are the update directions of the state and its duality respectively; F(z k ,ν k ) is the residual, and its component dual residual The original residual The specific expression is as follows:

[0179]

[0180] Calculate the coefficient matrix ▽F(z k ,ν k ) and put it into equation (38), and we get the update equation of the k-th Newton iteration:

[0181]

[0182] Among them, γP ú diag(d) 2 P is the Hessian matrix of the implicit inequality constraint γφ(z); Δz and Δν are the state update direction and the dual update direction, respectively.

[0183] Therefore, the updated iteration state (z k+1 ,ν k+1 )satisfy:

[0184] z k+1 =z k +sΔz,ν k+1 =v k +sΔν (41)

[0185] Different from the traditional Newton method, the initial point is not feasible. The Newton method introduces a weight coefficient s∈(0,1] in the state update process to ensure that the updated result satisfies the inequality constraint Pz k+1 <h.

[0186] From the process of constructing the infeasible Newton method for the initial point, we can know that the initial value of the state variable z 0 Need to strictly satisfy the inequality constraint Pz0 <h; and the equality constraint Cz k =b is corrected by the dual residual The iterations are completed step by step until the constraints are fully satisfied, so the initial value z 0 There is no need to satisfy the equality constraints. Similarly, the initial value of the dual variable ν 0 It can also be chosen arbitrarily. Finally, as the residual continues to decrease, the iterative solution (z k ,ν k ) will gradually approach the optimal solution within the error threshold (z * ,ν * ).

[0187] The solution of the iterative update direction (40) is the most computationally complex part of the Newton method with an infeasible initial point. Usually, solving the linear equation system (40) requires calculating the coefficient matrix ▽F(z k ,ν k ) ú DL decomposition and then solve its inverse matrix. ú The computational complexity of DL decomposition is high, which will significantly affect the computational efficiency of the initial point infeasible Newton method.

[0188] Compared with traditional explicit MPC, the efficiency of FMPC is mainly reflected in its ability to quickly solve high-dimensional linear equations (40). The coefficient matrix of equation (40) is a symmetric positive definite block diagonal matrix, which can be equivalently reduced to a low-dimensional matrix through Schur decomposition, thereby significantly reducing the dimension and computational complexity of the equation. The new equations after reduction still maintain a symmetric positive definite coefficient matrix, and further adopt the equation solution method combining the Clausky decomposition method with the matrix block elimination method. The computational efficiency is much higher than the traditional L ú The solution of the system of equations for matrix inversion after DL decomposition. Through the above simplification process, the dimension and solution complexity of the system of equations are greatly reduced, which significantly improves the overall solution efficiency.

[0189] First, this embodiment performs equivalent dimension reduction on the coefficient matrix of the equation group (40) by Schur decomposition, so that the dimension of the new equation group is the same as the submatrix Φ = 2H + γP ú diag(d) 2 The submatrix Φ is a block symmetric positive definite matrix, and its inverse matrix is:

[0190]

[0191] The original equation group (40) can be equivalently reduced in dimension and decoupled into equations only about the dual update direction Δν by Schur decomposition:

[0192] CΦ -1 C ú Δν=r p -CΦ-1 r d (43)

[0193] The coefficient matrix of the linear equation system (43) is Y = CΦ -1 C ú , the constant vector is β = -r p +CΦ -1 r d Therefore, it is easy to solve the root of the equation system Δν=-Y -1 β. The solution of the traditional linear equation system depends on the LDL of the coefficient matrix Y ú Decompose and then solve the inverse matrix Y -1 Obtain the solution to the equation; in contrast, using the Clausky LL ú After decomposing the coefficient matrix Y, the system of equations is solved by two eliminations, and the latter has more advantages in computational complexity.

[0194] After Clausky decomposition, the coefficient matrix Y can be decomposed into a lower triangular matrix L, satisfying the equation Y = LL ú The specific form of the matrix L is as follows:

[0195]

[0196] Each block sub-matrix of the matrix L can be calculated by recursive decomposition. The specific steps are as follows:

[0197] 1. Decompose the first element Y of the main diagonal 11 Solving for L 11 :

[0198] 2. Solve for the off-diagonal elements L 21 :

[0199] 3. Decomposition of the main diagonal elements To solve L 22 : Y 22 -L 21 L 21 =L 22 L 22 ;

[0200] 4. Follow the above process recursively until the entire matrix is ​​decomposed.

[0201] The above iterative calculation process can be summarized as the following formula:

[0202]

[0203] By performing Clausky decomposition on the coefficient matrix Y, the original equation YΔν=-β is equivalently transformed into the equation of two triangular coefficient matrices: L úΔν=y and Ly=-β. The triangular coefficient matrix equations obtained after decomposition can be efficiently solved for Δν using the matrix block elimination method. Then, Δν is substituted back into the original equation (40) to obtain the expression for Δz: Δz=Φ -1 (-r d -C ú Δν).

[0204] Implementation method 4: See Figures 6 to 11 This embodiment is described. This embodiment is to verify the simulation effect of a trajectory planning method for active observation of a non-cooperative target with high efficiency and strong security described in the above embodiment.

[0205] The non-cooperative target to be observed is in the GEO orbit, and the average angular velocity of the orbit is n = 7.277 × 10 -5 rad / s. The navigation information of the target determined by the ground orbit determination and navigation system includes the estimated state x 0 (t 0 )(mm / s) and its corresponding 98% confidence error covariance X 0 (t 0 )(m 2 m 2 / s 2 ), the initial state of our service spacecraft x(t 0 )(mm / s) and other scene information are summarized as follows Figure 4 shown.

[0206] The relative operation state of a spacecraft on orbit can be divided into two typical relative trajectories according to the solutions of the HCW equation: periodic solutions and drift solutions. This implementation method designs two representative simulation scenarios based on this classification. Specifically, the periodic solution corresponds to a periodic change in the relative orbit of the spacecraft, and the trajectory is an elliptical or circular closed curve, and its initial value must satisfy The drift solution corresponds to the spacecraft orbit gradually deviating from the initial state, and the trajectory appears as a spiral curve.

[0207] The range of the target navigation covariance is usually determined by the accuracy of the navigation system. For example, taking the GPS navigation system as a reference, its ground positioning accuracy is about 10 meters and its speed measurement accuracy is about 0.2m / s. The GPS navigation accuracy is affected by the target distance. The farther the distance, the lower the navigation accuracy. Usually, in the orbit determination mission of non-cooperative targets in space, the navigation positioning error is between 10 and 100m, and the speed measurement error is between 0.1 and 1m / s. For this article, the simulation scenario is designed on the GEO orbit, and the scenario scale is at the kilometer level, so the positioning error is about 50m and the speed measurement error is about 0.25m / s.

[0208] In addition, the envelope size of the satellite and its attached structures generally does not exceed 10 meters. This implementation assumes that the size of the target and service spacecraft is the upper limit value ρ P =ρ I =10m. There is no special range for the value of the fault tolerance distance ρ, which is selected as Figure 6 At the same time, it is assumed that the attitude angle of the non-cooperative target relative to the LVLH coordinate system is The rotation order is zyx.

[0209] The model predictive control used requires solving the discrete optimization problem (30), and setting the prediction time domain N = 50 and the discrete time step Δt = 20s. Discrete dynamic parameter A d ,B d Calculate according to the above formula (4). In the explicit model predictive control (32): the indicator matrix H∈R N(6+3)×N(6+3) The state and control weight sub-matrix Q = I 6×6 , R=I 3×3 ; Expected trajectory z p Calculate according to formula (31), where the expected observation trajectory x p is a curve that sequentially links each observation node and fits the safe collision avoidance ellipsoid; the inequality constraint coefficient matrix P, the state and control constraint submatrix F in h x ,f x Calculated according to formula (26), it is the result of the current iteration Function of u ,f u Calculated by formula (28), the control constraint u max =-u min =1m / s 2 The initial value of the weighted coefficient of the barrier function is γ = 1, the maximum number of iterations is 1000, and the termination error threshold is 1e-8. The preset time of the safe collision avoidance ellipsoid is selected as the time interval between the four observation nodes N / 5·Δt = 200s. The specific simulation results are as follows Figures 7 to 11 As shown; among them, Figure 7 and Fig. 9 The absolute observation trajectories under two typical scenarios are plotted respectively. Figure 8 and Fig.10 The corresponding relative observation trajectories are shown. Fig.11The non-cooperative target observation control sequence is demonstrated. The absolute observation trajectory is presented as a curve coiled in the tubular safety area drawn by the safe collision avoidance ellipsoid at different times; the corresponding relative observation trajectory is a closed curve surrounding the safe collision avoidance ellipsoid. Different safe collision avoidance ellipsoids will affect the shape of the absolute observation trajectory by affecting the "radius" of the tubular safety area; at the same time, they will directly affect the minimum circling radius of the relative observation trajectory. This implementation method is based on the trajectory planning strategy proposed by FMPC and can solve the observation trajectory that does not pass through the safe collision avoidance ellipsoid, so that the service spacecraft can achieve circling observation of non-cooperative targets in a covert, efficient and safe manner.

[0210] In summary, the present invention proposes a trajectory planning method for active observation of non-cooperative targets with high efficiency and strong security. The task time of observing non-cooperative targets is shortened from 1 day to 1000 seconds without entering the safety ellipsoid, which greatly shortens the task cycle while having good security.

[0211] The above description is only the implementation mode of the present invention and is not limited to the present invention. For those skilled in the art, the present invention may have various modifications and variations. Any modification, equivalent substitution, improvement, etc. made within the spirit and principle of the present invention shall be included in the scope of the claims of the present invention.

Claims

1. A highly secure and efficient trajectory planning method for active observation of non-cooperative targets, characterized in that: The method is: S1: Determine the target navigation estimation state and the target navigation error covariance matrix based on the navigation information output by the non-cooperative target spacecraft; S2: Determine the state gain reachable set of the non-cooperative target spacecraft according to the target navigation estimated state and the target navigation error covariance matrix; S3: Use the state gain reachable set of the non-cooperative target spacecraft to compensate the size of the structure of the non-cooperative target spacecraft and the service spacecraft, and determine the safe observation ellipsoid of the non-cooperative target spacecraft and the service spacecraft; S4: During active observation of non-cooperative targets, the service spacecraft must avoid entering the safe observation ellipsoid; S5: Select the diagonal plane of the non-cooperative target spacecraft as the planning plane of the observation trajectory, and use the four intersection points of the two diagonals in the planning plane and the safety observation ellipsoid as observation nodes, and form a reference closed trajectory that meets the observation task by sequentially connecting the four observation nodes; S6: Select the discretized observation trajectory sequence and its corresponding control sequence in the reference closed trajectory as the optimization variables, select the distance-energy mixed quadratic index as the optimization index, use the discretized inequality as the path constraint condition, and calculate the corresponding optimized trajectory sequence and the corresponding control sequence by using the trajectory planning solver.

2. According to claim 1, a trajectory planning method for active observation of non-cooperative targets with high security and high efficiency is characterized in that: The output navigation information is expressed in the form of navigation error ellipsoid; The navigation error ellipsoid is X0=ε(x0,p -1 X0), where p is the confidence, the center of the navigation error ellipsoid x0 represents the target navigation estimation state, and the shape of the navigation error ellipsoid X0 represents the target navigation error covariance matrix.

3. According to claim 2, a trajectory planning method for active observation of non-cooperative targets with high security and high efficiency is characterized in that: Substitute the target navigation estimated state x0 and the target navigation error covariance matrix X0 into the following differential equation: And solve it to get the state gain reachable set; Among them, A represents the coefficient matrix in the HCW dynamics equation, A ú Represents the transposed matrix of the A matrix.

4. According to claim 2, a trajectory planning method for active observation of non-cooperative targets with high security and high efficiency is characterized in that: The state gain reachable set is expressed as: Among them, [t0,t f ] is the exercise time range.

5. According to claim 3, a trajectory planning method for active observation of non-cooperative targets with high security and high efficiency is characterized in that: The safety observation ellipsoid is expressed as: Among them, e At represents the state transition matrix for time interval t, represents the transpose of the state transfer matrix, and I represents the 6×6 identity matrix.

6. According to claim 4, a trajectory planning method for active observation of non-cooperative targets with high security and high efficiency is characterized in that: The expression that the service spacecraft needs to avoid entering the safe observation ellipsoid is:

7. According to claim 1, a trajectory planning method for active observation of non-cooperative targets with high security and high efficiency is characterized in that: The applicable scenario of the trajectory planning method is: the non-cooperative target and the service spacecraft are located near the same near-circular orbit, and the navigation and guidance system measures that the non-cooperative target is located within the error ellipsoid ε(x0,X0) with a confidence level of p = 98%.

8. A highly secure and efficient trajectory planning system for active observation of non-cooperative targets, characterized in that: The system comprises a storage device for executing the method and steps described in claim 1.

9. A computer-readable storage medium, characterized in that: The computer-readable storage medium stores a computer program, which, when executed by a processor, executes a trajectory planning method for active observation of a non-cooperative target with high efficiency and strong security as described in any one of claims 1 to 7.

10. A computer device, characterized in that: The device includes a memory and a processor, wherein a computer program is stored in the memory. When the processor runs the computer program stored in the memory, the processor executes a trajectory planning method for active observation of a non-cooperative target with high efficiency and strong security as described in any one of claims 1 to 7.

Images