CN112009727A - Optimal low-thrust transfer sectional design method for translation point orbit - Google Patents

Optimal low-thrust transfer sectional design method for translation point orbit Download PDF

Info

Publication number
CN112009727A
CN112009727A CN202010850752.2A CN202010850752A CN112009727A CN 112009727 A CN112009727 A CN 112009727A CN 202010850752 A CN202010850752 A CN 202010850752A CN 112009727 A CN112009727 A CN 112009727A
Authority
CN
China
Prior art keywords
point
thrust
design
transfer
manifold
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202010850752.2A
Other languages
Chinese (zh)
Other versions
CN112009727B (en
Inventor
曾豪
彭坤
王平
侯振东
齐岳
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing Space Technology Research and Test Center
Original Assignee
Beijing Space Technology Research and Test Center
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beijing Space Technology Research and Test Center filed Critical Beijing Space Technology Research and Test Center
Priority to CN202010850752.2A priority Critical patent/CN112009727B/en
Publication of CN112009727A publication Critical patent/CN112009727A/en
Application granted granted Critical
Publication of CN112009727B publication Critical patent/CN112009727B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B64—AIRCRAFT; AVIATION; COSMONAUTICS
    • B64G—COSMONAUTICS; VEHICLES OR EQUIPMENT THEREFOR
    • B64G1/00—Cosmonautic vehicles
    • B64G1/22—Parts of, or equipment specially adapted for fitting in or to, cosmonautic vehicles
    • B64G1/24—Guiding or controlling apparatus, e.g. for attitude control
    • B64G1/242—Orbits and trajectories
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/04—Constraint-based CAD
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/06—Multi-objective optimisation, e.g. Pareto optimisation using simulated annealing [SA], ant colony algorithms or genetic algorithms [GA]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Mathematical Optimization (AREA)
  • General Engineering & Computer Science (AREA)
  • Mathematical Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Geometry (AREA)
  • Algebra (AREA)
  • Data Mining & Analysis (AREA)
  • Computational Mathematics (AREA)
  • Evolutionary Computation (AREA)
  • Computer Hardware Design (AREA)
  • Remote Sensing (AREA)
  • Combustion & Propulsion (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Chemical & Material Sciences (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Fluid Mechanics (AREA)
  • Operations Research (AREA)
  • Computing Systems (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Traffic Control Systems (AREA)

Abstract

The invention relates to an optimal low-thrust transfer sectional design method for a translation point orbit, which comprises the following steps of: a. establishing a dynamic equation of the aircraft under a centroid rotating coordinate system; b. determining the invariant manifold of the translation point track; c. global searching is carried out to determine initial values of design variables; d. and determining the optimal transfer track by local optimization. The method is based on a three-body system, global search is carried out by combining invariant manifold, the position of a proper splicing point and the initial value of a design variable are rapidly determined, then the optimal low-thrust transfer trajectory of the fuel is locally optimized and solved, the transfer design of the detector among system translation points is realized, the problems of low-thrust near-ground multi-circle divergence and the initial value can be effectively solved, the application range is wide, and the method can be applied to the transfer design of translation point tracks of the lunar L1 point, the solar L1 point and the L2 point.

Description

Optimal low-thrust transfer sectional design method for translation point orbit
Technical Field
The invention relates to an optimal low-thrust transfer sectional design method for a translation point orbit.
Background
The translational point is a space point where gravitational force and inertial force balance each other in the restrictive trisomy dynamics, and 7 translational points exist near the earth: the L1 and L2 points of the Sun-Earth system, and the L1-L5 translation points of the Earth-moon system. The periodic orbit near the translation point is the optimal foot drop point for space environment observation, and can be used as a transfer station for future lunar exploration personnel and material transportation, manned mars and asteroid exploration tasks. Meanwhile, the low-energy-consumption transfer track often exists in the translation point track task, and the detector can develop a detection task or a expansibility task based on the track, so that the application value can be effectively improved.
The design of the translation point orbit space transfer orbit is mainly divided into pulse control and small-thrust control. The pulse control transfer mode includes direct transfer, weak stable transfer and force-assisted transfer. Zeng adopts a layered differential correction strategy and an initial value polynomial, and better solves the problems of multiple constraints and initial value sensitivity in Halo pulse transfer track design (Zeng Hao, Zhang J ingrui. Qi Rui, Li Mingtao. Study of time-free transfer point with multiple constraints [ J ]. Journal of Guidance, Control, and dynamics, 2017,40(11):2752 and 2770). Parker researches weak and stable boundary transfer of Halo orbits of L2 points of aircraft flying earth moon by connecting space invariant manifold of earth-moon system and sun-earth system, but has the problems of long flying time and the like (Parker J S. families of low-energy lunar Halo transfers [ C ]. AAS/AIAA Spaceflight Dynamics Conference, 2016). Folta proposes an improved lunar gravity assisted triple pulse transfer scheme, gives a lunar gravity orientation and orbit selection strategy and quantitatively analyzes the influence of gravity constraint sets (Folta D C, Pavlak T A, Haapala A F, et al.
In recent years, some scholars have explored a track design method combining advanced low-thrust technology and a panning point task. Ozimek derives the first-order optimal conditions for inserting manifolds for the translational point orbital-to-variable-specific-impulse small-thrust transfer optimization problem, but the design does not set upper and lower limits for thrust and specific impulse (Ozimek M. Low-run traffic design and optimization of nuclear pole coverage [ D ]. Purde University: School of aeronautical and Astronautics, 2010). Zhang adopts Gauss-Lobatt matching method to solve the fuel optimal specific impulse low thrust transfer orbit from GTO to Earth moon L1 and L2 point Lyapunov orbit, and combines tangential thrust and optimal Control method to reduce the difficulty of complex dynamic orbit design (Zhang C, Topputo F, Zazzera F B, et al. Low-run minimum-fuel optimization in the circular estimated three-body protocol [ J ]. Journal of Guidant, Control, and dynamics.2015,38(8):1501 2509).
Compared with a pulse translational point task, the low-thrust transfer has the characteristic of high specific impulse, can greatly reduce the fuel consumption of the task, and is suitable for a space material transportation task with low time requirement, but the low-thrust translational point task has the problems of multi-celestial body strong gravitational field nonlinear characteristic, low near-earth-end multi-turn low-thrust transfer calculation efficiency, difficulty in convergence and the like. Therefore, how to quickly determine a proper splicing point, comprehensively consider the influence of an earth berthing orbit escape point and a manifold capture point, and how to effectively search and calculate feasible design initial values to reduce sensitivity is worthy of further research.
Disclosure of Invention
The invention aims to provide a translational point orbit optimal low-thrust transfer sectional design method capable of quickly determining an optimal track.
In order to achieve the aim, the invention provides a translational point orbit optimal low thrust transfer sectional design method, which comprises the following steps:
a. establishing a dynamic equation of the aircraft under a centroid rotating coordinate system;
b. determining the invariant manifold of the translation point track;
c. global searching is carried out to determine initial values of design variables;
d. and determining the optimal transfer track by local optimization.
According to one aspect of the present invention, in the step (a), the origin of the centroid rotating coordinate system is a centroid of a three-body system composed of the aircraft, the first main celestial body and the second main celestial body;
wherein, if the first main celestial body is the earth, the second main celestial body is the moon, and if the first main celestial body is the sun, the second main celestial body is the earth or the moon;
the X axis points to the second main celestial body from the first main celestial body, the Z axis is parallel to the angular momentum direction of the three-system, and the Y axis is determined by a right-hand rule;
the dynamic equations comprise an uncontrolled dynamic equation and a controlled dynamic equation;
the uncontrolled kinetic equation is as follows:
Figure BDA0002644639870000031
wherein the content of the first and second substances,
Figure BDA0002644639870000032
mu is the mass coefficient of the system,
Figure BDA0002644639870000033
the state quantity of the aircraft under the mass center rotation system,
Figure BDA0002644639870000034
and
Figure BDA0002644639870000035
the aircraft is respectively relative to the first main celestial body and the second main celestial bodyThe distance of the body;
the controlled dynamics equation is:
Figure BDA0002644639870000041
wherein the content of the first and second substances,
Figure BDA0002644639870000042
t is the magnitude of the thrust, uTIs a unit vector of the three axial directions of the thrust, and satisfies
Figure BDA0002644639870000043
P is the aircraft engine power, m is the aircraft mass, and v is the aircraft velocity.
According to an aspect of the present invention, in the step (b), the spatially invariant manifold is determined according to eigenvalues and eigenvectors of a matrix M, where the matrix M is:
M＝Φ(t0+T,t0)；
wherein the content of the first and second substances,
Figure BDA0002644639870000044
03and I3Respectively 3 × 3 order zero matrix and unit matrix;
the linear perturbation state quantities of the unstable manifold and the stable manifold of the translational point orbit at any time t are respectively:
Figure BDA0002644639870000045
wherein the superscripts "u" and "s" respectively represent unstable and stable manifolds in the invariant manifold,
Figure BDA0002644639870000046
is the eigenvector, perturbation d, of the matrix Mm200km in a daily terrestrial system and 50km in a terrestrial lunar system, the sign "±" corresponding to the invariant manifold in both directions.
According to an aspect of the present invention, in the step (c), the searching method includes:
s1, modeling the fuel optimization problem;
s2, determining and analyzing a design variable set;
s3, global search is carried out to determine the splicing point area and the initial value range of the design variable;
s4, repeating the step (S3) to screen the splicing point area and the initial value range of the design variable.
According to an aspect of the present invention, in the step (S1), a fuel optimum low thrust shift trajectory is first determined, and its optimum performance index maxJ can be expressed as:
maxJ＝k·mf；
wherein k is an arbitrary normal number, mfAnd determining the Hamilton function of the system for the residual mass of the aircraft at the moment of engine shutdown by combining the Pontryagin maximum principle and the controlled dynamics equation:
Figure BDA0002644639870000051
the co-state equation of the optimal problem is:
Figure BDA0002644639870000052
wherein λ m is a covariate.
According to one aspect of the present invention, the set of design variables in said step (S2) is:
Figure BDA0002644639870000053
during analysis, the optimal low thrust transfer track is divided into a first low thrust section, an uncontrolled gliding section and a second low thrust section, wherein the thrust direction of the first low thrust section is along the direction of a velocity vector;
wherein the content of the first and second substances,
Figure BDA0002644639870000054
and
Figure BDA0002644639870000055
the elevation angle, the elevation intersection point declination and the latitude argument are respectively used for describing the near-ground end mooring track under the centroid rotating coordinate system; tau and alpha are time variables, wherein tau is used for representing the state quantity of the acquisition point of the track of the translational point, and alpha represents the state quantity of the corresponding moment of the invariant manifold so as to represent the manifold state quantity of the splicing point;
Figure BDA0002644639870000056
for accompanying control parameters, for solving co-state vectors λ in said co-state equationrAnd λvAn initial value of (d); t is tF1And tF2The starting time lengths of the first small thrust section and the second small thrust section are respectively set; t is tcThe flight time of the uncontrolled taxiing section.
According to an aspect of the invention, in the step (S3), a global search is performed on the design variables based on a genetic algorithm, and the splice point region and the initial value range of the design variables are determined, wherein the optimized performance index is:
Figure BDA0002644639870000061
wherein | | | Δ r | | and | | | Δ v | | are the distance vector and velocity vector difference between the second small thrust segment end and the stable manifold respectively, Δ i is the end inclination angle difference, and the coefficient
Figure BDA0002644639870000062
Symbol
Figure BDA0002644639870000063
Is a rounding operation.
According to an aspect of the present invention, in the step (S4), J is obtained by repeating the step (S3)GAAggregate and select the smallest JGAThereby screening the splicing point region and the initial value range of the design variable, wherein the repetition time is 10 times.
According to an aspect of the present invention, in step (d), the local optimization algorithm is used to iteratively solve the initial values of the design variables to obtain an optimal transfer splice point position, the performance index of the optimal transfer splice point position is the same as the performance index set when the optimal low thrust transfer orbit is solved in step (S1), and the constraint equation is:
Figure BDA0002644639870000064
in the formula, rs(τf,αf) And vs(τf,αf) The end position and the speed of the invariant manifold are measured; r (t)F2) And v (t)F2) Transferring the position and the speed of the terminal of the orbit for the small thrust;
further, the partial derivative relation of the constraint vector with respect to the design variable is obtained as follows:
Figure BDA0002644639870000071
finally, variable D is subjected to Newton optimization algorithmsAnd carrying out iterative solution to obtain the optimal transfer track meeting the constraint equation.
According to one scheme of the invention, after an aircraft dynamics equation is analyzed, the invariant manifold of the earth berthing orbit and the translation point orbit is analyzed, the global optimization search is carried out on the thrust vector of the small thrust in the near-earth flight stage, and the splicing point positions of different orbit segments and the appropriate initial values of design variables are determined. And further adjusting the positions of the earth escape point and the invariant manifold capture point by combining a local optimization algorithm, and the magnitude and direction of the thrust in the transfer process, so that the optimal fuel transfer orbit meeting the task constraint can be quickly and effectively designed.
According to one scheme of the invention, in the local optimization step, the size of the parking track, the positions of the initial escape point and the target track capture point and the thrust vector are automatically searched and adjusted, so that the influence of uncertain feature point selection on fuel consumption is effectively avoided, and the design reliability and result optimality are improved.
According to one scheme of the invention, different low-thrust transfer schemes are determined according to different translation point track task requirements under different three-body systems, so that corresponding track transfer design tasks can be completed, and references can be provided for future detectors to perform tasks such as lunar exploration, space environment observation and cargo transportation.
Drawings
In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings needed in the embodiments will be briefly described below, and it is obvious that the drawings in the following description are only some embodiments of the present invention, and it is obvious for those skilled in the art to obtain other drawings without creative efforts.
FIG. 1 schematically represents a schematic view of a centroid rotating coordinate system of the design method of the present invention;
FIG. 2 is a schematic representation of the stable and unstable manifolds of the pan point trajectory for a design method according to an embodiment of the present invention;
FIG. 3 schematically represents a transfer trajectory (upper side) and design variables and strategies (lower side) schematic of the design method of an embodiment of the present invention;
FIG. 4 is a schematic diagram of the optimum low thrust transfer orbit of the Earth's moon L1 point halo orbit of the design method of one embodiment of the present invention (3D on the left and 2D on the right);
FIG. 5 is a graph of mass versus thrust for the transfer process in the embodiment of FIG. 4;
FIG. 6 is a schematic diagram of a DAY L1 Point halo trajectory optimal low thrust transfer trajectory (3D on the left and 2D on the right) illustrating a design method according to an embodiment of the present invention;
FIG. 7 is a graph of mass versus thrust for the transfer process in the embodiment of FIG. 6;
FIG. 8 is a schematic diagram of a DAY L2 Point halo trajectory optimal low thrust transfer trajectory (3D on the left and 2D on the right) illustrating a design method according to an embodiment of the present invention;
fig. 9 is a graph of mass versus thrust for the transfer process in the embodiment of fig. 8.
Detailed Description
In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings used in the embodiments will be briefly described below. It is obvious that the drawings in the following description are only some embodiments of the invention, and that for a person skilled in the art, other drawings can be derived from them without inventive effort.
In describing embodiments of the present invention, the terms "longitudinal," "lateral," "upper," "lower," "front," "rear," "left," "right," "vertical," "horizontal," "top," "bottom," "inner," "outer," and the like are used in an orientation or positional relationship that is based on the orientation or positional relationship shown in the associated drawings, which is for convenience and simplicity of description only, and does not indicate or imply that the referenced device or element must have a particular orientation, be constructed and operated in a particular orientation, and thus, the above-described terms should not be construed as limiting the present invention.
The present invention is described in detail below with reference to the drawings and the specific embodiments, which are not repeated herein, but the embodiments of the present invention are not limited to the following embodiments.
Referring to fig. 1, according to the optimal low thrust transfer sectional design method of the translational point orbit, firstly, a kinetic equation of an aircraft is established under a centroid rotating coordinate system. The centroid rotating coordinate system may also be referred to as a centroid converging coordinate system, as shown in particular in fig. 1. The origin B in the coordinate system is the mass center of the three-system, and the three-system is composed of an aircraft, a first main celestial body and a second main celestial body, wherein the first main celestial body and the second main celestial body do circular motion around the mass center. If the first main celestial body is the earth, the second main celestial body is the moon, and if the first main celestial body is the sun, the second main celestial body is the earth or the moon. Thus, a three-system may also be referred to as a Earth-moon or a Sun-Earth/moon three-system, and the coordinate system may also be referred to as an Earth-moon (or Sun-Earth/moon) centroid rotation coordinate system. In the coordinate system, the X-axis is directed from the first main celestial body to the second main celestial body, the Z-axis is parallel to the angular momentum direction of the three-system (which is vertically outward and therefore cannot be drawn in FIG. 1), and the Y-axis is determined by the right-hand rule. The determination method of these coordinate axes is common knowledge in the field of aerospace and is therefore not described in more detail (r in the coordinate system is the distance of the aircraft from the centroid). The dynamic equations of the aircraft comprise an uncontrolled dynamic equation and a controlled dynamic equation (or a model), wherein the controlled dynamic equation introduces a small thrust engine of the aircraft to provide the thrust angular velocity, so that the name is given. Under the centroid rotating coordinate system, the uncontrolled flight dynamics equation of an aircraft (also called a spacecraft) is as follows:
Figure BDA0002644639870000091
wherein the content of the first and second substances,
Figure BDA0002644639870000092
mu is the mass coefficient of the system,
Figure BDA0002644639870000093
the state quantity of the aircraft under the mass center rotation system,
Figure BDA0002644639870000094
and
Figure BDA0002644639870000095
respectively the distance of the aircraft relative to the first main celestial body and the second main celestial body.
And the controlled kinetic equation is:
Figure BDA0002644639870000096
wherein the content of the first and second substances,
Figure BDA0002644639870000097
t is the magnitude of the thrust, uTIs a unit vector of the three axial directions of the thrust, and satisfies
Figure BDA0002644639870000098
P is aircraft engine power. Where m is the flight-to-mass, and therefore T/m is the thrust acceleration provided by the low thrust engine. And v is the speed of the aircraft.
Through the steps, the uncontrolled and controlled dynamic equations of the aircraft are respectively established, and then the invariant manifold of the track of the translational point is determined. In this step, a space invariant manifold is first determined according to the eigenvalues and eigenvectors of a matrix M, where the matrix M is:
M＝Φ(t0+T,t0)；
wherein the content of the first and second substances,
Figure BDA0002644639870000101
03and I3Respectively 3 × 3 order zero matrix and unit matrix;
the linear perturbation state quantities of the unstable manifold (left equation) and the stable manifold (right equation) of the translational point orbit at any time t are respectively:
Figure BDA0002644639870000102
wherein the superscripts "u" and "s" respectively represent unstable and stable manifolds in the invariant manifold,
Figure BDA0002644639870000103
is the eigenvector, perturbation d, of the matrix Mm200km in a sun-earth system and 50km in a moon-earth system, and the perturbation d in the sun-earth system is small because of the small gravitational force of the moonmCan be considered the same as the daily-terrestrial system. The symbol "±" corresponds to the invariant manifolds in two directions, as shown in fig. 2, the left and right sides of the "cross section" are the stable manifold and the unstable manifold, respectively, the upper side is the internal manifold, and the lower side is the external manifold.
After the invariant manifold is determined, global search (or global solution) can be performed for initial values of the design variables. In the aerospace field, splice points are typically represented using state quantities or coordinates. And the design variables are understood to represent state quantities of the splice point. This step therefore searches for the initial values of the design variables while also determining the approximate area of the splice point. In a specific global search step, the fuel with the most problems is modeled firstly, so that the fuel-optimal small thrust track is solved. When modeling is solved, the optimal performance index maxJ is as follows:
maxJ＝k·mf；
wherein k is any normal number, mfThe remaining mass of the aircraft at the moment of engine shutdown.
Combining the Pontryagin maximum principle with the above controlled kinetic model equation, the Hamilton function of the system was determined as:
Figure BDA0002644639870000111
the co-state equation of the optimal problem is:
Figure BDA0002644639870000112
in the formula, λ m is a covariate.
After the model is built, a set of design variables needs to be determined and analyzed and solved. Wherein, the design variable set is as follows:
Figure BDA0002644639870000113
in order to solve the problem that the near-earth-end multi-circle low-thrust solution is not easy to converge, the optimal low-thrust transfer trajectory is divided into a first low-thrust section A, an uncontrolled gliding section D and a second low-thrust section B when the design variable set is analyzed and solved, wherein the first low-thrust section A is connected with an earth mooring track C, as shown in figure 3. In addition, the thrust direction of the first small thrust section A is along the direction of the velocity vector, so that the height of the orbit can be effectively and quickly raised, the influence of the earth strong gravitational field is reduced in a short time, and the convergence of the algorithm is improved.
In the set of design variables,
Figure BDA0002644639870000114
and
Figure BDA0002644639870000115
the three are used for describing a near-earth end parking track under a centroid rotation coordinate system, and the definition of the three is similar to that of a centroid inertia coordinate system.
And tau and alpha are respectively time variables, wherein tau is used for representing the state quantity of the capture point of the translational point track, and alpha is used for representing the state quantity of the corresponding moment of the invariant manifold and respectively solving the state quantities of the translational point track and the corresponding moment of the manifold so as to represent the manifold state quantity of the splicing point.
Figure BDA0002644639870000121
For accompanying control parameters, for solving the co-state vector lambda in the above-mentioned co-state equationrAnd λvInitial value of (d), tF1And tF2The starting time length t of the first small thrust section and the second small thrust section respectivelycThe flight time of the uncontrolled taxiing section.
After the design variable set is determined, a global optimization search strategy is required, so that the approximate area of the splicing point and the approximate range of the initial value of the design variable can be determined. The global search in the invention is mainly carried out based on a genetic algorithm, wherein the optimized performance indexes are as follows:
Figure BDA0002644639870000122
the determination condition can ensure that the space relative distance and the speed direction of the small thrust end point and the manifold end point are close. Coefficient of performance
Figure RE-GDA0002743293810000123
Symbol (sign)
Figure RE-GDA0002743293810000124
Is a rounding operation.
The global search performed in the previous step can obtain the rough area of the splicing point and the rough range of the initial value of the design variable. Therefore, a single global optimization may cause the search to have randomness, so the invention repeatedly performs the global search to avoid the defect. Specifically, the number of repeated searches is 10, and the J can be obtained finallyGAA collection of (a). Based on the principle of fuel economy, the smallest J should be selected from the setGAAnd screening out the final splicing point area and the initial value of the design variable. It should be noted that the initial value of the final design variable obtained in this step is also an approximate value, so the area of the splicing point is also an approximate area, which is only more accurate with respect to a single search. Therefore, the initial value of the design variable cannot guarantee the position continuity of the low-thrust trajectory and the stable manifold, and local optimization is needed to determine the accurate splice point position.
Specifically, the initial values of the design variables are further iteratively solved by adopting a local optimization algorithm, so that the position of the fuel optimal transfer splicing point is obtained. Because the position of the splicing point is positioned on the optimal transfer track, the setting of the performance index is the same as the optimal low-thrust transfer track, and the constraint equation is as follows:
Figure BDA0002644639870000131
in the formula, rs(τf,αf) And vs(τf,αf) The end position and the speed of the invariant manifold are measured; r (t)F2) And v (t)F2) The terminal position and the speed of the orbit are transferred for small thrust.
The positions and the speeds of the small thrust end and the stable manifold end can be made continuous through the above formula, and then the partial derivative relation of the constraint vector relative to the design variable is derived as follows:
Figure BDA0002644639870000132
finally, variable D is subjected to Newton optimization algorithmsAnd carrying out iterative solution to obtain the fuel optimal low-thrust transfer trajectory meeting the constraint equation.
In order to verify the feasibility of the method, the invention designs the optimal low-thrust transfer orbit of the Earth-moon L1 point Halo orbit (namely Halo orbit). Assuming that the initial mass of the aircraft is 1500kg, firstly, a Halo orbit with a normal amplitude Az of 10000km at a point L1 is selected as a target orbit, and by combining the design steps, the optimal small thrust trajectory is shown in fig. 4, and the corresponding mass change and thrust magnitude change are shown in fig. 5. Wherein, a translation point Halo track, a stable manifold, a small thrust transfer track (first and second are not distinguished) and an uncontrolled sliding track are respectively indicated in the figure. The analysis shows that: the aircraft is started from a GEO orbit low-thrust engine with the inclination angle of 4.7 degrees, the ascension intersection declination angle of 15.8 degrees and the latitude amplitude angle of 28.9 degrees, the thrust is 0.6N, the aircraft is started for 41.4 days along the speed direction, the uncontrolled gliding section flies for about 2.35 days, then the engine is started again for 38.72 days to enter the stable manifold, and the aircraft flies to the target Halo orbit after continuing to navigate along the stable manifold for 22.7 days. The whole process takes about 105.1 days, and the mass consumption is about 49 kg.
In order to verify the effectiveness of the proposed design strategy, a three-system-in-the-sun system L1 and a three-system-in-the-earth system L2 Halo track are selected to complete corresponding track design tasks, and the small thrust transfer track and the key parameter change under different working conditions are shown in FIGS. 6 to 9. Aiming at a transfer task of a day-ground system L1 point, the inclination angle and the ascension point right ascension of an earth mooring orbit are respectively 7.9 degrees and 308.9 degrees, the flight time of a first small thrust section A, an uncontrolled gliding section D and a second small thrust section B are respectively 56.3 days, 4.1 days and 35.8 days, the aircraft reaches a target orbit after going into a stable manifold and continuing to sail for about 188.8 days, the total flight time of the whole task is about 285.01 days, 77.05kg of propellant is consumed, and the maximum thrust of a small thrust arc section is less than 0.9N. Likewise, for the L2 point Halo orbital mission, the entire mission flies for a total of about 310.61 days, the propellant consumption is about 75.39kg, and the maximum thrust is less than 0.45N.
In conclusion, the invention is based on a circular restrictive three-system formed by an earth-moon-detector or a sun-earth-detector, global search is carried out by combining invariant manifold, the position of a proper splicing point and the initial value of a design variable are rapidly determined, then the optimal low-thrust transfer trajectory of fuel is locally optimized and solved, and the transfer design of the detector from the earth to the earth-moon system or the solar-earth system translation point is realized. The method can effectively solve the problems of low-thrust near-ground multi-circle divergence and initial value guessing, has wide application range, and can be suitable for the transfer design of the translation point track of the Earth-moon point L1, the Sun-Earth point L1 and the L2. The translation point orbit low thrust transfer strategy is suitable for a space task without time limitation, the propellant consumption of the whole task is low, and the first low thrust section ensures that the effectiveness and the feasibility of the design strategy can be effectively improved along the speed direction. The transfer characteristics of different translation point tracks of different three-system systems are researched and analyzed, and the method has reference value in the aspects of ground-moon space material transportation and expansibility tasks, strategy and parameter selection of a future continuous low-thrust aircraft.
The above description is only one embodiment of the present invention, and is not intended to limit the present invention, and it is apparent to those skilled in the art that various modifications and variations can be made in the present invention. Any modification, equivalent replacement, or improvement made within the spirit and principle of the present invention should be included in the protection scope of the present invention.

Claims (9)

1. A translational point orbit optimal low thrust transfer subsection design method comprises the following steps:
a. establishing a dynamic equation of the aircraft under a centroid rotating coordinate system;
b. determining the invariant manifold of the translation point track;
c. global searching is carried out to determine initial values of design variables;
d. and determining the optimal transfer track by local optimization.
2. The design method according to claim 1, wherein in the step (a), the origin of the centroid rotating coordinate system is a centroid of a three-body system consisting of the aircraft, the first main celestial body and the second main celestial body;
wherein, if the first main celestial body is the earth, the second main celestial body is the moon, and if the first main celestial body is the sun, the second main celestial body is the earth or the moon;
the X axis points to the second main celestial body from the first main celestial body, the Z axis is parallel to the angular momentum direction of the three-system, and the Y axis is determined by a right-hand rule;
the dynamic equations comprise an uncontrolled dynamic equation and a controlled dynamic equation;
the uncontrolled kinetic equation is as follows:
Figure FDA0002644639860000011
wherein the content of the first and second substances,
Figure FDA0002644639860000012
mu is the mass coefficient of the system,
Figure FDA0002644639860000013
the state quantity of the aircraft under the mass center rotation system,
Figure FDA0002644639860000014
and
Figure FDA0002644639860000015
the distances of the aircraft relative to the first main celestial body and the second main celestial body respectively;
the controlled dynamics equation is:
Figure FDA0002644639860000021
wherein the content of the first and second substances,
Figure FDA0002644639860000022
t is the magnitude of the thrust, uTIs a unit vector of the three axial directions of the thrust, and satisfies
Figure FDA0002644639860000023
P is the aircraft engine power, m is the aircraft mass, and v is the aircraft velocity.
3. The design method according to claim 2, wherein in step (b), the spatially invariant manifold is determined according to eigenvalues and eigenvectors of a matrix M, where M is:
M＝Φ(t0+T,t0)；
wherein the content of the first and second substances,
Figure FDA0002644639860000024
03and I3Respectively 3 × 3 order zero matrix and unit matrix;
the linear perturbation state quantities of the unstable manifold and the stable manifold of the translational point orbit at any time t are respectively:
Figure FDA0002644639860000025
wherein the superscripts "u" and "s" respectively represent unstable and stable manifolds in the invariant manifold,
Figure FDA0002644639860000026
is the eigenvector, perturbation d, of the matrix Mm200km in a daily terrestrial system and 50km in a terrestrial lunar system, the sign "±" corresponding to the invariant manifold in both directions.
4. The design method according to claim 3, wherein in the step (c), the search method comprises:
s1, modeling the fuel optimization problem;
s2, determining and analyzing a design variable set;
s3, global search is carried out to determine the splicing point area and the initial value range of the design variable;
s4, repeating the step (S3) to screen the splicing point area and the initial value range of the design variable.
5. The design method according to claim 4, characterized in that in said step (S1), a fuel optimum low thrust shift trajectory is first determined, the optimum performance index max J of which can be expressed as:
max J＝k·mf；
wherein k is an arbitrary normal number, mfAnd determining the Hamilton function of the system for the residual mass of the aircraft at the moment of engine shutdown by combining the Pontryagin maximum principle and the controlled dynamics equation:
Figure FDA0002644639860000031
the co-state equation of the optimal problem is:
Figure FDA0002644639860000032
wherein λ m is a covariate.
6. The design method according to claim 5, wherein the set of design variables in the step (S2) is:
Figure FDA0002644639860000033
during analysis, the optimal low thrust transfer track is divided into a first low thrust section, an uncontrolled gliding section and a second low thrust section, wherein the thrust direction of the first low thrust section is along the direction of a velocity vector;
wherein the content of the first and second substances,
Figure FDA0002644639860000034
and
Figure FDA0002644639860000035
the elevation angle, the elevation intersection point declination and the latitude argument are respectively used for describing the near-ground end mooring track under the centroid rotating coordinate system; tau and alpha are time variables, wherein tau is used for representing the state quantity of the acquisition point of the track of the translational point, and alpha represents the state quantity of the corresponding moment of the invariant manifold so as to represent the manifold state quantity of the splicing point;
Figure FDA0002644639860000041
for accompanying control parameters, for solving co-state vectors λ in said co-state equationrAnd λvAn initial value of (d); t is tF1And tF2The starting time lengths of the first small thrust section and the second small thrust section are respectively set; t is tcThe flight time of the uncontrolled taxiing section.
7. The design method according to claim 6, wherein in the step (S3), a global search is performed on the design variables based on a genetic algorithm to determine the splicing point region and the initial value range of the design variables, wherein the optimized performance index is:
Figure FDA0002644639860000042
wherein | | | Δ r | | and | | | Δ v | | are the distance vector and velocity vector difference between the second small thrust segment end and the stable manifold respectively, Δ i is the end inclination angle difference, and the coefficient
Figure FDA0002644639860000043
Symbol
Figure FDA0002644639860000044
Is a rounding operation.
8. The design method according to claim 7,in the step (S4), J is obtained by repeating the step (S3)GAAggregate and select the smallest JGAThereby screening the splicing point region and the initial value range of the design variable, wherein the repetition time is 10 times.
9. The design method according to claim 8, wherein in the step (d), the initial values of the design variables are iteratively solved by using a local optimization algorithm to obtain the optimal transfer splice point position, the performance index of the optimal transfer splice point position is the same as the performance index set when the optimal low-thrust transfer orbit is solved in the step (S1), and the constraint equation is as follows:
Figure FDA0002644639860000045
in the formula, rs(τf,αf) And vs(τf,αf) The end position and the speed of the invariant manifold are measured; r (t)F2) And v (t)F2) Transferring the position and the speed of the terminal of the orbit for the small thrust;
further, the partial derivative relation of the constraint vector with respect to the design variable is obtained as follows:
Figure FDA0002644639860000051
finally, variable D is subjected to Newton optimization algorithmsAnd carrying out iterative solution to obtain the optimal transfer track meeting the constraint equation.
CN202010850752.2A 2020-08-21 2020-08-21 Optimal low-thrust transfer sectional design method for translation point orbit Active CN112009727B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010850752.2A CN112009727B (en) 2020-08-21 2020-08-21 Optimal low-thrust transfer sectional design method for translation point orbit

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010850752.2A CN112009727B (en) 2020-08-21 2020-08-21 Optimal low-thrust transfer sectional design method for translation point orbit

Publications (2)

Publication Number Publication Date
CN112009727A true CN112009727A (en) 2020-12-01
CN112009727B CN112009727B (en) 2021-11-23

Family

ID=73505523

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010850752.2A Active CN112009727B (en) 2020-08-21 2020-08-21 Optimal low-thrust transfer sectional design method for translation point orbit

Country Status (1)

Country Link
CN (1) CN112009727B (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113968360A (en) * 2021-08-09 2022-01-25 中国空间技术研究院 Satellite autonomous electric propulsion orbit transfer method for stationary orbit satellite

Citations (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105574261A (en) * 2015-12-15 2016-05-11 北京理工大学 Method for designing earth-moon libration point transfer orbit via moon leveraging constraint
CN105912819A (en) * 2016-05-06 2016-08-31 北京理工大学 Quick design method of earth-moon L1 Lagrange point transfer orbit
CN106379555A (en) * 2016-09-05 2017-02-08 北京理工大学 Optimal orbital transfer method of low-earth-orbit satellite under limited thrust by taking J2 perturbation into consideration
CN107992682A (en) * 2017-12-05 2018-05-04 北京理工大学 A kind of optimal multiple-pulse transfer method of interplanetary multi-body system asteroid detection
CN108100307A (en) * 2017-12-05 2018-06-01 北京理工大学 One kind is for low energy small feature loss accurately detecting Orbital Transfer under Complex Constraints
CN108674695A (en) * 2018-07-19 2018-10-19 北京空间技术研制试验中心 Spacecraft reentry return trajectory planing method
CN109344449A (en) * 2018-09-07 2019-02-15 北京空间技术研制试验中心 The spacecraft moon ground transfer orbit Reverse Design
CN110046436A (en) * 2019-04-19 2019-07-23 北京空间技术研制试验中心 The autonomous observing and controlling covering analyzing method of spacecraft
CN110736470A (en) * 2019-11-06 2020-01-31 北京理工大学 method for hybrid search of continuous thrust tracks near small irregular celestial bodies
CN110733667A (en) * 2019-09-29 2020-01-31 北京空间技术研制试验中心 Design method for inter-orbit transfer of earth-moon translation point

Patent Citations (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105574261A (en) * 2015-12-15 2016-05-11 北京理工大学 Method for designing earth-moon libration point transfer orbit via moon leveraging constraint
CN105912819A (en) * 2016-05-06 2016-08-31 北京理工大学 Quick design method of earth-moon L1 Lagrange point transfer orbit
CN106379555A (en) * 2016-09-05 2017-02-08 北京理工大学 Optimal orbital transfer method of low-earth-orbit satellite under limited thrust by taking J2 perturbation into consideration
CN107992682A (en) * 2017-12-05 2018-05-04 北京理工大学 A kind of optimal multiple-pulse transfer method of interplanetary multi-body system asteroid detection
CN108100307A (en) * 2017-12-05 2018-06-01 北京理工大学 One kind is for low energy small feature loss accurately detecting Orbital Transfer under Complex Constraints
CN108674695A (en) * 2018-07-19 2018-10-19 北京空间技术研制试验中心 Spacecraft reentry return trajectory planing method
CN109344449A (en) * 2018-09-07 2019-02-15 北京空间技术研制试验中心 The spacecraft moon ground transfer orbit Reverse Design
CN110046436A (en) * 2019-04-19 2019-07-23 北京空间技术研制试验中心 The autonomous observing and controlling covering analyzing method of spacecraft
CN110733667A (en) * 2019-09-29 2020-01-31 北京空间技术研制试验中心 Design method for inter-orbit transfer of earth-moon translation point
CN110736470A (en) * 2019-11-06 2020-01-31 北京理工大学 method for hybrid search of continuous thrust tracks near small irregular celestial bodies

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113968360A (en) * 2021-08-09 2022-01-25 中国空间技术研究院 Satellite autonomous electric propulsion orbit transfer method for stationary orbit satellite
CN113968360B (en) * 2021-08-09 2024-03-15 中国空间技术研究院 Autonomous electric propulsion orbital transfer method on stationary orbit satellite

Also Published As

Publication number Publication date
CN112009727B (en) 2021-11-23

Similar Documents

Publication Publication Date Title
CN110733667B (en) Design method for inter-orbit transfer of earth-moon translation point
Singh et al. Eclipse-conscious transfer to lunar gateway using ephemeris-driven terminal coast arcs
CN105511493B (en) A kind of low rail constellation dispositions method based on martian atmosphere auxiliary
CN109911249A (en) The interspace transfer Finite Thrust of low thrust ratio aircraft is entered the orbit interative guidance method
CN109613822A (en) Control method is unfolded in spatial tether system based on Nonlinear Model Predictive Control
CN112009727B (en) Optimal low-thrust transfer sectional design method for translation point orbit
Ranieri et al. Indirect optimization of three-dimensional finite-burning interplanetary transfers including spiral dynamics
He et al. Solution set calculation of the Sun-perturbed optimal two-impulse trans-lunar orbits using continuation theory
Li et al. Indirect transfer to the Earth–Moon L1 libration point
Xie et al. Guidance, navigation, and control for spacecraft rendezvous and docking: theory and methods
Samsam et al. Multi-impulse shape-based trajectory optimization for target chasing in on-orbit servicing missions
Ploen et al. Dynamics of earth orbiting formations
Golan et al. Minimum fuel lunar trajectories for a low-thrust power-limited spacecraft
Tang et al. Low-thrust trajectory optimization of asteroid sample return mission with multiple revolutions and moon gravity assists
Ulybyshev Study of optimal transfers from L2 halo-orbits to lunar surface
CN110155370B (en) Transverse formation method based on solar sails
Chen et al. Trajectory design for the Moon departure libration point mission in full ephemeris model
Lee et al. Numerical analysis of relative orbit control strategy for CANYVAL-X mission
Atkins et al. Autonomous satellite formation assembly and reconfiguration with gravity fields
Nan et al. Global 4-D trajectory optimization for spacecraft
Cuevas et al. Dynamics, Guidance and Control for Autonomous Rendezvous and Docking in the Restricted Three Body Problem
Li et al. Earth-Phobos transfer with ballistic trajectory in the Sun-Mars system
Shirobokov et al. Parametric analysis of low-thrust lunar transfers with resonant encounters
Cai et al. Solar sailing trajectory optimization with planetary gravity assist
Lee et al. Poincaré filter for trajectory design in multi-body regime

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant