Satellite orbit prediction method based on physical information neural network, terminal and medium

Abstract

The application relates to the technical field of satellite orbit prediction, and discloses a satellite orbit prediction method based on a physical information neural network, a terminal and a medium. The method comprises the following steps: constructing a sample set based on satellite position data, sun position data and moon position data in the same time period; constructing a physical information neural network; taking historical satellite position data as input, the network predicts the satellite orbit position at a future time; training the physical information neural network by using the sample set, and performing a satellite orbit prediction task by using the trained physical information neural network; the loss function of the training comprises data loss, physical loss and smoothness loss; the data loss is used for measuring the difference between the network prediction output and the actual observation data; the physical loss is jointly calculated by the earth non-spherical perturbation, three-body gravity perturbation and solar radiation pressure physical perturbation suffered by the satellite; and the smoothness loss is used for constraining the smoothness of the predicted satellite orbit. The application improves the satellite orbit prediction accuracy and robustness.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of satellite orbit prediction, and particularly relates to a satellite orbit prediction method based on a physical information neural network, a computer terminal applying the method and a computer readable storage medium. BACKGROUND

[0002] Orbit prediction is a process of predicting the future position and related uncertainty of any given resident space object (RSO). An accurate satellite orbit prediction model can accurately predict the position of a satellite, which is crucial for preventing satellite collisions, improving satellite navigation and positioning accuracy, and enhancing the completion degree of satellite orbit spatial docking.

[0003] Existing orbit prediction methods can be mainly divided into two categories:

[0004] One category is the traditional physical method, which establishes a dynamic model according to the physical principles of the satellite itself and various perturbation forces received during the on-orbit period, and predicts the average state of the satellite at a future time through analytical or numerical solution. The numerical method predicts the position of the satellite at a future time by numerically solving the motion equation of the space object, and has high prediction accuracy, but the calculation cost and time cost are too large. The analytical method uses a closed solution to estimate the motion of the object to provide the position of the satellite at a future time, and has fast calculation but poor prediction accuracy. Due to the uncertainty of various perturbation forces, it is difficult to establish an accurate dynamic model for the traditional method, and the prediction accuracy of the commonly used SGP4 / SDP4 model will decrease with the increase of the prediction time.

[0005] The other category is the method based on artificial intelligence, which mainly has two ideas: one is to correct the prediction error through machine learning on the basis of the traditional dynamic prediction model to improve the accuracy of satellite orbit prediction, but this method will obviously decrease the error correction ability when the prediction time increases, and is still limited by the accuracy of the satellite dynamic model. The other is to use the ability of deep learning to fit high-order nonlinear models without constructing a complex dynamic model, and directly use a deep learning model for satellite orbit prediction. The existing deep learning method has two deficiencies in solving the orbit prediction problem: 1) satellite orbit prediction involves the coupling of multiple perturbation forces, and single orbit data is difficult to fully represent the physical relationship, which leads to a large deviation between the prediction result of the model and the actual orbit; 2) the existing method mainly focuses on modeling the time dependence of single coordinate component data in the orbit data, and ignores the dependence between different coordinate components, which limits the prediction ability of the model. SUMMARY

[0006] To solve the technical problems in the prior art, the present application provides a satellite orbit prediction method based on a physical information neural network, a terminal and a medium, which organically combines deep learning with classical orbit dynamics, introduces multiple physical constraints as additional loss terms, and optimizes them in cooperation with data-driven loss to ensure that the prediction results comply with the basic laws of orbit dynamics, thereby improving the satellite orbit prediction accuracy and robustness.

[0007] To achieve the above object, the present application provides the following technical solutions.

[0008] The present application discloses a satellite orbit prediction method based on a physical information neural network, comprising steps S1-S3.

[0009] S1. Based on satellite position data, sun position data and moon position data within the same time period, a sample set is constructed.

[0010] S2. A physical information neural network is constructed; the network takes satellite historical position data as input to predict the satellite orbit position at a future time.

[0011] S3. The physical information neural network is trained using the sample set, and the trained physical information neural network is used to perform a satellite orbit prediction task.

[0012] The training process of the physical information neural network includes: inputting a satellite position data sequence of a preset length in the sample set into a physical information neural network with a Crossformer model as the network backbone, extracting orbit multi-scale features using cross-variable attention and cross-time attention, and generating a prediction output sequence; calculating data loss, physical loss and smoothness loss based on the prediction output sequence to construct a joint loss, updating the network parameters of the Crossformer model using backpropagation for the joint loss, and iterating in batches until the training is completed; the data loss is used to measure the difference between the network prediction output and the actual observation data; the physical loss is jointly calculated by the earth's non-spherical perturbation, three-body gravitational perturbation and solar pressure physical perturbation; the smoothness loss is used to constrain the smoothness of the predicted satellite orbit.

[0013] As a further improvement of the above scheme, in step S3, the training method of the physical information neural network comprises the following specific steps, i.e. S31-S39.

[0014] S31. Batch selection of training samples from the sample set; each training sample is represented as ; wherein, represents a satellite position data sequence of a preset length, represents the corresponding label, represents The corresponding timestamp information, express The corresponding timestamp information, and These represent the position data of the sun and moon at the same time.

[0015] S32. Extract each training sample... The sequence is used as the encoder input of the Crossformer model, and the sequence is segmented into fixed-length segments using the Crossformer model's tile embedding; where, T The sequence length excluding the label length, for T The satellite position vector at each moment.

[0016] S33. Using the two-stage encoding structure of the Crossformer model, multi-scale features are extracted by alternating between cross-variable attention and cross-temporal attention; among them, cross-variable attention is used to capture the coupling features between the three components x, y, and z, and cross-temporal attention is used to capture the long- and short-term orbital evolution features.

[0017] S34. Generate future step orbital segments at the decoder end of the Crossformer model, and restore them to the predicted output sequence through de-embedding. ; T pred The length of the predicted output sequence is consistent with the label length.

[0018] S35. Calculate data loss :

[0019] ;

[0020] In the formula, and The first and second parts of the label and prediction output are respectively T + j The satellite position vector at each moment.

[0021] S36. Calculate physical losses :

[0022] ;

[0023] In the formula, Indicates the satellite in T + j The predicted total acceleration value at each moment; For the first T + j The total perturbation acceleration of the satellite at a given moment; ||·|| is the Euclidean norm of ·.

[0024] S37. Calculate the smoothness loss :

[0025] ;

[0026] where, is the jerk of the satellite at the T + j +1 time instant, ; is the fixed step size of the adjacent timestamps.

[0027] S38. Construct the joint loss :

[0028] ;

[0029] where, is the loss corresponding to the learnable noise variance; is the loss corresponding to the precision weight.

[0030] S39. Backpropagate the joint loss , calculate the gradients, and simultaneously update the network weights and uncertainty parameters of the Crossformer model according to the optimizer strategy, complete the first iteration, and return to step S31 for the next loop until all rounds of iterative training are completed.

[0031] As a further improvement of the above scheme, in step S36, for the total acceleration prediction value and total perturbation acceleration of the satellite at the T + j time instant, the calculation method includes the following specific steps, i.e., S361~S367.

[0032] S361. Denormalize the predicted output sequence .

[0033] S362. Approximate the total acceleration by second-order central difference:

[0034] .

[0035] S363. Calculate the acceleration of the satellite under the Earth's non-spherical perturbation:

[0036] ;

[0037] where, is the acceleration of the satellite under the Earth's non-spherical perturbation at the T + j time instant, is the geodetic flattening coefficient of the Earth, taking 1.08263e-3. Let be the Earth's gravitational constant, taken as 3.986e14. ; The radius of the Earth's equator is taken as 6.378e3 km; r The distance from the satellite to the Earth's center is obtained through the predicted satellite position vector; , and for T + j The satellite's three-dimensional coordinates at that moment.

[0038] S364. Calculate the acceleration of the satellite under the gravitational perturbation of the Sun as a third body:

[0039] ;

[0040] In the formula, For the first T + j The acceleration of the satellite under the gravitational perturbation of the Sun's third body at a given moment; The gravitational constant of the Sun as a third body; For the first T + j The solar position vector at each moment, through get; for The Euclidean norm.

[0041] S365. Calculate the acceleration of the satellite under the gravitational perturbation of the Moon as a third body:

[0042] ;

[0043] In the formula, For the first T + j The acceleration of the satellite under the gravitational perturbation of the third body on the Moon at a given moment; The third-body gravitational constant of the Moon; For the first T + j The lunar position vector at each instant, through Obtain.

[0044] S366. Calculate the satellite's acceleration under solar radiation pressure perturbation:

[0045] ;

[0046] In the formula, For the first T + j The acceleration of the satellite under solar radiation pressure perturbation at a given moment; It is the solar constant; is the albedo of the satellite; is the reflectance of the satellite; is the shadow factor, = 1 means the satellite is in the lighted region, = 0 means the satellite is in the Earth's shadow region.

[0047] S367. Combine the physical accelerations calculated in steps S363-S366 to form the total perturbation acceleration:

[0048] .

[0049] As a further improvement of the above scheme, in step S38, the precision weight is calculated by the formula:

[0050] ;

[0051] where exp(·) denotes the natural exponential function.

[0052] As a further improvement of the above scheme, step S1 includes the following specific steps, i.e., S11-S13.

[0053] S11. Define the position vector of a single satellite at the observation time as , and set the satellite position data set in the time period as ; similarly, set the sun position data set in the time period as , and the moon position data set as ; where N is the number of observations, i.e., the sequence length; and are the sun and moon position vectors at the observation time , respectively.

[0054] S12. According to the observation requirements, set the observation time sequence , and use ephemeris files or ground tracking systems to obtain the satellite position vector at each observation time point by point, and save it to ; use the same ephemeris files or astronomical calculation library to calculate the sun's center-of-mass position vector and the moon's center-of-mass position vector at each observation time point by point, and generate data sets and .

[0055] S13. Based on the satellite orbit coordinate data set, the sun position data set and the moon position data set, a plurality of samples are generated to form a sample set; wherein each sample includes a satellite position data sequence of a preset length, a label, timestamp information, and sun and moon position data at the same time.

[0056] As a further improvement of the above scheme, in step S1, the following preprocessing is further performed on the data in the sample set:

[0057] The position vectors in the satellite orbit coordinate data set, the sun position data set and the moon position data set are mapped to the same geocentric inertial coordinate system.

[0058] For each coordinate component in the satellite orbit coordinate data set, the sun position data set and the moon position data set, the mean and standard deviation are calculated respectively for normalization processing, thereby obtaining the normalized corresponding data set.

[0059] The satellite orbit coordinate data set is generated in the form of a sliding window to generate a plurality of sequences, the former part of each sequence is data and is used as the input of the physical information neural network, and the latter part is the label and is used for loss calculation with the predicted data of the network output; all samples are divided into a training set, a validation set and a test set according to a set proportion.

[0060] As a further improvement of the above scheme, in step S12, the satellite position data set, the sun position data set and the moon position data set are aligned according to a common index to ensure that there are , and at the same observation time.

[0061] The application also discloses a computer terminal, comprising a memory, a processor and a computer program stored in the memory and executable on the processor, wherein the processor executes the computer program to realize the steps of the satellite orbit prediction method based on the physical information neural network.

[0062] The application also discloses a computer readable storage medium having a computer program stored thereon, wherein the program is executed by a processor to realize the steps of the satellite orbit prediction method based on the physical information neural network.

[0063] Compared with the prior art, the application has the following beneficial effects:

[0064] 1. The satellite orbit prediction method disclosed by the application uses a physical information neural network, so that the model prediction result not only fits the historical trajectory, but also follows the law of orbit dynamics, thereby improving the robustness and interpretability to abnormal states (such as near-earth disturbance and thrust mutation).

[0065] 2. The application uses a physical information neural network with a Crossformer model structure, adopts a two-stage attention layer structure, and captures the coupling relationship between time sequence characteristics and three coordinate components, thereby improving the modeling capability of satellite orbit sequences.

[0066] 3. The application adopts an uncertainty weighting strategy to balance data loss and physical constraints, stabilizes the training process, and speeds up convergence. BRIEF DESCRIPTION OF DRAWINGS

[0067] Figure 1 A flowchart of the satellite orbit prediction method based on the physical information neural network in embodiment 1 of the application.

[0068] Figure 2 A three-dimensional display diagram of the prediction results of the first batch of the test set in embodiment 1 of the application.

[0069] Figure 3 A three-dimensional display diagram of the X-axis coordinate prediction results of the first batch of the test set in embodiment 1 of the application.

[0070] Figure 4 A three-dimensional display diagram of the Y-axis coordinate prediction results of the first batch of the test set in embodiment 1 of the application.

[0071] Figure 5 A three-dimensional display diagram of the Z-axis coordinate prediction results of the first batch of the test set in embodiment 1 of the application.

[0072] Figure 6 A prediction effect display diagram of the original Crossformer model in embodiment 1 of the application.

[0073] Figure 7 A prediction effect display diagram of the method proposed in embodiment 1 of the application.

[0074] Figure 8 A structural schematic diagram of a computer terminal in embodiment 2 of the application. DETAILED DESCRIPTION

[0075] The technical solutions in the embodiments of the application will be described clearly and completely below with reference to the drawings in the embodiments of the application. Obviously, the described embodiments are only part of the embodiments of the application, rather than all the embodiments of the application. Based on the embodiments in the application, all other embodiments obtained by those skilled in the art without creative labor fall within the protection scope of the application.

[0076] Embodiment 1

[0077] Please refer to Figure 1The embodiment provides a satellite orbit prediction method based on a physical information neural network, and comprises steps S1-S3.

[0078] S1. Constructing a sample set based on satellite position data, sun position data and moon position data in the same time period.

[0079] Step S1 comprises the following specific steps, namely S11-S13.

[0080] S11. Data definition:

[0081] The position vector of a single satellite at an observation time point is defined as , and the satellite position data set in the time period is set as ; similarly, the sun position data set in the time period is set as , and the moon position data set in the time period is set as ; wherein is the observation frequency, namely the sequence length; N and are the sun and moon position vectors at the observation time point .

[0082] The satellite position vector can be represented as , the sun position vector can be represented as , and the moon position vector can be represented as . It should be noted that the superscript "T" of the vector is the transpose symbol, which is different from the subsequent sequence length " T " in italics.

[0083] S12. Data acquisition:

[0084] According to the observation requirement, the observation time sequence is set, the satellite position vector at each observation time point can be obtained point by point by using a high-precision ephemeris file or a ground orbit measurement system, and is saved to ; the sun center position vector and the moon center position vector at each observation time point are calculated point by point by using the same ephemeris file or an astronomical calculation library, and the data sets and are generated.

[0085] In addition, the satellite position data set, the sun position data set and the moon position data set can also be aligned according to the public index, so that , and exist at the same observation time point .​

[0086] S13. Sample set construction:

[0087] A number of samples are generated based on the satellite orbit coordinate data set, the sun position data set, and the moon position data set, thereby constituting a sample set; each sample includes a satellite position data sequence of a preset length, a label, timestamp information, and sun and moon position data at the same time.

[0088] In step S1, the data in the sample set is also preprocessed as follows, i.e., steps 1-3.

[0089] 1. Coordinate system unification:

[0090] The position vectors in the satellite orbit coordinate data set, the sun position data set, and the moon position data set are mapped to the same geocentric inertial coordinate system.

[0091] The satellite and sun and moon position data are converted from their respective original reference systems (such as ITRF, ICRF, etc.) to a unified geocentric inertial system, ensuring physical consistency and eliminating the effects of the Earth's rotation. In this embodiment, the geocentric inertial coordinate system can use the J2000 coordinate system, which takes the Earth's mean equatorial plane and the vernal equinox on 2000-01-01T12:00:00 TT as the reference, with the Z axis pointing to the Earth's rotation axis at that time, the X axis pointing to the vernal equinox, and the Y axis determined by the right-hand rule.

[0092] 2. Feature normalization:

[0093] For each feature (coordinate component) in the satellite orbit coordinate data set, the sun position data set, and the moon position data set, the mean and standard deviation are calculated respectively for normalization processing, thereby obtaining the corresponding normalized data set.

[0094] Through normalization processing, the effects of different dimensions can be eliminated, and the training convergence can be accelerated. Taking the satellite orbit coordinate data set as an example, the sequence length is N , and then:

[0095] ;

[0096] In the formula, μ is the mean vector of the three-dimensional coordinate vector, σ is the standard deviation vector of the three-dimensional coordinate vector. For each sample point in the data set , the following transformation is performed:

[0097] ;

[0098] Referring to the above principles, the normalized data set can be obtained​ , and Similarly.

[0099] 3. The satellite orbit coordinate dataset is generated in the form of a sliding window to form multiple sequences, the former part of each sequence is data data and is used as the input of the physical information neural network, and the latter part is label label and is used for loss calculation with the predicted data pred of the network output; all samples are divided into a training set, a validation set and a test set according to a ratio of 6:2:2.

[0100] S2. Construct a physical information neural network; the network takes satellite historical position data as input to predict the satellite orbit position at a future time.

[0101] In this embodiment, the physical information neural network adopts a Crossformer model as the network backbone, and adds physical constraints as a new loss function on the basis of the original model. In some embodiments, a time series prediction model based on a Transformer architecture, a long short-term memory network (LSTM) long-term time series prediction model can also be used instead.

[0102] The core structure of Crossformer includes the following important components:

[0103] (1) Dimension-Segment (DSW) Embedding

[0104] Satellite orbit data usually contains multiple dimensions of measurement values. Through DSW embedding, the model can effectively integrate and structure these multi-dimensional information. This technology first divides the time series data into segments along each dimension, and segments according to a fixed time interval, so that each segment contains multiple measurement values, thereby ensuring the preservation of time information. Through linear projection and position embedding, each feature segment is embedded into a unified feature vector. Finally, the model generates a two-dimensional vector array H, where each embedded vector hid represents a one-dimensional slice of a time series.

[0105] (2) Two-stage attention layer

[0106] The two-stage attention layer (TSA) receives a two-dimensional array H as input, which can be the output of the DSW embedding or the result of the lower TSA layer. For each dimension, the model directly applies a multi-head self-attention (MSA) mechanism to capture the dependencies between different time steps within the same dimension. In the cross-dimension stage, the model introduces a routing mechanism that sets a small set of learnable vectors, called "routers", for each time step to aggregate information from all dimensions. These routers distribute the aggregated information to individual dimensions, effectively establishing a full connection between dimensions. Layer normalization and MLP are used to process the output of each layer in both the cross-time and cross-dimension stages.

[0107] (3) Hierarchical encoder-decoder (HED) structure

[0108] Based on the aforementioned two-stage attention network, different size time series segments are generated for the input data. Specifically, the sequence is divided into 2, 4, 8, etc. different number of time series segments from top to bottom, and the window length of the time series segments of each layer is different. The input of the model is initially composed of fine-grained time series segments, and as the number of layers increases, these time series segments are gradually aggregated into more coarse-grained representations. This structure allows the model to extract information at multiple granularities. In the decoding stage, the model uses the encoding information at different levels to make predictions, and the prediction results of each layer are then summed up to obtain the final prediction result.

[0109] S3. Training the physical information neural network using the sample set, and performing a satellite orbit prediction task using the trained physical information neural network; wherein the loss function for training the physical information neural network includes a data loss, a physical loss, and a smoothness loss; the data loss is used to measure the difference between the network prediction output and the actual observation data; the physical loss is jointly calculated by the earth's non-spherical perturbation, three-body gravitational perturbation, and solar radiation pressure physical perturbation on the satellite; the smoothness loss is used to constrain the smoothness of the predicted satellite orbit.

[0110] In this embodiment, the physical constraints are defined as follows:

[0111] The satellite orbit is a complex perturbed orbit, and high-precision orbit prediction is to calculate the spacecraft orbit according to the modeled various perturbation accelerations. The main perturbation force models include earth's non-spherical perturbation, tide, three-body gravitational perturbation, atmospheric drag perturbation, solar radiation pressure perturbation, and earth radiation pressure perturbation. The complete motion equation of the satellite in a given coordinate system can be written as:

[0112]

[0113] In the formula, is the vector pointing to the earth, ​is the velocity vector of the satellite, is the acceleration vector of the satellite, , is the gravitational constant multiplied by the mass of the Earth and the satellite, perturbation force is:

[0114] ;

[0115] In the formula, each term on the right side of the equal sign is, in turn, the Earth's non-spherical perturbation, three-body gravitational perturbation, atmospheric resistance perturbation, solar radiation pressure perturbation, tidal effect and other perturbation forces (including general relativity adjustment and Earth radiation pressure perturbation, etc.). For low-orbit spacecraft, the largest perturbation factor is the Earth's non-spherical gravity, followed by atmospheric resistance, and when the orbit position is higher, the deviation caused by solar radiation pressure is also large; for medium-orbit and high-orbit spacecraft, the perturbation factors that have a greater impact are the Earth's non-spherical gravity, the third body gravity of the sun and the moon, and the solar radiation pressure.

[0116] The present application selects the Earth's non-spherical perturbation (J2), the three-body gravitational perturbation (sun and moon), and the solar radiation pressure as the physical constraints to limit the satellite orbit predicted by the network within the range of physical laws, and introduces the jerk (derivative of acceleration) as a regularization term to constrain the smoothness of the predicted orbit.

[0117] Earth's non-spherical perturbation (J2):

[0118] ;

[0119] wherein, , , are the accelerations in the three coordinate components, respectively; is the Earth's dynamic flattening coefficient, taking 1.08263e-3; is the Earth's gravitational constant, taking 3.986e14 ; is the Earth's equatorial radius, taking 6.378e3 km; is the distance from the satellite to the center of the Earth, calculated by the predicted satellite orbit coordinates. Three-body gravitational perturbation (sun and moon):

[0120]

[0121] ; In the formula,

[0122] and are the third body gravitational constants corresponding to the sun and the moon, respectively, and the sun takes 1.327e20 , and the moon takes 1.215e20 .Take 4.903e12 , and is the vector of the third body pointing to the satellite, and is the vector of the third body pointing to the Earth, which is calculated by dynamically obtaining the positions of the sun and the moon through the timestamp of the predicted position.

[0123] Solar pressure perturbation (assuming uniform satellite surface reflection characteristics):

[0124] ;

[0125] where, is the solar constant, taking 4.56e-6, is the satellite surface specific mass ratio, is the reflection coefficient, which is determined according to different satellite parameters; is the vector of the sun pointing to the satellite, which is calculated by dynamically obtaining the position of the sun through the timestamp of the predicted position; is the shadow factor, used to determine whether the satellite is in the Earth's shadow (conical shadow approximation), in the light area, in the shadow area, the expression is:

[0126] ;

[0127] where, is the angle between the satellite pointing to the Earth vector and the sun pointing to the Earth vector; if indicates "if", and indicates "and", otherwise indicates "otherwise".

[0128] The loss function of the physical information neural network includes a data error term and a physical information error term, where the data error term is used to measure the difference between the network prediction output and the actual observation data, and the purpose is to enable the network to fit the data as much as possible; the physical information error term is unique to PINN, which considers whether the network prediction result satisfies the physical law, and the residual obtained by substituting the network prediction physical quantity into the corresponding physical law constitutes this part of the loss function, thereby ensuring physical consistency.

[0129] The data error term adopts mean square error MSE:

[0130] ;

[0131] where, is the sequence length, is the real position of the satellite at the th time, which is a three-dimensional column vector , is the model prediction of the The satellite's position at any given moment.

[0132] The physical information error term selects the above-mentioned perturbations of Earth's non-spherical shape, three-body gravitational force, and solar radiation pressure. The perturbation accelerations formed by these three factors are vector-superimposed to obtain the total perturbation acceleration:

[0133] ;

[0134] Furthermore, the actual acceleration predicted by the network in discrete time... Below, we use the second-order central difference approximation:

[0135] ;

[0136] In the formula, For the first The model predicts the orbital acceleration at each time point; To predict the orbital coordinates of the model over time A changing function; The fixed step size is defined for adjacent timestamps. The physical loss term is obtained by calculating the residual between the total perturbation acceleration calculated using physical constraints and the actual acceleration calculated from the network's predicted position:

[0137] ;

[0138] The smoothness constraint term is used to suppress abrupt changes in motion and constrains the smoothness of the predicted trajectory. The jerk (the derivative of acceleration) is used as the smoothness constraint in discrete time. Next, the first-order error is eliminated internally using the second-order central difference:

[0139] ;

[0140] in For the first The model predicts the jerk of the orbit at each moment. To predict the acceleration of the trajectory over time for the model A changing function.

[0141] The corresponding smoothness loss is:

[0142] ;

[0143] The combined three-phase loss is used, where the weight of each loss is determined using an uncertainty-weighted formula. , and The uncertainty parameters are learnable, allowing the network to adaptively assign optimal weights to each loss term based on its contribution to the loss, smoothing the training process without requiring manual curve adjustment. Specifically, for each loss term... Introduce a learnable noise variance When a certain loss Noise estimation When the gradient is large, it is divided by a larger value, automatically reducing its impact on parameter updates; conversely, it is strengthened. The final loss function is obtained as follows:

[0144] .

[0145] In step S3 of the present invention, the training method of the physical information neural network includes the following specific steps, namely S31 to S39.

[0146] S31. For a given epoch, multiple (batch_size) training samples are selected from the sample set in batches; each training sample is represented as... ;in, This represents a sequence of satellite position data of a preset length. express The corresponding tags express The corresponding timestamp information, express The corresponding timestamp information, and These represent the position data of the sun and moon at the same time.

[0147] S32. Extract each training sample... The sequence is used as the encoder input of the Crossformer model, and the sequence is segmented into fixed-length segments using the Crossformer model's tile embedding; where, T The sequence length excluding the label length, for T The satellite position vector at each moment.

[0148] S33. Using the two-stage encoding structure of the Crossformer model, multi-scale features are extracted by alternating between cross-variable attention and cross-temporal attention; among them, cross-variable attention is used to capture the coupling features between the three components x, y, and z, and cross-temporal attention is used to capture the long- and short-term orbital evolution features.

[0149] S34. Generate future step orbital segments at the decoder end of the Crossformer model, and restore them to the predicted output sequence through de-embedding. ; T pred The length of the predicted output sequence is consistent with the label length.

[0150] S35. Calculate data loss :

[0151] ;

[0152] wherein, and are the satellite position vectors at the T + j th time instant in the label and predicted output, respectively.

[0153] S36. Calculate the physical loss :

[0154] ;

[0155] wherein, denotes the total acceleration prediction of the satellite at the T + j th time instant; is the total perturbation acceleration of the satellite at the T + j th time instant.

[0156] In step S36, for the total acceleration prediction and the total perturbation acceleration of the satellite at the T + j th time instant, the calculation method comprises the following specific steps, i.e., S361-S367.

[0157] S361. Denormalize the predicted output sequence .

[0158] S362. Approximate the predicted total acceleration by second-order central difference:

[0159] .

[0160] S363. Calculate the acceleration of the satellite under the Earth's non-spherical perturbation:

[0161] ;

[0162] wherein, is the acceleration of the satellite under the Earth's non-spherical perturbation at the T + j th time instant; is the Earth's dynamic flattening coefficient, taking 1.08263e-3; is the Earth's gravitational constant, taking 3.986e14 ; is the Earth's equatorial radius, taking 6.378e3 km; r is the distance from the satellite to the Earth's center, obtained by the predicted satellite position vector; 、 and are Tj the satellite's three-dimensional coordinates at the time instant t.

[0163] S364. Calculate the satellite's acceleration under the third-body gravitational perturbation of the Sun:

[0164]

[0165] where is the satellite's acceleration under the third-body gravitational perturbation of the Sun at the time instant t; T j is the third-body gravitational constant of the Sun; is the position vector of the Sun at the time instant t, obtained by T j is the Euclidean norm of

[0166] S365. Calculate the satellite's acceleration under the third-body gravitational perturbation of the Moon:

[0167]

[0168] where is the satellite's acceleration under the third-body gravitational perturbation of the Moon at the time instant t; T j is the third-body gravitational constant of the Moon; is the position vector of the Moon at the time instant t, obtained by T j

[0169] S366. Calculate the satellite's acceleration under the solar radiation pressure perturbation:

[0170]

[0171] where is the satellite's acceleration under the solar radiation pressure perturbation at the time instant t; T j is the solar constant; is the satellite's areal density; is the satellite's albedo; is the shadow factor, = 1 if the satellite is in the illuminated region, = 0 if the satellite is in the Earth's shadow.

[0172] S367. Combine the various physical accelerations calculated in steps S363-S366 to form the total perturbation acceleration:​​​​​​​​​​​​​​​​​

[0173] .

[0174] S37. Calculate the smoothness loss :

[0175] ;

[0176] wherein, is the jerk of the satellite at the i-th time point, T + j +1; ; is the fixed step size of adjacent time stamps.

[0177] S38. Construct the joint loss :

[0178] ;

[0179] wherein, is the loss corresponding to the learnable noise variance; is the loss corresponding to the precision weight.

[0180] In step S38, the calculation formula of the precision weight is:

[0181] ;

[0182] wherein, exp(·) represents the natural exponential function.

[0183] S39. Backpropagate the joint loss , calculate the gradient, and simultaneously update the network weights and uncertainty parameters of the Crossformer model according to the optimizer strategy, complete the first iteration, and return to step S31 for the next cycle until all rounds of iterative training are completed.

[0184] In the experiment, after the training samples are iterated, the validation set iteration without updating the network parameters is performed, and the completion is the end of a round of training. According to the experimental setting, there are generally 10 to 40 rounds of training.

[0185] The embodiment also applies the above method to real prediction:

[0186] In this embodiment, the satellite selected is the Beidou medium orbit satellite BD-3M23, the observation sequence time is from January 1, 2024 to December 6, 2024, the observation interval is 5 minutes, and the obtained satellite position data is obtained from the precise ephemeris published by the Wuhan IGS data center. The positions of the sun and the moon are calculated by the astronomical calculation library Skyfield according to the observation sequence. The data set is normalized according to step S2, and the training set, the validation set and the test set are divided, and then the data is divided into sequences, the sequence seq length is set to 24x12x10, the label label length is 24x12x3, and the prediction pred length is 24x12x3, that is, one week of historical orbit data is used to predict the orbit position of the next three days, and the size of the sliding window is 24. The training round epoch=20, the batch size is 32, the hidden layer dimension is 512, the head is 8, the optimizer is Adam, and the initial value of the learning rate is set to 0.0004, which is halved with the training period.

[0187] The results of the orbit prediction for three consecutive days (3d and XYZ axes) on September 29, 2024 are shown in Figures 2 to 5 Figure 2 , in which the three coordinate axes are the x, y and z space coordinate axes in the J2000 coordinate system, with a unit of km; Figure 3 The X-axis coordinate prediction result display of the first batch of the test set (the vertical coordinate is the value of the satellite X coordinate in the J2000 coordinate system, and the horizontal coordinate is the sequence number, of which the last 864 points are prediction values; Figure 4 The Y-axis coordinate prediction result display of the first batch of the test set (the vertical coordinate is the value of the satellite Y coordinate in the J2000 coordinate system, and the horizontal coordinate is the sequence number, of which the last 864 points are prediction values); Figure 5 The Z-axis coordinate prediction display of the first batch of the test set (the vertical coordinate is the value of the satellite Z coordinate in the J2000 coordinate system, and the horizontal coordinate is the sequence number, of which the last 864 points are prediction values). It can be seen that the model shows high precision in time series prediction of satellite orbits, and the predicted trajectory (red line) is basically consistent with the true trajectory (blue line), indicating that the model can effectively fit the trend of satellite orbit changes. The mean square error (MSE) is 1.76e-4, and the mean absolute error (MAE) is 0.0109, which are very low, indicating that the model prediction error is small, with good precision and stability.

[0188] The prediction results of the first 50 data points after the start of the prediction data are compared in detail with the prediction results of the pure Crossformer model, and the results are shown in Figure 6 and Figure 7 Figure 6 ​​The prediction effect display for the original Crossformer model (under the same experimental conditions, the display data is the first 50 sequence points of the X-axis coordinate component of the orbit prediction value, the horizontal coordinate is the sequence number, and the vertical coordinate is the value of the X coordinate of the satellite in the J2000 coordinate system). Figure 7 The prediction effect display for the method proposed in the application (under the same experimental conditions, the display data is the first 50 sequence points of the X-axis coordinate component of the orbit prediction value, the horizontal coordinate is the sequence number, and the vertical coordinate is the value of the X coordinate of the satellite in the J2000 coordinate system). It can be seen that, under the same experimental parameter settings, the original Crossformer prediction result ( Figure 6 ) and the model prediction result after adding the physical constraint ( Figure 7 ) can be compared. It can be seen that the physical error term and the jerk smoothing constraint term help the model to better fit the orbit curve, reduce the sudden noise points, and better smooth the prediction curve.

[0189] Embodiment 2

[0190] The embodiment provides a computer terminal, including a memory, a processor and a computer program stored in the memory and executable on the processor, when the processor executes the computer program, the steps of the satellite orbit prediction method based on the physical information neural network are realized.

[0191] As shown in Figure 8 , the computer terminal provided in the embodiment includes at least one processor 101 and a memory 102 connected with the at least one processor 101, and the specific connection medium between the processor 101 and the memory 102 is not limited in the embodiment, Figure 8 and the processor 101 and the memory 102 are connected through a bus 100 in the embodiment. The bus 100 is represented by a thick line in the embodiment, Figure 8 and the connection mode between other components is only schematically described and is not limited. The bus 100 can be divided into an address bus, a data bus, a control bus and the like, and for convenience of representation, Figure 8 only one thick line is represented in the embodiment, but it does not mean that there is only one bus or one type of bus. Alternatively, the processor 101 can also be called a controller, and the name is not limited.

[0192] In the embodiment, the memory 102 stores instructions executable by the at least one processor 101, and the at least one processor 101 can execute the foregoing method by executing the instructions stored in the memory 102.

[0193] The processor 101 is the control center of the apparatus, and can connect all parts of the apparatus via various interfaces and lines. The processor 101 can monitor the apparatus as a whole by running or executing instructions stored in the memory 102 and calling data stored in the memory 102, and process data to realize various functions of the apparatus.

[0194] In a possible design, the processor 101 can include one or more processing units, and the processor 101 can integrate an application processor and a modem processor. The application processor can mainly process an operating system, a user interface, and an application program, and the modem processor can mainly process wireless communication. It can be understood that the modem processor can also not be integrated into the processor 101. In some embodiments, the processor 101 and the memory 102 can be implemented on the same chip, and in some embodiments, they can also be implemented on separate chips respectively.

[0195] The processor 101 can be a general-purpose processor, for example, a central processing unit (CPU), a digital signal processor, an application-specific integrated circuit, a field programmable gate array, or other programmable logic device, a discrete gate or transistor logic device, a discrete hardware component, and can implement or execute the methods, steps, and logic block diagrams disclosed in the embodiments. The general-purpose processor can be a microprocessor or any conventional processor. The steps of the satellite orbit prediction method based on a physical information neural network disclosed in Embodiment 1 can be directly embodied by a hardware processor for execution, or be executed by a combination of hardware and software modules in the processor 101.

[0196] The memory 102 is a non-volatile computer readable storage medium, and can be used to store non-volatile software programs, non-volatile computer executable programs, and modules. The memory 102 can include at least one type of storage medium, for example, can include a flash memory, a hard disk, a multimedia card, a card-type memory, a random access memory (RAM), a static random access memory (SRAM), a programmable read-only memory (PROM), a read-only memory (ROM), an electrically erasable programmable read-only memory (EEPROM), a programmable logic unit, and / or a register. The memory 102 can be any electronic, magnetic, optical, or other physical storage device that stores data including program code and / or data structures in the form of instructions, data structures, or other data. The memory 102 is not limited to the foregoing examples, but is intended to include any and all possible electronic, magnetic, optical, or other physical storage devices that store data including program code and / or data structures in the form of instructions, data structures, or other data. The memory 102 in the embodiments can also be a circuit or any other device capable of storing a program code and / or data.

[0197] The processor 101 is designed and programmed to implement the steps of the satellite orbit prediction method based on a physical information neural network as described in the foregoing embodiments, so that the chip can execute the steps of the satellite orbit prediction method based on a physical information neural network at runtime. Figure 1 The steps of the satellite orbit prediction method based on a physical information neural network are shown. How to design and program the processor 101 is a technology known to those skilled in the art, which will not be described here.

[0198] Embodiment 3

[0199] The embodiments provide a computer-readable storage medium having a computer program stored thereon, wherein the program is executed by a processor to implement the steps of the satellite orbit prediction method based on a physical information neural network as described in Embodiment 1.

[0200] The computer-readable storage medium can include a flash memory, a hard disk, a multimedia card, a card-type memory (e.g., SD or DX memory, etc.), a random access memory (RAM), a static random access memory (SRAM), a read-only memory (ROM), an electrically erasable programmable read-only memory (EEPROM), a programmable read-only memory (PROM), a magnetic memory, a magnetic disk, an optical disk, etc. In some embodiments, the storage medium can be an internal storage unit of a computer device, such as a hard disk or a memory of the computer device. In other embodiments, the storage medium can also be an external storage device of the computer device, such as a plug-in hard disk, a smart media card (SMC), a secure digital (SD) card, a flash card, etc. Of course, the storage medium can include both the internal storage unit and the external storage device of the computer device. In the embodiments, the memory is generally used to store an operating system and various application software installed on the computer device, etc. In addition, the memory can also be used to temporarily store various data that has been output or will be output.

[0201] The above description is only a preferred embodiment of the present application, but the protection scope of the present application is not limited thereto. Any person skilled in the art can make equivalent replacements or changes to the technical solution and the inventive concept of the present application within the technical scope disclosed by the present application, which should be covered by the protection scope of the present application.

Claims

1. A satellite orbit prediction method based on physical information neural networks, characterized in that, Including the following steps: S1. Construct a sample set based on satellite position data, solar position data, and lunar position data within the same time period; S2. Construct a physical information neural network; this network uses historical satellite position data as input to predict the satellite's orbital position at future moments; S3. Train the physical information neural network using the sample set, and use the trained physical information neural network to perform satellite orbit prediction tasks; The training process of the physical information neural network includes: inputting a satellite position data sequence of a preset length from the sample set into the physical information neural network with the Crossformer model as the network backbone; extracting multi-scale features of the orbit using cross-variable attention and cross-temporal attention to generate a predicted output sequence; calculating data loss, physical loss, and smoothness loss based on the predicted output sequence to construct a joint loss; updating the network parameters of the Crossformer model using backpropagation on the joint loss; and iterating in batches until training is complete; the data loss is used to measure the difference between the network's predicted output and the actual observed data; the physical loss is jointly calculated by the physical perturbations of the Earth's non-spherical perturbation, the three-body gravitational perturbation, and the solar radiation pressure perturbation experienced by the satellite; and the smoothness loss is used to constrain the smoothness of the predicted satellite orbit.

2. The satellite orbit prediction method based on a physical information neural network according to claim 1, characterized in that, In step S3, the training method for the physical information neural network includes the following specific steps: S31. Select training samples in batches from the sample set; each training sample is represented as... ;in, This represents a sequence of satellite position data of a preset length. express The corresponding tags express The corresponding timestamp information, express The corresponding timestamp information, and These represent the positions of the sun and moon at the same time; S32. Extract each training sample... The sequence is used as the encoder input of the Crossformer model, and the sequence is segmented into fixed-length segments using the Crossformer model's tile embedding; where, T The sequence length excluding the label length, for T Satellite position vector at each moment; S33. Using the two-stage encoding structure of the Crossformer model, multi-scale features are extracted by alternating between cross-variable attention and cross-temporal attention; among them, cross-variable attention is used to capture the coupling features between the three components x, y, and z, and cross-temporal attention is used to capture the long-term and short-term orbital evolution features. S34. Generate future step orbital segments at the decoder end of the Crossformer model, and restore them to the predicted output sequence through de-embedding. ; T pred The length of the predicted output sequence is consistent with the label length; S35. Calculate data loss : ； In the formula, and The first and second parts of the label and prediction output are respectively T + j Satellite position vector at each moment; S36. Calculate physical losses : ； In the formula, Indicates the satellite in T + j The predicted total acceleration value at each moment; For the first T + j The total satellite perturbation acceleration at each moment; for The Euclidean norm; S37. Calculate the smoothness loss : ； In the formula, For the first T + j The satellite's jerkness at time +1. ; A fixed step size between adjacent timestamps; S38. Constructing Joint Losses : ； In the formula, For loss The corresponding learnable noise variance; For loss The corresponding precision weights; S39. Joint Losses Perform backpropagation, calculate the gradient, update the network weights and uncertain parameters of the Crossformer model simultaneously according to the optimizer policy, complete the first iteration, and return to step S31 to start the next loop, until all rounds of iterative training are completed.

3. The satellite orbit prediction method based on a physical information neural network according to claim 2, characterized in that, In step S36, for the satellite in the... T + j The predicted total acceleration and total perturbation acceleration at each moment are calculated using the following specific steps: S361. The predicted output sequence Anti-normalization; S362. Predicting total acceleration using a second-order central difference approximation: ； S363. Calculate the satellite's acceleration under non-spherical perturbations of the Earth: ； In the formula, For the first T + j The satellite's acceleration at a given moment under nonspherical perturbations of the Earth; The flattening coefficient for geodynamics is taken as 1.08263e-3; Let be the Earth's gravitational constant, taken as 3.986e14. ; The radius of the Earth's equator is taken as 6.378e3 km; r The distance from the satellite to the Earth's center is obtained through the predicted satellite position vector; , and for T + j The satellite's three-dimensional coordinates at that moment; S364. Calculate the acceleration of the satellite under the gravitational perturbation of the Sun as a third body: ； In the formula, For the first T + j The acceleration of the satellite under the gravitational perturbation of the Sun's third body at a given moment; The gravitational constant of the Sun as a third body; For the first T + j The solar position vector at each moment, through get; S365. Calculate the acceleration of the satellite under the gravitational perturbation of the Moon as a third body: ； In the formula, For the first T + j The acceleration of the satellite under the gravitational perturbation of the third body on the Moon at a given moment; The third-body gravitational constant of the Moon; For the first T + j The lunar position vector at each instant, through Obtain; S366. Calculate the satellite's acceleration under solar radiation pressure perturbation: ； In the formula, For the first T + j The acceleration of the satellite under solar radiation pressure perturbation at a given moment; It is the solar constant; The surface mass ratio of the satellite; The reflection coefficient of the satellite; As the shadow factor, =1 indicates that the satellite is in the illuminated area. =0 indicates that the satellite is in the Earth's shadow area; S367. Combine the physical accelerations calculated in steps S363 to S366 to form the total perturbation acceleration: 。 4. The satellite orbit prediction method based on a physical information neural network according to claim 2, characterized in that, In step S38, precision weights The calculation formula is: ； In the formula, exp(·) represents the natural exponential function.

5. The satellite orbit prediction method based on a physical information neural network according to claim 1, characterized in that, Step S1 includes the following specific steps: S11. Define the observation time of a single satellite. The position vector is The satellite location dataset for the specified time period is then set as follows: Similarly, the dataset of the sun's position within a given time period is... The lunar position dataset is ;in, N The number of observations is the sequence length. and At the observation time The position vectors of the sun and moon; S12. Set the observation time series according to the observation requirements. Using ephemeris files or ground-based orbit determination systems, each observation time is acquired point by point. satellite position vector Save to ; Calculate each observation time point by point using the same ephemeris file or astronomical computing library. Sun's barycenter position vector and the position vector of the Moon's center of mass Generate dataset and ; S13. Generate several samples based on the satellite orbit coordinate dataset, the solar position dataset, and the lunar position dataset to form a sample set; each sample includes a satellite position data sequence of a preset length, labels, timestamp information, and solar and lunar position data at the same time.

6. The satellite orbit prediction method based on a physical information neural network according to claim 5, characterized in that, In step S1, the data in the sample set are also preprocessed as follows: Map the position vectors in the satellite orbit coordinate dataset, the solar position dataset, and the lunar position dataset to the same geocentric inertial coordinate system; For each coordinate component in the satellite orbit coordinate dataset, the solar position dataset, and the lunar position dataset, the mean and standard deviation are calculated for normalization, thus obtaining the corresponding normalized dataset. The satellite orbit coordinate dataset is used to generate multiple sequences in the form of a sliding window. The first part of each sequence is the data and is used as the input to the physical information neural network, and the second part is the label and is used to calculate the loss with the prediction data output by the network. All samples are divided into training set, validation set and test set according to a set ratio.

7. The satellite orbit prediction method based on a physical information neural network according to claim 5, characterized in that, In step S12, the satellite position dataset, solar position dataset, and lunar position dataset are also aligned according to a common index to ensure that they are observed at the same time. exist , and .

8. A computer terminal, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that, When the processor executes the computer program, it implements the steps of the satellite orbit prediction method based on a physical information neural network as described in any one of claims 1 to 7.

9. A computer-readable storage medium having a computer program stored thereon, characterized in that, When the program is executed by the processor, it implements the steps of the satellite orbit prediction method based on a physical information neural network as described in any one of claims 1 to 7.

Images