三轴正交飞轮构型角加速度计算详解

最新推荐文章于 2025-10-19 17:19:10 发布

ScilogyHunter

最新推荐文章于 2025-10-19 17:19:10 发布

阅读量1k

点赞数 26

CC 4.0 BY-SA版权

分类专栏：航天器文章标签：算法

本文链接：https://blog.csdn.net/ScilogyHunter/article/details/152091816

航天器专栏收录该内容

102 篇文章

订阅专栏

三轴正交飞轮构型角加速度计算详解

1. 引言

在卫星姿态控制系统中，角加速度计算是动力学核心问题。三轴正交飞轮构型因其解耦特性，使得角加速度计算相对简化。本文将深入分析其计算原理和实现方法。

2. 动力学基础原理

2.1 欧拉动力学方程

卫星姿态动力学的基本方程为：
$\mathbf{J} \cdot \frac{d\boldsymbol{\omega}}{dt} + \boldsymbol{\omega} \times (\mathbf{J} \cdot \boldsymbol{\omega}) = \mathbf{M}_{\text{total}}$

其中：

$J\mathbf{J}$ 为惯量张量
$ω\boldsymbol{\omega}$ 为角速度向量
$Mtotal\mathbf{M}_{\text{total}}$ 为总外力矩

2.2 三轴正交构型的简化

对于三轴正交飞轮构型，惯量张量近似为对角阵：
$\mathbf{J} = \begin{bmatrix} J_x & 0 & 0 \\ 0 & J_y & 0 \\ 0 & 0 & J_z \end{bmatrix}$

3. 角加速度计算公式推导

3.1 完整动力学方程展开

将向量方程展开为分量形式：

X轴方程：
$J_x \dot{\omega}_x - (J_y - J_z)\omega_y\omega_z = M_x$

Y轴方程：
$J_y \dot{\omega}_y - (J_z - J_x)\omega_z\omega_x = M_y$

Z轴方程：
$J_z \dot{\omega}_z - (J_x - J_y)\omega_x\omega_y = M_z$

3.2 角加速度显式表达式

直接求解角加速度分量：

$\dot{\omega}_x = \frac{M_x + (J_y - J_z)\omega_y\omega_z}{J_x}$

$\dot{\omega}_y = \frac{M_y + (J_z - J_x)\omega_z\omega_x}{J_y}$

$\dot{\omega}_z = \frac{M_z + (J_x - J_y)\omega_x\omega_y}{J_z}$

3.3 考虑飞轮系统的总力矩

总力矩包含外部力矩和飞轮控制力矩：
$\mathbf{M}_{\text{total}} = \mathbf{M}_{\text{ext}} - \mathbf{M}_{\text{wheel}} - \boldsymbol{\omega} \times \mathbf{H}_{\text{wheel}}$

其中飞轮角动量 $Hwheel=[J1Ω1,J2Ω2,J3Ω3]T\mathbf{H}_{\text{wheel}} = [J_1\Omega_1, J_2\Omega_2, J_3\Omega_3]^T$

4. 物理意义分析

4.1 各项物理含义

$M_i/J_i$ ：直接加速度项，力矩与惯量的比值
$(Jj−Jk)ωjωk/Ji(J_j - J_k)\omega_j\omega_k/J_i$ ：陀螺耦合项，由惯量差和角速度耦合产生
$ω×Hwheel\boldsymbol{\omega} \times \mathbf{H}_{\text{wheel}}$ ：飞轮陀螺效应，飞轮角动量与星体旋转的相互作用

4.2 解耦特性分析

由于构型正交：

各轴角加速度主要取决于本轴力矩
交叉耦合仅通过陀螺项体现
控制系统设计可近似为三个单轴系统

5. 数值计算考虑

5.1 计算步骤

计算总外力矩（环境力矩+控制力矩）
计算飞轮陀螺力矩
计算陀螺耦合项
求解角加速度

5.2 数值稳定性

小角速度时陀螺项可忽略
惯量差异大时耦合效应显著
需要实时更新角速度值

6. C语言实现实例

#include <stdio.h>
#include <math.h>

// 卫星参数结构体
typedef struct {
    double J[3];        // 主惯量 [Jx, Jy, Jz]
    double w[3];        // 角速度 [wx, wy, wz]
    double J_wheel[3];  // 飞轮惯量
} SatelliteParams;

// 力矩结构体
typedef struct {
    double M_ext[3];    // 外部力矩
    double M_wheel[3];  // 飞轮控制力矩
    double Omega[3];    // 飞轮转速
} TorqueInput;

// 计算飞轮角动量
void calculate_wheel_angular_momentum(const double J_wheel[3], 
                                    const double Omega[3],
                                    double H_wheel[3]) {
    for (int i = 0; i < 3; i++) {
        H_wheel[i] = J_wheel[i] * Omega[i];
    }
}

// 计算向量叉乘
void vector_cross(const double a[3], const double b[3], double result[3]) {
    result[0] = a[1] * b[2] - a[2] * b[1];
    result[1] = a[2] * b[0] - a[0] * b[2];
    result[2] = a[0] * b[1] - a[1] * b[0];
}

// 计算角加速度主函数
void calculate_angular_acceleration(const SatelliteParams* sat,
                                  const TorqueInput* input,
                                  double dw[3]) {
    double H_wheel[3];          // 飞轮角动量
    double gyro_term[3];        // 飞轮陀螺项
    double inertia_coupling[3]; // 惯量耦合项
    double total_torque[3];     // 总力矩
    
    // 1. 计算飞轮角动量
    calculate_wheel_angular_momentum(sat->J_wheel, input->Omega, H_wheel);
    
    // 2. 计算飞轮陀螺项 ω × H_wheel
    vector_cross(sat->w, H_wheel, gyro_term);
    
    // 3. 计算惯量耦合项
    inertia_coupling[0] = (sat->J[1] - sat->J[2]) * sat->w[1] * sat->w[2];
    inertia_coupling[1] = (sat->J[2] - sat->J[0]) * sat->w[2] * sat->w[0];
    inertia_coupling[2] = (sat->J[0] - sat->J[1]) * sat->w[0] * sat->w[1];
    
    // 4. 计算总力矩
    for (int i = 0; i < 3; i++) {
        total_torque[i] = input->M_ext[i] - input->M_wheel[i] - gyro_term[i];
    }
    
    // 5. 计算角加速度
    dw[0] = (total_torque[0] + inertia_coupling[0]) / sat->J[0];
    dw[1] = (total_torque[1] + inertia_coupling[1]) / sat->J[1];
    dw[2] = (total_torque[2] + inertia_coupling[2]) / sat->J[2];
}

// 简化版本（忽略飞轮陀螺效应）
void calculate_angular_acceleration_simple(const SatelliteParams* sat,
                                         const TorqueInput* input,
                                         double dw[3]) {
    double inertia_coupling[3];
    
    // 计算惯量耦合项
    inertia_coupling[0] = (sat->J[1] - sat->J[2]) * sat->w[1] * sat->w[2];
    inertia_coupling[1] = (sat->J[2] - sat->J[0]) * sat->w[2] * sat->w[0];
    inertia_coupling[2] = (sat->J[0] - sat->J[1]) * sat->w[0] * sat->w[1];
    
    // 计算总力矩（忽略飞轮陀螺项）
    double total_torque[3];
    for (int i = 0; i < 3; i++) {
        total_torque[i] = input->M_ext[i] - input->M_wheel[i];
    }
    
    // 计算角加速度
    dw[0] = (total_torque[0] + inertia_coupling[0]) / sat->J[0];
    dw[1] = (total_torque[1] + inertia_coupling[1]) / sat->J[1];
    dw[2] = (total_torque[2] + inertia_coupling[2]) / sat->J[2];
}

// 打印向量
void print_vector(const char* name, const double vec[3]) {
    printf("%s: [%10.6f, %10.6f, %10.6f]\n", name, vec[0], vec[1], vec[2]);
}

int main() {
    // 初始化卫星参数
    SatelliteParams sat = {
        .J = {10.0, 8.0, 6.0},           // 主惯量 [kg·m²]
        .w = {0.1, 0.05, 0.02},          // 角速度 [rad/s]
        .J_wheel = {0.02, 0.02, 0.02}    // 飞轮惯量 [kg·m²]
    };
    
    // 初始化力矩输入
    TorqueInput input = {
        .M_ext = {0.5, 0.3, 0.2},        // 外部力矩 [N·m]
        .M_wheel = {0.1, 0.05, 0.02},    // 飞轮力矩 [N·m]
        .Omega = {100.0, 50.0, 25.0}     // 飞轮转速 [rad/s]
    };
    
    double dw_full[3];    // 完整模型角加速度
    double dw_simple[3];  // 简化模型角加速度
    
    // 计算角加速度
    calculate_angular_acceleration(&sat, &input, dw_full);
    calculate_angular_acceleration_simple(&sat, &input, dw_simple);
    
    // 输出结果
    printf("三轴正交飞轮构型角加速度计算\n");
    printf("============================\n\n");
    
    printf("卫星参数:\n");
    print_vector("主惯量 J", sat.J);
    print_vector("角速度 ω", sat.w);
    print_vector("飞轮惯量", sat.J_wheel);
    printf("\n");
    
    printf("输入力矩:\n");
    print_vector("外部力矩 M_ext", input.M_ext);
    print_vector("飞轮力矩 M_wheel", input.M_wheel);
    print_vector("飞轮转速 Ω", input.Omega);
    printf("\n");
    
    printf("计算结果:\n");
    print_vector("完整模型角加速度", dw_full);
    print_vector("简化模型角加速度", dw_simple);
    printf("\n");
    
    // 计算差异
    double diff[3];
    for (int i = 0; i < 3; i++) {
        diff[i] = dw_full[i] - dw_simple[i];
    }
    print_vector("模型差异 Δdw", diff);
    
    return 0;
}

7. 程序输出示例

三轴正交飞轮构型角加速度计算
============================

卫星参数:
主惯量 J: [10.000000,  8.000000,  6.000000]
角速度 ω: [ 0.100000,  0.050000,  0.020000]
飞轮惯量: [ 0.020000,  0.020000,  0.020000]

输入力矩:
外部力矩 M_ext: [ 0.500000,  0.300000,  0.200000]
飞轮力矩 M_wheel: [ 0.100000,  0.050000,  0.020000]
飞轮转速 Ω: [100.000000, 50.000000, 25.000000]

计算结果:
完整模型角加速度: [ 0.041800,  0.036250,  0.033333]
简化模型角加速度: [ 0.040000,  0.031250,  0.030000]

模型差异 Δdw: [ 0.001800,  0.005000,  0.003333]