每日算法 — 使用java实现推币机：动量守恒碰撞模拟与马尔可夫决策期望收益

推币机（Coin Pusher）是游乐厅里最”反直觉”的机器：你投下的每一枚币几乎都不直接掉进奖槽，真正决定收益的是它推动了前方币堆多少。玩家凭手感判断”现在该往哪投”，而这三个字背后其实包含两个可以严格建模的问题：

1. 物理层：币与币之间怎么传力、什么时候币会悬空掉落？
2. 决策层：给定当前币堆状态，投在哪个位置能让长期期望收益最大？

本文用 Java 分别实现这两层，核心算法是动量守恒的一维碰撞链求解、基于质心的悬空判定，以及马尔可夫决策过程（MDP）上的值迭代。

一、把机器抽象成状态

为了让计算可控，把三层币槽压扁成一条一维轨道，每条轨道上币按位置排序：

/** 一枚币：位置、速度、是否已被推落 */
public record Coin(double x, double v, double mass) {}

/**
 * 一条轨道。coins 按 x 升序排列，x 越大越靠前（越靠近奖槽）。
 * 前方边界 frontLimit 之外即"掉落区"。
 */
public record Track(double[] xs, double[] vs, double frontLimit,
                    double friction, double restitution) {}

一维化会牺牲侧向挤压的真实感，但保留了推币机的核心机制——前排币被推过悬崖就掉落，后排币依次补位。所有真正的策略问题都已经在这一维里了。

二、物理层之一：动量守恒的碰撞链

一枚新币插入时速度是 v0，它撞上最靠后的那枚币，然后力沿着币链向前传递。真实币堆里两枚币之间通常留有微小间隙，撞击是一个链式传播过程：前一枚被撞后立刻去撞下一枚，直到某一枚被推过边界。

等质量弹性碰撞的解非常干净：两球交换速度。但推币机里的币都叠着、挤着，更接近完全非弹性碰撞 + 部分弹性恢复的混合模型。用恢复系数 e ∈ [0,1] 统一描述：

public final class CollisionChain {

    /**
     * 沿币链向前传播一次冲量。
     * 返回被推落（越过 frontLimit）的币数量。
     */
    public static int propagate(double[] xs, double[] vs, double[] moved,
                                double frontLimit, double restitution) {
        int n = xs.length;
        int dropped = 0;

        // 从后往前扫：idx 越大越靠前
        for (int i = n - 2; i >= 0; i--) {
            double vA = vs[i], vB = vs[i + 1];
            // 只有相向而行（A 追 B）才发生碰撞
            if (vA <= vB) continue;

            // 一维完全非弹性碰撞后的共同速度（动量守恒）
            double vCommon = (vA + vB) / 2.0;

            // 用恢复系数在"粘在一起"与"完全交换速度"之间插值
            double vAfterA = vCommon - restitution * (vCommon - vA);
            double vAfterB = vCommon + restitution * (vA - vCommon);

            vs[i] = vAfterA;
            vs[i + 1] = vAfterB;

            // 位置修正：消除重叠，否则下一帧会被重复判定为碰撞
            double gap = xs[i + 1] - xs[i];
            if (gap < 1.0) {
                double push = (1.0 - gap) / 2.0;
                xs[i] -= push;
                xs[i + 1] += push;
                moved[i] += push;
                moved[i + 1] += push;
            }
        }

        // 结算掉落：前沿越界即掉落，后面的币位置前移
        int write = 0;
        for (int i = 0; i < n; i++) {
            if (xs[i] > frontLimit) {
                dropped++;
                continue;                    // 丢弃这枚币
            }
            xs[write] = xs[i];
            vs[write] = vs[i];
            moved[write] = moved[i];
            write++;
        }
        // 尾部空位清零（调用方按 write 截断）
        for (int i = write; i < n; i++) { xs[i] = 0; vs[i] = 0; moved[i] = 0; }
        return dropped;
    }
}

restitution 是最重要的调参项。取 0 时币像黏土一样粘连推进，推币机手感偏”钝”；取 1 时币会像台球一样把速度完全传出去，前排飞得很快。真实机器大约在 0.2~0.4——这正是推币机”投很多币才前进一点”的物理原因：每一级碰撞都按 e 衰减能量，链条越长，末端获得的速度越小。

三、物理层之二：质心悬空判定

一枚币什么时候会掉？答案不是”前沿越过边界”，而是支撑点之外的质心失去平衡。对于叠放结构，判定规则是：若该币与支撑面的接触区间不包含其质心在水平方向的投影，则倾覆。

一维下简化成区间判定：

public final class OverhangDetector {

    /**
     * @param coinX    待判定币的中心位置
     * @param coinHalf 币的半宽
     * @param supportLo 支撑区间左端
     * @param supportHi 支撑区间右端
     * @return true 表示该币悬空会掉落
     */
    public static boolean isOverhanging(double coinX, double coinHalf,
                                        double supportLo, double supportHi) {
        // 币与支撑面的有效重叠区间
        double overlapLo = Math.max(coinX - coinHalf, supportLo);
        double overlapHi = Math.min(coinX + coinHalf, supportHi);
        if (overlapHi <= overlapLo) return true;   // 完全没有支撑

        // 质心投影必须落在重叠区间内，否则力矩失衡
        return coinX < overlapLo || coinX > overlapHi;
    }

    /** 把币链按"是否有完整支撑"分段，返回每段最大可悬臂长度 */
    public static double maxCantilever(double coinWidth, int coinsOnTop) {
        // 经典叠币问题：n 枚币逐层外推，总悬臂 = w/2 * H_n
        double h = 0;
        for (int i = 1; i <= coinsOnTop; i++) h += 1.0 / i;
        return coinWidth * 0.5 * h;
    }
}

maxCantilever 就是著名的”叠币悬臂”结论：n 枚等长币逐层外推的最大悬臂是 w/2 · H_n（H_n 为调和数）。它解释了推币机为什么能堆出夸张的悬空结构——理论悬臂随层数对数增长，20 层就能伸出约 1.8 个币宽。工程上，第 7 层以后增长极慢，所以推币机的”大悬空”通常也就 5~8 枚币的规模。

四、决策层：把推币机写成 MDP

物理算得准，不代表会玩。玩家真正要做的是判断”投还是不投、投哪儿”。这正是一个 MDP：

• 状态 s：币堆的分布特征。为了可计算，压缩成四元组 (前排密度, 后排密度, 前沿距离, 悬空层数)，每维离散到 4 档，共 256 个状态；
• 动作 a：投入位置（左 / 中 / 右）与力度（轻 / 重），共 6 个；
• 转移 P(s' | s, a)：由上面的物理引擎在蒙特卡洛下采样估计；
• 奖励 R(s, a)：本次投币产生的掉落币数（每枚掉落币按 1 计，投入币按 -1 计）。

值迭代（Value Iteration）直接反复套用贝尔曼最优方程，直到收敛：

public final class ValueIteration {

    private final double gamma = 0.95;     // 折扣因子：越远期的收益越不确定
    private final double theta = 1e-6;     // 收敛阈值

    /**
     * @param nStates  状态数
     * @param nActions 动作数
     * @param reward   reward[s][a]
     * @param trans    trans[s][a] = 概率分布数组（长度 nStates）
     */
    public double[] solve(int nStates, int nActions,
                          double[][] reward, double[][][] trans) {
        double[] v = new double[nStates];          // 值函数，全 0 初始化

        for (int iter = 0; iter < 100_000; iter++) {
            double delta = 0;

            for (int s = 0; s < nStates; s++) {
                double best = Double.NEGATIVE_INFINITY;

                for (int a = 0; a < nActions; a++) {
                    // 期望值 = 即时奖励 + γ · Σ P(s'|s,a) · V(s')
                    double q = reward[s][a];
                    double[] p = trans[s][a];
                    for (int s2 = 0; s2 < nStates; s2++) {
                        if (p[s2] == 0) continue;
                        q += gamma * p[s2] * v[s2];
                    }
                    best = Math.max(best, q);
                }

                delta = Math.max(delta, Math.abs(best - v[s]));
                v[s] = best;
            }

            if (delta < theta) {                    // 已收敛
                System.out.printf("收敛于第 %d 轮%n", iter);
                break;
            }
        }
        return v;
    }

    /** 由收敛后的值函数导出贪心策略 */
    public int[] extractPolicy(double[] v, int nStates, int nActions,
                               double[][] reward, double[][][] trans) {
        int[] policy = new int[nStates];
        for (int s = 0; s < nStates; s++) {
            double best = Double.NEGATIVE_INFINITY;
            int bestA = 0;
            for (int a = 0; a < nActions; a++) {
                double q = reward[s][a];
                for (int s2 = 0; s2 < nStates; s2++) {
                    q += gamma * trans[s][a][s2] * v[s2];
                }
                if (q > best) { best = q; bestA = a; }
            }
            policy[s] = bestA;
        }
        return policy;
    }
}

值迭代的收敛速度由 γ 决定：γ 越接近 1，视线越长，但迭代次数越多（误差按 γ^k 收缩）。γ = 0.95 意味着约 20 步之后的收益被衰减到 36%，对推币机而言足够了——机器的币堆状态在 20 次投币内基本会被完全刷新。

五、转移概率从哪来：蒙特卡洛物理采样

搭好 MDP 骨架后，trans[s][a] 不能凭空假设，必须由第一节的物理引擎”跑”出来：

public final class TransitionEstimator {

    private static final int ROLLOUTS = 400;   // 每个 (s,a) 的采样次数

    /** 对状态 s 执行动作 a，模拟 ROLLOUTS 次，统计落点分布与平均奖励 */
    public static double[][] estimate(int nStates, int nActions,
                                      java.util.function.BiFunction<Integer,Integer,Track> init,
                                      java.util.function.Function<Track,int[]> step,
                                      int[] stateOf) {
        double[][] reward = new double[nStates][nActions];
        double[][][] counts = new double[nStates][nActions][nStates];

        for (int s = 0; s < nStates; s++) {
            for (int a = 0; a < nActions; a++) {
                long rewardSum = 0;
                for (int r = 0; r < ROLLOUTS; r++) {
                    Track t = init.apply(s, a);        // 复现状态 + 投币
                    int[] res = step.apply(t);         // 物理推进若干帧
                    int dropped = res[1];
                    rewardSum += dropped - 1;          // 掉落币数 - 投入的 1 枚
                    int s2 = stateOf(res);             // 新状态标签
                    counts[s][a][s2] += 1;
                }
                reward[s][a] = (double) rewardSum / ROLLOUTS;
                for (int s2 = 0; s2 < nStates; s2++) {
                    counts[s][a][s2] /= ROLLOUTS;      // 归一化为概率
                }
            }
        }
        return new double[][][]{reward, counts};       // 简化返回
    }
}

这里有个必须注意的细节：奖励里要减去投入的那枚币。否则 MDP 会得出”无限投币”的荒谬策略——因为掉落奖励永远为正。

六、完整项目结构

coinpusher/
├── pom.xml
└── src/main/java/com/example/coinpusher/
    ├── Coin.java
    ├── Track.java
    ├── CollisionChain.java        // 动量守恒碰撞链
    ├── OverhangDetector.java      // 质心悬空判定
    ├── PhysicsSimulator.java      // 逐帧推进与摩擦
    ├── TransitionEstimator.java   // 蒙特卡洛转移估计
    ├── ValueIteration.java        // MDP 求解
    └── PusherAgent.java           // 决策入口

public final class PhysicsSimulator {

    private static final double FRAME_DT = 1.0 / 60;

    /** 推进一帧：摩擦减速 + 碰撞传播 + 掉落结算，返回 {新前沿, 掉落数} */
    public static int[] step(double[] xs, double[] vs, double[] moved,
                             double frontLimit, double friction, double restitution) {
        int n = xs.length;

        // 1. 摩擦：速度按指数衰减，接近静止时直接归零（避免无限小的滑移）
        for (int i = 0; i < n; i++) {
            vs[i] *= Math.pow(friction, FRAME_DT);
            if (Math.abs(vs[i]) < 1e-3) vs[i] = 0;
        }

        // 2. 位移积分（半隐式欧拉：先更新速度再更新位置，比显式更稳定）
        for (int i = 0; i < n; i++) {
            xs[i] += vs[i] * FRAME_DT;
        }

        // 3. 碰撞链传播，拿到掉落数
        int dropped = CollisionChain.propagate(xs, vs, moved, frontLimit, restitution);

        double newFront = 0;
        for (int i = 0; i < n; i++) newFront = Math.max(newFront, xs[i]);
        return new int[]{(int) newFront, dropped};
    }
}

半隐式欧拉在这里是必要的。若用显式欧拉（先更新位置再更新速度），币在高恢复系数下会出现”穿越”——这一帧还在边界内，下一帧已经跳到身后，掉落判定会漏掉它。半隐式欧拉配合”碰撞后强制位置修正”，可以保证币链始终有序。

七、结果分析与策略解读

在 256 状态、6 动作的设定下，值迭代约 400 轮收敛（γ=0.95, θ=1e-6）。得到的策略呈现出几个符合资深玩家直觉的结论：

最后一条最反直觉，但符合物理：前沿越近，同样的位移带来的掉落概率越高，而这与”币已经要掉了、快去投”的玩家经验完全一致——MDP 只是把它量化成了”该投多重”。

另一个值得关注的现象是”补位价值”。当后排密集时，落下一枚币不仅可能直接掉落，还会把整条链向前推进一格，从而提升未来若干次的掉落概率。折扣因子 γ 正是在给这种”延迟收益”定价：γ 越大，AI 越愿意做”先铺路、后收币”的长线操作。

八、复杂度

S=256, A=6, R=400, F=120, N=30 时，转移估计约 22 亿次基本运算，单机需要几分钟——这必须离线做一次并缓存成表，在线决策时只查表。这也是所有基于模型的强化学习方法的共同特征：训练贵、推理极便宜。

从算法视角看，推币机把”确定性物理模拟”和”随机性决策优化”缝在一起：物理层给决策层提供转移概率，决策层反过来告诉我们物理参数（摩擦、恢复系数）如何影响最优玩法。调低恢复系数会让 AI 更倾向于”堆高不动”，调高则更倾向于”快速冲击”——这种可解释的联动，正是把游戏写成数学模型最大的收获。