Satellite flexible attitude dynamics modeling method considering whole satellite centroid variation

Abstract

The invention provides a satellite flexible attitude dynamics modeling method considering the change of the mass center of a whole satellite, and the model is used for simulating the on-orbit flexible attitude behavior of the satellite. The method comprises the following steps: establishing an inertial coordinate system Oi, an orbital coordinate system Oo, a star coordinate system Ob and a sailboard coordinate system Of; establishing an attitude motion equation Eq.1 of a whole satellite system, an equation Eq.2 of flexible vibration of the sailboard and an attitude motion equation Eq.3 of the sailboard relative to a central body, wherein the Eq.1, Eq.2 and Eq.3 consider the influence of the mass center change of the whole satellite and the flexible vibration of the large sailboard; and calculating a direction cosine matrix RoI of the centroid orbital coordinate system, a direction cosine matrix RsI of the star coordinate system, an attitude conversion matrix Rso of the star coordinate system relative to the centroid orbital coordinate system and an Euler angle, and taking the Euler angle as attitude control input. According to the satellite attitude dynamics calculation method considering the satellite centroid change and the sailboard elastic deformation in real time, the problem that the instantaneous centroid change cannot be simulated in an existing attitude dynamics algorithm is solved.

Description

Technical Field

[0001] The present invention relates to the technical field of satellite simulation, in particular to a satellite flexible attitude dynamics modeling method, and in particular to a satellite flexible attitude dynamics modeling method that takes into account the change of the mass center of the entire satellite. Background Art

[0002] With the continuous advancement of aerospace technology and the expansion of satellite applications, the energy supply demand of modern satellite systems is increasing. The use of large solar panels is not only an inevitable technological trend but also a necessary means to meet the high power consumption and long life requirements of satellites. However, as the panel area increases, its mass also increases significantly, thereby increasing the proportion of large panels in the total satellite mass. This change in mass distribution not only alters the satellite's inertial parameters but also complicates the dynamic characteristics of the entire satellite. In particular, during the process of attitude adjustment and sun tracking by large panels, the position of the center of mass of the entire satellite changes dynamically, which places higher demands on the satellite's attitude control. In addition, due to their large size and complex structure, large panels are highly susceptible to various disturbances in space, causing flexible vibrations. These flexible vibrations not only affect the stability and energy conversion efficiency of the panels themselves but also, through structural coupling, significantly affect the attitude motion of the entire satellite.

[0003] Therefore, to accurately describe and predict the attitude behavior of satellites during on-orbit operation and improve the accuracy and stability of satellite attitude control, a satellite flexible attitude dynamics modeling method that fully considers the influence of the satellite's center of mass fluctuations and the flexible vibration of large sail panels is required. Existing attitude dynamics modeling methods, whether based on vector mechanics such as the Newton-Euler method, analytical mechanics such as the Lagrange and Hamilton methods, or methods combining the characteristics of vector mechanics and analytical mechanics such as the Kane method, all require the establishment of a coordinate system during the modeling process. Based on the selected coordinate system, multibody system modeling methods can be further categorized as hybrid coordinate methods, rotational coordinate methods, and inertial coordinate system methods. Furthermore, the description of the deformation of flexible structures is also crucial to the accuracy of system dynamics modeling. Currently, the main description methods include distributed parameter methods, discrete parameter methods, and modal coordinate methods. In summary, current research on flexible multibody system dynamics modeling methods is relatively mature.

[0004] However, the current dynamic attitude equations based on the reference frame with the center of mass of the entire satellite as the origin all assume that the origin does not change during the attitude motion. But in fact, the center of mass of the system changes dynamically. Especially in satellite systems where the mass of the sailboard accounts for a large proportion, the impact of this center of mass change cannot be ignored. Therefore, there is an urgent need for a satellite attitude dynamics calculation method that can consider the center of mass change while taking into account the flexibility of the sailboard. Summary of the Invention

[0005] To address some or all of the issues in the prior art, the present invention provides a satellite flexible attitude dynamics modeling method that considers the changes in the satellite's center of mass. This method uses Newton's method to establish the attitude mechanics equations for the entire satellite system, Lagrangian methods to establish the equations for the flexible vibration of the sailboard itself, hybrid coordinates to describe the coupling between the elastic deformation of the sailboard and the attitude motion of the central body, and finite element methods to discretize the flexible deformation of the sailboard. This satellite attitude dynamics calculation method considers the changes in the satellite's center of mass and the elastic deformation of the sailboard in real time, overcoming the problem of existing attitude dynamics algorithms that cannot simulate instantaneous center of mass changes.

[0006] The present invention provides a satellite flexible attitude dynamics modeling method taking into account the change of the center of mass of the entire satellite, the method comprising the following steps:

[0007] S1: Establish inertial (global) coordinate system O i , orbital (reference) coordinate system O o , Star (satellite) coordinate system O b , sailboard (body) coordinate system O f ;

[0008] S2: Establish the attitude motion equation of the entire satellite system, Eq.1, the equation for the flexible vibration of the sailboard itself, Eq.2, and the attitude motion equation of the sailboard relative to the central body, Eq.3. Eq.1, Eq.2, and Eq.3 take into account the influence of the center of mass change of the entire satellite and the flexible vibration of the large sailboard;

[0009] S3: Start calculating satellite attitude dynamics, with position vector As the state variable, solve the orbit determination equation to obtain the satellite's velocity at any time and location Then calculate the direction cosine matrix R of the center of mass orbit coordinate system oI ;

[0010] S4: The generalized angular velocity ω of the central body s , the generalized coordinate η of the flexible vibration of the sailboard Ai and its derivatives constitute the state variables And solve the attitude equation to get the satellite's ω s , then ω s After converting to attitude velocity, integrate to get attitude quaternion, and then calculate the direction cosine matrix R of the stellar coordinate system sI ;

[0011] S5: Calculate the attitude transformation matrix R of the star coordinate system relative to the center of mass orbit coordinate system so , calculate the Euler angle according to the conversion formula from direction cosine matrix to Euler angle, and use the two-variable inverse tangent function arctan2 to expand the value range, and the obtained Euler angle is used as the attitude control input;

[0012] S6: After each step iterative calculation, determine whether the simulation end time has been reached. If not, return to S3 to enter the next iteration. If the simulation end time has been reached, terminate the entire calculation process, complete the satellite flexible attitude dynamics modeling, and generate a satellite flexible attitude dynamics model; and

[0013] S7: simulating the flexible attitude behavior of the satellite on orbit using the satellite flexible attitude dynamics model.

[0014] Furthermore, the direction cosine matrix R of the mass center orbit coordinate system oI It is used to describe the directional relationship between the center-of-mass orbit coordinate system and other coordinate systems, which is of great significance for the subsequent analysis of the satellite's attitude and motion on orbit; the direction cosine matrix R of the stellar coordinate system sI It is used to describe the directional relationship between the celestial coordinate system and other coordinate systems, which is crucial for analyzing the attitude of the satellite itself; the attitude transformation matrix R of the celestial coordinate system relative to the center of mass orbit coordinate system is so Describes the rotation relationship between two coordinate systems.

[0015] Furthermore, in step S1, the orbital coordinate system O o -X o Y o Z o The origin is located at the center of mass of the satellite, X o The Y axis is along the tangent direction of the satellite's orbit and points to the direction of the satellite's movement. o The axis points in the direction of the normal of the track, Z o Given by the right-hand rule, perpendicular to the orbital plane; the stellar coordinate system O b -X b Y b Z b The origin of the orbital coordinate system O o Origin coincides, X b Axis, Y b Axis and Z b The axis is fixed on the entire satellite; the sailboard coordinate system O f -X f Y f Z f The origin is located at the connection between the sailboard and the satellite body.

[0016] Furthermore, in step S2, Eq.1, Eq.2 and Eq.3 are respectively:

[0017]

[0018] Where I s is the inertia tensor of the entire star system in the star coordinate system; ωs and is the generalized angular velocity and angular acceleration of the central body in the stellar coordinate system; C s is the static moment matrix of the entire star system relative to the zero-position center of mass of the system in the star coordinate system; r b0 、 and h is the position vector of the instantaneous center of mass of the entire star system pointing to the zero center of mass in the star coordinate system and its first and second order derivatives (velocity, acceleration); wj and is the angular momentum and control torque of momentum wheel j in the stellar coordinate system; J Ais and J sAi is the mutual coupling inertia tensor between the entire satellite system and the sailboard in the satellite and sailboard coordinate systems; F Ait is the coupling coefficient matrix of the sailboard vibration and its translation in the sailboard coordinate system; F Air is the coupling coefficient matrix of the sailboard vibration and its rotation in the sailboard coordinate system; F Ais is the coupling coefficient matrix of the sailboard vibration to the rotation of the entire star system in the star coordinate system; F Aits is the coupling coefficient matrix of the sailboard vibration and its translation in the stellar coordinate system; η Ai 、 and is the generalized coordinates of the flexible vibration of the sailboard and its first and second order derivatives; Ω Ai is the windsurfing board modal / vibration shape frequency matrix; ξ Ai is the windsurfing board modal / mode damping ratio coefficient / matrix; I Ai is the inertia tensor of the sailboard relative to the mounting point in the sailboard coordinate system; C Ais is the static moment matrix of the sailboard in the stellar coordinate system; and is the driving angular velocity and driving angular acceleration of the sailboard in the sailboard coordinate system; g s 、g Ai and g BAi is the generalized force exerted by the entire star, the sail and the central body on the sail.

[0019] Furthermore, the Eq.1 includes an inertia term, a center of mass motion-related term, and a sailboard-related term. The inertia term is Shows I s With ω s The inertial effect related to its derivative reflects the inertial characteristics of the satellite during rotation; the center of mass motion related term is This part considers the influence of center of mass motion on satellite attitude, which reflects the consideration of instantaneous mass in the present invention.

[0020]

[0021] Furthermore, Eq. 2 includes vibration-related terms and coupling terms between the sailboard and the entire satellite. The vibration-related terms are: It is a typical vibration equation, which reflects the vibration characteristics of the sailboard itself; the coupling term between the sailboard and the whole satellite is This item takes into account the coupling effects of the angular velocity change of the center body, the driving angular velocity change of the sailboard, and the acceleration of the center of mass on the flexible vibration of the sailboard.

[0022] Furthermore, the Eq. 3 includes the terms related to inertia and center of mass motion and the windsurfing board's own motion. The terms related to inertia and center of mass motion are: This term takes into account the influence of the angular velocity of the center body, the driving angular velocity of the sailboard and the center of mass motion on the attitude of the sailboard relative to the center body; the sailboard's own motion term is This item takes into account the influence of the sailboard's own rotation and flexible vibration on its attitude relative to the center body.

[0023] Furthermore, step S3 marks the start of the satellite attitude dynamics calculation process.

[0024] Furthermore, in step S3, the satellite at any time and Obtained by:

[0025] The orbit determination equation is based on Newton's second law. Without loss of generality, a simplified satellite orbit motion equation is given here that only considers the influence of J2 perturbation (higher precision orbit determination equations can consider more perturbation terms):

[0026]

[0027] Where: μ = Gm e =398600.44km 3 / s 2 is the Earth's gravitational constant; R e =6371km is the average equatorial radius; J2 = 1082.63×10 -6 is the perturbation coefficient with harmonic terms;

[0028] remember The orbit determination equation can be simplified as The initial position is given at the beginning of the simulation and speed By using numerical integration method (as long as the integration step is small enough, various numerical integration methods can guarantee convergence), the orbit determination equation is integrated once to obtain the satellite's orbit at any time. Performing quadratic integration on the orbit determination equation yields In this way, the motion state of the satellite at different times can be obtained.

[0029] Furthermore, the numerical integration method includes the Euler method, the Runge-Kutta method, etc.

[0030] Furthermore, the R oI Get it as follows:

[0031] Direction cosine matrix (transpose) of the center-of-mass orbital coordinate system

[0032] In the formula

[0033] Furthermore, in step S4, solving the posture equation is mainly divided into three sub-steps:

[0034] ① Solve the dynamic equations using generalized angular velocity coordinates to analyze the influence of the satellite's motion and force on its attitude;

[0035] ② Update the kinematic equations using Euler quaternion coordinates, using the mathematical tool quaternion to describe and update the satellite's attitude changes;

[0036] ③ Update the control equation with Euler angle parameter coordinates to convert the attitude information into parameters that are easy for the control system to process.

[0037] Furthermore, in step S4, the ω s Obtained by:

[0038] Since the sailboard is fully controlled, its motion law is known, namely ω Ai 、 It is known that the posture motion equations Eq.1 and Eq.2 can be simplified as The matrix equation is as follows:

[0039]

[0040] Given the initial value p(0), the above matrix equation is solved by numerical integration method (such as explicit Runge-Kutta method or implicit backward Euler method) to obtain

[0041] Furthermore, the left side of the matrix equation is the derivative of the state variable p, and the right side of the equation is composed of multiple matrix multiplications and additions, including I s 、F Aits , center of mass motion related items, g Ai , sailboard-related items, etc., which reflect the comprehensive effect of multiple factors in satellite attitude dynamics.

[0042] Furthermore, the ω s Converted to attitude velocity:

[0043]

[0044] Furthermore, the attitude quaternion is [q4 q1 q2 q3] T .

[0045] Furthermore, according to the attitude quaternion [q4 q1 q2 q3] T And the following formula, to obtain the R sI :

[0046]

[0047] Furthermore, in step S5, the R so Get it as follows:

[0048]

[0049] Furthermore, the conversion formula from the direction cosine matrix to Euler angles is:

[0050]

[0051] Furthermore, considering that the range of arcsin and arctan is It is not enough to cover the entire sphere. We can use the bivariate inverse tangent function arctan2 to replace the inverse tangent function and expand its range to [-π,π]. so The corresponding Euler angles.

[0052] Furthermore, in step S7, the specific steps of simulating the flexible attitude behavior of the satellite on orbit using the satellite flexible attitude dynamics model include:

[0053] 1) Satellite overall structural parameter calculation: Through 3D modeling and mass estimation, the approximate overall mass, static moment, moment of inertia, and the position coordinates of the center of mass in the satellite coordinate system are measured;

[0054] 2) Pre-calculation of satellite flexibility parameters: Through 3D modeling, meshing, boundary constraint application, and finite element analysis, the natural frequencies and vibration modes of the flexible sailboards, as well as their coupling coefficients to the translational and rotational motions of the entire satellite, are measured.

[0055] 3) Setting the initial motion state of the satellite: by obtaining the real telemetry data on the satellite, selecting the appropriate time point, and calibrating the initial motion state of the satellite;

[0056] 4) Satellite space environment parameter setting: by obtaining the real telemetry data on the satellite, calibrate the satellite space environment state at the start of the simulation;

[0057] 5) Select an appropriate simulation step size (10ms) and synchronously update the orbit position and attitude position using the orbit determination equation and attitude motion equation;

[0058] 6) Based on the updated orbital position and attitude position, the time-varying structural parameters in the satellite flexible dynamics are recalculated and the next iterative calculation is performed until the end.

[0059] Furthermore, in step 3), the real telemetry data on the satellite includes: orbital position, body attitude quaternion / Euler angle, body attitude angular velocity, sailboard angle and angular velocity, etc.

[0060] Furthermore, in step 4), the second on-board real telemetry data includes: ephemeris data, atmospheric data, optical pressure data, geomagnetic data, etc.

[0061] Furthermore, in step 5), the numerical integration method involved in the orbit determination equation and the attitude motion equation is a 4th-order fixed-step Runge-Kutta method.

[0062] Furthermore, in step 6), the time-varying structural parameters include the moment of inertia of the entire satellite, the position of the center of mass of the entire satellite, etc.

[0063] This paper establishes four reference coordinate systems: an inertial coordinate system, an orbital coordinate system, a stellar coordinate system, and a sailboard coordinate system. Based on force analysis of the satellite and the sailboard, the paper establishes the attitude motion equations for the entire satellite system (Eq. 1), the equation for the flexible vibration of the sailboard itself (Eq. 2), and the equation for the attitude motion of the sailboard relative to the central body (Eq. 3). These equations contain numerous parameters related to the satellite and the sailboard, such as the inertia tensor, static moment matrix, and coupling coefficient matrix. This paper proposes a satellite flexible attitude dynamics modeling method that considers the center of mass change of the entire satellite. This model is used to simulate the flexible attitude behavior of a satellite on orbit. Comparisons have shown that the attitude accuracy calculated using this model, when compared to actual telemetry, has a maximum error of no more than 13.7%. This paper provides important theoretical support and technical assurance for achieving precise satellite attitude control and optimizing satellite design, overcoming the inability of existing attitude dynamics algorithms to simulate instantaneous center of mass changes.

[0064] The present invention has at least the following beneficial effects: 1) by establishing the attitude motion equation of the entire satellite system, the flexible vibration equation of the sailboard itself, and the attitude motion equation of the sailboard relative to the central body, factors such as the instantaneous center of mass change of the satellite and the flexible vibration of the sailboard are fully considered, and the attitude dynamic characteristics of the satellite in space can be described more accurately. It is more in line with the actual operation conditions than the traditional modeling method, and the generated model is more accurate. For example, the model can be closer to the actual satellite attitude dynamic parameters; 2) multiple reference coordinate systems are constructed, and relevant equations and calculation processes are derived based on them, so that the analysis of the motion and interaction relationship of various parts of the satellite is more accurate, providing a more reliable model basis for satellite attitude control and orbit calculation, so that the generated model can improve the control accuracy of the satellite in orbit operation, attitude adjustment, etc. during simulation. BRIEF DESCRIPTION OF THE DRAWINGS

[0065] To further illustrate the above and other advantages and features of various embodiments of the present invention, a more detailed description of various embodiments of the present invention will be presented with reference to the accompanying drawings. It will be understood that these drawings depict only typical embodiments of the present invention and are not to be considered as limiting the scope thereof. In the drawings, for clarity, identical or corresponding components will be represented by the same or similar reference numerals.

[0066] Figure 1 shows a schematic diagram of the reference coordinate system in the present invention;

[0067] Figure 2 The figure shows the schematic diagram of the force analysis of the whole satellite and sailboard;

[0068] Figure 3 The diagram shows the implementation principle of the present invention. DETAILED DESCRIPTION

[0069] In the following description, the present invention is described with reference to various embodiments. However, those skilled in the art will recognize that the various embodiments can be implemented without one or more of the specific details or with other alternative and / or additional methods, materials, or components. In other cases, well-known structures, materials, or operations are not shown or described in detail to avoid obscuring the inventive aspects of the present invention. Similarly, for the purpose of explanation, specific quantities, materials, and configurations are described to provide a comprehensive understanding of the embodiments of the present invention. However, the present invention is not limited to these specific details. In addition, it should be understood that the various embodiments shown in the drawings are illustrative representations and are not necessarily drawn to scale.

[0070] In this specification, reference to "one embodiment" or "the embodiment" means that a particular feature, structure, or characteristic described in connection with the embodiment is included in at least one embodiment of the present invention. The appearances of the phrase "in one embodiment" in various places in this specification are not necessarily all referring to the same embodiment.

[0071] It should be noted that the embodiments of the present invention describe the process steps in a specific order. However, this is only for the purpose of illustrating the specific embodiment and does not limit the order of the steps. On the contrary, in different embodiments of the present invention, the order of the steps can be adjusted according to the process.

[0072] The following embodiment establishes four reference coordinate systems: inertial coordinate system, orbital coordinate system, celestial coordinate system, and sailboard coordinate system ( Figure 1 ), according to the force analysis of the whole satellite and sailboard ( Figure 2 ), establish the attitude motion equation of the entire satellite system Eq.1, the equation for the flexible vibration of the sailboard itself Eq.2, and the attitude motion equation of the sailboard relative to the central body Eq.3. The equations contain many parameters related to the satellite and the sailboard, such as the inertia tensor, static moment matrix, coupling coefficient matrix, etc. Eq.1, Eq.2 and Eq.3 are respectively:

[0073]

[0074] Where I s is the inertia tensor of the entire star system in the star coordinate system; ω s and is the generalized angular velocity and angular acceleration of the central body in the stellar coordinate system; C s is the static moment matrix of the entire star system relative to the zero-position center of mass of the system in the star coordinate system; r b0 、 and h is the position vector of the instantaneous center of mass of the entire star system pointing to the zero center of mass in the star coordinate system and its first and second order derivatives (velocity, acceleration); wj and is the angular momentum and control torque of momentum wheel j in the stellar coordinate system; J Ais and J sAi is the mutual coupling inertia tensor between the entire satellite system and the sailboard in the satellite and sailboard coordinate systems; F Ait is the coupling coefficient matrix of the sailboard vibration and its translation in the sailboard coordinate system; F Air is the coupling coefficient matrix of the sailboard vibration and its rotation in the sailboard coordinate system; F Ais is the coupling coefficient matrix of the sailboard vibration to the rotation of the entire star system in the star coordinate system; F Aits is the coupling coefficient matrix of the sailboard vibration and its translation in the stellar coordinate system; η Ai 、 and is the generalized coordinates of the flexible vibration of the sailboard and its first and second order derivatives; Ω Ai is the windsurfing board modal / vibration shape frequency matrix; ξ Ai is the windsurfing board modal / mode damping ratio coefficient / matrix; I Aiis the inertia tensor of the sailboard relative to the mounting point in the sailboard coordinate system; C Ais is the static moment matrix of the sailboard in the stellar coordinate system; and is the driving angular velocity and driving angular acceleration of the sailboard in the sailboard coordinate system; g s 、g Ai and g BAi is the generalized force exerted by the entire star, the sail and the central body on the sail.

[0075] Eq.1 includes inertia terms, center of mass motion related terms, and sailboard related terms. The inertia term is Shows I s With ω s The inertial effect related to its derivative reflects the inertial characteristics of the satellite during rotation; the center of mass motion related item is This part considers the influence of the center of mass motion on the satellite attitude, which reflects the characteristic of considering the instantaneous center of mass change in the present invention;

[0076] The influence of the sailboard's flexible vibration on the entire satellite's attitude is considered. Eq.2 includes the vibration-related term and the coupling term between the sailboard and the entire satellite. The vibration-related term is It is a typical vibration equation, which reflects the vibration characteristics of the sailboard itself; the coupling term between the sailboard and the whole satellite is This term takes into account the coupling effect of the angular velocity change of the center body, the driving angular velocity change of the sailboard, and the acceleration of the center of mass on the flexible vibration of the sailboard. Eq.3 includes the terms related to the inertia and center of mass motion and the sailboard's own motion. The terms related to the inertia and center of mass motion are: This term takes into account the influence of the angular velocity of the center body, the driving angular velocity of the sailboard, and the center of mass motion on the attitude of the sailboard relative to the center body; the sailboard's own motion term is This item takes into account the influence of the sailboard's own rotation and flexible vibration on its attitude relative to the center body.

[0077] The technical solution of the present invention will be further described below in conjunction with the accompanying drawings of the embodiments.

[0078] A satellite flexible attitude dynamics modeling method considering the change of the center of mass of the entire satellite includes the following steps:

[0079] Establish inertial coordinate system O i , orbital coordinate system O o , stellar coordinate system O b , windsurfing coordinate system O f , Figure 1 The reference coordinate system diagram of the present invention is shown in FIG. A is a dotted line representing the trajectory of the satellite in orbit. o -X o Yo Z o The origin is located at the center of mass of the satellite, X o The Y axis is along the tangent direction of the satellite's orbit and points to the direction of the satellite's movement. o The axis points in the direction of the normal of the track, Z o Given by the right-hand rule, perpendicular to the orbital plane; the stellar coordinate system O b -X b Y b Z b The origin and orbital coordinate system O o Origin coincides, X b Axis, Y b Axis and Z b The axis is fixed on the entire satellite; the sailboard coordinate system O f -X f Y f Z f The origin of the inertial coordinate system O is located at the connection between the sailboard and the satellite body. i Provides an absolute reference base for the entire system, used to describe the absolute position and movement of the satellite in space; orbital coordinate system O o The coordinate axis direction is related to the orbital motion of the satellite and is used to describe the motion state of the satellite in orbit; the stellar coordinate system O b The origin of the satellite will rotate with the rotation of the satellite and is used to describe the satellite's own attitude and internal motion; the sailboard coordinate system O f Used to describe the motion and attitude of the sailboard relative to the entire star, Figure 1 The velocity V of the sailboard is also marked. The relative position and direction relationship between these coordinate systems helps to establish the satellite attitude dynamics model and analyze the motion of the satellite body, sailboard and the interaction between them.

[0080] Figure 2 The figure shows the force analysis diagram of the whole satellite and sailboard. Figure 2 The rectangle on the left represents the satellite's central body, and the polygon on the right represents the satellite's sail. Figure 2 Middle: B refers to the centrosome, A i Refers to the No. i windsurfing board, For B in O b The external force For B in O b The external torque at A i In O fi The external force A i In O fi The external torque at For B and A i In O fi The internal force of For B and A i In O fi Based on the force analysis of the entire satellite and the sailboard, the attitude motion equation of the entire satellite system, Eq.1, the equation for the flexible vibration of the sailboard itself, Eq.2, and the attitude motion equation of the sailboard relative to the central body, Eq.3, are established. Eq.1, Eq.2, and Eq.3 take into account the influence of the center of mass change of the entire satellite and the flexible vibration of the large sailboard;

[0081] Start calculating satellite attitude dynamics, with position vector As the state variable, the orbit determination equation is:

[0082]

[0083] Abbreviated as Where: μ = Gm e =398600.44km 3 / s 2 is the Earth's gravitational constant R e =6371km is the average equatorial radius; J2 = 1082.63×10 -6 is the perturbation coefficient with harmonic terms;

[0084]

[0085] At the beginning of the simulation, the initial position of the satellite needs to be given and initial velocity Then, a reasonable numerical integration method is used (as long as the integration step is small enough, various numerical integration methods can guarantee convergence, such as the Euler method and the Runge-Kutta method) to perform a single and double integration on the orbital motion equation. The single integration can obtain the velocity of the satellite at any time. The second integral can be used to obtain the position of the satellite at any time In this way, the motion state of the satellite at different moments can be calculated;

[0086] Then calculate the direction cosine matrix R of the center of mass orbit coordinate system oI , R oI It is used to describe the directional relationship between the center-of-mass orbital coordinate system and other coordinate systems, which is of great significance for the subsequent analysis of the satellite's attitude and motion on orbit. The direction cosine matrix (transpose) of the center-of-mass orbital coordinate system is expressed as:

[0087]

[0088] Where,

[0089] Then calculate the direction cosine matrix R of the stellar coordinate system sI , R sIIt is used to describe the directional relationship between the stellar coordinate system and other coordinate systems, which is crucial for analyzing the attitude of the satellite itself. Since the sailboard is fully controlled, its motion law is known, that is, ω Ai 、 It is known that the generalized angular velocity ω of the central body s , the generalized coordinate η of the flexible vibration of the sailboard Ai and its derivatives constitute the state variables The posture motion equations Eq.1 and Eq.2 are simplified as The matrix equation is as follows (the left side of the matrix equation is the derivative of the state variable p, and the right side of the equation is composed of multiple matrix multiplications and additions, including I s 、F Aits , center of mass motion related items, g Ai , sailboard-related items, etc., reflecting the comprehensive effect of multiple factors in satellite attitude dynamics. ):

[0090]

[0091] Given the state variable p(0), solving the above matrix equation by a reasonable numerical integration method (such as the explicit Runge-Kutta method or the implicit backward Euler method) yields Then the obtained ω s Converted to attitude velocity:

[0092]

[0093] Integrating the above attitude velocity, we get the attitude quaternion [q4 q1 q2 q3] T , according to the attitude quaternion [q4q1 q2 q3] T And the following formula, find R sI :

[0094]

[0095] According to the direction cosine matrix R of the mass center orbit coordinate system oI And the direction cosine matrix R of the stellar coordinate system sI , R so Describes the rotation relationship between the two coordinate systems, and the attitude transformation matrix R of the stellar coordinate system relative to the center of mass orbit coordinate system is calculated by the following formula so :

[0096]

[0097] The Euler angle is calculated according to the conversion formula from direction cosine matrix to Euler angle:

[0098]

[0099] Considering that the range of arcs in and arctan is It is not enough to cover the entire sphere. We can use the bivariate inverse tangent function arctan2 to replace the inverse tangent function and expand its range to [-π,π]. so The corresponding Euler angles, where ψ is the yaw angle, is the pitch angle, θ is the roll angle, and these Euler angles can intuitively represent the attitude of the satellite and can be used as attitude control input to control the attitude adjustment of the satellite.

[0100] like Figure 3 As shown in the figure, it is a schematic diagram of the implementation principle of the present invention. After each step-by-step iterative calculation, it is determined whether the simulation end time has been reached. If not (i.e., entering a loop), it returns to the step of "solving the orbit determination equation", continues the loop calculation, and enters the next iteration. If the simulation end time has been reached (i.e., not entering a loop), the entire calculation process is terminated, the satellite flexible attitude dynamics modeling is completed to generate a satellite flexible attitude dynamics model, and the satellite flexible attitude dynamics model is used to simulate the flexible attitude behavior of the satellite in orbit.

[0101] The following uses a real satellite model as an example to illustrate the implementation process of the above principles:

[0102] 1) Satellite overall structural parameter calculation: Through 3D modeling and mass estimation, the satellite's approximate overall mass, static moment, moment of inertia, and the position coordinates of the center of mass in the satellite coordinate system are measured;

[0103] 2) Pre-calculation of satellite flexibility parameters: Through 3D modeling, meshing, boundary constraint application, and finite element analysis, the natural frequencies and vibration modes of the flexible sailboards, as well as their coupling coefficients to the translational and rotational motions of the entire satellite, are measured.

[0104] 3) Setting the initial motion state of the satellite: By obtaining the real telemetry data on the satellite, selecting the appropriate time point, and calibrating the initial motion state of the satellite, the real telemetry data on the satellite includes: orbital position, body attitude quaternion / Euler angle, body attitude angular velocity, sailboard rotation angle and angular velocity, etc.

[0105] 4) Space environment parameter setting: The satellite space environment state at the start of the simulation is calibrated by obtaining the second real telemetry data on the satellite. The second real telemetry data on the satellite includes: ephemeris data, atmospheric data, optical pressure data, geomagnetic data, etc.

[0106] 5) Select an appropriate simulation step size (10ms) and synchronously update the orbit position and attitude position using the aforementioned orbit determination equations and attitude equations (the core calculation steps are as described above, and the numerical integration method used is the 4th-order fixed-step Runge-Kutta method);

[0107] 6) Based on the updated values, recalculate the time-varying structural parameters in the satellite's flexible dynamics, such as the moment of inertia of the entire satellite and the position of the center of mass of the entire satellite, and perform the next iterative calculation until the end.

[0108] Comparisons have shown that the maximum error between the attitude accuracy calculated using this method and actual telemetry measurements is no more than 13.7%. This demonstrates that this method has a certain degree of accuracy and reliability in measuring satellite attitude. While some error exists, it is within an acceptable range and can provide a valuable reference for satellite attitude control and related research.

[0109] The present invention establishes four reference coordinate systems, namely the inertial coordinate system, the orbital coordinate system, the stellar coordinate system, and the sailboard coordinate system. Based on the force analysis of the entire satellite and the sailboard, the attitude dynamics equations of the entire satellite and the sailboard are established. The equations contain many parameters related to the satellite and the sailboard, such as the inertia tensor, the static moment matrix, the coupling coefficient matrix, etc., and proposes a satellite flexible attitude dynamics modeling method that takes into account the change of the center of mass of the entire satellite. This method will provide important theoretical support and technical guarantee for achieving precise control of satellite attitude and optimizing satellite design, and overcome the problem that the existing attitude dynamics algorithm cannot simulate instantaneous center of mass changes.

[0110] Although various embodiments of the present invention have been described above, it should be understood that they are presented by way of example only and not limitation. It will be apparent to those skilled in the relevant art that various combinations, modifications, and variations may be made thereto without departing from the spirit and scope of the present invention. Therefore, the breadth and scope of the present invention disclosed herein should not be limited by the exemplary embodiments disclosed above, but should be defined solely in accordance with the appended claims and their equivalents.

Claims

1. A satellite flexible attitude dynamics modeling method considering the change of the center of mass of the entire satellite, characterized by: The method comprises the following steps: S1: Establish inertial coordinate system O i , orbital coordinate system O o , stellar coordinate system O b , windsurfing coordinate system O f ; S2: Establish the attitude motion equation of the entire satellite system, Eq.1, the equation for the flexible vibration of the sailboard itself, Eq.2, and the attitude motion equation of the sailboard relative to the central body, Eq.

3. Eq.1, Eq.2, and Eq.3 take into account the influence of the center of mass change of the entire satellite and the flexible vibration of the large sailboard; S3: Start calculating satellite attitude dynamics, with position vector As the state variable, solve the orbit determination equation to obtain the satellite's velocity at any time and location Then calculate the direction cosine matrix R of the center of mass orbit coordinate system oI ; S4: The generalized angular velocity ω of the central body s , the generalized coordinate η of the flexible vibration of the sailboard Ai and its derivatives constitute the state variables And solve the attitude equation to get the satellite's ω s , then ω s After converting to attitude velocity, integrate to get attitude quaternion, and then calculate the direction cosine matrix R of the stellar coordinate system sI ; S5: Calculate the attitude transformation matrix R of the star coordinate system relative to the center of mass orbit coordinate system so , calculate the Euler angle according to the conversion formula from direction cosine matrix to Euler angle, and use the two-variable inverse tangent function arctan2 to expand the value range, and the obtained Euler angle is used as the attitude control input; S6: After each step iterative calculation, determine whether the simulation end time has been reached. If not, return to S3 to enter the next iteration. If the simulation end time has been reached, terminate the entire calculation process, complete the satellite flexible attitude dynamics modeling, and generate a satellite flexible attitude dynamics model; and S7: simulating the flexible attitude behavior of the satellite on orbit using the satellite flexible attitude dynamics model.

2. The satellite flexible attitude dynamics modeling method considering the change of the center of mass of the entire satellite according to claim 1 is characterized in that: In step S1, the orbital coordinate system O o -X o Y o Z o The origin is located at the center of mass of the satellite, X o The Y axis is along the tangent direction of the satellite's orbit and points to the direction of the satellite's movement. o The axis points in the direction of the normal of the track, Z o Given by the right-hand rule, perpendicular to the orbital plane; the stellar coordinate system O b -X b Y b Z b The origin of the orbital coordinate system O o Origin coincides, X b Axis, Y b Axis and Z b The axis is fixed on the entire satellite; the sailboard coordinate system O f -X f Y f Z f The origin is located at the connection between the sailboard and the satellite body.

3. The satellite flexible attitude dynamics modeling method considering the change of the center of mass of the entire satellite according to claim 1 is characterized in that: In step S2, Eq.1, Eq.2 and Eq.3 are respectively: Where I s is the inertia tensor of the entire star system in the star coordinate system; ω s and is the generalized angular velocity and angular acceleration of the central body in the stellar coordinate system; C s is the static moment matrix of the entire star system relative to the zero-position center of mass of the system in the star coordinate system; r b0 、 and h is the position vector of the instantaneous center of mass of the entire star system pointing to the zero center of mass in the star coordinate system and its first and second order derivatives; wj and is the angular momentum and control torque of momentum wheel j in the stellar coordinate system; J Ais and J sAi is the mutual coupling inertia tensor between the entire star system and the sailboard in the star and sailboard coordinate systems; F Ait is the coupling coefficient matrix of the sailboard vibration and its translation in the sailboard coordinate system; F Air is the coupling coefficient matrix of the sailboard vibration to its rotation in the sailboard coordinate system; F Ais is the coupling coefficient matrix of the sailboard vibration to the rotation of the entire star system in the star coordinate system; F Aits is the coupling coefficient matrix of the sailboard vibration and its translation in the stellar coordinate system; η Ai 、 and is the generalized coordinates of the flexible vibration of the sailboard and its first and second order derivatives; Ω Ai is the windsurfing board modal / shape frequency matrix; ξ Ai is the windsurfing board modal / mode damping ratio coefficient / matrix; I Ai is the inertia tensor of the sailboard relative to the mounting point in the sailboard coordinate system; C Ais is the static moment matrix of the sailboard in the stellar coordinate system; and is the driving angular velocity and driving angular acceleration of the sailboard in the sailboard coordinate system; g s 、g Ai and g BAi is the generalized force exerted by the entire star, the sail and the central body on the sail.

4. The satellite flexible attitude dynamics modeling method considering the change of the center of mass of the entire satellite according to claim 1 is characterized in that: In step S3, the satellite is and Obtained by: The orbit determination equation is established according to Newton's second law, and the orbit determination equation is abbreviated as The initial position is given at the beginning of the simulation and speed By using numerical integration method, the orbit determination equation is integrated once to obtain the satellite's orbit at any time. Performing quadratic integration on the orbit determination equation yields 5. The satellite flexible attitude dynamics modeling method considering the change of the center of mass of the entire satellite according to claim 4 is characterized in that: The R oI Get it as follows: In the formula 6. The satellite flexible attitude dynamics modeling method considering the change of the center of mass of the entire satellite according to claim 1 is characterized in that: In step S4, the ω s Obtained by: Since the board is fully controlled, ω Ai 、 It is known that the posture motion equations Eq.1 and Eq.2 can be simplified as The matrix equation is as follows: Given the initial value p(0), the above matrix equation is solved by numerical integration method to obtain 7. The satellite flexible attitude dynamics modeling method considering the change of the center of mass of the entire satellite according to claim 6 is characterized in that: The ω s Converted to attitude velocity:

8. The satellite flexible attitude dynamics modeling method considering the change of the center of mass of the entire satellite according to claim 7 is characterized in that: The attitude quaternion is [q4 q1 q2 q3] T ; According to the attitude quaternion [q4 q1 q2 q3] T And the following formula, to obtain the R sI :

9. The satellite flexible attitude dynamics modeling method considering the change of the center of mass of the entire satellite according to claim 8 is characterized in that: In step S5, the R so Get it as follows:

10. The satellite flexible attitude dynamics modeling method considering the change of the center of mass of the entire satellite according to claim 9 is characterized in that: In step S7, the specific steps of simulating the flexible attitude behavior of the satellite on orbit using the satellite flexible attitude dynamics model include: 1) Satellite overall structural parameter calculation: Through 3D modeling and mass estimation, the approximate overall mass, static moment, moment of inertia, and the position coordinates of the center of mass in the satellite coordinate system are measured; 2) Pre-calculation of satellite flexibility parameters: Through 3D modeling, meshing, boundary constraint application, and finite element analysis, the natural frequencies and vibration modes of the flexible sailboards, as well as their coupling coefficients to the translational and rotational motions of the entire satellite, are measured. 3) Setting the initial motion state of the satellite: By acquiring real telemetry data on the satellite, selecting a suitable time point, and calibrating the initial motion state of the satellite; the real telemetry data on the satellite includes: orbital position, body attitude quaternion / Euler angle, body attitude angular velocity, and sailboard rotation angle and angular velocity; 4) Satellite space environment parameter setting: The satellite space environment state at the start of the simulation is calibrated by obtaining the second real telemetry data on the satellite; the second real telemetry data on the satellite includes: ephemeris data, atmospheric data, light pressure data, and geomagnetic data; 5) Selecting a simulation step size of 10ms, and synchronously updating the orbit position and attitude position using the orbit determination equation and attitude motion equation; 6) Based on the updated orbital position and attitude position, the time-varying structural parameters in the satellite flexible dynamics are recalculated and the next iterative calculation is performed until the end.

Images