新型太空电梯——半摆卫星太空电梯的设计与验证

最新推荐文章于 2025-06-11 18:04:45 发布

原创最新推荐文章于 2025-06-11 18:04:45 发布

· 1.1k 阅读

12 ·

版权

文章标签：

#航天 #数学 #仿真 #科技 #验证

新想法专栏收录该内容

1 篇文章

订阅专栏

一、同步轨道太空电梯方案简介

1.1 简介

当一个航天器运行在地球赤道上方的同步轨道时，其运行周期与地球的自转周期相同。利用这一性质，如果在该轨道处垂下一条绳缆到地表，则该绳缆会与地面保持相对静止，绳缆上可以通过类似电梯的设备运输物品到太空，其成本将远低于当前的火箭。

据估计火箭运载成本约7500 美元/kg，而太空电梯预计成本仅约200
美元/kg，相差约三十多倍。

传统太空电梯主要结构如图1.1、1.2所示，其中卡门线（Kármán line）通常被认为是大气与太空的分界线，卡门线以下被认为是各国领空，其上则为太空。

在这里插入图片描述

图1.1 太空电梯示意图

在这里插入图片描述

图1.2 卡门线

在这里插入图片描述

图1.3 a 《流浪地球2》中的太空电梯
在这里插入图片描述

图1.3 b 《流浪地球2》中的太空电梯

1.2 技术细节

（1）绳缆

为保证电梯绳缆各个截面处所受应力相同，绳缆各截面的面积将不断变化，下窄上宽，且并非为简单的线性变化。定义绳缆上截面与下截面的比为锥度比，通过对绳缆所受重力与离心力的综合分析，可得到该锥度比的计算公式：
在这里插入图片描述

（1-1）

公式中各符号含义如下表所示：

符号	物理意义	单位
( A(x) )	距离地心 ( x ) 处的横截面积	m²
( A(x_0) )	地球表面（( x = x_0 )）的初始横截面积	m²
( \rho )	材料密度	g/cm³
( \sigma_{\text{max}} )	材料抗拉强度	GPa (GN/m²)
( G )	万有引力常数	N·m²/kg²
( M )	地球质量	kg
( \Omega )	地球自转角速度	rad/s
( x )	太空电梯某点距地心的距离	m
( x_0 )	地球平均半径	m

公式中物理常数和参考值如下表所示：

符号	物理意义	参考值	单位
( G )	万有引力常数	( 6.674 \times 10^{-11} )	N·m²/kg²
( M )	地球质量	( 5.972 \times 10^{24} )	kg
( \Omega )	地球自转角速度	( 7.292 \times 10^{-5} )	rad/s
( x_0 )	地球平均半径	( 6.371 \times 10^6 )	m

当确定绳的长度 ( x ) 时，会发现式（1-1）左侧的锥度比只与 ( \rho / \sigma_{\text{max}} ) 项有关，即抗拉强度与材料密度的比值，该比值又称为材料的比强度。

例如作为当前唯一有潜力成为同步轨道太空电梯材料的碳纳米管，其抗拉强度约80-100GPa，密度约1.3
g/cm³，比强度约61.5-76.9
GPa/(g/cm³)。利用"附件1------Matlab计算锥度比"中的matlab代码进行计算，得其对应锥度比约2.2，在可承受范围之内。

然而由于技术问题，当前能够制造碳纳米管的最大长度只有0.5米，与3.6万千米的要求相差甚大，且由于宏观缺陷等问题，实际制造的该材料绳缆的性能远低于理论，仅10GPa左右。下表列出了可实际制造的比强度最高的几种材料的参数：

材料名	抗拉强度（GPa）	密度（g/cm³）	比强度（GPa/(g/cm³)）	备注
碳纳米管	80-100	1.3	61.5-76.9	宏观连续制备困难
碳纤维T1100G	7	1.8	3.89	剖析日本东丽S、H、G型高强碳纤维典型特性
PBO纤维（Zylon AS）	5.8	1.54	3.77	pbo纤维_百度百科
超高分子量聚乙烯（Dyneema SK99）	3.5	0.97	3.61	解密高强T1100G碳纤维匹配树脂的典型特征
芳纶纤维（Kevlar KM2）	3.6	1.44	2.5	—

一般认为锥度比为个位数级别时太空电梯才能够实际制造。附录1根据式（1-1）计算各材料与高度下的锥度比，计算得从地面到同步轨道处的各绳缆的锥度比为：

碳纳米管：2.2
其它材料：最小也达10⁵级别。

平衡锤

另一方面，为了平衡下方绳缆的重力，需要一运行轨道直径大于同步轨道的平衡锤与太空电梯相连。

二、方案1简介

2.1 初步设想——单半摆卫星

在这里插入图片描述

图2.1 单半摆卫星模型

如图2.1设地球静止不动，由其赤道上方垂直固连一杆( AO ) ，( AO ) 又与另一杆 ( AB ) 连接，设两杆为刚体，且无质量，铰接处光滑。杆 ( AB ) 以角速度 ( \omega ) (( rad/s ))绕 ( A ) 点逆时针旋转。

现有一卫星 ( C )，质量 ( m_{卫} ) (kg)，为 ( \angle BAO ) 如图所示为直角时，其位于轨道上的 ( H )点(H、I点沿轴AO对称)。此时卫星的速度为 ( v_{卫}(t_H) ) = ( v_H = v_0 )。

且卫星到达(I)点时杆杆(AB)上的(B)点也与(I)点重合；且卫星(C)到达(I)点时自动与(B)上固连，到(H)点时自动与( B ) 点分离。

在该过程中，只要半摆卫星质量足够，则杆(AB)运动过程中产生的离心力的平均合力会与杆(AB)、(AO)自身的重力相等。

不过上述模型受力并不均匀，一个周期中(A)点在不同时刻所受的力的大小、方向都不相同，解决这个问题可以使用多半摆卫星。

2.2 完善设想——多半摆卫星

在这里插入图片描述

图2.2 多半摆卫星模型

进一步的如图2.2，摆臂不再只是一个杆，而是如同摩天轮一样，由多个杆A1、AB2、AB3…组成，每个Bn点都可以挂载半摆卫星以产生向上离心力。则当半摆卫星足够多时，A点所受的合力近似均匀，持续受到垂直向上的拉力。

同时，由于A点受力均匀，AO、ABn等可以不再是理论中的刚体，而可以换作实际可行的柔性绳缆。

另一方面，对于C点连接处，可以使用一根轴上多个轴承的结构K，AB绳缆固定于轴K中间位置，ABn等摆臂各安装于两侧。为使两侧受力平衡，可两端可以同步运行一个“摩天轮”。

另一方面，上述假设地球不动，相当于太空电梯位于南北两极，于赤道等处是否可行留待验证。

三、卫星轨道仿真验证

下图中各部分含义：

蓝色球体地球
橙色线网大气与太空分界线（200km高度）
亮蓝线条卫星运行轨道
红色点卫星起始位置
黄色点卫星实时位置

在这里插入图片描述

图3.1 正常卫星运行的轨道
在这里插入图片描述

在这里插入图片描述

图3.2 半摆卫星低速时实际运行轨迹
在这里插入图片描述

在这里插入图片描述

图3.3 半摆卫星高速时实际运行轨迹

然而，上述椭圆部分的轨道在其它无外力作用下不能实现的，由开普勒第一定律得，卫星绕地球运行的轨道是一个以地球为焦点之一的椭圆。

经附录2程序对卫星运行轨迹进行仿真，会发现该卫星离开H点且不与地球相撞的轨迹就当是这样的，如图3.3所示。卫星在近地点D后，由于速度过大，会一骑绝尘地远离地球，设想中的轨迹在该条件下不能实现。

四、方案2简介

在这里插入图片描述

图4.1 双半摆卫星协作模型

参考“2.2”小节，基于方案1进行优化，采用两个半摆卫星太空电梯进行协作，半摆卫星在两个太空电梯间循环传递。经数学仿真，该对称式的卫星运行轨道是可以实现的，进一步要对缆绳材料的性能进行验证。

五、缆绳材料性能验证

先来看一下能连接至同步轨道处的绳缆其锥度比是怎么变化的，使用附录1中的程序进行计算：

在这里插入图片描述

图5.1 碳纳米管锥度比变化
在这里插入图片描述

图5.2 碳纤维T1100G锥度比变化

如图5.1、5.2，发现各材料锥度比附高度的升高而增大，且由于临近地表处的重力加速度要比远地处更大，而远地处向心力又大于近地处，其增大的速度也不是线性的。

经计算，碳纤维T1100G材料绳缆锥度比为2.2左右时，其离地高度仅约350km，如图5.3所示。而卫星运行高度一般大于200km，而方案2要求绳长度最好长于一个地球的直径，即约6371km。随绳缆长度减小，半摆卫星轨迹中的有效作用轨迹将越来越短，350Km的高度几乎已经贴近地面，不认为具有可行性。

在这里插入图片描述

图5.3 碳纤维T1100G绳缆高度分析

六、项目总结

经验证方案1、方案2分别由于卫星运行轨迹不可实现、现有材料可制造强缆高度过小而不具有工程上的可实现性。

然而，以另一种方式绕开当前太空电梯当前绳缆材料方面的困难，仍是值得思考的。本项目的主体思路是以多卫星的作用代替同步轨道太空电梯中的“单卫星”作用，那还有没有其它更好、更简单的思路呢，太空电梯一定要建设到同步轨道那么高吗？可不可以通过减小材料密度方式提高比强度？近地的大气层是否可以利用？能否与其它低成本太空运输方案结合？能否利用其它的作用，电、磁、光、机械波、真空、其它天体……
未来的太空电梯也许与今天普遍设想的方案截然不同。

附录

附录1 Matlab锥度比计算

注：程序主体由Deepseek生成

a1_Rope_compute.m

close all;
clear all;
clear;

% --- 定义常数和材料参数（使用新单位） ---
G = 6.674e-11; % 万有引力常数 [N·m²/kg²]
M = 5.972e24; % 地球质量 [kg]
Omega = 7.292e-5; % 地球自转角速度 [rad/s]
x0 = 6.371e6; % 地球半径 [m]
% x_sync = 4.2164e7; % 同步轨道半径[m]
% x_sync = x0*3; % 地球直径高度时 [m]
x_sync = x0+350e3; % 地球直径高度时 [m]

% % 碳纳米管
% rho_gcm3 = 1.3; % 材料密度 [g/cm³] 
% sigma_max_GPa = 80; % 抗拉强度 [GPa] 

% 碳纤维T1100G
rho_gcm3 = 1.8; % 材料密度 [g/cm³] 
sigma_max_GPa = 7; % 抗拉强度 [GPa] 


% --- 定义计算范围（从地球表面到同步轨道）---
x = linspace(x0, x_sync, 100); % 生成 100 个点，范围：x0 到 x_sync 米


% --- 调用函数计算锥度比 ---
A_ratio = calculate_cone_ratio(x, rho_gcm3, sigma_max_GPa, G, M, Omega, x0);

% --- 可视化结果 ---
plot(x, A_ratio);
xlabel('距离地心的距离 x [m]');
ylabel('A(x)/A(x_0)');
xlim([0 x_sync]); % 限制绘图区域 
ylim([0 A_ratio(100)]);

title('太空电梯锥度比随高度的变化（单位：g/cm³, GPa）');
grid on;

a1_Rope_compute.m

function A_ratio = calculate_cone_ratio(x, rho_gcm3, sigma_max_GPa, G, M, Omega, x0)
    % 计算太空电梯横截面积比 A(x)/A(x0)
    %
    % 输入参数：
    %   x             : 距离地心的距离（单位：米，m），标量或数组
    %   rho_gcm3      : 材料密度（单位：克每立方厘米，g/cm³）
    %   sigma_max_GPa : 材料抗拉强度（单位：吉帕斯卡，GPa）
    %   G             : 万有引力常数（单位：牛·米²/千克²，N·m²/kg²）
    %   M             : 地球质量（单位：千克，kg）
    %   Omega         : 地球自转角速度（单位：弧度每秒，rad/s）
    %   x0            : 地球平均半径（单位：米，m）
    %
    % 输出参数：
    %   A_ratio       : 横截面积比 A(x)/A(x0)，无量纲
    %
    % 公式：
    %   A_ratio = exp[ (ρ/σ_max) * (GM(1/x0 - 1/x) + 0.5Ω²(x0² - x²)) ]

    % --- 单位转换 ---
    rho = rho_gcm3 * 1000;          % 将 g/cm³ 转换为 kg/m³ (1 g/cm³ = 1000 kg/m³)
    sigma_max = sigma_max_GPa * 1e9; % 将 GPa 转换为 Pa (1 GPa = 1e9 Pa)

    % --- 引力项计算（单位：m²/s²） ---
    GM_term = G * M * (1/x0 - 1./x); 
    
    % --- 离心力项计算（单位：m²/s²） ---
    Omega_term = 0.5 * Omega^2 * (x0^2 - x.^2); 
    
    % --- 综合指数项（无量纲） ---
    exponent = (rho / sigma_max) * (GM_term + Omega_term); 
    
    % --- 计算面积比（无量纲） ---
    A_ratio = exp(exponent); 
end

附录2 Matlab卫星运行轨道仿真

注：程序主体由Deepseek生成

main.m

close all;
clear;

satellite_trajectory();

satellite_trajectory.m

function satellite_trajectory()
    % --- 地球参数与物理常数 ---
    Re = 6371e3;                % 地球半径 [m]
    mu = 3.986e14;              % 地球标准引力参数 [m³/s²]
    atm_boundary = Re + 200e3;  % 太空分界线高度（地表+200km）[m]
    V_I = 7.9e3;                % 第一宇宙速度[m/s]
    
    % --- 初始条件（示例：近地圆轨道）---
%     r0 = [Re + 400e3, 0, 0];    % 初始位置矢量 [m]（400km高度赤道轨道）
%     v0 = [0, sqrt(mu/norm(r0)), 0]; % 初始速度矢量 [m/s]（圆轨道速度）
%     r0 = [-1.2*Re, 2.5*Re, 0];    % 初始位置矢量 [m]（400km高度赤道轨道）
%     v0 = [-0.3*V_I, -0.8*V_I, 0]; % 初始速度矢量 [m/s]（圆轨道速度）

%     % 标准运行轨道
%     r0 = [0, 2.5*Re, 0];    % 初始位置矢量 [m]（400km高度赤道轨道）
%     v0 = [ -0.8*V_I, 0, 0]; % 初始速度矢量 [m/s]（圆轨道速度）

%      % 低速半摆卫星
%     r0 = [2*Re, 2.5*Re, 0];    % 初始位置矢量 [m]（400km高度赤道轨道）
%     v0 = [ -0.5*V_I, 0, 0];    % 初始速度矢量 [m/s]（圆轨道速度）

     % 高速半摆卫星
    r0 = [2*Re, 2.5*Re, 0];    % 初始位置矢量 [m]（400km高度赤道轨道）
    v0 = [ -0.8*V_I, 0, 0];    % 初始速度矢量 [m/s]（圆轨道速度）
    
    % --- 仿真参数 ---
%     T = 2*pi*sqrt(norm(r0)^3/mu) * 3; % 轨道周期 [s]
*     T = 2*pi*sqrt(norm(r0)^3/mu) * 0.5; % 轨道周期 [s]    %高速半摆卫星
    tspan = [0, 1.5*T];         % 仿真时间（1.5个轨道周期）
    
    % --- 求解运动方程（ODE45）---
    options = odeset('RelTol',1e-9,'AbsTol',1e-9);
    [t, y] = ode45(@orbital_dynamics, tspan, [r0, v0], options, mu);
    
    % --- 可视化 ---
    figure('Color','k','Position',[100 100 1200 800])
    axes('Color','k', 'XColor','w','YColor','w','ZColor','w')
    hold on; grid on; axis equal; view(3)
    xlabel('X [m]'); ylabel('Y [m]'); zlabel('Z [m]')
    title('绕地卫星轨迹仿真','Color','w')
    
    % 绘制地球
    [X,Y,Z] = sphere(50);
    surf(Re*X, Re*Y, Re*Z, 'FaceColor',[0.2 0.5 0.8], 'EdgeColor','none',...
        'FaceAlpha',0.8, 'DisplayName','地球');
    
    % 绘制太空分界线（200km高度）
    surf(atm_boundary*X, atm_boundary*Y, atm_boundary*Z,...
        'FaceColor','none', 'EdgeColor',[1 0.5 0], 'LineWidth',1.5,...
        'DisplayName','太空分界线（200km）');
    
    % 标注地球直径
    plot3([-Re, Re], [0,0], [0,0], 'r--', 'LineWidth',1.5,...
        'DisplayName','地球直径');
    text(Re,0,0, sprintf('直径: %.0f km',2*Re/1e3),...
        'Color','r', 'HorizontalAlignment','left')
    
    % 绘制卫星轨迹
    plot3(y(:,1), y(:,2), y(:,3), 'c', 'LineWidth',1.5,...
        'DisplayName','卫星轨迹');
    scatter3(y(1,1), y(1,2), y(1,3), 100, 'r', 'filled',...
        'DisplayName','初始位置');
    
    % 动画演示
    sat = scatter3(y(1,1), y(1,2), y(1,3), 100, 'y', 'filled',...
        'DisplayName','卫星实时位置');
    for k = 1:10:length(t)
        sat.XData = y(k,1);
        sat.YData = y(k,2);
        sat.ZData = y(k,3);
        drawnow
    end
    
    legend('show','TextColor','w')
end

% --- 轨道动力学方程 ---
function dydt = orbital_dynamics(~, y, mu)
    r = y(1:3);       % 位置矢量
    v = y(4:6);       % 速度矢量
    r_norm = norm(r); % 地心距
    
    % 加速度计算（仅考虑地球引力）
    a = -mu * r / r_norm^3;
    
    dydt = [v; a];    % 状态导数
end