3. 太阳翼对日定向PID控制器设计

任务定义

本案例要求设计一个用于卫星太阳翼对日定向控制的PID控制器。被控对象为太阳翼驱动机构（SADA），其输入为驱动电机指令角速度，输出为太阳翼法线相对太阳方向的夹角（即对日定向误差角）。控制目标为：在轨道光照角变化及机构摩擦扰动下，使定向误差角收敛至±0.1°以内，且超调量小于15%，稳态无振荡。

被控对象模型

太阳翼驱动机构可等效为一阶惯性环节加积分环节（电机转速环闭环后），传递函数为：

G(s) = K / [s(τs + 1)]

其中K=0.5 (rad/s)/V为增益，τ=0.2s为机电时间常数。控制量u(t)为驱动电压（V），被控量为太阳翼转角θ(t)（rad），扰动d(t)为轨道光照角变化率（rad/s），等效为输入端的角速度扰动。

PID控制器设计

采用位置式PID控制器，控制律为：

u(t) = Kp·e(t) + Ki·∫e(t)dt + Kd·de(t)/dt

其中e(t) = θ_ref(t) − θ(t)为定向误差角。根据被控对象参数，采用零极点配置法整定：取期望闭环阻尼比ξ=0.8，自然频率ωn=3 rad/s，得到PID参数：Kp=3.6，Ki=2.0，Kd=1.2。控制器输出限幅为±10V，积分器设置抗饱和限幅（输出超过限幅时冻结积分累加）。

代码实现

以下Python代码实现PID控制器与卫星太阳翼对日定向闭环仿真，包含轨道光照角变化输入、机构摩擦扰动、控制器限幅与抗饱和逻辑：


import numpy as np
import matplotlib.pyplot as plt

# 被控对象参数
K = 0.5      # 增益 (rad/s)/V
tau = 0.2    # 时间常数 s
dt = 0.01    # 仿真步长 s
T_total = 30 # 仿真时长 s
N = int(T_total / dt)
t = np.linspace(0, T_total, N)

# PID参数
Kp, Ki, Kd = 3.6, 2.0, 1.2
u_max = 10.0   # 控制量限幅 V
integral = 0.0
prev_error = 0.0

# 状态变量: theta为太阳翼转角, omega为角速度
theta = 0.0
omega = 0.0

# 输入信号: 轨道光照角变化(目标角速度扰动)
# 模拟卫星入影出影过程: 0-10s光照角匀速变化, 10-20s恒定, 20-30s反向变化
def light_angle_rate(t_val):
    if t_val < 10:
        return 0.02   # rad/s, 正向变化
    elif t_val < 20:
        return 0.0    # 恒定
    else:
        return -0.02  # 反向变化

# 摩擦扰动: 库仑+粘性摩擦等效为输入扰动
def friction_torque(omega_val):
    return 0.05 * np.sign(omega_val) + 0.1 * omega_val

# 数据记录
theta_log = np.zeros(N)
error_log = np.zeros(N)
u_log = np.zeros(N)
theta_ref_log = np.zeros(N)

# 仿真主循环
for i in range(N):
    # 目标角度: 光照角积分
    if i == 0:
        theta_ref = 0.0
    else:
        theta_ref = theta_ref_log[i-1] + light_angle_rate(t[i-1]) * dt
    theta_ref_log[i] = theta_ref
    
    # 当前误差
    error = theta_ref - theta
    error_log[i] = error
    
    # PID控制律(含抗饱和)
    integral += error * dt
    u_pid = Kp * error + Ki * integral + Kd * (error - prev_error) / dt
    
    # 限幅
    u = np.clip(u_pid, -u_max, u_max)
    
    # 抗饱和: 若输出饱和则冻结积分
    if np.abs(u_pid) > u_max:
        integral -= error * dt  # 撤销本次积分累加
    
    u_log[i] = u
    prev_error = error
    
    # 被控对象更新(一阶惯性+积分)
    # 等效输入: 电压-摩擦扰动
    d_omega = (K * u - friction_torque(omega) - omega) / tau
    omega += d_omega * dt
    theta += omega * dt
    
    theta_log[i] = theta

# 计算性能指标
error_steady = np.abs(error_log[-500:])  # 最后5s稳态误差
max_error = np.max(np.abs(error_log))
overshoot = 0.0
# 计算超调: 第一次到达目标后的最大偏差
for i in range(1, N):
    if error_log[i-1] > 0 and error_log[i] <= 0:
        overshoot = np.max(np.abs(error_log[i:]))
        break

print(f"最大误差角: {max_error:.4f} rad = {np.degrees(max_error):.2f}°")
print(f"稳态误差(最后5s最大): {np.max(error_steady):.6f} rad = {np.degrees(np.max(error_steady)):.4f}°")
print(f"超调量: {overshoot:.4f} rad = {np.degrees(overshoot):.2f}°")

# 绘图
plt.figure(figsize=(10, 6))
plt.subplot(3, 1, 1)
plt.plot(t, np.degrees(theta_ref_log), 'b--', label='目标角度')
plt.plot(t, np.degrees(theta_log), 'r-', label='实际角度')
plt.ylabel('角度 (°)')
plt.legend()
plt.grid(True)

plt.subplot(3, 1, 2)
plt.plot(t, np.degrees(error_log), 'g-')
plt.ylabel('误差角 (°)')
plt.axhline(y=0.1, color='k', linestyle=':', label='±0.1° 边界')
plt.axhline(y=-0.1, color='k', linestyle=':')
plt.legend()
plt.grid(True)

plt.subplot(3, 1, 3)
plt.plot(t, u_log, 'm-')
plt.ylabel('控制电压 (V)')
plt.xlabel('时间 (s)')
plt.grid(True)
plt.tight_layout()
plt.savefig('pid_solar_tracking.png', dpi=100)
plt.show()

验证方法

运行上述代码，观察输出结果并判定是否满足指标：

跟踪精度验证：检查误差角曲线，在10s和20s光照角变化方向切换后，误差角应在2s内收敛至±0.1°以内，且稳态误差不超限；
稳定性验证：观察整个30s仿真过程，误差角曲线不应出现持续振荡或发散；
超调验证：在10s和20s输入突变后，误差角反向超调量应小于0.26°（对应15%的输入阶跃幅值）；
控制量检查：控制电压始终在±10V限幅范围内，且未出现频繁饱和切换。

判定标准：若打印输出的稳态误差最大值小于0.0017 rad（0.1°），且超调量小于0.26 rad的15%（即0.04 rad），则控制器设计满足要求。若超调过大，可适当增大Kd或减小Kp；若稳态误差不满足，则增大Ki。

运行与结果判定

在Python环境中运行上述代码（需安装numpy和matplotlib），程序将输出最大误差角、稳态误差和超调量三个指标，并生成误差角随时间变化的曲线图。判定标准为：

稳态误差（最后5s内最大值）≤ 0.0017 rad（即0.1°）；
超调量 ≤ 0.04 rad（即输入阶跃幅值的15%）；
误差角曲线在10s和20s扰动切换后2s内进入±0.1°误差带且不振荡。

若三项指标均满足，则PID控制器设计合格。若超调超标，增大Kd至1.5~2.0；若稳态误差不达标，增大Ki至2.5~3.0后重新仿真验证。

📚 课程目录