如何使用ReachabilityAnalysis.jl提取微分方程解的x值？

我在Julia中使用ReachabilityAnalysis.jl包求解带区间初始条件的微分方程，代码如下：

prob = @ivp(x' = f(x), x(0) ∈ X0, dim=1)
sol = solve(prob, alg=TMJets21a(abstol=1e-18), T=15)

当前得到的结果是TaylorModelReachSet类型，形式如下：

TaylorModelReachSet
  55.0 + 5.0 x₁ + ( 44.9801344648411 - 4.990972005107259 x₁ - 0.0024613652164294903 x₁²) t + ( - 22.449459189998198 + 2.4688376478960143 x₁ + 0.010727058992613947 x₁²) t² + ( 7.40324329161496 - 0.757126031600233 x₁ - 0.03224951488485455 x₁²) t³ + ( - 1.702781535410468 + 0.04388098980804231 x₁ + 0.08253929478277094 x₁²) t⁴ + ( 0.07895091288058795 + 0.2882499345694954 x₁ - 0.1908707483937028 x₁²) t⁵ + ( 0.43218402721400556 - 0.6203143736740161 x₁ + 0.4098322353914475 x₁²) t⁶ + ( - 0.7971949696667339 + 1.141845422435601 x₁ - 0.8306607029702969 x₁²) t⁷ + ( 1.2840090159806252 - 2.0106344691518556 x₁ + 1.6070995456266788 x₁²) t⁸ + [-2.21465e-07, 2.51483e-07]
[0, 0.00263234]

我希望将其转换为类似x = [[1,2], [4,8], ..., [23,45]]格式的区间数组，每个子数组代表对应时间步的x值区间。

解决方法 #

要将TaylorModelReachSet转换为目标格式的区间数组，可按以下步骤操作：

1. 转换为区间集合 #

使用overapproximate函数将泰勒模型形式的可达集近似为区间类型，这是保留原解区间信息的关键步骤：

interval_sol = overapproximate(sol, Interval)

2. 提取区间上下界生成数组 #

通过迭代转换后的区间集合，提取每个区间的上下界并组织为二维数组：

x_intervals = [[low(set), high(set)] for set in interval_sol]

3. （可选）获取对应时间区间 #

如果需要同时得到每个x区间对应的时间范围，可提取每个ReachSet的时间边界：

time_ranges = [[t_start(set), t_end(set)] for set in interval_sol]

完整示例代码 #

假设f(x)为具体的微分方程（如f(x) = -x），完整可运行代码如下：

using ReachabilityAnalysis

# 定义微分方程与初始区间
f(x) = -x
X0 = Interval(50.0, 60.0)  # 示例初始区间
prob = @ivp(x' = f(x), x(0) ∈ X0, dim=1)

# 求解微分方程
sol = solve(prob, alg=TMJets21a(abstol=1e-18), T=15)

# 转换为目标格式的区间数组
interval_sol = overapproximate(sol, Interval)
x_intervals = [[low(set), high(set)] for set in interval_sol]

# 输出结果
println(x_intervals)

说明 #

• overapproximate(sol, Interval)会为每个泰勒模型可达集生成最紧的区间包络，完整保留原解的区间约束。
• 生成的x_intervals中，每个子数组[low, high]对应解在某一时间段内的x值范围。

内容的提问来源于stack exchange，提问作者Kanishk Sharma