阅读说明：本技术 一种微型电推进推力大小的在轨标定方法 (On-orbit calibration method for magnitude of micro electric propulsion thrust ) 是由 王菲 张众正 李明翔 吴彤 郭晓华 牟邵君 于 2021-06-26 设计创作，主要内容包括：本发明公开一种微型电推进推力大小的在轨标定方法,包括：计算升轨和降轨过程中半长轴的变化；理论分析升轨和降轨过程中地球形状摄动和大气阻力摄动的影响；联立升轨和降轨过程中计算得到的半长轴变化量公式,求得推力加速度；进而根据推力加速度和卫星质量得到电推推力。本方案不需要通过姿态或角速度变化转换获得推力大小,也就是说推力器不需要偏心安装产生偏心力矩改变姿态或角速度,消除了与微型推力器推力大小量级较为接近的地球形状摄动、大气阻力摄动对推力大小标定的影响,并且可以实现对单个过质心安装推力器的推力大小在轨标定。(The invention discloses an on-orbit calibration method for the magnitude of micro electric propulsion thrust, which comprises the following steps: calculating the change of the semi-long axis in the rail lifting and lowering processes; theoretically analyzing the influence of earth spherical perturbation and atmospheric resistance perturbation in the rail ascending and descending processes; calculating a semi-major axis variable quantity formula obtained by calculation in the process of simultaneous rail lifting and rail lowering to obtain thrust acceleration; and then the electric thrust is obtained according to the thrust acceleration and the satellite mass. The thrust size does not need to be obtained through attitude or angular speed change conversion, namely the thruster does not need to be eccentrically installed to generate eccentric torque to change the attitude or the angular speed, the influence of spherical perturbation and atmospheric resistance perturbation which are relatively close to the magnitude order of the thrust of the miniature thruster on the calibration of the thrust size is eliminated, and the on-orbit calibration of the thrust size of a single thruster passing through the mass center can be realized.)

1. An on-orbit calibration method for the magnitude of micro electric propulsion thrust is characterized by comprising the following steps:

step 1, calculating semimajor axis change in a rail lifting process;

a_tis the thrust acceleration to which the satellite is subjected during the orbit raising process, a_d1For acceleration of rail resistance, a_e1For the perturbed acceleration of the earth's spherical shape during the rail-lifting process, t₁Total thrust application time for the rail lifting process, n₁Is the total orbit number of the satellite in the orbit raising process, T is the average orbit period in the orbit raising process, Delta a is the variation of the orbit semi-major axis in the orbit raising process,denotes the average orbit half-major axis, μ 398600km³/s²Is an earth gravity parameter;

step 2, calculating the change of the semimajor axis in the rail descending process;

a_tthrust acceleration to which the satellite is subjected during the course of falling into orbit, a_d2For falling rail resistance acceleration, a_e2For the acceleration of the earth's spherical perturbation during the rail lowering process, t₂For the total thrust application time in the rail lowering process, n₂The total number of the running orbits of the satellite in the orbit descending process is T, the average orbit period in the orbit descending process is T, and delta a is the variation of the orbit semi-major axis in the orbit descending process;

step 3, theoretically analyzing the influence of earth spherical perturbation and atmospheric resistance perturbation in the rail ascending and descending processes;

1) for a low earth orbit satellite, the earth gravity perturbation mainly considers J₂Determining that the J2 perturbation does not cause the change of the semi-major axis by integrating one circle along the satellite orbit;

2) during the process of raising and lowering the same orbit height, the atmospheric resistance suffered by the satellite is equal in magnitude and same in direction;

step 4, further according to the theoretical analysis of the step 3, equations (3) and (5) are combined to obtain the magnitude of the thrust acceleration;

step 5, calculating the electric thrust F according to the thrust acceleration and the satellite mass_TThe following are:

wherein m is the satellite mass.

2. The on-orbit calibration method for the magnitude of the micro electric propulsion thrust according to claim 1, characterized in that: in the step 3, the specific theoretical analysis process is as follows:

(1) the Gaussian perturbation equation of the orbit root can know that perturbation forces influencing the semi-major axis of the orbit are radial force and tangential force, and the following is shown:

wherein a is the semi-major axis of the track, e is the eccentricity, theta is the true paraxial point angle,is the average velocity, F_r、F_t、F_nRadial, tangential and normal external forces respectively applied to the satellite for the near-earth orbitSatellite, earth gravity perturbation mainly considering J₂Influence of item perturbation, J₂The components of the perturbation force in the radial direction, the tangential direction and the normal direction are as follows:

in the formula, J₂＝1.08263×10^-3Is J₂Coefficient of term, R_e6378.14km is the radius of the earth, r is a (1-e)/(1+ ecos theta) is the center distance of the earth, u is omega + theta is the latitude amplitude angle, the perturbation force is substituted into a Gaussian equation, and the Gaussian equation is integrated for one circle along the orbit, so that J is found₂Item perturbation does not cause variation in the semi-major axis, so J is every complete revolution during long run₂Item perturbation all self-counteracts, incomplete rounds of J₂The item perturbation influence is ignored;

(2) because thin high-rise atmosphere exists on the near-earth orbit, the satellite runs under atmospheric resistance for a long time in orbit, the orbit height is gradually attenuated, and the perturbation of the atmospheric resistance and the aerodynamic coefficient C of the satellite surface_dThe atmospheric density ρ, the windward area S and the satellite velocity v are related, and the equation is as follows:

because the space atmosphere model and the satellite surface aerodynamic coefficient can not be accurately estimated, the atmospheric resistance received by the satellite can not be accurately calculated, but the atmospheric resistance received by the satellite is considered to be equal in size and same in direction in the process of raising and lowering the same orbit height, and therefore the influence of the atmospheric resistance is counteracted through the simultaneous calculation of the subtraction of the lifting orbit thrust equation.

Technical Field

The invention belongs to the field of on-orbit calibration of electric propulsion thrust, and particularly relates to an on-orbit calibration method of micro electric propulsion thrust.

Background

The micro electric propulsion has the advantages of small volume, light weight, large specific impulse, high total impulse and the like, and is more and more widely applied in the field of commercial aerospace in recent years. However, the general thrust of the micro electric propulsion product is tens or hundreds of micro-newtons, and the temperature and pressure conditions consistent with the space environment are difficult to establish on the ground due to the limitation of the ground environment; and the accurate calibration of the thrust of the miniature electric propeller is difficult to realize due to the limitation of ground measurement means. In order to ensure the effectiveness of the in-orbit operation of the miniature electric thruster, the thrust of the miniature electric thruster needs to be calibrated in the in-orbit mode.

The existing thrust calibration methods comprise two methods, one is an orbit calibration method, and the orbit calibration method calculates to obtain the satellite speed increment according to the orbit change parameters; and calculating to obtain the thrust of the propulsion system according to the satellite speed increment and the electric propulsion working time. The method does not consider the influence of space perturbation force such as earth spherical perturbation and atmospheric resistance perturbation on thrust calibration. The other method is an attitude calibration method, in the method, a thruster which is not beyond the centroid is started to generate a control moment, and the thrust is calculated through the change of the satellite attitude and the angular speed. But the method is only suitable for installing a plurality of thrusters, each thruster is not installed with an excessive satellite mass center, and the on-orbit calibration of the thrust of a single thruster installed with the excessive mass center cannot be met.

Disclosure of Invention

Aiming at the defects of the existing micro electric thruster thrust magnitude on-orbit calibration technology, the invention provides an on-orbit calibration method of the micro electric thruster thrust magnitude, namely a lifting rail combined calibration method, which considers the influences of earth spherical perturbation and atmospheric resistance perturbation and can realize the on-orbit calibration of the thrust magnitude of a single mass center passing mounting thruster.

The invention is realized by adopting the following technical scheme: an on-orbit calibration method for the magnitude of micro electric propulsion thrust comprises the following steps:

step 1, calculating semimajor axis change in a rail lifting process;

a_tis the thrust acceleration to which the satellite is subjected during the orbit raising process, a_d1For acceleration of rail resistance, a_e1For the perturbed acceleration of the earth's spherical shape during the rail-lifting process, t₁Total thrust application time for the rail lifting process, n₁The total orbit number of the satellite in the orbit raising process is shown, T is the average orbit period in the orbit raising process, delta a is the variation of the orbit semimajor axis in the orbit raising process, a represents the average orbit semimajor axis, and mu is 398600km³/s²Is an earth gravity parameter;

step 2, calculating the change of the semimajor axis in the rail descending process;

a_tthrust acceleration to which the satellite is subjected during the course of falling into orbit, a_d2For falling rail resistance acceleration, a_e2For the acceleration of the earth's spherical perturbation during the rail lowering process, t₂For the total thrust application time in the rail lowering process, n₂The total number of the running orbits of the satellite in the orbit descending process is T, the average orbit period in the orbit descending process is T, and delta a is the variation of the orbit semi-major axis in the orbit descending process;

step 3, theoretically analyzing the influence of earth spherical perturbation and atmospheric resistance perturbation in the rail ascending and descending processes;

1) for a low earth orbit satellite, the earth gravity perturbation mainly considers J₂Determining that the J2 perturbation does not cause the change of the semi-major axis by integrating one circle along the satellite orbit;

2) during the process of raising and lowering the same orbit height, the atmospheric resistance suffered by the satellite is equal in magnitude and same in direction;

step 4, further according to the theoretical analysis of the step 3, equations (3) and (5) are combined to obtain the magnitude of the thrust acceleration;

step 5, calculating the electric thrust F according to the thrust acceleration and the satellite mass_TThe following are:

wherein m is the satellite mass.

Further, in step 3, a specific theoretical analysis process is as follows:

(1) the Gaussian perturbation equation of the orbit root can know that perturbation forces influencing the semi-major axis of the orbit are radial force and tangential force, and the following is shown:

wherein a is the semi-major axis of the track, e is the eccentricity, theta is the true paraxial point angle,is the average velocity, F_r、F_t、F_nThe external forces of the satellite in the radial direction, the tangential direction and the normal direction are respectively considered, and J is mainly considered for the perturbation of the earth gravity of the satellite in the near-earth orbit₂Influence of item perturbation, J₂The components of the perturbation force in the radial direction, the tangential direction and the normal direction are as follows:

in the formula, J₂＝1.08263×10^-3Is J₂Coefficient of term R_e6378.14km is the radius of the earth, r is a (1-e)/(1+ ecos theta) is the center distance of the earth, u is omega + theta is the latitude amplitude angle, the perturbation force is substituted into a Gaussian equation, and the Gaussian equation is integrated for one circle along the orbit, so that J is found₂Item perturbation does not cause variation in the semi-major axis, so J is every complete revolution during long run₂Item perturbation all self-counteracts, incomplete rounds of J₂The item perturbation influence is ignored;

(2) because thin high-rise atmosphere exists on the near-earth orbit, the satellite runs under atmospheric resistance for a long time in orbit, the orbit height is gradually attenuated, and the perturbation of the atmospheric resistance and the aerodynamic coefficient C of the satellite surface_dThe atmospheric density ρ, the windward area S and the satellite velocity v are related, and the equation is as follows:

because the space atmosphere model and the satellite surface aerodynamic coefficient can not be accurately estimated, the atmospheric resistance received by the satellite can not be accurately calculated, but the atmospheric resistance received by the satellite is considered to be equal in size and same in direction in the process of raising and lowering the same orbit height, and therefore the influence of the atmospheric resistance is counteracted through the simultaneous calculation of the subtraction of the lifting orbit thrust equation.

Compared with the prior art, the invention has the advantages and positive effects that:

the calibration method provided by the scheme eliminates the influence of spherical perturbation and atmospheric resistance perturbation on the calibration of the thrust magnitude, wherein the magnitude order of the thrust magnitude of the micro thruster is relatively close to that of the spherical perturbation and the atmospheric resistance perturbation; and the on-orbit calibration of the thrust magnitude of a single through-center-of-mass mounting thruster can be realized.

Drawings

FIG. 1 is a schematic diagram of a combined calibration process of a lifting rail for micro electric propulsion thrust in an embodiment of the present invention;

FIG. 2 is a schematic diagram of a variation curve of a semi-major axis of a satellite orbit with time when a micro electric propulsion device is used to complete orbit lifting according to an embodiment of the present invention.

Detailed Description

In order to make the above objects, features and advantages of the present invention more clearly understood, the present invention will be further described with reference to the accompanying drawings and examples. In the following description, numerous specific details are set forth in order to provide a thorough understanding of the present invention, however, the present invention may be practiced in other ways than those described herein, and thus, the present invention is not limited to the specific embodiments disclosed below.

The method for calibrating the magnitude of the micro electric propulsion thrust in the on-orbit requires a satellite to perform 2 stages of operation of ascending and descending the orbit, eliminates the influence of perturbation force by selecting the part with the consistent height of the orbit of the ascending and descending orbit, and realizes the accurate on-orbit calibration of the magnitude of the micro electric propulsion thrust.

The on-orbit calibration method for the magnitude of the micro electric propulsion thrust provided by the invention obtains the following calculation formula for the magnitude of the micro electric propulsion thrust:

in the formula, F_TRepresenting micro-electric propulsion thrust, m representing satellite mass, n₁Representing the total orbit number n of the satellite in the orbit raising stage₂Represents the total orbit number t of the satellite in the orbit reduction stage₁Representing the total thrust application time, t, in the rail lifting phase₂The total thrust application time in the rail lowering stage is shown, Δ a represents the variation of the semi-major axis of the track, a represents the average semi-major axis of the track, and μ is 398600km³/s²Is the gravity parameter.

Specifically, as shown in fig. 1, a schematic diagram of a combined calibration flow of a lifting rail for micro electric propulsion thrust is shown, the objective of the invention is to accurately calibrate the magnitude of the micro electric propulsion thrust, and consider the influences of spherical perturbation and atmospheric resistance perturbation in the calibration process, and simultaneously realize the on-rail calibration of the magnitude of the thrust of a single over-center-of-mass mounted thruster, and the specific implementation steps of the combined calibration of the lifting rail for micro electric propulsion thrust are as follows:

step 1, calculating the semimajor axis change in the rail lifting process according to external force and time;

suppose the thrust acceleration of the satellite during the orbit raising process is a_tAcceleration of resistance a_d1The earth's spherical perturbation acceleration is a_e1. Recording the total thrust application time t during the rail lifting process₁Total number of orbits n of satellite₁The average track period T and the variation Δ a of the track semimajor axis are as follows:

a_t·t₁+a_e1·n₁·T-a_d1·n₁·T＝Δv (2)

in the formula (I), the compound is shown in the specification,the velocity increment generated when the satellite is subjected to external force in the process of orbit rising.

Thus obtaining:

step 2, calculating the change of the semimajor axis in the rail descending process according to the external force and time;

assuming that the thrust acceleration received by the satellite in the process of orbit reduction is a_tAcceleration of resistance a_d2The earth's spherical perturbation acceleration is a_e2. Recording the total application time t of the thrust during the rail descending process₂Total number of orbits n of satellite₂The average track period T and the variation Δ a of the semi-major axis of the track are

a_t·t₂+a_e2·n₂·T+a_d2·n₂·T＝Δv (4)

Step 3, theoretically analyzing the influence of earth spherical perturbation and atmospheric resistance perturbation in the rail ascending and descending processes;

the Gaussian perturbation equation of the orbit root can know that perturbation forces influencing the semi-major axis of the orbit are radial force and tangential force, and the following is shown:

wherein a is the semi-major axis of the track, e is the eccentricity, theta is the true paraxial point angle,is the average velocity, F_r、F_t、F_nRadial, tangential and normal external forces to which the satellite is subjected are provided. For a low earth orbit satellite, the earth gravity perturbation mainly considers J₂Influence of item perturbation, J₂The components of the perturbation force in the radial direction, the tangential direction and the normal direction are as follows:

in the formula, J₂＝1.082×63^-3Is J₂Coefficient of term, μ 398600km³/s²As a parameter of Earth's gravity, R_e6378.14km is the earth radius, r a (1-e)/(1+ ecos θ) is the earth center distance, and u ω + θ is the latitude argument. The perturbation force is substituted into the Gaussian equation, and the J can be found by integrating one circle along the track₂Item perturbation does not cause variation in the semi-major axis. Thus, J for each complete revolution during long run₂The item perturbation can be automatically counteracted, and J of incomplete circle₂Item perturbation effects are negligible in small quantities.

Due to the thin high atmosphere on the low earth orbit, the satellite runs under atmospheric resistance for a long time in orbit, and the orbit height gradually attenuates. Atmospheric drag perturbation and satellite surface aerodynamic coefficient C_dThe atmospheric density ρ, the windward area S and the satellite velocity v are related, and the equation is as follows:

because the space atmosphere model and the satellite surface aerodynamic coefficient can not be accurately estimated, the atmospheric resistance of the satellite can not be accurately calculated. However, the atmospheric resistance experienced by the satellite during the same orbital altitude elevation and lowering may be considered equal and equal. Therefore, the influence of the atmospheric resistance can be offset by the subtraction simultaneous calculation of the lifting rail thrust equation.

Step 4, simultaneous equations (3) and (5) are used for solving the magnitude of the thrust acceleration;

and 5, calculating the electric thrust force according to the thrust acceleration and the satellite mass, wherein the electric thrust force comprises the following steps:

from the above formula, it can be seen that: the method does not need to obtain the thrust magnitude through attitude or angular speed change conversion, namely the thruster does not need to be eccentrically installed to generate eccentric torque to change the attitude or the angular speed, so that the method can realize the on-orbit calibration of the thrust magnitude of a single thruster passing through the center of mass.

The above method was verified by simulation as follows:

and setting the satellite mass to be 38kg and the electric thrust to be 500 mu N, and establishing a high-precision space environment ground simulation model. Firstly, the average semi-major axis of the initial orbit of the satellite is set to be 6878.5km, the thruster is opened for 150 tracks, each track is opened for 40 minutes, and the average semi-major axis of the orbit which can be reached by the satellite is 6886.75km (see the curve ascending segment of fig. 2).

Then, the thruster is started at 146 tracks, each track is started for 40 minutes, and the simulation finds that the satellite can return to the initial track (see the curve descending segment of fig. 2) with the average semi-major axis of the track being 6878.5km, and the thrust of the micro electric propulsion is 498 mu N through the calculation of the formula (1), and is basically consistent with the setting. The analysis can show that the proposed method can effectively calibrate the magnitude of the micro electric propulsion thrust.

The above description is only a preferred embodiment of the present invention, and not intended to limit the present invention in other forms, and any person skilled in the art may apply the above modifications or changes to the equivalent embodiments with equivalent changes, without departing from the technical spirit of the present invention, and any simple modification, equivalent change and change made to the above embodiments according to the technical spirit of the present invention still belong to the protection scope of the technical spirit of the present invention.