Development of method for selecting motion control laws of space optical system on based diffractive membranes during transfer into geostationary orbit G.P. Anshakov1, V.V. Salmin1, A.S. Chetverikov1, K.V. Peresypkin1, I.S. Tkachenko1 1 Samara National Research University, 34 Moskovskoe Shosse, 443086, Samara, Russia Abstract The article presents a formulation of a problem of trajectory optimization using low-thrust engines for an optical space system based on diffractive membranes. A methodology, where first stage nominal trajectories and control programs are selected and then corrected at the long- range guidance, has been developed for solving the problem of optimizing the trajectories of a flight to a geostationary orbit. At the final stage, algorithms for terminal control are formed, which allows to deliver a cosmic optical system based on diffraction membranes to a given point in the geostationary orbit. The end result is acquisition of Pareto-optimal solutions in the coordinates "characteristic speed-duration of the flight", where each point of the set of solutions has a corresponding a measure of accuracy of payload delivery to a geostationary orbit at a given set of coordinates. Keywords: cosmic optical system; diffraction membrane; interorbital flight; geostationary orbit; multicriteria problem; terminal control algorithms; low-thrust engine 1. Introduction Modern remote Earth sensing satellites function on low orbits. This imposes several constraints on their performance. For example, it is impossible to perform continuous monitoring of an object from a single satellite. In order to do so the satellite must be delivered to the geostationary orbit. Since the geostationary orbit is quite high (about 36000 km), obtaining high-quality images requires complex optical systems. Remote geostationary Earth sensing satellites fitted with traditional refractive payloads capable of obtaining high- quality images are inevitably heavy which complicates their delivery to the orbit. To solve this problem, the US Agency for Advanced Defense Research and Development of the US DARPA is developing a project of a modern space telescope with a membrane diffractive optical system. The project was named MOIRE or Membrane Optical Imager for Real-Time Exploitation [1-3]. The spacecraft with a membrane diffractive optical system will have a much smaller mass compared to a spacecraft with a refractive optical system. However, such a spacecraft would have large dimensions - the diameter of the lens is of the order of 10-20 m, the distance from the lens to the spacecraft body is of the order of 50 to 70 m. The paper [4] presents a design of the load-bearing structure for an Earth remote sensing satellite with diffractive optics (Figure 1). Fig. 1. Load-bearing structure for an Earth remote sensing satellite with diffractive optics [4]. When such a structure is being propelled by a chemical booster block, significant overloads can occur, which can lead to undesirable structural changes of the diffractive observation system. In this case, it is preferable to use low-thrust electric propulsion engines, that are creating accelerations of the order of 0.5.1.0 mm/s2, in order to bring the cosmic optical system on the basis of a large diffraction membrane to the geostationary orbit (GSO). The duration of transportation from low Earth orbit to GSO would be in the range from 100 to 200 days. When optimizing ballistic schemes for such flights, it is necessary to seek a compromise between the mass of the payload and the duration of the flight - the main efficiency criteria [5]. The problem of ensuring the required accuracy of delivery is also very important. 3rd International conference “Information Technology and Nanotechnology 2017” 35 Computer Modeling / G.P. Anshakov et al. Since spacecraft of this type have significant dimensions and, accordingly, mass-inertial characteristics, the process of motion control is significantly hampered. This calls for a solution of the problem of joint optimization of trajectory and angular movements. 2. Statement of the Problem Let us consider a ballistic scheme for a transfer of a spacecraft from initial circular orbit to operational (geostationary) orbit, with insertion into a given station and return of the space tug to the initial orbit. The optimization problem in a general statement is shaped as the problem of tied choice of design parameters, p  P , ballistic parameters, b B and set of functions u(t,x), x(t) out of allowable multitude D, ensuring the transfer of a spacecraft from the initial state x (t 0 )  x0 to a finite multitude x(t f )  X f with the maximum value of the optimization criterion. Relative payload mass (relation of the payload mass Мpl to launch mass of the spacecraft М0) is adopted as the main optimality criterion  : М pl   max . М0 Another criterion, no less important, is the overall transfer time ТΣ, which is a sum of the powered flight time ТP and the time of unpowered legs ТUP. The latter result from the need to make navigational observation and to "phase" precision targeting of an EP-powered spacecraft. Let us represent the launch weight of a spacecraft as the sum of the masses of its functional elements: power plant (consisting of the reactor, energy transformer and radiator); thruster unit (consisting of cruising EP thrusters and controlling thrusters with their controls); propellant supply for the direct and return transfer, including additional expenses for control; propellant storage and supply system; payload; the body of the spacecraft: М 0  М PP  М EP  М P1  М P 2  М PSS  М PL  М B , (1) where М0 is the launch mass of the spacecraft; МP1, МP2 is the propellant mass for the direct and the return transfer respectively; МPL is payload mass; МPP is the power plant mass; МEP is electric propulsion mass; МPSS is the propellant storage and supply system; МB is the spacecraft body mass. The criterion µ may be represented as [6]:  1  Az1  z 2  z1 z 2    B    max   p, b, u (t ), T , М 0   max   , (2) u (t )U u (t )U  1  Az 2  pP pP с   PP  c  where A  1   PSS     EP  ; αPP - is the relative mass of the power plant; γEP - is the relative mass of electric T  2  P  propulsion; μB is the relative mass of the body; γPSS is the relation of the propellant storage and supply system mass to propellant  W f1  mass; Р is the thrust of the cruising and controlling thrusters; с is the thruster exhaust velocity; z1  1  exp   ;  c   Wf2  z 2  1  exp    ; Wf1, Wf2 – is the final characteristic velocity expense for direct and return transfer.  c   In this statement, the optimization problem is traditionally divided into two parts: dynamic (selection of optimal trajectories and control laws) and parametrical (selection of optimal design parameters). Thus, two main optimality criteria are set: the relative payload mass μ or, considering the expression (2), the final characteristic velocity expense Wf, as well as the total transfer time TΣ. Then the optimization of the ballistic schemes, trajectories, and control modes as the ultimate goal is reduced to building a Pareto set in coordinates Wf - TΣ. In addition, every point in the Pareto set must correspond to a measure of meeting the final boundary conditions Φ, for example, in the following expression: Φ  x T f x f . Here x f is the vector of deviation of the final values of the state vector from required values; Λ is a positively determined weighted coefficient matrix. The general mathematical model of motion includes the equation of the motion of the center of mass in equinoctial elements, dynamic equations of angular motion, kinematic equations in quaternion form. 3. The expansion principle in solution of a dynamic problem The optimal control problem in this statement is extremely difficult, as the state of the controlled object is characterized by a set of variables, some of which are "slow-changing" (trajectorial motion) and some are relatively "fast-changing" (angular motion); besides, the mathematical model of motion includes a perturbation vector. It is apparently problematic to approach this problem with the classical methods of optimal control (the maximum principle and dynamic programming). Therefore, we propose an approach to the solution of the dynamic problem, which is based on the principle of the extension – narrowing of the class of admissible states and controls [7]. 3rd International conference “Information Technology and Nanotechnology 2017” 36 Computer Modeling / G.P. Anshakov et al. The control vector u will include only the angles of the thrust orientation vector, and other components (controlling torques МХ, МY, MZ) will be subject to constraints. Let us impose constraints also on the state vector x at the final point of the trajectory. In this manner we obtain the following mathematical model: x*  f x* ,u*  , x* t 0   x*0 (3) where x*  h, e X , eY , i X , iY , F T , u*   , T , with extended functional T  L  V f  adt  min (4) 0 and constraints x* t f   x* f  h f , e xf , e yf , i xf , i yf , F f T  X f , max М x (t )  М max_x ; max М y (t )  М max_ y ; max М z (t )  М max_z . (5) t t t Here h, e x , e y , i x , i y , F - are the equinoctial elements;  , - are the thrust orientation vector angles; а – is reactive acceleration; МMAX_Х, МMAX_Y, MMAX_Z are the maximum possible controlling torques that can be created by controlling thrusters of the spacecraft. The conditions of (5) are checked as the result of a numerical modeling of the basic system of equations, with optimal control found by solving the problem for the simplified mathematical model. Assume that the spacecraft is moving along a near circular orbit, so the eccentricity can be taken as zero. Assume also that the radial component of the acceleration is zero. Then the thrust direction is determined by the angle ψ between the transversal and the thrust vector, and the projections of the reactive acceleration to the axis of the orbital system of coordinates are: P P aT  cos , aS  0, aW  sin (6) M M Let us also exclude from the differential components the equations that describe the changes in the longitude of the ascending node and the perigee argument. Fig. 2. Position of the spacecraft orbit and thrust vector control. The system (3), taking into consideration (6), is transformed to the equations of a near-circular motion of a small-thrust spacecraft. By averaging the "slow-changing" variables r (average radius of a near-circular orbit, equivalent to semi-major axis) and i (orbit inclination) along the "fast-changing" variable u (polar angle of the argument of latitude) we can obtain the "asymptotic" model of motion [8]. Thrust vector control for transfers between non-coplanar orbits requires the sign of the аW to change twice per revolution. The thrust angle orientation control is set in the following way:  (W , u )   m (W ) signcos u  , (7) where ψm is the amplitude of oscillations of the angle ψ, W is current characteristic velocity and u is the argument of latitude. Analytical solutions, obtained within the "asymptotic" model, describe a "universal" trajectory of a transfer between non- coplanar orbits, that does not depend on the spacecraft's design parameters (Figure 3). 4. Method of selecting motion control laws of the space optical system on based diffractive membranes during transfer to given point on the geostationary orbit The method includes an algorithm of developing nominal programs for thrust vector control, an algorithm of EP thrust magnitude adjustment, algorithm of terminal control development using motion models in discrete setting, numerical algorithm of building a set of Pareto-optimal solutions. The goal of the control at the stage of lifting to GEO is narrowing the area G to an allowable area Ga. The problem will be solved consequentially, in two stages: - development of an algorithm for moving the final state deflection vector ΔХf into an area G’, where one or more of the components of the vector ΔХf satisfy the required precision (e.g. deflection by inclination ∆if); 3rd International conference “Information Technology and Nanotechnology 2017” 37 Computer Modeling / G.P. Anshakov et al. - developing the control laws and algorithms for narrowing the area G’ to area Ga, where all components of the vector ΔХf satisfy the required precision conditions. Fig. 3. Phase trajectory of a transfer to GEO. 4.1 Goal of control at the long-range guidance stage If the deflections are considerable, or it is impossible to build in advance a precise model of perturbing accelerations, which may change significantly during the transfer, it is advisable to use multistage control algorithms with end state prognosis and identification of perturbations. Since sufficiently precise models of perturbations from the Earth's gravity field, solar and lunar perturbations are presently available, the parameter to be adjusted is the magnitude of reactive acceleration. Thrust magnitude has been chosen as the adjusted parameter. Adjustment of actual thrust magnitude will be carried out according to the expression:  Ti  P i  Pnom 1  , (8)  T calc     where Tcalc is the oscillating orbit time at the measured unpowered leg; Ti is the deviation of oscillating orbit time from calculated value. An example of modeling the trajectorial motion of an EP-powered spacecraft during a transfer to GEO with set values of initial orbit (А0, i0, e0), exhaust speed с and initial acceleration а0 with adjustment of the thrust magnitude depending on the number of corrections N is represented in Table 1. Here ∆P is the systematic error in thrust magnitude; tup is the length of unpowered legs, where the low thrust propulsion unit is, as a rule, shut down to ensure more precise orbit parameter calculation and consequent adjustment of thrust magnitude. Table 1. Results of modeling the motion of a EP OTV during transfer to GEO with adjustment of thrust magnitude (a0 = 0,0006 m/s2, c = 25 km/s, tup = 0,5 day, i0 = 51,60, e0 = 0, A0 = 7171 km, ΔP = -1,5% (Рr = 11,82 N)) N if, deg ef Wf, km/s Thrust and length of the i-th powered A f  A f  Areal f , km flight leg i Pi , N Tai , days 1 0 12 63,278 -464,1 -0,72 0,002 7,670 1 11,665 66,079 2 0 12 63,278 -12,8 -0,02 0,004 7,614 1 11,665 33,039 2 11,789 32,222 4.2 Development of control laws and algorithms at the final stage of GEO transfer The goal of the terminal control. The goal of the control is to move the end state deviation vector ΔХК from area G’ to area Ga. Assume that the correction maneuver is performed with the help of transversal low thrust, which creates the transversal acceleration аТ. This problem is set as a terminal control problem with functional Φ  х Tf х f  min . (9)  Here Λ is the fixed coefficient matrix; x f  T f ,  f , e f T .   The control is structured as a sequence of powered and unpowered leg lengths u   1 ... i , t up1 ...t upi T .  3rd International conference “Information Technology and Nanotechnology 2017” 38 Computer Modeling / G.P. Anshakov et al. Deviation of the semi-major axis ΔA is equivalent to the deviation of the orbit time ΔT = T – TS. Here the orbit time on GEO equals star day TS = 86164,09 s. In addition, the position of the OTV on GEO is described by longitude λ, which deviates from the required value of the station longitude λS by Δλ = λ – λS. We shall solve this problem by increasing the sophistication of the motion model. 1. Neglecting the eccentricity value due to its insignificance, we obtain a near optimal analytical solution of the problem by Hamilton-Jacobi-Bellman's formalism. The characteristic velocity budget for the correction maneuver with аТ = const is determined by the relation: T0 W f  . (10) TS  T0 3(TS  T0 )3 2 E The total time of the maneuver is determined by the equation  0 T   . (11)  2     Т  Т   E   S 0  Here τ and tup are the lengths of the powered and unpowered flight legs respectively. Therefore, the characteristic velocity budget and time required for execution of the maneuver do not depend on design parameters (transversal acceleration аТ), and are determined only by the initial deflections ΔТ0 and Δλ0. Near optimal control by orbit time and longitude is performed in one step. Here by step we mean a sequence of a powered and an unpowered legs. 2. The problem of terminal control is solved with the help of multistage control algorithm with control parameters adjustment. Let the control law be set by a sequence of powered legs that is taken as decreasing and is defined by the expression [12]:   i  1 b   i  а  1    , (12)   n   where i, n are the number of adjustment and total number of adjustments respectively; а and b are the parameters that describe the law of decreasing the lengths of the powered legs. Then the problem of determining the optimal control law is reduced to a bi-parametric optimization problem, which is stated in the following way: for set initial values of orbital elements, acceleration аТ, number of adjustments n, lengths of unpowered flight legs tup one should find such parameters а and b that ensure the minimum of the functional (9). Control parameters а and b are found as the result of minimizing the functional of the view (9) and at that for better precision the dependency of the functional from parameter а is approximated by the least squares methods. When the control is adjusted (during motion modeling accounting for perturbations) at every unpowered leg the number of adjustment steps n is also corrected. A series of calculations of control laws for a transfer of and EP-powered spacecraft to a given station by orbit time and longitude were carried out [9]. The characteristic velocity W budget, depending on the initial deviation by orbit time (∆Т0 = 300…1000 s) is on the order of 4 to 14 m/s. 3. A discrete model for flat motion of a spacecraft under small transversal acceleration (with changes in orbit time, eccentricity and station longitude) was developed in [11], and an analytical solution for the problem of the search for the optimal control (lengths of powered and unpowered legs), minimizing final orbit time, eccentricity and station longitude errors, was obtained. In [23] the authors propose an approximate method for solving the problem based on the three-step control algorithm for circulation period, eccentricity and longitude of the point of standing. Orbit correction is carried out using low-thrust electric rocket engine that produces acceleration in the transversal direction. Discrete model of plane motion of geostationary satellites under the influence of small transversal acceleration is presented in as [11]: T3  T (k ) T (k  1)  T (k )  3aТ T3  T (k ) 3 τ( k ), 2π  μ  2π  Δλ(k  1)  Δλ(k )    ω3   tп (k )  τ( k ) , T  3  Δ T (k)  1/3  Т  ΔT (k )  e(k  1)  Δe(k )  2  aT  3  τ( k ) ,  2π  μ  where k = 0,…, N – 1, tП(k) is determined by formulas Т τ T tп  0 (1  2m)   0 0 , for аТ > 0, 2 2 2π 3rd International conference “Information Technology and Nanotechnology 2017” 39 Computer Modeling / G.P. Anshakov et al. τ T tп  Т 0 (1  m)   0 0 , for аТ < 0. 2 2π Functional I  X KT X K  min . Where aT – transversal acceleration,  0 - true anomaly angle before correction; ∆ХК = {∆TК, ∆λК, ∆eК}Т – the final state vector, where ΔТК = ТК – ТЗ, ΔеК = еК – еGEO, ΔλК = λК – λР; ТК, еК, λК – values of the orbital period, eccentricity and longitude of the satellite’s standing point on the orbit at the end of the correction maneuver; ТЗ – circulation period of the spacecraft in geostationary orbit, equal to a star day ТЗ = 86164,09 с; еGEO –eccentricity of the geostationary orbit; λР – longitude of the working point of standing of a satellite; ΔТ0 = Т0 – ТЗ, Δе0 = е0 – еGEO, Δλ0 = λ0 – λР, where Т0, е0, λ0 – values of the orbital period, eccentricity and longitude of the point of standing on the orbit before the spacecraft correction maneuver. Based on the proposed control structure an analytical solution is obtained for  0 , 1 ,  2 , t П1 , t П 2 [11] e0 1  TС  τ0  , τ1  1  1  , T3  T0  T  TС  1/3   3  T3  TС   2  aт 3 aт  3  2π  μ  2π  μ  1  TС  TС τ2  1  1    1 , (13)  T3  TС   1/3  3  T3  T  С   3T3  TС  a т    2π  μ  T3  T0 where TС  T0  3aт T3  T0 3 0 ; 2π  μ  1   1 S 2  tп2  1  3S TС   m   τ1     , for аТ > 0, (14)   2   2(1  3S ) 2  3S    1 S 2  tп2  1  3S TС 1  m  τ1     , for аТ < 0;   2(1  3S ) 2  3S  TС where m Z , S  1  1  ; 3TС  2π  tп1  λ С  λ B    ω3  , (15)  T3  TС       2π   2π  where λC  2π     T  ΔT  ω3   τ1    ω3   tп2  τ 2 ,  3   T3  TС    T3  TС   1  3aт 3 С τ1  2ππ       2π  λ B  λ 0    ω3   τ 0 .  T3  T0  The presented algorithm has shown high accuracy in modeling of the correction of orbit for a surveillance satellite, fitted with an electric propultion engine. For example, for ΔТ0 = 1000 s, e0 = 0,005, Δλ0 = 0,087 rad the final orbit parameters deviations were: orbital period ΔТK = 1,3 s, standing point longitude ΔλK = 0,150, eccentricity ΔeK = 1×10-4. Durations of the active and passive sections were τ0 = 7758 s, τ1 = 1997 s, , τ2 = 1998 s, tп1 = 260200 s ≈ 3 days, tп2 = 40170 s ≈ 0,46 days. Figure 4 and 5 shows an example simulation of the orbit control of geostationary spacecraft by using low thrust of EP. Fig. 4. Simulation the orbit correction for geostationary spacecraft using low thrust of EP (а0 = 0,001 м/s2, ΔТ0 = 1000 s, e0 = 0,005, Δλ0 = 0,087 rad). 3rd International conference “Information Technology and Nanotechnology 2017” 40 Computer Modeling / G.P. Anshakov et al. Fig. 5. The eccentricity change of geostationary spacecraft orbit in simulation the orbit correction using low-thrust of EP (а0 = 0,001 м/s2, ΔТ0 = 1000 s, e0 = 0,005, Δλ0 = 0,087 rad). 4.3 Building a set of Pareto-optimal solutions of a dynamic problem Suggested control algorithms make it possible to build a set of Pareto-optimal solutions of the dynamic optimization problem of a space optical system on based diffractive membranes transfer to GEO using low thrust of EP with insertion into a given station. The sequence of building the set of Pareto-optimal solutions is represented in Table 2. Table 2. Building Pareto-optimal set № Control algorithms for insertion of powered EP spacecraft to GEO Final error margin 1 1. Control program (7) with EP thrust adjustment algorithm (8) and correction of control program without Φ1 precision GEO formation stage. 2 1. Control program (7) with EP thrust adjustment algorithm (8) and correction of control program. Φ2 2. A near-optimal control by extended set used on the final stage. 3 1. Optimal control program (7) with EP thrust adjustment algorithm (8) and correction of control program. Φ3 2. Control algorithm (12) used on the final stage. 4 1. Optimal control program (7) with EP thrust adjustment algorithm (8) and correction of control program. Φ4 2. Three-stage control algorithm of terminal control (13) – (15) used on the final stage. Figure 6 shows a sample calculation of a set of Pareto-optimal solutions for a transfer to GEO with systematic thrust error ΔP = - 2,5%. Fig. 6. Sample set of Pareto-optimal solutions (ΔP = - 2,5%, i0 = 51,60, r0 = 7171 km, c = 25 km/s). Conclusion A method of solving the dynamic optimization problem of a transfer to a given station on geostationary orbit was developed, including: an algorithm for obtaining nominal thrust control vector programs; algorithm for EP thrust magnitude adjustment on the basis of actual orbit time measurement at the long-range targeting stage; algorithm for obtaining terminal control, using discrete motion models; algorithm for obtaining a set of Pareto-optimal solutions. 3rd International conference “Information Technology and Nanotechnology 2017” 41 Computer Modeling / G.P. Anshakov et al. Acknowledgements This work is supported by the Ministry of Education and Science of the Russian Federation in the framework of The Federal purpose-oriented program "Research and development on priority directions of development of scientific-technological complex of Russia for 2014-2020" (agreement № 14.578.21.0229, unique identificator of the project RFMEFI57817X0229). References [1] Early J, Hyde R, Baron R. Twenty meter space telescope based on diffractive Fresnel lens. Proceedings of SPIE - The International Society for Optical Engineering 2004; 5166: 148–156. [2] Atcheson P, Stewart C, Domber J, Whiteaker K, Cole J, Spuhler P, Seltzer A, Smith L. MOIRE - Initial demonstration of a transmis-sive diffractive membrane optic for large lightweight optical telescopes. Proceedings of SPIE - The International Society for Optical Engineering 2012; 8442: 844221. [3] Atcheson P, Domber J, Whiteaker K, Britten JA, Dixit SN, Farmer B. MOIRE - Ground demonstration of a large aperture diffractive transmissive telescope. Proceedings of SPIE – The International Society for Optical Engineering 2014; 9143: 91431W. [4] Salmin VV, Karpeev SV, Peresypkin KV, Chetverikov AS, Tkachenko IS. Feasibility study and modeling of com-ponents for an informational space system based on a large diffractive membrane. CEUR Workshop Proceedings 2016; 1638: 132–148. [5] Grodzovsky GL, Ivanov YH, Tokarev VV. Space Flight Mechanics (Optimization Problems). Moscow: Nauka, 1975; 704 p. [in Russian] [6] Salmin VV, Ishkov SA. Trajectory and Parameter Optimization of Inter-Orbital Transfer Vehicles with Low-Thrust Propulsion Systems. Kosmicheskie Issledovaniya 1989; 27(1): 42–53. [in Russian] [7] Gurman VI. Extension Principle in Control Problems. Moscow: Nauka, 1985; 284 p. [in Russian] [8] Lebedev VN. Calculations of Low Thrust Spacecraft Motion. Moscow: CS AS USSR Press, 1968; 108 p. [in Russian] [9] Chernyavsky GM, Bartenev BA, Malyshev VA. Controlling the orbit of a geostationary satellite. Moscow: Mashinostroyenie, 1984; 144 p. [in Russian] [10] Salmin VV, Chetverikov AS. Selecting control laws for trajectorial and angular motion of a nuclear powered electric propulsion spacecraft during non-coplanar inter-orbital transfers. Bulletin of Samara Science Center of RAS 2013; 15(6). [in Russian] [11] Salmin VV, Chetverikov AS. Controling the flat orbit parameters of a geostationary spacecraft with low thrust propulsion. Bulletin of Samara State Aerospace University 2015; 14(4): 92–101. [in Russian] 3rd International conference “Information Technology and Nanotechnology 2017” 42