Ultra-agile earth observation satellite full-precision model scheduling method and system

Abstract

The application discloses a hypersensitive satellite scheduling method, comprising the following steps: step one, accurately calculating adjacent task attitude switching time through Euler angles and rotation matrix; step two, searching task insertion time based on adjacent task attitude switching satisfaction deviation value, and realizing local search algorithm of hypersensitive satellite scheduling; step three, realizing optimal path exchange and reservation crossover and mutation operators based on discrete integer coding and minimum graph decoding, and further realizing global search evolution algorithm of hypersensitive satellite scheduling; and step four, realizing efficient solution of hypersensitive satellite full-precision model scheduling based on the combination of step two and step three.

Description

Technical Field

[0001] This invention relates to the field of satellite resource allocation and scheduling technology, and in particular to a method and system for scheduling ultra-agile Earth observation satellites with full-precision models. Background Technology

[0002] The agile satellite mission planning system is the control center for resolving conflict constraints, optimizing scheduling and operation, meeting user preferences, scientifically allocating telemetry resources, and maximizing satellite resource utilization and mission fulfillment rates.

[0003] Agile Earth observation satellites can be categorized into semi-agile and ultra-agile satellites based on their attitude angle maneuverability. Semi-agile satellites primarily adjust their pitch within a certain range, while their side-slip angle remains almost fixed and their yaw angle approaches 0°. In contrast, ultra-agile satellites offer significantly greater freedom in attitude adjustment, with unrestricted variations in their three attitude angles, further enhancing their observation capabilities. With the advancement of aerospace and satellite technologies, agile Earth observation satellites are increasingly becoming a focus of research and development for various countries.

[0004] The flexible and maneuverable hardware performance of ultra-agile satellites manifests itself in two ways: firstly, more and longer visible time windows, increasing the chances of mission fulfillment; secondly, it determines that the time-dependent attitude switching time is non-linear, complicating the scheduling problem. The intelligent mission planning of ultra-agile satellites aims to maximize the sum of the benefits of all scheduled missions. The solution process requires selecting the specific missions to be scheduled and determining the start time of the observation mission. Therefore, the intelligent mission planning problem of ultra-agile satellites is a mixed integer programming problem.

[0005] Currently, agile satellite mission planning mainly focuses on research into semi-agile satellite mission planning that satisfies specific properties (FIFO property). It is primarily implemented using iterative heuristic algorithms, mainly based on preprocessing by determining the upper and lower bounds of the observation start time using a bisection method, and using permutation encoding based on the upper and lower bounds of the start time for mission insertion to construct and improve the quality of solutions. The aim is to transform the start time within the search time window into a time determined by the upper and lower bounds based on the mission insertion order, so that integer solutions within the same sequence that fall between the upper and lower bound solutions correspond to a single scheduling scheme. This transforms the integer programming problem into a problem of searching for the insertion order, thus reducing the search space. However, this simplified model comes at the cost of scheduling accuracy: firstly, limiting the range of attitude angle changes reduces the time window; secondly, dividing the problem into subproblems on each orbit and limiting the observation mission to at most one visible time window on any orbital orbit reduces the computational scale, but loses some time windows at the boundary points between orbital orbits (or adjacent orbital orbits); finally, simplifying the attitude switching time and approximating it as a linear piecewise function of the sum of three angle changes to satisfy the FIFO property increases the switching time and reduces the opportunity for mission insertion.

[0006] For full-precision ultra-agile satellite mission planning without model simplification, existing encoding / decoding methods based on the FIFO property (bisection method for determining the upper and lower bounds of start time) and iterative heuristic search algorithms are unusable. To address this issue, a technique based on minimum graph optimization to find the optimal observation time series has been proposed. This technique is based on a continuous encoding and minimum graph model, enabling continuous evolutionary algorithms to solve the ultra-agile satellite mission planning problem. However, the above methods still cannot perform more accurate and effective local searches of scheduling schemes, resulting in insufficient optimization capabilities. Furthermore, the mapping from discrete integers to continuous encoding not only increases the search space but also weakens the interpretability of related operators and the preservation of constraint satisfaction. In addition, the mutual conversion of encodings is required, leading to low solution generation quality and high time overhead. Therefore, this invention proposes a full-precision model scheduling method and system for ultra-agile Earth observation satellites to solve the problems existing in the prior art. Summary of the Invention

[0007] To address the aforementioned problems, the present invention aims to propose a scheduling method and system for ultra-agile Earth observation satellites with full-precision model. This method and system solves the scheduling optimization problem of ultra-agile Earth observation satellites under a full-precision model without any simplification of the scheduling model.

[0008] To achieve the objectives of this invention, the invention is implemented through the following technical solution: a method for scheduling a full-precision model of ultra-agile Earth observation satellites, comprising the following steps:

[0009] Step 1: Accurately calculate the attitude switching time between adjacent tasks using Euler angles and rotation matrices;

[0010] Step 2: Search for the task insertion time based on the attitude switching satisfaction deviation value of adjacent tasks on both sides, and simultaneously implement the local search algorithm for hypersensitive satellite scheduling;

[0011] Step 3: Based on discrete integer encoding and minimum graph decoding, implement crossover and mutation operators for optimal path exchange and retention, and further implement a global search evolutionary algorithm for ultrasensitive satellite scheduling;

[0012] Step 4: Based on the meme algorithm combining Step 2 and Step 3, achieve high-efficiency solution for scheduling of the hypersensitive satellite full-precision model.

[0013] In step one, the attitude angle of the staring target is converted into a unit vector pointing from the satellite to the target. The specific method is as follows:

[0014] At a certain moment t, the satellite j The attitude angle of staring (pointing) at target j is represented by the vector Θ. j (t j )=(αj (t j ),β j (t j ),γ j (t j )) indicates that α j (t j ) is the lateral swing angle, β j (t j ) is the pitch angle, γ j (t j The yaw angle is 0. At this time, the coordinates of target j in the VVLH satellite orbital coordinate system are C. j (t j )=(x j (t j ),y j (t j ),z j (t j )), then Θ j (t j ) and C j (t j The conversion relationship between ) is shown in equation (1).

[0015]

[0016] The satellite's center of mass points to the unit vector U of target j. j (t j ) =(u x j (t j ),u y j (t j ),u z j (t j According to equation (1), u z j (t j ), u x j (t j ) and u y j (t j The following can be calculated using formulas (2)-(4):

[0017]

[0018]

[0019] The further calculations after the transformation involve, in step one, solving for the angles of rotation around the three degrees of freedom (i.e., Euler angles) and the attitude adjustment time when observing attitude adjustments between adjacent tasks. The specific method is as follows:

[0020] Given any two adjacent observation tasks p and s, where task p is the preceding task and task s is the succeeding task, the unit vectors pointing to the adjacent targets when the satellite is staring at a certain moment are denoted as U. p (t p ) and U s (t s ), can be determined by observing the attitude angle Θ p (t p ) and Θ s (t s The observation attitude angle is calculated according to the methods corresponding to equations (1)-(4), and is given and known in the form of a calculation example based on satellite parameters. The rotation axis K can be obtained from U. p (t p ) and U s (t s The cross product of vector U and vector U is calculated as follows: k = K / ||K||. p (t p ) to vector U s (t s The rotation angle φ = arccos(U) p (t p )∙ U s (t s Then the rotation matrix R can be calculated by the Rodrigues rotation formula, as shown in equation (5), where I is the identity matrix and [k] is the cross product matrix of k, given by equation (6).

[0021]

[0022]

[0023] The satellite attitude adjustment defaults to rotating in the order of Z, Y, X axes, with rotation angles denoted as θ. z , θ y , θ x The corresponding rotation matrix is ​​R. z , R y , R x Calculate according to formulas (7)-(9).

[0024]

[0025] Then we can construct equation (10):

[0026]

[0027] Among them, s x =sinθ x c x =cosθ x s y =sinθ y c y =cosθ y s z =sinθ z c z =cosθ z By observing the equation, we can obtain: tanθ x =-R 23 / R 33 sinθ y = R 13 tanθ z =-R 12 / R 11 If a gimbal lock occurs, then the rotation order of X, Y, Z is used to construct equation (11), which yields: tanθ x =R 23 / R 33 sinθ y = -R 31 tanθ z =R 21 / R 11 .

[0028]

[0029] If gimbal lock still occurs, then revert to the rotation order of Z, Y, X and process it separately, directly taking θ. x =0, sinθ y =R 13 tanθ z =R 21 / R 22。

[0030] Note: During the above rotation process, if the rotation axis and the vector to be rotated are parallel, no further judgment will be made, and the rotation angle will not be reset to 0, in order to improve calculation efficiency.

[0031] After the above rotations, the rotation time of the satellite around the three degrees of freedom is denoted as TR. x (t p ,t s ), TR y (t p ,t s ) and TR z (t p ,t sThen, the attitude adjustment time from adjacent observation task p to task s is calculated according to the following formula (12):

[0032]

[0033] Among them, t p and t s The observation times represented correspond to the starting and ending attitudes of the satellite attitude adjustment, respectively. Therefore, TR... x (t p ,t s )= TR x (Θ p (t p ),Θ s (t s ))= TR x (U p (t p ),U s (t s )) = TR x (θ x Given that the three degrees of freedom of the satellite have the same rotational performance, with a maximum angular velocity of ω and an angular acceleration of a, and using trapezoidal velocity planning, the time for attitude adjustment via X-axis rotation is calculated using equation (13) as follows. The calculation methods for the other two degrees of freedom are the same:

[0034]

[0035] In step two, when inserting an unscheduled observation task into a scheduled feasible sequence, for the case where there are tasks on both sides of the insertion position, the search range for the start time is determined by the intersection of the interval formed by the start times of the tasks on both sides of the insertion position and its own time window interval. Using the midpoint of the interval determined by the start time search range as the initial value, it is checked whether the switching time constraint is satisfied on both sides. If both sides are satisfied, this time is taken as the start time of the task to be inserted; if neither is satisfied, there is no suitable time for the task to be inserted in the local region, and the search is terminated. Otherwise, the search for the start time of the task to be inserted is performed. The specific steps of the search method are as follows:

[0036] S1: Initial value t of the observation start time for task j to be inserted. j 0 Let t be the midpoint of the search interval, and let t be the lower bound of the search interval. j min The maximum value of the start time of the selected observation time window and the start time of the predecessor task p, i.e., t j min =max(tws j ,t p ); the upper bound t of the search intervalj max The minimum value of the latest start time of the selected observation time window and the start time of the subsequent task s is t. j max =min(twe j -d j ,t s ), where d j The observation duration for the observation task j to be inserted;

[0037] S2: Determine the currently determined start time t j If the search has been performed before, execute S3; otherwise, calculate the time deviation lt between the task to be inserted j and the predecessor task p and successor task s that satisfy the attitude adjustment constraint according to equations (14) and (15), respectively. p,j and rt j,s If lt p,j and rt j,s If all values ​​are less than 0, then execute S3, lt. p,j and lt j,s If all values ​​are greater than 0, then execute S4; if none of the above conditions are met, then execute S5.

[0038]

[0039] S3: End search;

[0040] S4: with t j Insert task j into the scheduled feasible sequence as the observation start time;

[0041] S5: If lt p,j Less than 0 and rt j,s If it is greater than 0, then use equation (16) to calculate t. j Increase; otherwise, use equation (17) to increase t j Decrease; for t j Perform boundary constraint checks; if t j Less than the lower bound t of the search interval j min , then t j =t j min , t j Greater than the upper bound of the search interval t j max , then t j =t j max Finally, execute S2;

[0042]

[0043] In equations (16) and (17), u is a random decimal number between the interval [0.75-1.25].

[0044] A further implementation involves using a greedy insertion strategy in step two, which means arranging all tasks to be inserted in descending order of the observed task rewards, and then trying to insert them one by one in this order to construct and improve the solution.

[0045] The efficient global search evolutionary algorithm for solving the ultra-agile Earth observation satellite mission planning under the full-precision model in step three includes population initialization based on discrete integer encoding, individual selection based on minimum graph optimization decoding, crossover and mutation operators, and the specific method is as follows:

[0046] The solution to the satellite mission planning problem is represented using integer encoding, and the individual pop in the population... i Let pop represent a solution to the scheduling problem, and let vector pop be a solution to the scheduling problem. i =( pop i,1 ,…, pop i,j ,…, pop i,n ), where pop i,j =(vtw j i ,t j i ) represents an individual pop i Observational information for mission j, vtw j i Pop represents an individual. i The time window selected for performing observation task j, t j i Pop represents an individual. i The start time of observation task j must satisfy the time window constraint, as expressed by equation (18):

[0047]

[0048] In equation (18), tws j i and TWE j i Pop represents an individual in the population. i The time window selected for task j is vtw j i The start and end times, d j This indicates the working time of task j.

[0049] The initialization of integer encoding involves, in step three, specifically popping individuals from the population. i Its observation task j is based on the random selection of the observation window vtw. ji Then randomly select the observation window VTW j i Start time t j i Perform minimum graph optimization decoding on all individuals in the initialized population, and then pop the individuals. i The generated optimal feasible sequence is s i This indicates that its objective function value fs i Calculated from equation (19), s i Any two adjacent tasks p and s in the equation satisfy the attitude switching time constraint given by equation (20).

[0050]

[0051] Where, p j This represents the gains achieved by observation task j.

[0052]

[0053] Where, d p For the observation duration of task p, trans p,s (t p +d p ,t s Let be the attitude adjustment time from when the satellite completes mission p to when it begins observing mission s. If missions p and s satisfy the constraint of equation (20), they are called reachable, denoted as p→s; otherwise, they are called unreachable.

[0054] The individual selection process specifically involves sorting individuals in the population according to their objective function values ​​from largest to smallest, and selecting the first... % of individuals participate in cross-operations as elite individuals.

[0055] The crossover operator specifically selects the current individual for pop. i The optimal feasible sequence s i Pop with elite individuals el The optimal feasible sequence s of (el≠i) el To perform a crossover operation, first select observation sets from the two crossover sequences where both task and time window indices are the same. Then, randomly select an observation from these observation sets as the crossover position for the two crossover sequences. Finally, swap all tasks after the crossover position of the two crossover sequences to generate two new sequences s. i ’ and s el ’ The feasibility of the new sequence is repaired using a given strategy; finally, by comparison, the new sequence with the larger objective function value is selected as the final sequence generated in this crossover and denoted as cs. i,el If sequence csi,el Superior to s i Then cs i,el Replace s i and using sequence cs i,el Replace and update pop with observation information from the mission. i The observation information corresponding to the task, pop i The observation information for other missions remains unchanged;

[0056] After the crossover operation, the new sequence is used to perform feasibility repairs for two illegal cases: the task near the crossover position does not meet the attitude switching time constraint and repeated observations exist on both sides of the crossover position.

[0057] For illegal cases where tasks near the intersection do not meet the switching time constraints, the observation time of the corresponding task at the intersection is first adjusted. If no feasible observation time can be found within a limited number of iterations, the task is deleted to relax the time. For illegal cases of repeated observations on both sides of the intersection, time windows and tasks with larger time intervals are retained first.

[0058] In step three, for the invalid case where tasks near the intersection do not meet the switching time constraint, the specific adjustment method is as follows:

[0059] S1: Check whether the observation tasks at the intersection and their successor observation tasks meet the switching time constraints. If they do, there is no need to adjust the sequence; otherwise, execute S2.

[0060] S2: Adjust the start time of observations at the intersection. If a valid observation time can be found to satisfy the switching time constraint using the method in step two, then update the start time of observations at the intersection and end the adjustment; otherwise, execute S3.

[0061] S3: Delete the subsequent observations at the intersection position, and re-check whether the tasks on both sides of the deleted position meet the switching time constraint. If they do, end the adjustment; otherwise, execute step S2.

[0062] The specific adjustment steps for the illegal duplicate observations existing on both sides of the intersection are as follows:

[0063] S1: Iterate through each group of repeated observation tasks. If the time windows corresponding to the repeated observation tasks are the same, execute S2; otherwise, execute S3.

[0064] S2: Randomly select a duplicate observation as the observation to be deleted, and execute S4;

[0065] S3: Select observations with smaller time window intervals as observations to be deleted, and execute S4;

[0066] S4: If the switching time constraint can be met after the observation to be deleted is deleted, the observation to be deleted is deleted directly and the adjustment ends. Otherwise, it is determined whether the switching time constraint of the observations on both sides is met after the other observation is deleted. If it is met, the observation is deleted and the adjustment ends. Otherwise, S5 is executed.

[0067] S5: At this point, for a task with repeated observations, when any repeated observation is deleted, a sequence that does not meet the switching time constraint will be generated. Using the adjustment method for the first illegal case mentioned above, the start time of the adjacent task at the position where the constraint is not met is adjusted or the observation is deleted, and a new legal sequence can be formed. Compare the new sequences formed after deleting two repeated observations respectively, select the better new sequence and end the adjustment.

[0068] A further implementation involves the mutation operator in step three employing the idea of ​​prioritizing observations on the optimal sequence, performing mutation operations on each individual in the population. The specific adjustment method is as follows: for each individual popped... i The observation pop below i,j First, determine whether task j is popped by individual. i The optimal sequence s after minimum graph decoding i If the value is above the threshold, the mutation probability is set to a smaller probability PM; otherwise, the mutation probability is set to probability 1-PM. A random number between [0,1] is generated. If this number is less than the preset mutation probability, task j is mutated; otherwise, the observation information of task j remains unchanged. The observation pop... i,j The mutation operation is to randomly select a time window (vtw). j i Then randomly select the observation window VTW j i Start time t j i .

[0069] A high-precision satellite scheduling system based on a full-precision model includes an attitude switching time calculation module, a global search module, and a local search module. The global search module employs an evolutionary algorithm based on discrete integer encoding and minimum graph decoding techniques. The local search module inserts tasks and constructs solutions based on the attitude switching satisfaction deviation values ​​of adjacent tasks. The full-precision model solving module completes model task planning, solving, and scheduling based on a memetic framework combining global search and greedy insertion local search algorithms.

[0070] The beneficial effects of this invention are as follows: By constructing and solving a global search algorithm, and using the operation operators of optimal path preservation and crossover strategies, this invention provides an efficient method for solving the global search of ultra-agile Earth observation satellites. Furthermore, by utilizing the satisfaction deviation value of attitude switching between adjacent tasks on both sides, it efficiently searches for task insertion time, thus solving the local search problem of ultra-agile satellite scheduling. At the same time, the meme algorithm constructed based on this can efficiently solve the ultra-agile satellite scheduling problem of full-precision model, providing a more efficient scheduling method for ultra-agile satellite full-schedule model scheduling. Attached Figure Description

[0071] Figure 1 This is a flowchart of a method according to an embodiment of the present invention.

[0072] Figure 2 This is a framework diagram of the meme algorithm in an embodiment of the present invention.

[0073] Figure 3 This is a schematic diagram of the attitude switching of a hypersensitive satellite according to an embodiment of the present invention.

[0074] Figure 4 This is an example diagram of discrete integer encoding for the satellite 3 mission, which is an application example of the present invention.

[0075] Figure 5 This is a flowchart illustrating the cross-operation process in an application example of the present invention.

[0076] Figure 6 This is a diagram illustrating an application example of the cross-operation of the present invention. Detailed Implementation

[0077] To enhance understanding of the present invention, the present invention will be further described in detail below with reference to embodiments. These embodiments are only used to explain the present invention and do not constitute a limitation on the scope of protection of the present invention.

[0078] The problem concerning satellite mission planning in this invention is described as follows:

[0079] The user provides a set of observation tasks, T = {j1, j2, ..., j...} n}, where n represents the number of pre-scheduled tasks, and the working time required to execute observation task j is d. j The reward obtained by completing observation task j is p. j Within the scheduling timeframe, the set of selectable observation periods (i.e., visible time windows) for observation task j is represented by VTW. j It indicates that |VTW j | represents the number of visible time windows for observation task j; if task j is selected and observations are performed, then the decision variable y j =1, otherwise y j =0; When performing observation task j, only a specific time window (vtw) is selected. j (vtw)j ∈VTW j At a certain integer time t within ) j An observation is performed as a decision variable, satisfying a time window constraint, as shown in Equation 1, where the selected specific time window is vtw. j =[tws j ,twe j ], tws j The start time of this time window, twe j This is the end time of the time window.

[0080]

[0081] Figure 3 This is a schematic diagram of attitude switching for a hypersensitive satellite. The attitude switching time between any two adjacent observation tasks p and s is denoted by trans. p,s (t p +d p ,t s ) represents the attitude switching constraint; adjacent observation tasks that satisfy the constraint described in equation (2) are called reachable, denoted as p→s, otherwise they are called unreachable.

[0082]

[0083] Where trans in equation (2) p,s (t p +d p ,t s Calculate as follows:

[0084] At a certain moment t, the satellite j The attitude angle of staring (pointing) at target j is represented by the vector Θ. j (t j )=(α j (t j ),β j (t j ),γ j (t j )) indicates that α j (t j ) is the lateral swing angle, β j (t j ) is the pitch angle, γ j (t j The yaw angle is 0; at this time, the coordinates of target j in the VVLH satellite orbit coordinate system are C. j (t j )=(x j (t j ),y j (tj ),z j (t j )), then Θ j (t j ) and C j (t j The conversion relationship between ) is shown in equation (3).

[0085]

[0086] The satellite's center of mass points to the unit vector U of target j. j (t j ) =(u x j (t j ),u y j (t j ),u z j (t j According to equation (1), u z j (t j ), u x j (t j ) and u y j (t j The following can be calculated using formulas (4)-(6):

[0087]

[0088] Given any two adjacent observation tasks p and s, where task p is the preceding task and task s is the succeeding task, the unit vectors pointing to the adjacent targets when the satellite is staring at a certain moment are denoted as U. p (t p ) and U s (t s ), can be determined by observing the attitude angle Θ p (t p ) and Θ s (t s ) Calculated according to equations (4)-(6), the observation attitude angle is given and known in the form of a calculation example based on satellite parameters; the rotation axis K can be obtained from U p (t p ) and U s (t s The cross product of vector U is calculated, then the normalized rotation axis k = K / ||K||; vector U p (t p ) to vector U s (t s The rotation angle φ = arccos(U)p (t p )∙ U s (t s Then the rotation matrix R can be calculated by the Rodrigues rotation formula, as shown in equation (7), where I is the identity matrix and [k] is the cross product matrix of k, given by equation (8).

[0089]

[0090]

[0091] The satellite attitude adjustment defaults to rotating in the order of Z, Y, X axes, with rotation angles denoted as θ. z , θ y , θ x The corresponding rotation matrix is ​​R. z , R y , R x Calculate according to formulas (9)-(11).

[0092]

[0093] Then we can construct equation (12):

[0094]

[0095] Among them, s x =sinθ x c x =cosθ x s y =sinθ y c y =cosθ y s z =sinθ z c z =cosθ z By observing the equation, we can obtain: tanθ x =-R 23 / R 33 sinθ y = R 13 tanθ z =-R 12 / R 11 If a gimbal lock occurs, then the rotation order of X, Y, Z is used to construct equation (13), which yields: tanθ x =R 23 / R 33 sinθ y = -R 31 tanθz =R 21 / R 11 .

[0096]

[0097] If gimbal lock still occurs, then revert to the rotation order of Z, Y, X and process it separately, directly taking θ. x =0, sinθ y =R 13 tanθ z =R 21 / R 22。

[0098] Note: During the above rotation process, if the rotation axis and the vector to be rotated are parallel, no further judgment will be made, and the rotation angle will not be reset to 0, in order to improve calculation efficiency.

[0099] After the above rotations, the rotation time of the satellite around the three degrees of freedom is denoted as TR. x (t p ,t s ), TR y (t p ,t s ) and TR z (t p ,t s Then, the attitude adjustment time from adjacent observation task p to task s is calculated according to the following formula (14):

[0100]

[0101] Among them, t p and t s The observation times represented correspond to the starting and ending attitudes of the satellite attitude adjustment, respectively. Therefore, TR... x (t p ,t s )= TR x (Θ p (t p ),Θ s (t s ))= TR x (U p (t p ),U s (t s )) = TR x (θ x The satellite's three degrees of freedom are constrained to have the same rotational performance, with a maximum angular velocity of ω and an angular acceleration of a. Using trapezoidal velocity planning, the time for attitude adjustment via X-axis rotation is calculated using equation (15) as follows. The calculation methods for the other two degrees of freedom are the same:

[0102]

[0103] The feasible observation time series that satisfies the constraints of equations (1)-(2) above is a feasible scheduling scheme, denoted as s; the objective function of the feasible scheduling scheme, f(s), is the total revenue of the scheduled task, calculated by equation (16); the optimization objective is to maximize f(s).

[0104]

[0105] Example

[0106] This embodiment provides a method for scheduling ultra-agile Earth observation satellite full-precision models. This method is a memetic evolution algorithm that includes global search and local search, and its algorithm framework is as follows: Figure 2 As shown, individuals in the population continuously generate new superior individuals through crossover and mutation. Then, by introducing local search, partial solutions are optimized, thereby continuously improving the convergence speed and the quality of the solutions. Considering that the proposed optimal path crossover and mutation operators can enable a large number of ordinary individuals in the population to learn the solution structure of elite individuals, the local search in the meme algorithm is an improvement on the solution construction by performing greedy insertion on a subset of elite individuals (the top b% of individuals with larger objective functions). This design achieves a trade-off between computational cost and obtaining high-quality solutions.

[0107] Step 1: Implement local search.

[0108] The greedy insertion algorithm constructs a system that sorts all tasks to be inserted in descending order of observation task reward, and then attempts to insert them one by one in this order. When inserting a specific observation task, the algorithm first sorts all time windows of the task according to a given strategy (e.g., sorting by time window interval in descending order). Then, it sequentially traverses these windows to find feasible insertion positions and sorts these feasible insertion positions according to a certain strategy (e.g., random sorting). Finally, it traverses these feasible insertion positions and performs the task insertion operation. If the insertion is successful, the insertion of the task ends; otherwise, it continues to traverse other insertion positions to insert the task. If all insertion positions under a given time window are not feasible, it continues to traverse other time windows until the task is successfully inserted. If insertion fails, the task to be inserted is discarded.

[0109] This algorithm can quickly construct some tasks onto the scheduling sequence, improving the overall solution quality of the population in a short time. The key is how to insert unscheduled tasks into the already scheduled sequence to achieve an improved solution. Existing methods utilize the FIFO property, using a binary search or preprocessing method to obtain the upper and lower bounds of the start time of the task to be inserted. Simultaneously, the upper and lower bounds of the start times of other tasks in the new scheduled sequence are updated, transforming the problem of finding the task start time into finding the upper and lower bounds of the start time. This converts the integer programming problem into a sorting problem, reducing the search space and complexity. In ultra-agile satellite mission planning, the mission observation start time does not satisfy the strong FIFO constraint property. The feasible and infeasible regions for task insertion alternate within the time window and cannot be divided by a single upper (lower) bound time. Therefore, a binary search cannot be used for searching, although the constraint satisfaction still exists locally.

[0110] For the aforementioned insertion task, when inserting an unscheduled observation task into a scheduled feasible time series, if there are tasks on both sides of the insertion position, the search range for the start time is determined by the intersection of the interval formed by the start times of the tasks on both sides of the insertion position and its own time window interval. Using the midpoint of the interval determined by the start time search range as the initial value, it is checked whether the switching time constraint is satisfied on both sides. If both sides are satisfied, this time is taken as the start time of the task to be inserted; if neither is satisfied, there is no suitable time for the task to be inserted in the local region, and the search is terminated. Otherwise, the search for the start time of the task to be inserted is performed, with the specific steps as follows:

[0111] S1: Initial value t of the observation start time for task j to be inserted. j 0 Let t be the midpoint of the search interval, and let t be the lower bound of the search interval. j min The maximum value of the start time of the selected observation time window and the start time of the predecessor task p, i.e., t j min =max(tws j ,t p ); the upper bound t of the search interval j max The minimum value of the latest start time of the selected observation time window and the start time of the subsequent task s is t. j max =min(twe j -d j ,t s ), where d j The observation duration is the time required for the observation task j to be inserted.

[0112] S2: Determine the currently determined start time t jIf the search has been performed before, execute S3; otherwise, calculate the time deviation lt between the task to be inserted j and the predecessor task p and the successor task s that satisfy the attitude adjustment constraint according to equations (17) and (18), respectively. p,j and rt j,s If lt p,j and rt j,s If all values ​​are less than 0, then execute S3, lt. p,j and lt j,s If all values ​​are greater than 0, then execute S4; if none of the above conditions are met, then execute S5.

[0113]

[0114] S3: End search;

[0115] S4: with t j Insert task j into the scheduled feasible sequence as the observation start time;

[0116] S5: If lt p,j Less than 0 and rt j,s If it is greater than 0, then use equation (19) to calculate t. j Increase; otherwise, use equation (20) to increase t j Decrease; for t j Perform boundary constraint checks; if t j Less than the lower bound t of the search interval j min , then t j =t j min , t j Greater than the upper bound of the search interval t j max , then t j =t j max Finally, execute S2;

[0117]

[0118] In equations (19) and (20), u is a random decimal number between the interval [0.75-1.25].

[0119] Step 2: Implement global search.

[0120] The global search includes population initialization based on discrete integer encoding, individual selection based on minimum graph optimization decoding, and crossover and mutation operators, as detailed below:

[0121] The solution to the satellite mission planning problem is represented using integer encoding, and the individual pop in the population... i Let pop represent a solution to the scheduling problem, and let vector pop be a solution to the scheduling problem. i =( popi,1 ,…, pop i,j ,…, pop i,n ), where pop i,j =(vtw j i ,t j i ) represents an individual pop i Observational information for mission j, vtw j i Pop represents an individual. i The time window selected for performing observation task j, t j i Pop represents an individual. i The start time of observation task j must satisfy the time window constraint, as expressed in Equation 1 above:

[0122] Population initialization specifically involves popping individuals from the population. i Its observation task j is based on the random selection of the observation window vtw. j i Then randomly select the observation window VTW j i Start time t j i Perform minimum graph optimization decoding on all individuals in the initialized population, and then pop the individuals. i The generated optimal feasible sequence is s i This indicates that its objective function value fs i Calculated from equation (21), s i Any two adjacent tasks p and s in the equation satisfy the attitude switching time constraint given by equation (2) above.

[0123]

[0124] Where, p j This represents the gains achieved by observation task j.

[0125] Individual selection involves ranking individuals in the population according to their objective function values ​​from largest to smallest, with the top individuals selected first. % of individuals participate in cross-operations as elite individuals.

[0126] The crossover operator specifically selects the current individual for pop. i The optimal feasible sequence s i Pop with elite individuals el The optimal feasible sequence s of (el≠i) elTo perform a crossover operation, first select observation sets from the two crossover sequences where both task and time window indices are the same. Then, randomly select an observation from these observation sets as the crossover position for the two crossover sequences. Finally, swap all tasks after the crossover position of the two crossover sequences to generate two new sequences s. i ’ and s el ’ The feasibility of the new sequence is repaired using a given strategy; finally, by comparison, the new sequence with the larger objective function value is selected as the final sequence generated in this crossover and denoted as cs. i,el If sequence cs i,el Superior to s i Then cs i,el Replace s i and using sequence cs i,el Replace and update pop with observation information from the mission. i The observation information corresponding to the task, pop i The observation information for other missions remains unchanged.

[0127] After the crossover operation, the new sequence is used to perform feasibility repairs for two illegal cases: the task near the crossover position does not meet the attitude switching time constraint and repeated observations exist on both sides of the crossover position.

[0128] For illegal cases where tasks near the intersection do not meet the switching time constraints, the observation time of the corresponding task at the intersection is first adjusted. If no feasible observation time can be found within a limited number of iterations, the task is deleted to relax the time. For illegal cases of repeated observations on both sides of the intersection, time windows and tasks with larger time intervals are retained first.

[0129] For invalid cases where tasks near the intersection do not meet the switching time constraints, the specific adjustment steps are as follows:

[0130] S1: Check whether the observation tasks at the intersection and their successor observation tasks meet the switching time constraints. If they do, there is no need to adjust the sequence; otherwise, execute S2.

[0131] S2: Adjust the start time of observations at the intersection. If the insertion method in the local search can find a legal observation time that satisfies the switching time constraint, then update the start time of observations at the intersection and end the adjustment; otherwise, execute S3.

[0132] S3: Delete the subsequent observations at the intersection position, and re-check whether the tasks on both sides of the deleted position meet the switching time constraint. If they do, end the adjustment; otherwise, execute step S2.

[0133] The specific adjustment steps for the illegal duplicate observations existing on both sides of the intersection are as follows:

[0134] S1: Iterate through each group of repeated observation tasks. If the time windows corresponding to the repeated observation tasks are the same, execute S2; otherwise, execute S3.

[0135] S2: Randomly select a duplicate observation as the observation to be deleted, and execute S4;

[0136] S3: Select observations with smaller time window intervals as observations to be deleted, and execute S4;

[0137] S4: If the switching time constraint can be met after the observation to be deleted is deleted, the observation to be deleted is deleted directly and the adjustment ends. Otherwise, it is determined whether the switching time constraint of the observations on both sides is met after the other observation is deleted. If it is met, the observation is deleted and the adjustment ends. Otherwise, S5 is executed.

[0138] S5: At this point, for tasks with repeated observations, deleting any repeated observation will generate a sequence that does not meet the switching time constraint. Using the adjustment method described in the first illegal case above, the start time of adjacent tasks at the position where the constraint is not met can be adjusted or observations can be deleted to form a new legal sequence. Compare the new sequences formed after deleting two repeated observations respectively, select the better new sequence, and end the adjustment.

[0139] The mutation operator prioritizes observations on the optimal sequence and performs mutation on each individual in the population. The specific adjustment method is as follows: For each individual popped... i The observation pop below i,j First, determine whether task j is popped by individual. i The optimal sequence s after minimum graph decoding i If the value is above the threshold, the mutation probability is set to a smaller probability PM; otherwise, the mutation probability is set to probability 1-PM. Generate a random number between [0,1]. If this number is less than the preset mutation probability, mutate task j; otherwise, the observation information for task j remains unchanged. Pop the observation. i,j The mutation operation is to randomly select a time window (vtw). j i Then randomly select the observation window VTW j i Start time t j i .

[0140] Application examples

[0141] according to Figure 3-6 As shown, this embodiment provides examples of each part of the ultra-agile Earth observation satellite full-precision model scheduling method.

[0142] Example 1: Calculation of discrete coding and attitude switching time.

[0143] This example uses a self-built full-precision ultra-agile satellite observation dataset to construct an ultra-agile satellite with the same orbital elements as AS-01. Points within mainland China (3°N to 53°N, 74°E to 133°E) are randomly selected as observation targets. No constraints are imposed on the satellite's rotation attitude angles, and the orbit is not further subdivided. Information from three observation tasks is extracted from the dataset and used for discrete integer encoding of the example. The known data information from the satellite observations is shown in Table 1 below.

[0144] A time window for the task is randomly selected, and a moment within that time window is randomly chosen as the start time of the observation task. The discrete integer codes for the three observation tasks are as follows:

[0145] pop i,1 =([3410,3585],3500)

[0146] pop i,2 =([2414,3300],3000)

[0147] pop i,3 =([2014,2884],2500)

[0148] pop i,1 To illustrate the specific meaning of the encoding, we first select time 3500 of the first time window to execute task 1. [3410, 3585] represents the interval formed by the start and end times of the time window.

[0149] Based on the task start time directly given by the discrete encoding, the three attitude angles required during the task execution period can be obtained by looking up the table in the example set, as shown in Table 2 below.

[0150] Next, the switching time between adjacent tasks needs to be calculated to further determine whether the adjacent tasks meet the attitude transition constraints. First, the attitude angles are converted into unit vectors of the staring target using equations (3)-(6), and the conversion results are shown in Table 3.

[0151] Based on the start time of the observation tasks, the observation order was determined to be j3, j2, j1. Taking the unit vector U3(2524) at time 2524 after completing task 3 and rotating to the unit vector U2(3000) at time 3000 when task 2 begins as an example, the calculation of the switching time is explained. Table 3 shows:

[0152] U3(2524) = [-0.45359,-0.43964,0.775226]

[0153] U2(3000) = [-0.73908, 0.137262, 0.659481]

[0154] The rotation axis and rotation angle from U3(2524) to U2(3000) are as follows:

[0155] Rotation axis k = [-0.64128403 -0.44304123 -0.62647366]

[0156] Rotation angle φ=38.17

[0157] From equations (7)-(8), we can obtain the Rodriguez rotation as follows:

[0158]

[0159] First, rotate along the Z, Y, X axes in that order. According to equations (9)-(12), we can obtain the following equations.

[0160]

[0161] According to tanθ x =-R 23 / R 33 sinθ y = R 13 tanθ z =-R 12 / R 11 The rotation angle can be obtained as follows:

[0162] θ x = -27.64, θ y =-10.83, θ z =-27.13, according to Equation 14-15, taking the maximum angular velocity w=3 and angular acceleration a=0.3, the switching time trans can be obtained. 3,2 (t3+d3,t2)=50.24

[0163] Following the same calculation process as above, we can obtain: trans 2,1 (t2+d2,t1)=132.03

[0164] Finally, check the attitude switching constraints:

[0165] t3+d3+ trans 3,2 (t3+d3,t2)=2574.2< t2

[0166] t2+d2+ trans2,1 (t2+d2,t1)=3158.03< t1

[0167] In summary, the discrete integer solutions in the above examples are feasible solutions, and the resulting feasible sequences and scheduling schemes are shown in the appendix to the specification. Figure 4 As shown.

[0168] Example 2: Task insertion based on deviation values ​​satisfied by constraints on both sides.

[0169] The task insertion method is summarized as follows: For a time window of a task to be inserted, firstly, a suitable insertion position is found in the feasible scheduling sequence based on the time window information. Then, all insertion positions are traversed, and within the time range determined by the insertion position and the time window, a task start time that satisfies the switching time constraint of adjacent tasks on both sides is found. Among these, efficiently and quickly finding the start time of the observation task is crucial. A specific example is given below to illustrate this method.

[0170] A time window interval vtw for task j to be inserted j =[18504,18817], working time d j =20, the predecessor task started at 18762 seconds and the successor task started at 18831 seconds.

[0171] First, determine the valid search range for the task j to be inserted: the lower bound t of the search range. j min =max(18054,18762)=18762, upper bound t j max =min(18837-20,18831)=18817, then the initial value of the start time t j 0 = (18762+18817) / 2=18789

[0172] Next, the constraint satisfaction on both sides is calculated according to equations (17) and (18). At this point, lt p,j =-10.41<0, rt j,s =9.87>0, which means that the left constraint is not satisfied while the right constraint is satisfied. Therefore, the next step is to increase t according to equation (19). j Up to 18800. Recheck the constraints on both sides, and calculate lt. p,j =1.53>0, rt j,s =-0.66<0, which means that the left-side constraint is satisfied while the right-side constraint is not satisfied. Therefore, the next step is to reduce t according to equation (20). j Up to 18799. Re-examine the constraints on both sides, and calculate lt. p,j =0.71>0, rt j,s=0.34>0, all adjacent tasks satisfy the switching time constraint, stop searching, and 18799 seconds is the start time of the task to be inserted.

[0173] Example 3: Global search and related operators.

[0174] Instruction manual attached Figure 5 The operation flow of the crossover operator is given, including steps such as elite selection, sequence crossover, sequence repair, and individual update.

[0175] Instruction manual attached Figure 6 An application example of the crossover operator is given. Circles represent observation tasks on the optimal sequence after individual decoding. The numbers inside the circles are the task and time window numbers (separated by commas), respectively. The number above the circle represents the start time of the task. The arrows between the circles represent the reachability between two observation tasks (black arrows indicate that the predecessor task is reachable from the successor task, and red arrows indicate that reachability is not yet known). The cyan circles in the figure represent observation tasks at the crossover positions.

[0176] From the instruction manual Figure 6 It is evident that this crossover operator is more problem-specific and interpretable, and can more effectively promote communication between individuals, enabling them to obtain node information from elite individuals that form the optimal observation sequence, thus greatly increasing the population evolution speed. Secondly, the new individuals formed after the execution of this crossover operator directly satisfy the time window constraint. Finally, selecting observations corresponding to a common time window as the optimal sequence crossover position ensures that the observation sequences on both sides of the crossover position remain feasible, requiring only feasibility checks and repairs for the vicinity of the crossover position and duplicate observation positions.

[0177] The foregoing has shown and described the basic principles, main features, and advantages of the present invention. Those skilled in the art should understand that the present invention is not limited to the above embodiments. The embodiments and descriptions in the specification are merely illustrative of the principles of the invention. Various changes and modifications can be made to the invention without departing from its spirit and scope, and all such changes and modifications fall within the scope of the present invention as claimed. The scope of protection of the present invention is defined by the appended claims and their equivalents.

Claims

1. A method for scheduling a full-precision model of ultra-agile Earth observation satellites, characterized in that, Includes the following steps: Step 1: Accurately calculate the attitude switching time between adjacent tasks using Euler angles and rotation matrices; Step 2: Search for the task insertion time based on the attitude switching satisfaction deviation value of adjacent tasks on both sides, and simultaneously implement the local search algorithm for hypersensitive satellite scheduling; This also includes inserting unscheduled observation tasks into the already scheduled feasible sequence, as detailed below: For cases where there are tasks on both sides of the insertion position, the search range for the start time is determined by the intersection of the interval formed by the start times of the tasks on both sides of the insertion position and its own time window interval. The midpoint of the interval determined by the search range for the start time is used as the initial value. It is checked whether the switching time constraint is met on both sides. If both sides are satisfied, the time is taken as the start time of the task to be inserted. If neither is satisfied, there is no suitable time in the local area for the task to be inserted, and the search is terminated. Otherwise, the search for the start time is performed. The search for the start time is as follows: S1: Task to be inserted j Initial value of observation start time t j 0 The midpoint of the search interval, and the lower bound of the search interval. t j min The start time of the selected observation window and the predecessor mission. p The maximum value of the start time, i.e. t j min = max ( tws j ,t p ); Upper bound of the search interval t j max The latest start time of the selected observation window and the subsequent task. s The minimum start time, i.e. t j max = min ( twe j - d j ,t s ),in d j For the observation mission to be inserted j The duration of observation work; S2: Determine the currently determined start time. t j If the search has been performed before, execute S3; otherwise, calculate the task to be inserted according to equations (14) and (15). j With the precursor mission p and subsequent tasks s Time deviation that satisfies attitude adjustment constraints lt p,j and rt j,s ;like lt p,j and rt j,s If all values ​​are less than 0, then execute S3. lt p,j and lt j,s If all values ​​are greater than 0, then execute S4; if none of the above conditions are met, then execute S5. Equation (14) Equation (15) S3: End search; S4: with t j The mission will be used as the start time of observation. j Insert a scheduled feasible sequence; S5: If lt p,j Less than 0 and rt j,s If it is greater than 0, then use equation (16) to... t j Increase; otherwise, use equation (17) to increase t j Decrease; t j Perform boundary constraint checks, if t j Less than the lower bound of the search interval t j min ,but t j =t j min , t j Greater than the upper bound of the search interval t j max ,but t j =t j max Finally, execute S2; Equation (16) Equation (17) In equations (16) and (17) u All are random decimals within the range [0.75-1.25]; Step 3: Based on discrete integer encoding and minimum graph decoding, implement the crossover and mutation operators for optimal path exchange and retention, and further implement the global search evolutionary algorithm for hypersensitive satellite scheduling; It also includes an efficient global search evolutionary algorithm for solving the ultra-agile Earth observation satellite mission planning under the full-precision model, the specific method of which is as follows: The crossover operator specifically selects the current individual. pop i Optimal feasible sequence s i With elite individuals pop el ( el≠i The optimal feasible sequence of ) s el To perform a crossover operation, first select observation sets from the two crossover sequences that have the same task and time window indices. Then, randomly select an observation from the observation set as the crossover position of the two crossover sequences. Finally, swap all tasks after the crossover position of the two crossover sequences to generate two new sequences. s i ’ and s el ’ The feasibility of the new sequence is repaired using a given strategy; finally, by comparison, the new sequence with the larger objective function value is selected as the final sequence generated in this crossover and denoted as [the sequence is not specified in the original text]. cs i,el , if the sequence cs i,el Superior s i ,but cs i,el replace s i and using sequence cs i,el Replace and update the observation information of the mission. pop i The observation information corresponding to the task in the middle, pop i The observation information for other missions remains unchanged; After the crossover operation, the new sequence is used to perform feasibility repairs for two illegal cases: tasks near the crossover position do not meet the attitude switching time constraint and repeated observations exist on both sides of the crossover position. The specific methods are as follows: For illegal cases where tasks near the intersection do not meet the switching time constraints, the observation time of the corresponding task at the intersection is first adjusted. If no feasible observation time can be found within a limited number of time iterations, the task is deleted to relax the time. For illegal cases of repeated observations on both sides of the intersection, time windows and tasks with larger time intervals are retained first. For the invalid case where tasks near the intersection do not meet the switching time constraint, the specific adjustment steps are as follows: S1. Check whether the observation tasks at the intersection and their successor observation tasks meet the switching time constraints. If they do, there is no need to adjust the sequence; otherwise, execute S2. S2. Adjust the start time of observations at the intersection. If a valid observation time can be found to satisfy the switching time constraint using the method in step two, update the start time of observations at the intersection and end the adjustment; otherwise, execute S3. S3. Delete the subsequent observations at the intersection position, and re-check whether the tasks on both sides of the deleted position meet the switching time constraint. If they do, end the adjustment; otherwise, execute step S2. The specific adjustment steps for the illegal duplicate observations existing on both sides of the intersection are as follows: S1. Iterate through each group of repeated observation tasks. If the time windows corresponding to the repeated observation tasks are the same, then execute S2; otherwise, execute S3. S2. Randomly select a duplicate observation as the observation to be deleted, and execute S4; S3. Select observations with smaller time window intervals as observations to be deleted, and execute S4. S4. If the switching time constraint can be met by the observations on both sides after the observation to be deleted is deleted, then the observation to be deleted is deleted directly and the adjustment ends. Otherwise, it is determined whether the switching time constraint of the observations on both sides is met after the other observation is deleted. If it is met, then the observation is deleted and the adjustment ends. Otherwise, S5 is executed. S5. At this point, for a task with repeated observations, when any repeated observation is deleted, a sequence that does not meet the switching time constraint will be generated. The adjustment method for the first illegal case mentioned above is used to adjust the start time or delete the observations of the adjacent tasks at the position where the constraint is not met, forming a new legal sequence. Compare the new sequences formed after deleting two repeated observations respectively, select the better new sequence and end the adjustment. Step 4: Based on the meme algorithm combining Step 2 and Step 3, efficiently solve the scheduling of the hypersensitive satellite full-precision model.

2. The ultra-agile Earth observation satellite full-precision model scheduling method according to claim 1, characterized in that: Step one further includes: converting the attitude angle of the staring target into a unit vector pointing from the satellite to the target, the specific method of which is as follows: At a certain moment, the satellite t j Stare at the target j The attitude angle is represented by a vector. Θ j ( t j )=( α j ( t j ), β j ( t j ), γ j ( t j )) indicates that among them α j ( t j () is the lateral swing angle. β j ( t j () represents the pitch angle. γ j ( t j The yaw angle is 0; at this time, the target j Coordinates in the VVLH satellite orbital coordinate system C j ( t j )=( x j ( t j ), y j ( t j ), z j ( t j )),but Θ j ( t j )and C j ( t j The conversion relationship between ) is shown in equation (1); Equation (1) Satellite's center of mass is pointing towards the target j unit vector U j ( t j ) =( u x j ( t j ), u y j ( t j ), u z j ( t j According to equation (1), u z j ( t j ), u x j ( t j )and u y j ( t j The following can be calculated using formulas (2)-(4): Equation (2) Equation (3) Equation (4).

3. The ultra-agile Earth observation satellite full-precision model scheduling method according to claim 2, characterized in that: In step one, when observing attitude adjustment between adjacent tasks, the angles of rotation around the three degrees of freedom, i.e., Euler angles, and the attitude adjustment time are calculated. The specific method is as follows: Given any adjacent observation task p and tasks s Among them, the task p As a prerequisite task, the task s For subsequent missions, when the satellite is staring at a certain moment, the unit vectors pointing to adjacent targets are denoted as follows: U p ( t p )and U s ( t s (This can be determined by observing the attitude angle) Θ p ( t p )and Θ s ( t s The observation attitude angle is calculated according to the method corresponding to equations (1)-(4), and is given and known in the form of a calculation example based on the satellite parameters; Rotation axis K can be U p ( t p )and U s ( t s The cross product calculation of ) then the normalized rotation axis k = K / || K ||; vector U p ( t p ) to vector U s ( t s rotation angle φ = arccos ( U p ( t p )∙ U s ( t s Then the rotation matrix R It can be calculated using the Rodriguez rotation formula, as shown in equation (5), where I For the identity matrix, [ k ]for k The cross product matrix is ​​given by equation (6); Equation (5) Equation (6) Satellite attitude adjustment defaults to Z , Y , X The axes are rotated sequentially, with rotation angles of 10° and 20° respectively. θ z , θ y , θ x The corresponding rotation matrix is R z , R y , R x Calculate according to formulas (7)-(9); Equation (7) Equation (8) Equation (9) Then construct equation (10): Equation (10) in, s x = sinθ x , c x = cosθ x , s y = sinθ y , c y = cosθ y , s z = sinθ z , c z = cosθ z By observing the equation, we can obtain: tanθ x =- R 23 / R 33 , sinθ y = R 13 , tanθ z =- R 12 / R 11 If a gimbal lock appears, replace it with... X , Y , Z By constructing equation (11) based on the rotation sequence, we can obtain: tanθ x = R 23 / R 33 , sinθ y = -R 31 , tanθ z = R 21 / R 11 ; Equation (11) If gimbal lock still occurs, then return to... Z , Y , X If the rotation order is handled separately, then directly take... θ x = 0, sinθ y = R 13 , tanθ z = R 21 / R 22 ; During the above rotation process, if the rotation axis and the vector to be rotated are parallel, no further judgment is made, and the rotation angle is not reset to 0, in order to improve computational efficiency; After the above rotations, the rotation times of the satellite around the three degrees of freedom are denoted as follows: TR x ( t p , t s ), TR y ( t p , t s )and TR z ( t p , t s ), then adjacent observation tasks p To the mission s The attitude adjustment time is calculated according to the following formula (12): Equation (12) in, t p and t s The observation times indicated correspond to the starting and ending attitudes of the satellite attitude adjustment, respectively. TR x ( t p , t s )= TR x ( Θ p ( t p ), Θ s ( t s ))= TR x ( U p ( t p ), U s ( t s )) = TR x ( θ x The satellite is constrained to have identical rotational properties across its three degrees of freedom, with a maximum angular velocity of [missing value]. ω Angular acceleration is a And using trapezoidal velocity programming, then through X The time for attitude adjustment by axis rotation is calculated using equation (13) as follows; the calculation methods for the other two degrees of freedom are the same. Equation (13).

4. The ultra-agile Earth observation satellite full-precision model scheduling method according to claim 1, characterized in that: In step two, the local search strategy is greedy insertion, which means that all tasks to be inserted are arranged in descending order of the observation task reward, and then the insertion is attempted one by one in this order.

5. The ultra-agile Earth observation satellite full-precision model scheduling method according to claim 1, characterized in that: In step three, discrete integer encoding represents the solution to the satellite mission planning problem using integer encoding, and the individuals in the population... pop i Let represent a solution to the scheduling problem, and let be a vector. pop i =( pop i,1 ,…, pop i,j ,…, pop i,n ),in pop i,j =( vtw j i ,t j i ) represents an individual pop i Medium task j Observational information, vtw j i Represents an individual pop i Perform observation tasks j The selected time window, t j i Represents an individual pop i Perform observation tasks j The start observation time must satisfy the time window constraint, as expressed by equation (18): Equation (18) In formula (18) tws j i and twe j i Individuals in a population pop i Execute the task j Selected time window vtw j i The start and end times, d j Indicates task j Working hours.

6. The ultra-agile Earth observation satellite full-precision model scheduling method according to claim 1, characterized in that: Step three also includes: population initialization, the specific method of which is as follows: initializing the individuals in the population. pop i Its observation mission j First, randomly select the observation window. vtw j i Then randomly select the observation window. vtw j i Start time within t j i Perform; after performing minimum graph optimization decoding on all individuals in the initialized population, the individual... pop i The generated optimal feasible sequence is used s i This indicates that its objective function value fs i Calculated by equation (19), s i Any two adjacent tasks p and s , satisfying the attitude switching time constraint given by equation (20); Equation (19) in, p j Indicates observation task j The gains obtained; Equation (20) in, d p For the task p The duration of observation work, trans p,s ( t p + d p , t s To complete the mission for the satellite p The attitude after the start of the observation mission s The time required to adjust the posture; if the task p and s The constraint that satisfies equation (20) is called reachable, denoted as p → s Otherwise it is considered unreachable; Individual selection is specifically based on ranking the individuals in the population from largest to smallest objective function value. Before selection... b % of individuals participate in cross-operations as elite individuals.

7. The ultra-agile Earth observation satellite full-precision model scheduling method according to claim 1, characterized in that: The mutation operator in step three adopts the idea of ​​prioritizing the observations on the optimal sequence, and performs the mutation operation on each individual in the population. The specific adjustment method is as follows: For each individual... pop i The following observation pop i,j First, determine the task. j Whether in an individual pop i The optimal sequence after minimum graph decoding s i If above, then the mutation probability is set to probability. PM Otherwise, the mutation probability is set to probability 1- PM Generate a random number between [0,1]. If this number is less than the preset mutation probability, then the task... j Mutation; otherwise, task j The observation information remains unchanged; Observation pop i,j The mutation operation involves randomly selecting a time window. vtw j i Then randomly select the observation window. vtw j i Start time within t j i .

8. A super-agile satellite full-precision model scheduling system, wherein the system executes the super-agile Earth observation satellite full-precision model scheduling method according to any one of claims 1-7, characterized in that: The system includes a full-precision model attitude switching time calculation module, a global search module, and a local search module. The global search module is an evolutionary algorithm based on discrete integer encoding and minimum graph decoding techniques. The local search module inserts tasks and constructs solutions based on the attitude switching satisfaction deviation values ​​of adjacent tasks on both sides. The full-precision model attitude switching time calculation module completes model task planning, solving, and scheduling based on a meme algorithm framework that combines global search and greedy insertion local search algorithms.

Images