用进化模拟退火算法求解Lambert方程

(1
−
cos
∆f
)
rK1
K V2
=
np r1r2 sin ∆f
(
rK2
−
rK1 )
+
r1 P
(1
−
cos
∆f
) rK2
而 Newton 迭代方法如下
给定方程 F(a) = ∆t − ∆t(a) = 0 给定初值 a0 = amin = s/ 2 然后给定迭代公式
an +1 = an − F (an ) / F / (an ), n = 0,1,2,....
进行迭代 直到满足收敛精度
在求得 a 以后 我们便可以通径的公式求得
P
=
4 a ( s − r1 )( s − r2 ) sin c2
2 α
± 2
β
其中对于 Lambert 定K理中K 1 3 用 + 号 2 4 用
-
号
从而得到V1, V2K
在用进化模拟退火算法求得 a 以 K
后 仍用上述公式算得V1, V2
第 19 卷第 2 期 2007 年 1 月
系 统 仿 真 学 报© Journal of System Simulation
Vol. 19 No. 2 Jan., 2007
用进化-模拟退火算法求解 Lambert 方程
王 石 1,2 文援兰 2 戴金海 2
(1.国防科技大学 ATR 实验室, 湖南 长沙 410073; 2.国防科技大学航天与材料工程学院, 湖南 长沙 410073)
几何函数展开方法并不收敛 因为椭圆的偏心率超过了极限
值 EA-SA 和 Newton 迭代法都具有较高的收敛精度 但
• 451 •
பைடு நூலகம்
万方数据
第 19 卷第 2 期 2007 年 1 月
=
1 2
(1 −
λ
−
xη ),η
=
y
−
λx
y=
1 − λ2 (1 − x 2 )
利用超几何级数迭代方法 可得 x 的值 进而可求出转移轨
道的长半轴和半通径
a
=
am 1− x2
P
=
r1 r2 a mη 2
sin 2
∆f 2
进而可求得
K V1 =
np r1r2 sin ∆ f
(
rK2
−
rK1 )
+
r2 P
step5 选择出前面最好的若干个个体组成 P(t+1)且 P(t+1)
和 P(t)的群体大小一样
step6 若不满足终止要件 则 t ←t +1转向 step3 继续进
化 若满足终止条件 则输出当前最优值 算法结束
2 算法的仿真实现及其和经典方法的对比
为了说明算法的收敛精度 我们首先仿真一个 真实
Evolutionary Algorithm, EA 和模拟退火算法[2] Simulated Algorithm, SA 来求解 Lambert 方程 进化算法具有全局优 化的能力 而模拟退火算法具有较强的局部探索能力 二者 相结合求解 Lambert 飞行时间问题 具有精度高 对模型依 赖少 操作容易等优点
是模拟自然界生物进化的一种算法 是一种智能搜索算法
其在工程中得到了广泛的应用 具体可参见[2-3] EA-SA
算法是把 SA 算法嵌入到 EA 中 因为 EA 算法的局部搜索
能力不强 而 SA 恰恰相反 因此两种算法的配合可达到更
好的效果 EA-SA 算法如下
step1 进化代数的初始化 t ← 0
轨道 然后用超几何函数展开 Newton 迭代法和 EA-SA 算
法进行精度比较 说明 EA-SA 的优越性
超几何函数展开的算法如下
令 x2 = 1 − am a 则 Lambert 定理可表示 如下的超几何级数形式[3-4]
∆t =
a
2 m
η
4 3
η
3
F
(3,1,
5 2
,
S1
)
+
4λη
其中 S1
退火的机理而建立起的一种优化算法 一般设置能量函数 当达到冷却状态时能量函数最小[2] 模拟退火算法具有较强
的局部搜索能力 简单的模拟退火算法如下(其中 xi 为迭代 序列)
step1 设置初始温度 T0 冷却温度 εT k = 0 to ←T0 xi ←x0 step2 若满足终止循环条件 则 step3 否则 从 x 的邻域 中任选一 x′
step2 随机产生初始群体 P(t)
step3 实行遗传算子操作
1 交叉 对 P′(t) 进行两两随机交叉 P(t) → P′(t) 2 变异 对每个个体进行正态变异操作 P(t) → P′′(t) 令 P′′′(t) = P(t) ∪ P′(t ) ∪ P′′(t )
step4 对 P ′′′( t ) 模拟退火算子操作 P (4) (t ) ← SA [P ′′′(t )]
摘 要 Lambert 方程在轨道拦截和初始轨道确定起着重要作用 求解 Lambert 方程的传统算法主
要有 Newton 迭代方法和超几何级数展开方法等 但这些算法都有一定的局限性 如有可能出现迭
代收敛过慢 级数展开收敛性问题 采用进化-模拟退火算法 EA-SA 算法求解 Lambert 方程
其中进化算法具有全局搜索能力 而模拟退火具有局部锁搜索能力 该方法克服了某些情况下梯
数展开 Newton 迭代方法进行计算比较 令 ∆a =| ac − a |表
示计算的长半轴和 真实 长半轴之差 即长半轴的计算精
度
而 ∆ri =| rGic − rGi |,
GH ∆Vi =|Vic − Vi |,
i = 1,2 分别表示 ti 时刻
的位置和速度的精度
从表 2 表 3 和表 4 可以看出 对于大椭圆轨道情形 超
计算 ∆ = E(x′) − E(x)
若∆ < 0 x ← x′
否则 当 exp( −∆ / t k ) > rand (0,1)
则 x ← x′ 重复 step2
step3
tk +1
=
T0 1 + tk
,k
←
k
+1
则 回到 step2
若 t k +1 < ε T
则终止 否
在我们讨论的问题中 E(x) = ∆t − ∆t(x) 2 进化算法(EA)
引言
已知卫星在空间的两点位置矢量以及卫星飞过两点的 时间 则可以确定卫星的轨道要素 主要的方法有高斯方法 f 和 g 级数方法 真近点角迭代方法 偏心率迭代法等 它 们都有一定的适用范围 还有一种重要的方法 就是求解 Lambert 方程 可用 Newton 迭代方法加以求解 但有可能 出现在某些点迭代收敛过慢 也可用超几何级数展开 可以 获得高精度 速度快等优点 但对模型的依赖性过强 也就 是说 需要过多的模型知识 而且用连分数展开时 级数收 敛性也是必须注意的一个问题 本文提出用进化算法[1]
模与仿真 计算机实时控制与实时系统 并行工程与产品数据管理
OE 图 1 Lambert 定理示意图
则 Lambert 飞行时间定理为[3]
• 450 •
万方数据
第 19 卷第 2 期 2007 年 1 月
王 石, 等 用进化-模拟退火算法求解 Lambert 方程
Vol. 19 No. 2 Jan., 2007
无论上述那一种情况 总可以写成如下的形式 ∆t − ∆t(a) = 0
1
其中 ∆t(a) 为上述(1) (2) (3)和(4)等式子的右边除以 µ 2 引入能量函数 E(Χ) = ∆t − ∆t(Χ) 2 很显然满足上述
Lambert 定理的 a 的能量为零
模拟退火算法 Simulated Annealing, SA 是基于金属
Solution of Lambert Theorem By Evolutionary Algorithm and Simulated Annealing
WANG Shi1,2, WEN Yuan-lan2, DAI Jin-hai2
(1.ATR State Key Lab of NUDT, ChangSha 410073, China; 2.College of Aerospace and material Engineering of NUDT, ChangSha 410073, China)
1 进化-模拟退火算法的设计
如图 1 所示 设地心为 OE 航天器由给定点 rK1 沿椭圆
轨道飞至给定点 rK2 设两点的地心距和夹角分别为 r1 r2 和
∆f
令c =
rK2 − rK1 ,
s=
r1 + r2 + C 2
即 S 为∆OE12 的半周长 设
椭圆的半长轴为 a 偏心率为 e 半通径为 P 令
a min = s / 2
K
=1− c/s
=
r1 + r2 r1 + r2
−c +c
<1
1
α
=
2
arcsin
a min a
2
1
β = 2 arcsin ( K a min / a ) 2
再令
a m = π
1
β m = 2 arcsin( K 2 )
∆tm
=
S2 2µ
π
− (βm
− sin 2
我们仿真一个 真实 的二体轨道 分大椭圆轨道 小
偏心率轨道 近圆轨道等 3 种情况讨论 其轨道参数见表 1
假定初始时到另外一时刻的转移时间为
∆t = t2 − t1 = T / 6 其中 T 表示轨道周期 即 T = 2π a 3 / µ
式中 µ = 3.986005 × 1014 为引力常数
下面我们对这 3 种情况分别用 EA-SA 算法 超几何函
βm )
G
1
V1
c
rG 1
rG 2
其中 µ 为引力常数
2
收稿日期 2005-10-28
修回日期 2006-11-23
作者简介 王石(1972-), 男, 河南濮阳人, 博士后, 研究方向为飞行器系统
建模与仿真; 文援兰(1965-), 男 ,湖南人, 副教授, 研究方向为轨道力学;
戴金海(1943-), 男 , 河北昌黎人, 教授, 博导, 研究方向为飞行器系统建
度下降法有时收敛过慢和超几何级数不收敛的缺点 并具有通用性 便于操作和理解 通过仿真计
算对比表明 EA-SA 具有普适性 而且精度优于其它两种算法
关键词 Lambert 定理 进化算法 模拟退火算法 仿真分析
中图分类号 V412
文献标识码 A
文章编号 1004-731X (2007) 02-0450-03
假设 ∆t 为从 1 点到 2 点的飞行时间 则 (1) ∆f ≤ π ∆t ≤ ∆tm
1
µ 2 ∆ t = [ a 3 / 2 (α − sin α ) − ( β − sin β ) ]
(2) ∆ f ≤ π , ∆ t > ∆ t m
1
µ 2 ∆ t = [ a 3 / 2 2π − (α − sin α ) − ( β − sin β ) ]
Abstract: Lambert theorem plays an important role in orbital interception and initial orbit determination. The classical methods, such as Newton iteration, hyper-geometry expansion and et al, are always valid in limited conditions. Evolutionary Algorithm and Simulated Annealing (EA-SA) was proposesed to solving Lambert theorem, which EA has global search performance and SA has local search performance. The algorithm overcomes the defects of low convergent speed of Newton iteration method in some cases and it is easy to be operated and understood. Finally, the computation verifies that EA-SA is universal adaptive and it is more precise than the classical methods. Key words: Lambert theorem; Evolutionary Algorithm (EA); Simulated Annealing (SA); Simulation Analysis
(3) ∆ f ≥ π ∆ t ≤ ∆ t m
1
µ 2 ∆ t = [ a 3 / 2 (α − sin α ) + ( β − sin β ) ]
(4) ∆f ≥ π , ∆ t > ∆tm
1
µ 2 ∆ t = [ a 3 / 2 2π − (α − sin α ) + ( β − sin β )]