Skip to content

Plant Physics: Constitutive Force Algebra

> Aligned with PCR Master Blueprint v1.0 — see Blueprint §2.4, §4.1, §4.3. > 职责:在 BodyRWS Monad 内,根据 BodyEnv(已被上游 prober 装配好的局部 Reader)与外部传入的设备激励(EngineEffect/DynInFrame),产出 Forces 贡献。是 monadic 算子族,不是纯函数族。 > C-Distillation noteBodyRWS<Forces> 是一层 monad 包装;蒸馏阶段会被 inline 成 void compute_*(BodyEnv*, Forces* out) 风格的 C 函数。Monad 不会留到 RTOS。


1. 在双层 RWS 中的位置

mermaid
graph LR
    SimProbe[sim::probe::*<br/>看不见的上游] -->|装配 BodyEnv| Reader
    Reader[BodyEnv<br/>= WorldEnv& + BodyAsset& + traj/aero/mass Ctxs]
    AvIn[avionics 上游<br/>std::vector EngineEffect<br/>DynInFrame] -->|外部参数| Phys
    Reader -->|body_ask 取| Phys
    Phys[plant::physics::compute_*_contribution<br/>BodyRWS Forces]
    Phys -->|Monoid +| Total[Forces 累加]
    Total --> ODE[dynamics::ode::rk4_step]

两条输入通道

  1. 隐式 Reader 通道BodyEnv 通过 body_ask() 在 monad 内取。算子不直接接 BodyEnv 参数。
  2. 显式参数通道:来自 avionics 的当周期激励(EngineEffect[]FinDeflection[]DynInFrame),作为外部参数传入。

physics 函数本身不感知 prober 在哪、不感知 WorldEnv 何时被装配。它只面对 BodyEnv 这个干净的 Reader。


2. 实际函数签名(与代码对齐)

cpp
// src/plant/physics/Thrust.h
namespace plant {
dynamics::BodyRWS<dynamics::Forces>
compute_thrust_contribution(const std::vector<EngineEffect>& initial_effects);
}

// src/plant/physics/Drag.h
namespace plant {
dynamics::BodyRWS<dynamics::Forces>
compute_aero_contribution(const DynInFrame& dyn_in);
}

// src/plant/physics/Gravity.h
namespace plant {
dynamics::BodyRWS<dynamics::Forces>
compute_gravity_contribution();
}

签名约定:

  • 返回类型一律 dynamics::BodyRWS<dynamics::Forces>
  • BodyEnv& 参数(通过 body_ask 在 monad 内取)
  • RocketBody& 参数(RWS 的 State 部分隐式持有)
  • 只接 avionics 上游的"外部激励"(vector / frame)

2.1 BodyRWS<A> 是什么

cpp
// 在 src/dynamics/DynMonad.h
template<typename A>
using BodyRWS = monad::RWS<BodyEnv, DynLog, plant::RocketBody, A>;
//                          ^Reader  ^Writer  ^State            ^Result

三件事被同时携带:

  • Reader = BodyEnv:通过 body_ask()
  • Writer = DynLog:通过 body_tell(log) 累加(Monoid)
  • State = RocketBody:通过 body_modify([](body){…}) 演化
  • Result = Forces:通过 body_pure(value) 包装

physics 算子最常用的两个原语:

  • body_ask():取 Reader(BodyEnv)
  • body_pure(value):把裸值升入 monad

需要时还可用:

  • body_tell(log):在 monad 内写日志,不破坏纯度
  • body_modify(λ):演化 RocketBody.engines[i] 的 EngineForceState 等设备态

3. 标准实现模式

所有 compute_*_contribution 都遵循同一骨架:

cpp
dynamics::BodyRWS<dynamics::Forces> compute_X_contribution(/* 外部参数 */) {
    return dynamics::body_ask() >>= [/* 捕获参数 */](const dynamics::BodyEnv& env) {
        dynamics::Forces f = dynamics::Forces::zero();

        // 1. 从 env 取已 probe 好的局部真值
        //    env.traj.*       客观位置/速度真值(TrajectoryCtx)
        //    env.aero.*       流体力学真值(AeroCtx)
        //    env.mass.*       质量瞬态真值(MassPropsCtx)
        //    env.asset.*      本物体的本构资产(BodyAsset 引用)
        //    env.world.*      全局常量与场(WorldEnv 引用)

        // 2. 结合外部参数计算 Forces 各分量
        f.total_force_b   = /* ... */;
        f.total_moment_b  = /* ... */;
        f.gravity_lic     = /* ... */;   // 重力专走 LIC 通道
        f.d_mass_dt       = /* ... */;   // 推力相关:质量流

        return dynamics::body_pure(f);
    };
}

关键点body_ask() >>= λ 把"取 Reader + 在 lambda 里用"合成一个新的 BodyRWS,不会破坏 monadic 性质。


4. Forces 实际结构

cpp
// src/dynamics/state/Forces.h
struct Forces {
    Vec3 total_force_b;     // 体坐标合力
    Vec3 total_moment_b;    // 体坐标合力矩
    Vec3 gravity_lic;       // LIC 系重力(特殊通道,不并入 total_force_b)

    double d_mass_dt;       // 总质量流率
    double d_fuel_dt;       // 燃料质量流率
    double d_oxidizer_dt;   // 氧化剂质量流率

    static Forces zero();
    Forces operator+(const Forces& other) const;   // Monoid
};

注意两点偏离教科书的设计:

  1. 重力走独立通道 gravity_lic:因为重力是 LIC 系矢量,与体坐标分离表达更方便积分器后续做 LIC→BODY 转换。不是把重力旋转到体坐标后加进 total_force_b
  2. 质量流是 Forces 的一部分d_mass_dt / d_fuel_dt / d_oxidizer_dt 在 monoid 加法下也自然相加。Forces 不仅承载"力",还承载"质量演化"——这是把动量与质量方程合一处理。

5. 实例:推力

来自 src/plant/physics/Thrust.cpp

cpp
dynamics::BodyRWS<dynamics::Forces>
compute_thrust_contribution(const std::vector<EngineEffect>& initial_effects)
{
    return dynamics::body_ask() >>= [initial_effects](const dynamics::BodyEnv& env) {
        dynamics::Forces total = dynamics::Forces::zero();
        double p_amb = env.aero.static_pressure;          // 来自 AeroCtx

        for (const auto& eff : initial_effects) {
            // 1. 从 BodyAsset 查 nozzle 几何
            double area = 0.65;
            Vec3 mount_pos;
            for (const auto& asset_eng : env.asset.engines) {
                if (asset_eng.slot_id == (int)eff.engine_id) {
                    area      = asset_eng.nozzle_area_m2;
                    mount_pos = asset_eng.mount_pos_body;
                    break;
                }
            }

            // 2. 背压补偿:real = vacuum − p_amb · Sa
            double real_mag = eff.vacuum_thrust - p_amb * area;

            // 3. 喷管偏摆 → 体坐标力
            Vec3 f_body = thrust_in_nozzle_frame(real_mag,
                              eff.nozzle_yaw.rad(),
                              eff.nozzle_pitch.rad());

            // 4. 累加力 + 力矩 + 质量流
            total.total_force_b  += f_body;
            total.total_moment_b += mount_pos.cross(f_body);
            total.d_mass_dt      += (eff.mass_flow_fuel + eff.mass_flow_ox);
            total.d_fuel_dt      += eff.mass_flow_fuel;
            total.d_oxidizer_dt  += eff.mass_flow_ox;
        }

        return dynamics::body_pure(total);
    };
}

注意:

  • 外部参数 initial_effects 来自 avionics 的 DynInFrame.engine_effects
  • Reader 通道 env.aero.static_pressure(来自 prober)+ env.asset.engines[](来自 BodyAsset 配置)
  • 这里没有直接访问 env.world.atmosphere——降维原则成立

未来如果要把 EngineForceState 持久化进 RocketBody.engines[i],可以把 body_pure(total) 替换为 body_modify(...) >> body_pure(total)——monad 留好了通道。


6. 实例:气动力

来自 src/plant/physics/Drag.cpp

cpp
dynamics::BodyRWS<dynamics::Forces>
compute_aero_contribution(const DynInFrame& dyn_in)
{
    (void)dyn_in;   // 当前未消费 fin_deflections,签名预留
    return dynamics::body_ask() >>= [](const dynamics::BodyEnv& env) {
        dynamics::Forces f = dynamics::Forces::zero();

        double mach     = env.aero.mach_number;      // AeroCtx
        double q        = env.aero.dynamic_pressure; // AeroCtx
        double ref_area = env.asset.aero.reference_area;  // BodyAsset

        // 上升/下降段判定(用 TrajectoryCtx 的 ECF 位置)
        bool is_ascending = (env.traj.ecf_pos.magnitude() < 6400000.0);

        double cd = calc_drag_coeff(mach, is_ascending);
        double drag_mag = q * ref_area * cd;

        // 来流方向已在 prober 里转到 BODY 系(AeroCtx.flow_in_body)
        Vec3 flow_dir = Vec3(env.aero.flow_in_body.x(),
                             env.aero.flow_in_body.y(),
                             env.aero.flow_in_body.z()).normalized();

        f.total_force_b = flow_dir * drag_mag;
        return dynamics::body_pure(f);
    };
}

注意:

  • AeroCtx.flow_in_body 已是 BODY 系矢量 —— prober 完成了 LIC→BODY 旋转,physics 不再做坐标变换
  • 攻角/侧滑/总攻角/气动滚转 4 个角 都在 AeroCtx.angles,物理算子可直接消费
  • (void)dyn_in; 是惯例:签名预留 fin_deflections 通道,当前实现先吃 BodyEnv 即可

7. 实例:重力

来自 src/plant/physics/Gravity.cpp

cpp
dynamics::BodyRWS<dynamics::Forces> compute_gravity_contribution() {
    return dynamics::body_ask() >>= [](const dynamics::BodyEnv& env) {
        dynamics::Forces f = dynamics::Forces::zero();

        double mu    = env.world.constants.mu;        // 全局常量
        double r_mag = env.traj.ecf_pos.magnitude();  // TrajectoryCtx

        Vec3 pos_ecf(env.traj.ecf_pos.x(),
                     env.traj.ecf_pos.y(),
                     env.traj.ecf_pos.z());
        Vec3 g_ecf = pos_ecf.normalized() * (-mu / (r_mag * r_mag));

        // 走专用通道 gravity_lic(后续由 integrator 旋转到所需 frame)
        f.gravity_lic = g_ecf;
        return dynamics::body_pure(f);
    };
}

注意:

  • 无外部参数 —— 重力完全由 (位置 + 重力常数) 决定
  • Forces.gravity_lic 通道 —— 不并入 total_force_b
  • ECF 到 LIC 的精确旋转目前为 TODO(代码注释里有说明),未来通过 env.world.frame 提供

8. Forces Monoid 与累加

物理算子族产出的多个 BodyRWS<Forces> 由调用方组合:

cpp
// 调用方(属 simulation/pipeline/body_tick.cpp)
auto thrust  = plant::compute_thrust_contribution(dyn_in.engine_effects);
auto aero    = plant::compute_aero_contribution(dyn_in);
auto gravity = plant::compute_gravity_contribution();

// 在 monad 内累加(Forces 的 Monoid + 满足结合律)
auto total = (thrust  >>= [](Forces f1) {
return (aero    >>= [f1](Forces f2) {
return (gravity >>= [f1, f2](Forces f3) {
    return dynamics::body_pure(f1 + f2 + f3);
}); }); });

或用 fluent operator(代码中尚未实现,但 dynamics/AGENTS.md 提到惯用约定):

cpp
auto total = Forces::zero()
             >> compute_thrust_contribution(...)
             >> compute_aero_contribution(...)
             >> compute_gravity_contribution();

详见 05_Dynamics_Core/Forces_Monoid.md


9. 设计原则

原则说明
Monadic 返回所有 compute_*_contribution 返回 BodyRWS<Forces>,不返回裸 Forces
隐式 Reader通过 body_ask()BodyEnv,不在签名里出现
显式激励参数avionics 上游的 EngineEffect[] / DynInFrame 作为外部参数
不感知 proberphysics 看到的就是已经 ctx 化的 BodyEnv;prober 在哪、何时执行不可见
Forces 含 mass flow推力贡献也产出 d_mass_dt,与力一起在 Monoid 加法下传播
重力走独立通道LIC 系重力在 Forces.gravity_lic,不污染 total_force_b
坐标变换前置LIC→BODY 等坐标旋转应在 prober 内完成;physics 算子尽量在单一 frame 内工作

10. 加新力源的步骤

例如想加"太阳辐射压力":

  1. src/plant/physics/SolarRadiation.h/cpp 新建:
    cpp
    dynamics::BodyRWS<dynamics::Forces> compute_solar_contribution();
  2. 若需要新的 Ctx 字段(如太阳方向矢量),不直接读 environment:而是让 sim::probe 加一个相应的 ctx 输出(如 SolarCtx solar 加进 BodyEnv),由 prober 装配
  3. 若需要新的全局场(太阳星历),在 environment/SolarField,由 prober 读取
  4. 在调用方(sim::body_tick)的 Forces 累加链里加上 compute_solar_contribution()

plant/physics/ 不应直接 #include <environment/...>。新场必须先经 prober 落到 BodyEnv,再被 physics 消费。


11. 反模式

反模式为什么不行
compute_thrust_contribution(BodyEnv&, ...) 显式接 BodyEnv破坏 monadic 约定;BodyEnv 应通过 body_ask
返回裸 Forces 而非 BodyRWS<Forces>失去 Writer/State 通道,未来加日志/演化设备态时被迫改签名
在 lambda 里 env.world.atmosphere.pressure_at(...)越过降维:应消费 env.aero.static_pressure(prober 已算好)
在 lambda 里做 LIC→BODY 旋转坐标变换应在 prober;physics 用 AeroCtx.flow_in_body 即可
把重力加到 total_force_b违反 Forces 结构约定;重力专走 gravity_lic 通道
在 physics 算子里调用 body_modify 改 spatial越权;spatial 演化只能由 dynamics_core 的积分器执行
在 lambda 里持有静态缓存或全局 mutable破坏纯度,多线程 body 并行执行时不安全

12. 测试策略

每个 compute_*_contribution 适合 L1 bench 测试:

cpp
TEST(Thrust, VacuumThrust) {
    // 构造 fixture: WorldEnv + BodyAsset + 一组 EngineEffect
    auto env = make_fixture_body_env(...);
    std::vector<EngineEffect> effects = { /* 单台 RD180 真空标称 */ };

    auto rws = compute_thrust_contribution(effects);
    auto [forces, log, body_out] = rws.run(env, /* RocketBody fixture */);

    EXPECT_NEAR(forces.total_force_b.x(), 4'000'000.0, 1e3);  // 4 MN
    EXPECT_GT(forces.d_mass_dt, 0.0);
}

prober 的测试在另一个 bench:给定 fixture WorldState + RocketBody,验证 TrajectoryCtx / AeroCtx / MassPropsCtx 三个产物字段是否符合预期。两条测试链互不耦合。

详见 08_Cross_Cutting/Testing_Framework.md


13. Cross References

  • BodyRWS / body_ask / body_pure 等 monad 原语 → 01_Foundation/Monad_Toolkit.md
  • BodyEnv 字段全集与 Ctx 三分(TrajectoryCtx/AeroCtx/MassPropsCtx)→ 06_Simulation/Body_World_Tick.md §3(待写)
  • prober 函数签名与装配责任 → 06_Simulation/Body_World_Tick.md §4(待写)
  • 资产 spec → 02_Physical_World/Plant_Model_Assets.md
  • Forces Monoid 完整代数 → 05_Dynamics_Core/Forces_Monoid.md(待写)
  • 完整推力流五阶段 → Blueprint §4.1