1. Introduction Japan Aerospace Exploration Agency (JAXA) is currently planning the next generation Earthmagnetosphere observation mission called gSCOPEh (cross-Scale Coupling in Plasma universE, Fig.1)1. Its predecessor gGEOTAILh, launched in 1992, played a big role for understanding the ion behavior in the Earthfs magnetotail and triggered interests toward further microscopic scale (electrons scale) phenomena. Following this successful work of GEOTAIL, SCOPE aims at observing the Earthfs magnetotail where the ions and electrons interact with each other, with 5 satellites flying in formation. To fully resolve the time-domain behavior and spatial distribution of the magnetospheric phenomena, a simultaneous observation by spatially distributed electro-magnetic instruments is essential. From the formation design point of view, there are several requirements to be satisfied, which include not only the optimality of observation, but also constraints stemming from the operationality, the physical layout of the onboard components and so on. In the SCOPE mission, the fundamental requirements for formation design come from (1) observation optimality, (2) relative orbit determination accuracy, Fig.1: SCOPE Mission Image 2 and (3) inter-satellite communication link capability constrained by antennae layout of each satellite. The (2) requirement is needed because in the SCOPE mission, the inter-satellite radiometric ranging and clock synchronization via inter-satellite communications are employed. (3) is needed to comply with the actual antennae directivity limitations As the first step of designing the SCOPE formation, one wants to obtain all the possible answers which satisfy (1), (2) and (3) at the same time. Here what is needed is surveying a whole solution space, not focusing on the highly accurate local minimum, but rather on approximated answers which represent all the possible types of solutions. This paper proposes a formation design strategy which surveys a whole solution space. This method employs the simulated annealing (SA) algorithm to obtain the globally optimum solution. In this paper, the method is validated through analyzing the well-known optimum formation in low-earth orbit, and then applied to SCOPE formation design. Finally, based on the solution types obtained by SA, a heuristics to design the optimum SCOPE formation is introduced. 2. Formation Optimality Index Before constructing the design method, this section introduces several indices used in the following formation design. 2.1. Observation Optimality (FDOP-O) Consider the multi-point simultaneous observations by spatially distributed satellites fleet. The observation objective is to obtain the spatial differential coefficient of the observation field. The observation of i-th satellite, which is located at the relative position vector ( , , )T x y z = r can be written as; ( ) Y f = r (1) Suppose the field to be observed is linear, i.e.; ( ) , , T f f f f x y z ? ? =? ? ? ? r (2) then from the observation of i-th and j-th satellites, the difference of observations ij j i Z Y Y = ? can be expressed as; ( )T T ij j i ij Z f f = ? = r r r (3) Eq.(3) implies that the quantity f can be calculated from the observation differences. To obtain 3-dimensional f value, at least 3 different ij Z are required, which means at least 4 satellites are needed to obtain f . Practically f is determined by the least-square method as follows; ( )( ) 1 T ij ij ij ij f Z ? = r r r (4) where rij is the relative position vector from i-th to j-th satellite. If the field is supposed to be isotropic, the error covariance for equation (4) can be calculated from the following equation; ( )1 2 T ij ij ? = P rr (5) where 2 is the error variance of the observation. Therefore the estimation quality can be estimated from the following quantity; ( )1 Tr T ij ij ? ? ? ? ? ? ? r r (6) which represents the dilution of precision for the observation quantity. (6) implies that the well-known tetrahedron formation is optimum, and in addition, the larger formation is better. By eliminating the size effect from (6), we finally obtain the formation dilution of precision for observation as follows; ( )1 2 ij max(r )Tr T ij ij FDOPO ? ? ? = ? ? ? ? r r (7) Here the index is normalized by the maximum relative distance among all the two satellites combinations. 2.2. Formation Geometry Determination Optimality (FDOP-G) Consider the quality of formation geometry determination from the relative ranging information. The state vector to be estimated is the formation shape, 3 which is written as; [ ]T 1 2 2 3 3 3 1 1 1 n n n x x y x y z x y z ? ? ? ? = x ?? (8) where the 0th satellite is supposed to be at the origin, and the 1st satellite is supposed to be located on the x-axis, the 2nd satellite is on the xy-plane. Hence the dimension of estimation state vector becomes 3(n-2) for n-satellite system. From the relation between ranging output and relative position vector ij i j r ? = ? r r , the following observation equation is obtained; ij ij ij i j i j r r r ? ? ? = + r r r r (9) Eq.(9) for all the satellites combinations are combined to form the following matrix relation; _ _ = y H x (10) where 01 02 0 1 12 T 1 1 3 2 2 1 n n n n n n r r r r r r r ? ? ? ? ? ? ? = ? ? ? ? ? ? ? ? ? ? y ?? ?? ?? (11) is the observation vector. The formation geometry estimation is now translated to minimize ( )( ) T _ _ _ _ 1 2 J= ? ? y H x y H x hence the dilution of precision for geometry determination is derived as follows; ( )1 _ _ Tr /( 1) T FDOPG n ? ? ? = ? ? ? ? ? H H (12) In (12), division by n-1 is introduced to eliminate the satellite number dependency. 2.3. Formation Dynamics Determination Optimality (FDOP-D) Consider the quality of relative motion (or formation dynamics) determination from the relative ranging history information. The state vector to be estimated is the relative position vector (with reference to the inertial space), which can be written as; [ ]T 1 1 1 2 2 2 1 1 1 n n n x y z x y z x y z ? ? ? = x ?? (13) Different from (8), x includes the information not only about the shape of the formation, but also the orientation of the formation with reference to the inertial space. Therefore the dimension of estimation state vector becomes 3(n-1) for n-satellite system. Denoting the transition matrix of x as ( ) 2 1 , r t t , the observation equation is written as r = y H x . Then the formation dynamics estimation is translated to minimize ( ) ( ) 1 ( ) ( ) ( ,0) (0) 2 ( ) ( ) ( ,0) (0) T r r J t t t t t t dt = ? ? y H x y H x i ?? hence the dilution of precision for dynamics determination is derived as follows; Fig.2: Samples Distribution in E-P Plane Fig.3: Correlations Between FDOPO, FDOPG and E,P 4 1 T 1 Tr ( ,0) ( ) ( ) ( ,0) /( 1) T FDOPD t t t t dt n T ? ? ? ? ? = ? ? ? ? ? ? ? ? ? ? ? H H ?? (14) where T is the orbital period. Again division by n-1 in (14) is introduced to eliminate the satellite number dependency. 2.4. Basic Characteristics of the Indices The proposed indices FDOPO, FDOPG are compared to the indices used in ESAfs CLUSTER mission2. In CLUSTER, a generic evaluation of formation shape is done by abbreviating the formation shape to the pseudo-ellipsoid, and extracting the characteristic values gelongationh E and gplanarityh P. These values are derived in summary as follows. First, define the tensor R; T ij ij = R rr (15) The principle axes of the pseudo-ellipsoid are given by the eigenvectors R(n) of R. If we suppose (1) (2) (3) R R R ? ? (16) then, their square roots represent respectively the major, middle and minor semiaxis of the pseudo-ellipsoid; (1) (2) (3) , , a R b R c R = = = (17) The elongation E and planarity P are defined using (17) as follows; 1 ( / ) 1 ( / ) E b a P c b = ? = ? (18) E and P take the value between 0 and 1. The sphere corresponds to (E,P)=(0,0). For 4-satellite formation, (E,P)=(0,0) corresponds to the tetrahedron formation. To survey the relation between E,P and FDOPO, FDOPG, 1000 4-satellite formation samples are generated as is shown in Fig.2. Each satellite in each sample is generated randomly. The correlations between FDOPO, FDOPG and E,P are shown in Fig.3. Seen from these figures, the optimality of FDOPO and FDOPG corresponds to (E,P)=(0,0). It can also be said that the optimality of FDOPO also guarantees the FDOPG optimality. 3. Global Optimization Using Simulated Annealing Algorithm 3.1. Simulated Annealing Optimization To survey the whole solution space to find the optimum formation shapes, the simulated annealing (SA) optimization is employed. Some characteristics of SA are described below; ? SA updates the solution by perturbing the previous solution randomly. ? At first SA searches around the whole solution space, rather than insatiably pursuing the optimality. As the optimization step goes by, updates which increase the optimality becomes preferred more by analogy of a cooling process of steel gannealingh. ? SA is suitable for searching the discrete space, but it is also useful for obtaining the approximated optimum of the multimodal continuous space. ? The results of SA optimization are remained to be refined, by investigating the proximity of the SA solution in detail (by local minimum search). The overall algorithm of SA is illustrated in Fig.4. In the formation optimization, the SA process becomes as follows; (1) The initial solution, which is a set of relative state vectors xi=(xi, yi, zi, vxi, vyi, vzi)T at a certain orbital phase (i.e. at apogee), is generated randomly. (2) Chose 1 satellite (i-th satellite), and regenerate the state vector xfi randomly. (3) If the objective function J(xf) is smaller than J(x), apply x:=xf, otherwise stochastically apply x:=xf according to the process temperature T. Fig.4: Simulated Annealing Algorithm 5 (4) Terminate the algorithm when J(x) is sufficiently converged, otherwise update T by T:=RT and then go to (2). 3.2. Validation: Cart-Wheel Motion To validate the SA algorithm works properly, a cart-wheel motion in LEO is tested. The fleet of satellites in a circular orbit, having relative motion plane 60 degree inclined from the orbital plane, satisfying a certain relation between relative position and velocity, would hold their relative geometry, while the whole formation is rotating one cycle per orbital period. It is called a gcart-wheel motionh, and it can be solved via Hillfs equations easily. If we chose the objective function so as just to keep the geometry of the formation, and as a result of that, if a cart-wheel motion is obtained by SA, it can be said that the SA works properly for formation optimization. The objective function is designed so as to keep the relative distances between all the satellites are kept to L for whole through the orbit; ( )2 2 1 ( ) ij i j n J t L dt C T = ? r ?? (19) To calculate the relative orbital motion over one orbital period, the solution of the following relative orbital motion, which can be found in many papers3 and described briefly in appendix, is used; ( ) T 3 2 3 2 i i i i i R R ? ? + ~ + ~ + ~ ~ =? ? ? ? ? ? RR r ? r ? r ? ? r 1 r ?? ???? ?? (20) where R is the radius, ? the orbital angular velocity, the gravity constant, respectively. Because the relative motion must be periodic in the practical formation flight, the following additional condition must hold; 3 1/ 2 3 / 2 (0) (2 ) (0) (1 ) (1 ) y e x a e e + = ? + ? ?? or (21) 3 3/ 2 1/ 2 ( ) (2 ) ( ) (1 ) (1 ) y e x a e e ? = ? + ? ?? where the depedent variable for x,y is the true anomaly . The result of SA optimization is shown in Fig.5, where the conditions are chosen as follows; ? 3 satellites in 200km circular LEO, L=1.0 ? search space: (x,y,z)=[-1:1], (vx,vy,vz)=[-1e3:+1e3]. ? SA setting: T0=100, R=0.9 ? Number of trials: 100 In Fig.5, x is the radial direction, z the outer-plane direction, and y completes the right-hand system. The solid lines indicate the analytical solution of the cart-wheel motion. As the numerical solution distribution well aligns the analytical one, it can be said that the SA algorithm works properly to generate the optimum formation. 4. High Elliptic Formation Optimizations 4.1. FDOP-O Minimization Problem The primary objective of SCOPE mission is to perform multi-point simultaneous observations around apogee, therefore consider the following objective function; 210 150 ( ) J FDOPO d = ?? ?? (22) where the reference orbit is chosen to be a high elliptic orbit for SCOPE((RE+3000km)~30RE). As is shown in section 2.4, the minimization of FDOPO automatically optimizes FDOPG as well. The result of SA optimization is shown in Fig.6, where the conditions are chosen as follows; ? 4 satellites in (RE+3000km)~30RE, sat0 at origin. Fig.5: Cart-Wheel Motion Generated by SA 6 ? Search space: (x,y,z)=[-1:1], (vx,vy,vz)= [-1e-5:1e-5]. Periodicity condition (21) must hold. ? SA setting: T0=100, R=0.9 ? Number of trials: 100 The coordinate system for Fig.6 is the same as that of Fig.5. The linear approximated analysis based on Eq.(20) and Eq.(A1) imply that there are 4 geometrically equivalent formations for a given FDOPO and FDOPG history, which are (x,y,z)=(+,+,+), (+,+,-), (-,-,+) and (-,-,-). This sign selectivity has been removed and unified to (+,+,+) to clarify the pure formation shape dependency on FDOPO. Fig.6 indicates that, of course it shapes a tetrahedron formation, but in addition to that, there is a plane symmetry about yz-plane, and viewed from x-axis, 3 satellites (sat1,2,3) forms a cocentric circle. In summary, the FDOPO optimum formation is illustrated as Fig.7. 4.2. SCOPE Formation Optimization In SCOPE mission, FDOPO is not an only thing to be optimized. Because of the communication performance limitation of the member satellite, there are following two constraints in real operation; (1) The relative distance between all two satellites must be within 500km to assure the inter-satellite communications. (2) Each satellite is equipped with the limited directivity antennae for inter-satellite communication. All the satellites have orbit-normal spin axes, and have the antenna link coverage of }60deg. Therefore the relative location of each satellite is limited to be within }60deg elevation with reference to the orbital plane. Taking above (1) and (2) into account, the objective function for SCOPE formation optimization is written as follow; ( ) ( ) 210 150 1 2 2 , 2 ( ) asin max ( ) min ( ) ij i j ij ij ij ij J FDOPO d z k d x y r k r = ? ? ? ? ? + ?? ? ? ? +? ? ? ? + ? ?? ?? ?? (23) where the two constraints are implemented as penalty functions with weighting parameters k1 and k2. The second term restricts the maximum relative elevation between all the satellites combinations, and the third term tries to minimize the variation of the relative distance over one orbital period. Fig.6: FDOPO Minimum Formation Generated by SA (Left:position, Right:velocity distributions) reference orbit (around apogee) orbital plane top view side view Fig.7: Illustration of FDOPO Minimum Formation 7 Fig.8 is the result of SA optimization under the following conditions; ? 3 satellites in 200km circular LEO, L=1.0 ? k1=1.0, k2=0.01, ?=40deg ? search space: (x,y,z)=[-1:1], (vx,vy,vz)=[-1e3:+1e3]. Periodicity condition (21) must hold. ? SA setting: T0=100, R=0.9 ? Number of trials: 100 Fig.8 indicates that the optimum formation for SCOPE is not so much different in shape as that of Fig.6, but the velocity distribution is narrowed. It is understood as the effect of the two constraints, restricting the relative motion of the satellites in a whole orbital period. 5. Heuristics for Optimum SCOPE Formation Based on the results obtained in section 4, the heuristics of designing the optimum SCOPE formation can be constructed as follow; (1) Periodicity condition All the satellites must have the same orbital period. (2) Symmetry condition Sat1 and Sat2 are on the yz-plane symmetric positions at apogee. Sat3 on the yz-plane at apogee. (3) Cocentric condition The yz-plane projection of Sat1,Sat2 and Sat3 positions are cocentric at apogee. If we approximate the condition (3) to the coplanar condition, in which Sat0, Sat1 and Sat3 are on the same orbital plane, the conditions for each satellite at apogee can be listed as follows; ? Sat0: all 0 (on the reference orbit) ? Sat1: / , 0, / tan( /2) z y x p z v y x l = = = = ?? ? Sat2: / , 0, / tan( /2) z y x p z v y x l = = = =? ?? ? Sat3: 0, / tan( ) x x y z z y l = = = = = ?? ?? ?? where 3 3/ 2 1/ 2 (2 ) (1 ) (1 ) e p a e e ? = ? + ? and l is the size of the formation, the inplane angle (Sat2-0-1), the outer-plane angle (elevation of Sat3 with reference to the orbital plane), a the semimajor axis. Substituting above conditions to Eq.(A9), we get the following formation state vector for each satellite at apogee; [Sat1][Sat2] ( ) ( ) ( ) ( ) ( ) ( ) cos / 2 , sin / 2 , 0 ( )sin /2 , 1 1 2 ( )cos /2 , 0 1 x y z x l y l z e v l e e e v l v e = =} = ? = ? ? + + ? = = ? ?? ?? ? (24) [Sat3] 0, cos , sin 0, 0, 0 x y z x y l z l v v v = = = = = = (25) Fig.8: SCOPE Optimum Formation Generated by SA (Left:position, Right:velocity distributions) 8 Originally there are 6(n-1)=18 degree-of-freedom design parameters for 4-satellite formation, but now from (24) and (25), the number of design parameters is reduced to only 4, that are ,, and l. Because we know, from section 2, that the tetrahedron is the optimum shape, the unit size formation can be assumed as =60deg, =53.5deg, l=1.0. Therefore only remaining parameter now is . The simulation examples are shown in Fig.9 and Fig.10 with different . FDOPD is 5.0 for Fig.9 and 8.0 for Fig.10, respectively. Seen from these figure, the outer-plane relative elevations are kept always below 60deg, and the FDOPO and FDOPG become small enough around apogee. It is also important to note that, approximately corresponds to the perigee/apogee relative distance ratio. This perigee/apogee relative distance ratio can be approached to 1.0 at the cost of increasing FDOPO, FDOPG and FDOPD. 6. Conclusions The new formation design strategy using simulated annealing (SA) optimization was discussed. It was shown that SA is useful to survey a whole solution space of optimum formation, taking into account realistic constraints. It was also revealed that this method was not only applicable for circular orbit, but also for high-elliptic orbit formation flying. Fig.9: Heuristic Design of SCOPE Formation (=2.0) Fig.10 Heuristic Design of SCOPE Formation (=5.0) 9 The developed algorithm was applied to SCOPE formation design. A distinctive and useful heuristics was found by investigating SA results, showing the effectiveness of the proposed design process. Reference 1. Fujimoto, M., Saito,Y., Tsuda,Y., Shinohara,I., Kojima,H., gThe scientific targets of the SCOPE missionh, COSPAR04-A-01348, D2.3/E3.3/ PSW2-0016-04, COSPAR 2004 2. Paschmann, G.., Daly, P. W., Analysis Methods for Multi-Spacecraft Data., SR-001, ISSI Scientific Reports, 1998 3. How, J. P., Inalhan, G., gRelative Dynamics & Control of Spacecraft Formations in Eccentric Orbitsh, AIAA Paper, AIAA-2000-4443, 2000 10 Appendix The analytical solution for Eq.(20) is derived as follows3; ( ) ( ) 2 1 2 3 2 1 2 3 4 5 6 cos ( ) sin 2 sin ( ) cos 1 cos ( ) 1 cos 2 1 cos ( ) 1 1 1 sin 1 cos 1 cos sin cos ( ) 1 cos 1 cos ( ) c e x e d e H d d e y e d e e H d d d e e z d d e e x e ? ? ? ? = + ? + ? ? ? ? ? ? ? ? ? ? ? + ? ? ? ? = + + + ? ? ? ? ? ? ? ? ? ? ? ? + + + ? ? ? ? ? ? + + ? ? ? ? ? ? ? ? ? ? = + ? ? ? ? + + ? ? ? ? = ( ) ( ) ( ) 2 2 2 1 2 2 3 3 2 2 1 2 2 3 2 sin 2 sin cos os 2 cos ( ) 2 sin ( ) 1 cos 1 cos sin ( ) sin 2 ( ) 2 sin ( ) 2 cos ( ) cos sin sin cos 1 cos 1 cos 1 cos e e d e H e H d e e d y e d eH e H e H d e e d e e e ? ? ? ? + + + ? ? ? ? ? ? ? + + ? ? + ? ? ? ? ? ? = ? + ? + ? ? ? ? ? ? ? ? ? ? + + + + + ? ? + + ? ? ( ) ( ) ( ) 4 2 2 5 6 2 2 cos sin sin sin cos ( ) 1 cos 1 cos 1 cos 1 cos d e e z d d e e e e ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? = + + ? + + + ? ? ? ? + + ? ? ? ? (A1) where ( ) ( ) 5 / 2 2 2 3 ( ) 1 1 sin sin cos 2 2 eE e H e e E E E ? ? ? =? + ? + + ? ? ? ? (A2) (E)f is -differentiation of (E), e the eccentricity and E the eccentric anomaly. -differentiation and t-differentiation is related by the following equation; 2 3 2 3 / 2 (1 cos ) (1 ) e a e + = ? ?? (A3) (A1) indicates that the inplane motion and outer-plane motion are independent of each other, and the relative distance history is equivalent for the sign combinations of (x,y,z)=(+,+,+), (+,+,-), (-,-,+), (-,-,-). The periodicity condition (21) is equivalent to d2=0. When this periodicity holds, the state vector at perigee can be expressed using d1,c,d6 as follows; 3 1 4 6 1 3 5 ( ) 1 ( ) 1 1 1 ( ) 1 ( ) 1 ( ) 1 1 1 ( ) 1 x d y ed d e z d e x ed y d e z d e =? ? ? = + + ? ? ? ? ? ? + ? ? ? ? =? ? + ? ? =? ? ? ? ? ? = + ? ? + ? ? ? ? =? ? + ? ? (A4) In the same way, the state vector at apogee can be expressed as follows; 3 1 4 6 1 3 5 ( ) 1 ( ) 1 1 1 ( ) 1 ( ) 1 ( ) 1 1 1 ( ) 1 x d y e d d e z d e x ed y d e z d e = ? ? = ? + ? ? ? ? ? ? ? ? ? ? ? =?? ? ? ? ? = ? ? ? ? ? ? ? = ? + ? ? ? ? ? ? ? = ?? ? ? ? ? (A5) If d1 is selected so that the inplane relative distance is preserved between apogee and perigee, it is derived from (A4) and (A5) that; [ ] [ ] * * 1 4 1 4 * 1 4 2 1 1 1 1 1 1 1 1 e d d e d d e e d d e ? ? ? ? + + = ? + ? ? ? ? + ? ? ? ? ? ? = ? (A6) Let denote the offset ratio of d1 from d1 *.; * 1 1 d d = (A7) Substituting (A6),(A7) into (A5), we get; 3 * 1 * 1 3 (0) 1 (0) 1 1 (0) 1 (0) 1 1 x d y d e e x ed y d e = ? ? ? = ? + ? ? ? + ? ? = ? ? ? = + ? ? + ? ? (A8) and 3 * 1 * 1 3 ( ) 1 ( ) 1 1 ( ) 1 ( ) 1 1 x d y d e e x ed y d e = ? ? = ? + ? ? + ? ? ? = ? ? ? = ? + ? ? ? ? ? (A9) If <1, the formation size at apogee becomes larger than the size at perigee, while if >1, the size at apogee becomes smaller than the size at apogee. When =1, the inplane distance at apogee is equal to that of perigee.