Design concept for In-Drag Sail with individually controllable elements

Aerospace Science and Technology 89 (2019) 382–391

Contents lists available at ScienceDirect

Aerospace Science and Technology www.elsevier.com/locate/aescte

Design concept for In-Drag Sail with individually controllable elements Qingyu Qu, Ming Xu ∗ , Tong Luo School of Astronautics, Beihang University, Beijing 100191, China

a r t i c l e

i n f o

Article history: Received 28 June 2018 Received in revised form 3 February 2019 Accepted 9 April 2019 Available online 12 April 2019 Keywords: In-drag sail Pneumatic expansion Decentralized layout Attitude control De-orbiting

a b s t r a c t The lifetime of the existing attitude and orbit control systems (AOCS) generally depends on the amount of the fuel carried on a spacecraft. This paper proposes a novel design concept for a low Earth orbit (LEO) spacecraft named In-Drag Sail with Individually Controllable Elements (IDSICE), which can release the fuel’s limitation to its attitude control system. IDSICE is designed to have four different working modes, which are Compression and Storage mode (CS mode), Orbital Flight mode (OF mode), Attitude Control mode (AC mode), and De-Orbiting mode (DO mode). By switching between different modes, it can achieve different missions. When IDSICE is in AC mode, the aerodynamic torques can be applied to completed the attitude control phase. When IDSICE is in DO mode, its area-to-mass ratio is enlarged so that the orbit altitude can decrease rapidly to achieve the de-orbiting phase at the end of the life. The decentralized layout is applied in IDSICE and most components have redundant backups, which allows a higher fault tolerance compared with the conventional layout. Finally, some necessary simulations are made, which verified the feasibility of the attitude control and de-orbiting methods proposed based on IDSICE. © 2019 Elsevier Masson SAS. All rights reserved.

1. Introduction In general, a spacecraft’s attitude has to be controlled and stabilized for a variety of reasons. Attitude control can be implemented by many methods, such as thrusters, spin stabilization, momentum wheels, control moment gyros, solar sails, gravity-gradient stabilization, and magnetic torquers. Thrusters are the most common actuators and the main limitation is the fuel consumption. Spin stabilization is mainly applied in the missions where the orientation of a spacecraft’s primary axis does not need to change dramatically in the entire lifetime. Momentum wheels are rotors driven by motors, which are made to spin opposite to the required direction. Control moment gyros are rotors spinning at constant speed, which are mounted on gimbals to provide attitude control. Because of a quite large mass and cost, control moment gyros are usually applied in the large spacecraft. Solar sails can be used to achieve attitude control in the long-time missions so that large amounts of fuel can be saved. Gravity-gradient stabilization is usually used to stabilize the attitude of a LEO spacecraft with one axis much longer than the other two axes, making the long axis point to the Earth centroid. Magnetic torquers work only where there is a magnetic field against which to react, and only two-axis control is available at any given time with this method. Among these

*

Corresponding author. E-mail address: [email protected] (M. Xu).

https://doi.org/10.1016/j.ast.2019.04.016 1270-9638/© 2019 Elsevier Masson SAS. All rights reserved.

attitude control methods, solar sails attracted increasing interests of scholars owing to its negligible fuel consumption. Besides, it was found that the solar pressure and the drag pressure are approximately equal at an altitude of about 800 km. Therefore, some scholars proposed the concept of in-drag sails which are similar to the traditional solar sails, orbiting in the LEO below an altitude of 800 km [1]. The first solar sail IKAROS inspired engineers to study a novel attitude control method. In the IKAROS mission, the attitude control was accomplished by varying the reflectance of 80 liquid crystal panels embedded in the sail [2]. Similarly, Luo et al. proposed another attitude control method by varying the angles between the normal vector of the sails and the direction of solar radiation based on a new configuration of the solar sail with individually controllable elements [3]. Based on Luo’s works, Yao et al. constructed an improved imperfectly reflecting sail model to generate the trajectory around an artificial Lagrangian point for idea solar sails [4]. Steyn and Lappas analyzed a CubeSat mission with a deployable solar sail of 5 m × 5 m in a sun-synchronous LEO to demonstrate solar sailing using active attitude stabilization of the sail panel [5]. Fu and Eke presented further investigation of the idea of performing orientation change and attitude control of solar sail spacecraft by moving one or more of the attachment points of the sail material [6]. In Liu’s paper [7], the pitch dynamics is addressed by considering the torques by the center-of-mass and center-of-pressure offset, the gravity gradient, the internal damping and the control vane. By checking the Poincare surface of section and the power

Q. Qu et al. / Aerospace Science and Technology 89 (2019) 382–391

spectral density, the validity of the Melnikov method is numerically verified. Usually, the solar sails are designed to be an expansion mechanism. Rong et al. [8] built the finite element model of solar sail support tube for the issue of solar sail spacecraft support tube inflatable deployment. And they proved that solar sail can be deployed smoothly and effectively as long as appropriate restraint controls were applied to the support tube. Some other papers of Liu also obtained some useful conclusions. Lin et al. [9] utilized moment of momentum theorem to deduce rigid-flexible coupling dynamics equations including the attitude dynamics and vibration equations for super-flexible solar sail spacecraft with control vanes, and studied the yaw axis earth-pointing problems of solar sail with Bang-Bang control. Considering that the aerodynamic drag cannot be ignored in LEO, the attitude control for the LEO spacecraft could be accomplished by changing the aerodynamic force and torque. Leonard et al. proposed an orbital formation keeping method with differential drag and examined its feasibility [1]. Kumar et al. investigated an orbital formation keeping controller that combines the efficacies of linear control and time optimal bang-bang control methodologies [10]. Bevilacqua et al. proposed a novel two-phase hybrid controller to optimize propellant consumption during multiple spacecraft rendezvous maneuvers in LEO, and this controller exploited generated differentials in aerodynamic drag on each involved chaser spacecraft to effect a propellant-free trajectory near to the target spacecraft [11]. Varma et al. developed a control methodology based on sliding mode control to maintain formation flying in LEO with differential aerodynamic drag, and analyzed the stability of such a formation in the presence of external perturbations [12]. Psiaki designed a novel configuration for the LEO spacecraft whose aerodynamic design resembles a badminton shuttlecock, and used aerodynamic and magnetic torques to develop a new three-axis attitude stabilization system [13]. In this paper, a new configuration for the LEO spacecraft named In-Drag Sail with Individually Controllable Elements (IDSICE) is proposed. Based on this configuration, a novel attitude control method is given by applying aerodynamic force and torque. Currently, there are more and more potentially hazardous population of objects in Earth orbits. Thus, space debris has been considered as a serious problem for operational space missions. In the past decades, a large number of space debris capturing and removal methods have been proposed, especially aimed at the space debris in LEO [14]. The existing space debris in LEO is mainly caused by the antisatellite weapons tests and accidental collisions, and the major source of new debris could be in-orbit collisions, between both intact spacecraft and existing debris [15]. The mitigation guidelines for LEO implies a disposal by re-entry in atmosphere [14]. It is extremely effective for de-orbiting in LEO by means of aerodynamic drag, especially at an altitude below 1000 kilometers. The spacecraft in this region will decay naturally, and the decay time depends on the orbit altitude and the area-to-mass ratio. Some scholars have been concentrating on the design of deorbiting in LEO by drag augmentation. Steyn et al. described the development and simulation results of a full attitude and orbit control system (AOCS) for the DeorbitSail mission, and they proved that the expected de-orbiting time would be less than 50 days by maximizing the drag force [16]. Yoo et al. presented the design results of the attitude control and de-orbiting process for CNUSAIL-1, and the feasibility and performance of the design are verified in high-fidelity nonlinear simulations [17]. Visagie et al. presented a model to evaluate the collision risk of a decaying object under drag conditions, and showed that using drag augmentation from a deployable drag sail to de-orbit a satellite in LEO will lead to a reduction in collision risk [18]. Lappas et al. demonstrated the concept of solar sailing and end-of-life de-orbiting using the sail membrane as a drag-sail in CubeSail mission [19]. Alexandra and

383

David proposed the passively stable pyramid sail and provided its overview by describing the design of the ESPA-class system, followed by discussion of the prototype system for the aerodynamic deorbit experiment [20]. In this paper, a de-orbiting method is proposed based on the concept of in-drag sail. By varying the status of the sails, the amplitude of drag force can be increase rapidly so that IDSICE will decay in a limited time. In this paper, IDSICE is proposed to improve the existing AOCS of LEO spacecraft. There are mainly four advantages compared with the existing spacecraft. Firstly, IDSICE is designed to have four different working modes. By switching between different modes, it can achieve different missions. Therefore, it allows the CubeSat to adapt to more working scenarios. Secondly, when IDSICE is in Attitude Control mode, the aerodynamic torques is used for attitude control. Compared with the spacecraft with conventional attitude control system, the fuel consumption in attitude control will be less because the aerodynamic torques can serve as the control torques and the individually controllable elements require a tiny amount of electrical energy, instead of fuel. Thirdly, when IDSICE is in De-Orbiting mode, its area-to-mass ratio is enlarged so that the orbit altitude can decrease rapidly to achieve the de-orbiting process at the end of the life, which is quite different from the existing CubeSat. Thus, it is supposed to be a feasible scheme to install this kind of de-orbiting device in order to reduce the space debris effectively. Fourthly, the decentralized layout is applied in IDSICE. Most components have redundant backups, and it allows a higher fault tolerance compared with the conventional layout. To verify the feasibilities of the attitude control and de-orbiting methods, some simulations are carried out in the end of this paper. 2. Atmospheric pressure model In general, continuum mechanics is applied to describe physical properties of fluids, and it is highly accurate when the characteristic length of object are much greater than that of inter-atomic distance. However, when the characteristic length of object is approximately equal to molecular mean free path, the shear and normal stress cannot be expressed by macroscopic quantities, which means that continuum mechanics fails. Therefore, rarefied gas dynamics is supposed to be applied to deal with this kind of problem. Free molecular flow describes the fluid dynamics of gas where the molecular mean free path is larger than the characteristic length, and it is a kind of idealized model in rarefied gas dynamics [21]. Considering that the molecular mean free path is about 1 × 106 m which is far larger than the size of spacecraft when the orbit altitude is about 300–1000 km, the atmospheric pressure including shear and normal pressure can be obtained by applying the free molecular flow model. Based on momentum theorem, the aerodynamic force can be computed by analyzing the momentum of free molecular flow. Let Q i = Q i τ τ + Q in n represent the momentum of incident flow in one unit of area, and Q r = Q r τ τ + Q rn n represent the momentum of reflex flow in one unit of area, where τ and n represent the tangential and normal direction respectively. In order to describe the effect of spacecraft surface on momentum, two parameters are defined as

⎧ Q in − Q rn ⎪ ⎪ ⎨ σn = Q in − Q w Q ⎪ iτ − Q rτ ⎪ ⎩ στ =

(1)

Q iτ

where Q w is the momentum caused by the molecular thermal motion, which is related to the temperature of the spacecraft surface.

384

Q. Qu et al. / Aerospace Science and Technology 89 (2019) 382–391

Thus, the shear pressure p τ and normal pressure pn can be yielded that



pn = (2 − σn ) Q in + σn Q w p τ = στ Q i τ

(2)

Based on the Maxwell–Boltzmann distribution, the number of the molecules whose velocity equals c is

f =

ρ m

3

1

(2π R T )− 2 e − 2R T |c −c T |

2

(3)

where ρ is the atmospheric density, m is the mass of a molecule, R is √ the gas constant, T is the atmospheric temperature, and c T = 2R T . Then the incident velocity v of the molecules whose velocity equals c can be obtained, as well as the corresponding maximum probable velocity u, shown as



v = c − v0 u = cT − v0

(4)

Then, it can be solved from Eq. (4) that

c − cT = v − u

(5)

By substituting Eq. (5) into Eq. (3), the number of the molecules whose velocity equals c can be deduced as

f =

ρ m

3

1

(2π R T )− 2 e − 2R T | v −u|

2

(6)

The sheer component u τ and the sheer component un of u are formulized as



un = |u | sin α u τ = |u | cos α

(7)

where α is angle between the spacecraft surface and the direction of incident flow. Then, it is derived as

| v − u |2 = ( v n − u sin α )2 + ( v τ − u cos α )2 + v b2

+∞ mv n f vdΩ

(9)

−∞

which can be expanded as:

⎧  +∞ ⎪ Q = mv n2 f dΩ ⎪ ⎨ in −∞ +∞ Q i τ = −∞ mv n v τ f dΩ ⎪ ⎪ ⎩ Q =  +∞ mv 2 f dΩ w n w −∞

(10)

By substituting Eq. (10) into Eq. (2), the shear pressure p τ and normal pressure pn can be obtained as

⎧

⎫ ⎪ σn T w ⎪ ⎪ ⎪ −s2 sin2 α 2 − σn ⎪ ⎪ √ s sin α + ⎪e ⎪ ⎪ ⎪ ⎪ ⎪ 2 T π ⎪ ⎪

⎪ ⎪ ⎪ ⎪ sin α ) ⎨ +1 + erf (s ⎬  1 2 2 + s sin α × (2 − σn ) ⎪ 2s ⎪ ⎪ ⎪ ⎪ ⎪ (11) 2 ⎪ ⎪

⎪ ⎪  ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ σ π Tw n ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ (s sin α ) ⎩+ ⎭ ⎪ ⎪ 2  T ⎪ ⎪  ⎪ 2 ⎪

2 ρu 1 2 ⎪ ⎪ ⎩p τ = στ cos α √ e−s sin α + s sin α 1 + er f (s sin α ) ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ρ u2 ⎪ ⎪ ⎨pn = 2

2s

π

where T w is the temperature of the spacecraft surface, and x 2 er f (x) = √2 0 e − y dy. π

3. Operating modes of IDSICE To accomplish different working requirements during the entire lifetime, IDSICE is designed to have four different working modes, including CS, OF, AC and DO modes. 3.1. Compression and storage mode

(8)

where v n is the normal component of v , v τ is the sheer component of v, v b is the component of v along b = n × τ . Let dΩ denote a finite volume of elements, then we can represent f dΩ as the number of the molecules in dΩ . Thus, the number of the molecules in per unit time and per unit square is v n f dΩ , and the corresponding momentum is mv n f vdΩ . Therefore, the total momentum is

Qi=

Fig. 1. Structure of IDSICE in CS mode. (For interpretation of the colors in the figure(s), the reader is referred to the web version of this article.)

IDSICE is mainly constituted by a CubeSat, sails, individually controllable elements and an inflatable frame. In this paper, the CubeSat is designed as a 1U one, and the number of sails is determined as four, which can be changed to adapt to different working conditions. In fact, for the attitude control phase, there are three degrees of freedom. Therefore, three independent control variables should be introduced to ensure that the attitude control phase can be completed. In IDSICE, each individually controllable element can provide one control variable, and there are two individually controllable elements on each sail. Thus, the minimum number of sails is two theoretically. And there is no limit on the maximum number of sails. In order to make use of the rocket space as effective as possible, IDSICE in CS mode is designed as a folding structure. As shown in Fig. 1, component A (the white cube) represents the 1U CubeSat whose size is 10 cm × 10 cm × 10 cm, and all the payloads are put into it. Component B (the yellow board) represents the one piece of sail. When it is expanded (in OF/AC/DO mode), the size is determined as 2 m × 8 m. The mass and thickness can be ignored compared with other components of IDSICE. While it is folded (in CS mode), it is put on the side face of the CubeSat. Component C (the blue cube) represents individually controllable elements whose main functional part is a stepper motor, and its size is variable depending on different motor (for example, an OT-16GA030PA-115 stepper motor has a size of 29 mm × 15 mm × 12 mm). Component D (the red frame) represents the inflatable frame, and it is folded in CS mode. The expanding process of the inflatable frame is shown in Fig. 2. In general, an inflatable frame will be expanded to a circle one when there is no constrain. But in fact, we put four rigid constrains at the corners of the frame, which are shown to be the green components in Fig. 2. After ap-

Q. Qu et al. / Aerospace Science and Technology 89 (2019) 382–391

385

Fig. 2. Expanding process of the inflatable frame.

Fig. 3. Sail setting of IDSICE in OF mode.

plying these rigid constrains, the inflatable frame will maintain the square shape.

Fig. 4. Configuration of IDSICE in OF mode.

3.2. Orbital flight mode IDSICE in OF mode is supposed to work as a normal spacecraft. Therefore, the status of four sails is adjusted to guarantee that the aerodynamic force and torque can be ignored. When the angle between the sail surface and the incident flow is set to 0◦ , as shown in Fig. 3, the aerodynamic force and torque are approximately equal to zero. In this way, the IDSICE in OF mode has a configuration as shown in Fig. 4. 3.3. Attitude control mode In AC mode, the angle between the sail surface and the incident flow can be changed to obtained the required torques for attitude control. A body coordinate system S b is defined as follows. The origin is located at the center of IDSICE; x axis is along with normal direction of the frame plane, pointing to the front direction; y axis is located in the frame plane, pointing to the right direction; and z axis is determined by the right-hand rule. (See Fig. 5.) An orbital coordinate system S o is defined as follows. The origin is located at the center of IDSICE; x axis is located in the orbital plane, pointing to the front direction; z axis is located in the orbital plane, pointing to the center of Earth; and y axis is along with normal direction of the orbital plane, determined by the right-hand rule. S o can be rotated into S b by three basic coordinate rotations, which is shown as

Fig. 5. The body coordinate system of IDSICE.

R x (ϕ )

R y (θ )

R z (ψ)

S O −→ ◦ −→ ◦ −→ S b where ϕ , θ , and ψ are roll, pitch and yaw angles respectively. Then the coordinate transformation matrix is

 cos ψ cos θ

L b O = − sin ψ cos θ sin θ



sin ψ cos ϕ + cos ψ sin θ sin ϕ sin ψ sin ϕ − cos ψ sin θ cos ϕ cos ψ cos ϕ − sin ψ sin θ sin ϕ cos ψ sin ϕ + sin ψ sin θ cos ϕ − cos θ sin ϕ cos θ cos ϕ

(12)

386

Q. Qu et al. / Aerospace Science and Technology 89 (2019) 382–391

Thus the aerodynamic force at the point z is

⎡

⎛

d F bi

⎤

⎡

⎤⎞

cos δzi − sin δzi = ⎝− pn ⎣ sin δzi ⎦ + p τ ⎣ cos δzi ⎦⎠ wdz, 0 0

i = 1, . . . , n

(16)

The subscript b represents the body coordinate system, and the subscript i represents the label of different piece of sail. Then the aerodynamic force on one piece of sail can be obtained by an integration, which is

⎧  2l ⎪ ⎪ ⎪ F xbi = − l (− pn cos δzi − p τ sin δzi ) wdz ⎪ 2 ⎨  2l (− pn sin δzi + p τ cos δzi ) wdz F ybi = ⎪ ⎪ − 2l ⎪ ⎪ ⎩ F zbi = 0

(17)

Fig. 6. The local coordinate system defined on one piece of sail.

Supposing that the mass of CubeSat is m1 , the mass of each individually controllable element is m2 , the mass of each sail is m3 (in the simulation, m3 is set as zero), the mass of the inflatable frame is m4 , the length and width of each sail are noted as l and w respectively, the length of the CubeSat is a, and the number of sails is n. Then, the inertial matrix of IDSICE can be obtained as

⎡

J =⎣

⎤

Jx

⎦

Jy

⎧  n/2 ⎪

2 nm2l2 ⎪ m1 a2  m3 ⎪ ⎪ Jx = + m2 + a + (2k − 1) w + ⎪ ⎪ 6 2 2 ⎪ ⎪ k =1 ⎪ ⎪ ⎪  

nm m ⎪ 3 4 ⎪ ⎪ + l2 + w 2 + l2 + (a + nw )2 ⎪ ⎪ 12 12 ⎪ ⎪ ⎨ m1 a2 nm2l2 nm3l2 m4 l 2 Jy = + + + ⎪ 6 2 12 12 ⎪ ⎪ ⎪  n/2 ⎪ 2  ⎪

2 m a m 1 3 ⎪ ⎪ + m2 + a + (2k − 1) w Jz = ⎪ ⎪ ⎪ 6 2 ⎪ k =1 ⎪ ⎪ ⎪ ⎪ 2 ⎪ nm3 w m4 ⎪ ⎩ + + (a + nw )2 12

2

i =1

12

(δ1 − δ2 ) z l

l

l

2

2

− ≤z≤

(14)

Thus the normal vector of the sail at the point z can be expressed as {n}b = [cos δzi sin δzi 0] T in S b . Supposing the unit vector along with the incident flow is r 1 , and it is clear that {r 1 } O = [−1 0 0] T . Therefore, the angle between the normal direction of sail and the incident flow can be yielded that



η = arccos

−2

⎧  n T n   ⎪ ⎪ ⎪ F xbi F ybi 0 ⎪ F = ⎪ ⎨ b i =1 i =1  n n n ⎪    ⎪ ⎪ ⎪ ⎪ M M M zbi M = xbi ybi ⎩ b

As shown in Fig. 6, a one-dimensional local coordinate system is defined on one piece of sail. The origin is located at the center of the given piece of sail; the axis is along with longitudinal direction of the sail. Supposing the rotation angle of the motor located on the top of the sail is δ1 , while that of the motor located at the bottom is δ2 . A linear torsional deformation is performed on the sail so that the angle between the sail surface and the incident flow is no longer consistent along the axis of sail. Therefore, the rotation angle of the sail at the point z is

+

(18)

Then the aerodynamic force and torque on the entire spacecraft can be obtained as

where

δzi =

⎧ l ⎪ ⎪ M xbi = 2 l (− pn sin δzi + p τ cos δzi ) w zdz ⎪ ⎪ −2 ⎪ ⎪ ⎨  2l (− pn cos δzi − p τ sin δzi ) w zdz M ybi = − 2l ⎪ ⎪ ⎪ ⎪ l ⎪ ⎪ ⎩ M zbi = 2 l (− pn cos δzi − p τ sin δzi ) w ydz

(13)

Jz

δ1 + δ2

Similarly, the aerodynamic torque on one piece of sail can be obtained as well, which is

{n}b · {r 1 }b {n}b  · {r 1 }b 

(19)

i =1

Besides, the gravity gradient torque should be considered as well. Let dm denote a finite mass of elements. The distance vector from the centroid of the spacecraft to dm is denoted by r, and {r }b = [xb yb zb ] T . Then the distance vector from the center of the Earth to dm is R = R 0 + r, where R 0 is the distance vector from the center of Earth to the centroid of the spacecraft. Then the gravity gradient torque with respect to the centroid of the spacecraft in S b is

Mg = −

3μ R 30

 R0 ×

r

R0 · r

m

R 20

dm

(20)

Based on the definition of S o , the unit vector along with the z axis of S o is

ko = −

Ro

(21)

Ro

According to the definition of the inertia tensor, which is

!

I=

  m



r 2 U − r · r dm

ko × U · ko = 0

(22)

The gravity gradient torque can be adapted further, which is



(15)

i =1

T

Mg =

3μ r3

ko × ( I · ko )

(23)

Q. Qu et al. / Aerospace Science and Technology 89 (2019) 382–391

387

Thus the gravity gradient torque in S b is

⎧ 3μ ⎪ ⎪ M gxb = 3 ( J z − J y ) ⎪ ⎪ ⎪ R0 ⎪ ⎪ ⎪ ⎪ ⎪ × (cos ψ sin ϕ + sin ψ sin θ cos ϕ ) cos θ cos ϕ ⎪ ⎪ ⎪ ⎪ 3 μ ⎪ ⎨ M g yb = ( J x − J z) 3 R

(24)

0 ⎪ ⎪ ⎪ × (sin ψ sin ϕ − cos ψ sin θ cos ϕ ) cos θ cos ϕ ⎪ ⎪ ⎪ ⎪ 3μ ⎪ ⎪ ⎪ M gzb = 3 ( J y − J x )(sin ψ sin ϕ − cos ψ sin θ cos ϕ ) ⎪ ⎪ R0 ⎪ ⎪ ⎪ ⎩ × (cos ψ sin ϕ + sin ψ sin θ cos ϕ )

Based on the small-angle approximation, the gravity gradient torque is approximately as following:

⎧ 3μ ⎪ ⎪ M gxb = 3 ( J z − J y )ϕ ⎪ ⎪ R0 ⎪ ⎨ 3μ

(25)

M g yb = 3 ( J z − J x )θ ⎪ ⎪ ⎪ R0 ⎪ ⎪ ⎩ M gzb = 0

Fig. 7. Configuration of IDSICE in DO mode.

A spacecraft’s absolute angular velocity is a vector sum of its orbital coordinate system’s angular velocity and the angular velocity relative to orbital coordinate system, namely

ω = ωorb + ωr

(26)

˙ + θ˙ + ψ˙ , ωorb = where ωr = ϕ H represents the spacecraft’s moment of momentum, and r represents the orbital radius. Based on the small-angle approximation, the attitude dynamical equation can be described as H , r2

A ζ¨ + B ζ˙ + C ζ = M b

(27)

where ζ = [ϕ θ ψ] , and the coefficient matrixes are given by T

⎡

A=⎣

⎡ C =⎣

⎦

Jy Jz

⎡ B =⎣

⎤

Jx

⎤

0 0

( Jx − J y +

J z ) rH2

0 −( J x − J y + J z ) rH2 ⎦ 0 0 0 0

( J y − J z )( rH2 )2 − 3μ3 ( J z − J y )

0

0

0

− 3μ3 ( J z − J x )

0

R0

0

R0

0

( Jx −

Fig. 8. History of the coefficient of drag C d with respect to altitude.

⎤

where ρ is the atmospheric density, A is the area presented to the incident flow, v is the velocity of spacecraft, and C d is the coefficient of drag, which is given by

⎦

J y )( rH2 )2

Therefore, the attitude control can be accomplished by solving Eq. (19) and Eq. (27).

Cd =

F drag =

2

ρ AC d v



2 − σn

√

s2

π

s+

σn 2



Tw



2

2

T

Then, the variation of C d with respect to altitude can be shown in Fig. 8. It can be seen that the variation range of C d is quite small. Considering that the IDSICE’s altitude is about 300–800 km, the value of C d is set to 2.2. 4. Numerical simulations 4.1. Attitude control numerical simulation In this section, a simulation is achieved to verify the feasibility of the attitude control method based on IDSICE. Eq. (27) can be modified into

d (28)



T

s

A thin atmosphere will have an influence on the orbits of LEO spacecraft, which can cause them to spiral towards Earth. Some scholars have proved that a larger area-to-mass ratio will result in a shorter decay time, which means that the expanded configuration can be applied in the de-orbiting mission. As shown in Fig. 7, the area-to-mass ratio of IDSICE in DO mode turns to be quite large, which can lead to an enormous aerodynamic force. And the direction of the disturbance force is exactly opposite to the direction of spacecraft’s velocity. The expression of the corresponding aerodynamic force can be derived based on Eq. (11), which is 2

2

  1 + er f (s) 1 σn s π T w 2 ( + + 2 − σ ) + s n 2

3.4. De-orbiting mode

1

e −s

dt



ζ ζ˙





=

O − A −1 C

I − A −1 B



ζ ζ˙





+

0 A −1 M b



(29)

388

Q. Qu et al. / Aerospace Science and Technology 89 (2019) 382–391

Table 1 Simulation parameters in AC mode. Parameters

Values

Moment of inertia (kg m2 ) Initial Euler angles (◦ ) Initial Euler angle velocities (◦ /s) Orbital altitude (km)

ϕ=θ =ψ =1 ϕ˙ = θ˙ = ψ˙ = 0.001

J x = 1.93; J y = 1.55; J z = 0.39

500

A linear-quadratic regulator can be employed" to stabilize the # atO I titude. Let X e = [ζ ζ˙ ] T , U = [ O M ] T , D = − A −1 C − A −1 B , E =

"

O A −1

J=

#

. Then the quadratic cost function is defined as

1

t f



2

X eT (t ) Q (t ) X e (t ) + U T (t ) R (t )U (t ) dt

(30)

t0

where Q (t ) is a nonnegative weighted matrix, and R (t ) is a positive weighted matrix. In this simulation, they are chosen as

⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ Q =⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ R =⎣

180

⎤

2

π



180

2

π



180

2

π 0 0

1 × 10

⎤

12

1 × 1012

Fig. 9. Time histories of roll angle, pitch angle and yaw angle.

⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥, ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ 0

⎦ 1 × 1012

The feedback control law that minimizes the value of the cost function J is

U (t ) = − G (t ) · X e (t )

(31)

where G (t ) is given by

G (t ) = R −1 (t ) · E T · K (t )

(32)

And K(t) is found by solving the continuous time Riccati differential equation, i.e.

Fig. 10. Time history of control torques.

˙ (t ) = − K (t ) · D − D T · K (t ) + K (t ) · E · R −1 (t ) · E T · K (t ) K − Q (t )

(33)

with the boundary condition K (t f ) = 0. When the simulation parameters are set to the values in Table 1, we can obtain following simulation results. The variation of the roll angle, pitch angle and yaw angle is shown in Fig. 9. It can be shown that the Euler angles will converge to zero within 1500 s, which proves that the attitude control method proposed in this paper is effective. The variation of control torques is shown in Fig. 10. The variation of each individually controllable element’s rotation angle is shown in Fig. 11 (the labels of each individually controllable element in Fig. 11 are corresponding one by one to those in Fig. 7). It has been mentioned that the fuel consumption in attitude control will be less compared with the spacecraft with conventional attitude control system. In this section, the motion in yaw direction is taken as an example to show the comparison in fuel consumption. For the proposed method in this paper, the aerodynamic torques serve as the control torques and only some of

Table 2 The related technical parameters of the chosen steeper motor. Parameters

Values

Size (mm) Rated power (W) Rated revolution per minute Rated torque (mN m)

29 × 15 × 12 0.36 34 23.52

the individually controllable elements need to be controlled. The related technical parameters of the chosen stepper motor in the individually controllable elements is shown in Table 2. When the yaw angle is required to turn from 1◦ to 0◦ , the energy consumption can be calculated to be 0.6353 J, which is tiny enough to be covered by the energy from solar panels. And it does not require additional fuel consumption. However, when the conventional attitude control methods are adopted, it is necessary to require some additional fuel consumption.

Q. Qu et al. / Aerospace Science and Technology 89 (2019) 382–391

389

Fig. 11. Time histories of each individually controllable element’s rotation angle. Table 3 Simulation parameters in DO mode. Parameters

Values

Mass (kg) Effective area (m2 ) Apogee altitude (km) Perigee altitude (km) Inclination (◦ )

2 64 685 330 9

4.2. De-orbiting simulations In this section, a simulation is achieved to verify the feasibility of the de-orbiting method based on IDSICE. The orbital elements are determined to be the same as those of NanoSail-D, and the specific simulation parameters are shown in Table 3. The simulation result is shown in Fig. 12, which describes the variation of the semi-major axis. It can be concluded that the spacecraft will crash within 3 hours, which proves that the de-orbiting method is feasible. For the de-orbiting phase, the area-to-mass ratio of IDSICE is enlarged so that the orbital altitude can decrease rapidly to achieve the de-orbiting process at the end of the life, which is quite different from the existing CubeSat. Similar de-orbiting method is adapted by Willem et al. They [22] applied the deployed solar panels to the deOrbitSail FP7 mission. Considering that the deployed solar panels are rigid, the area-to-mass ratio cannot turn large enough. Thus the de-orbiting time is between 30 to 50 days when the initial orbit altitude is 650 km. However, flexible sails are adopted in the proposed design, which allows that the area-tomass ratio can be 32. And the spacecraft can crash within 3 hours. 5. Engineering tips In engineering applications, lots of practical issues need to be considered. In this part, we focus on the number of moles of the gas in the frame and the scale of the sails.

Fig. 12. Time history of semi-major axis of the orbit.

5.1. The number of moles of the gas in the frame A symmetric bending occurs to the two horizontal frames in Fig. 13a due to the atmospheric pressure. The maximum deflection is supposed to appear at the edge of the horizontal frame, and the maximum stress is supposed to appear at cross section A in Fig. 13a which is located at the middle. In order to maintain the shape of the expanded frame, it is essential to check and estimate whether the maximum deflection exceeds the allowable deflection value, which can be set as ε of the frame’s side length in this paper. In order to maintain the strength of the frame, it is necessary to check and estimate whether the maximum stress exceeds the allowable stress value [σ ], which is determined by the number of moles of the gas in the frame. As shown in Fig. 13b, the cross section of the frame is designed as a ring, where the thickness of the frame is represented by d and

390

Q. Qu et al. / Aerospace Science and Technology 89 (2019) 382–391

Fig. 13. The expanded frame and its cross section.

the inner radius of the frame is represented by r. Based on the beam bending theory, the maximum deflection can be expressed by

δmax =

5ql4

(34)

384( E 1 I 1 + E 2 I 2 )

where q is the uniform load caused by atmospheric pressure on p l the horizontal frame which is given by q = 2n , l is the side length of the expanded frame, E 1 is the elastic modulus of the gas in frame which is given by E 1 = γ nRV T , I 1 is the moment of inertia 4

of the inner circle which is given by I 1 = π4r , E 2 is the elastic modulus of the frame’s material which is determined by the material, I 2 is the moment of inertia of the outer ring which is given by I 2 = π4 [(r + d)4 − r 4 ], γ represents the specific heat capacity ratio of the gas, n represents the number of moles of the gas in frame, R represents the ideal gas constant, T represents the inner temperature of the frame which is supposed to be equal to the temperature of the gas, and V represents the capacity of the expanded frame for gas. Similarly, the maximum stress can be obtained based on the beam bending theory, which is expressed as

σmax =

ql2 (r + d) 8( E 1 I 1 + E 2 I 2 )

E2

(35)

The maximum deflection and the maximum stress must satisfy the following equations, which are



σmax ≤ [σ ] δmax ≤ ε · l

(36)

Thus the number of moles of the gas in frame can be obtained by solving Eq. (36) as

n≥

γV R T I1



5ql3 384ε

 − E2 I2

(37)

5.2. The scale of the sail of IDSICE One of the important performance indicators of the AOCS of IDSICE is the maximum aerodynamic torque that can be generated, and the amplitude of the maximum torque is depend on the scale of the sail. Therefore, by giving a proper scale of the sail, the required control torques can be obtained. The variation of the maximum torque’s amplitude with respect to the sail’s scale is shown

Fig. 14. The variation of the maximum torque’s amplitude with respect to the sail’s scale.

in Fig. 14. Besides, we also give the analytical expression by taking the amplitude of M z as an example. Clearly, when the two left sails are in the plane Oyz and the right two sails are perpendicular to the plane Oyz or vice versa, M z will reach the maximum value. And the expression of the maximum torque’s amplitude is

M max = wl(a + 2w )( pn − p τ )

(38)

where the expression of pn and p τ are given in Eq. (11). In this paper, the sail’s scale is set as 64 m2 . 6. Conclusion This paper proposed a design concept for IDSICE, which can improve the existing AOCS. In order to construct its attitude dynamic equation, we deduced atmospheric pressure model based on the rarefied gas dynamics in Section 2. IDSICE has four modes, which are CS mode, OF mode, AC mode, DO mode. And we introduced the details of the four modes in Section 3. Then, the attitude control numerical simulation and de-orbiting simulation were made in Section 4, which verified the feasibilities of the attitude control and de-orbiting methods proposed based on IDSICE. In Section 5,

Q. Qu et al. / Aerospace Science and Technology 89 (2019) 382–391

we computed the number of moles of the gas in the frame, which is necessary to maintain the strength of the inflatable frame; besides, we also obtained the variation of the maximum torque’s amplitude with respect to the sail’s scale, and we can choose proper sail’s scale according to the required torques. Conflict of interest statement None declared. Acknowledgements The research is supported by the National Natural Science Foundation of China (11772024 and 11432001), Aerospace Science Foundation by China Aerospace Science and Industry Corporation, the Fundamental Research Funds for the Central Universities (YWF-16-BJ-Y-10), and the Foundation of Key Laboratory of Spacecraft Design Optimization and Dynamic Simulation Technologies, Ministry of Education of China and the Innovation Practice Foundation of Beihang University for Postgraduate Students (YCSJ01-2018-08). References [1] C.L. Leonard, W.M. Hollister, E.V. Bergmann, Orbital formation keeping with differential drag, J. Guid. Control Dyn. 10 (10) (2012) 755–765. [2] Y. Tsuda, O. Mori, R. Funase, et al., Flight status of IKAROS deep space solar sail demonstrator, Acta Astronaut. 69 (9) (2011) 833–840. [3] T. Luo, M. Xu, Q. Qu, Design concept for a solar sail with individually controllable elements, J. Spacecr. Rockets 54 (27) (2017) 1–9. [4] C. Yao, M. Xu, T. Luo, Dynamics and control for nonideal solar sails around artificial Lagrangian points, J. Spacecr. Rockets 55 (3) (2018) 575–585. [5] W.H. Steyn, V. Lappas, CubeSat solar sail 3-axis stabilization using panel translation and magnetic torquing, Aerosp. Sci. Technol. 15 (6) (2011) 476–485. [6] B. Fu, F. Eke, Further investigation of the body torques on a square solar sail due to the displacement of the sail attachment points, Aerosp. Sci. Technol. 50 (2016) 281–294.

391

[7] J. Liu, L. Chen, N. Cui, Solar sail chaotic pitch dynamics and its control in earth orbits, Nonlinear Dyn. 90 (3) (2017) 1755–1770. [8] S. Rong, J. Liu, B. Zhao, N. Cui, Inflatable deployment dynamics simulation of solar sail spacecraft support tube, in: International Conference on Mechatronics and Automation, Hangchun, China, August 9–12, 2009, INSPEC-11177459. [9] J. Liu, S. Rong, J. Li, N. Cui, Attitude dynamics modeling and control of large flexible solar sail spacecraft, in: 3rd International Symposium on Systems and Control in Aeronautics and Astronautics, Harbin, China, June 8–10, 2010, INSPEC-11649611. [10] B.S. Kumar, A. Ng, A bang-bang control approach to maneuver spacecraft in a formation with differential drag, in: AIAA Guidance, Navigation and Control Conference and Exhibit, Honolulu, Hawaii, August 18-21, 2008, 2008-6469. [11] R. Bevilacqua, J.S. Hall, M. Romano, Multiple spacecraft rendezvous maneuvers by differential drag and low thrust engines, Celest. Mech. Dyn. Astron. 106 (1) (2010) 69–88. [12] S. Varma, K.D. Kumar, Multiple satellite formation flying using differential aerodynamic drag, J. Spacecr. Rockets 49 (2) (2015) 325–336. [13] M.L. Psiaki, Nanosatellite attitude stabilization using passive aerodynamics and active magnetic torquing, J. Guid. Control Dyn. 27 (3) (2004) 347–355. [14] Inter-Agency Space Debris Coordination Committee, IADC Space Debris Mitigation Guidelines, 2007. [15] J.C. Liou, N.L. Johnson, Instability of the present LEO satellite populations, Adv. Space Res. 41 (7) (2008) 1046–1053. [16] W.H. Steyn, H.W. Jordan, An active attitude control system for a drag sail satellite, Acta Astronaut. 128 (2016) 313–321. [17] Y. Yoo, S. Kim, J. Suk, J. Kim, Attitude control system design & verification for CNUSAIL-1 with solar/drag sail, Int. J. Aeronaut. Space Sci. 17 (4) (2016) 579–592. [18] L. Visagie, V. Lappas, S. Erb, Drag sails for space debris mitigation, Acta Astronaut. 109 (2015) 65–75. [19] V. Lappas, N. Adeli, L. Visagie, et al., CubeSail: a low cost cubesat based solar sail demonstration mission, Adv. Space Res. 48 (11) (2011) 1890–1901. [20] A.C. Long, D.A. Spencer, A scalable drag sail for the deorbit of small satellites, J. Small Satell. 7 (3) (2018) 773–788. [21] S.A. Schaaf, Recent progress in rarefied gas dynamics, J. Am. Rocket Soc. 30 (5) (2015) 443–447. [22] H.S. Willem, W.J. Hendrik, An active attitude control system for a drag sail satellite, Acta Astronaut. 128 (2016) 313–321.