use std::cmp::Ordering;

use crate::{
    N_INF, N00,
    arena::ScratchArena,
    compute_margin,
    drift::Drift,
    factor::{Factor, trunc::TruncFactor},
    gaussian::Gaussian,
    player::Player,
    tuple_gt, tuple_max,
};

#[derive(Debug)]
pub struct Game<'a, D: Drift> {
    teams: Vec<Vec<Player<D>>>,
    result: &'a [f64],
    weights: &'a [Vec<f64>],
    p_draw: f64,
    pub(crate) likelihoods: Vec<Vec<Gaussian>>,
    pub(crate) evidence: f64,
}

impl<'a, D: Drift> Game<'a, D> {
    pub fn new(
        teams: Vec<Vec<Player<D>>>,
        result: &'a [f64],
        weights: &'a [Vec<f64>],
        p_draw: f64,
        arena: &mut ScratchArena,
    ) -> Self {
        debug_assert!(
            result.len() == teams.len(),
            "result must have the same length as teams"
        );
        debug_assert!(
            weights
                .iter()
                .zip(teams.iter())
                .all(|(w, t)| w.len() == t.len()),
            "weights must have the same dimensions as teams"
        );
        debug_assert!(
            (0.0..1.0).contains(&p_draw),
            "draw probability must be >= 0.0 and < 1.0"
        );
        debug_assert!(
            p_draw > 0.0 || {
                let mut r = result.to_vec();
                r.sort_unstable_by(|a, b| a.partial_cmp(b).unwrap());
                r.windows(2).all(|w| w[0] != w[1])
            },
            "draw must be > 0.0 if there are teams with draw"
        );

        let mut this = Self {
            teams,
            result,
            weights,
            p_draw,
            likelihoods: Vec::new(),
            evidence: 0.0,
        };

        this.likelihoods(arena);
        this
    }

    fn likelihoods(&mut self, arena: &mut ScratchArena) {
        arena.reset();

        let n_teams = self.teams.len();

        // Sort teams by result descending; reuse arena.sort_buf to avoid allocation.
        arena.sort_buf.extend(0..n_teams);
        arena.sort_buf.sort_by(|&i, &j| {
            self.result[j]
                .partial_cmp(&self.result[i])
                .unwrap_or(Ordering::Equal)
        });

        // Team performance priors (TeamSumFactor logic inlined).
        let team_prior: Vec<Gaussian> = arena
            .sort_buf
            .iter()
            .map(|&t| {
                self.teams[t]
                    .iter()
                    .zip(self.weights[t].iter())
                    .fold(N00, |p, (player, &w)| p + (player.performance() * w))
            })
            .collect();

        let n_diffs = n_teams.saturating_sub(1);

        // One TruncFactor per adjacent sorted-team pair; each owns a diff VarId.
        let mut trunc: Vec<TruncFactor> = (0..n_diffs)
            .map(|i| {
                let tie = self.result[arena.sort_buf[i]] == self.result[arena.sort_buf[i + 1]];
                let margin = if self.p_draw == 0.0 {
                    0.0
                } else {
                    let a: f64 = self.teams[arena.sort_buf[i]]
                        .iter()
                        .map(|p| p.beta.powi(2))
                        .sum();
                    let b: f64 = self.teams[arena.sort_buf[i + 1]]
                        .iter()
                        .map(|p| p.beta.powi(2))
                        .sum();
                    compute_margin(self.p_draw, (a + b).sqrt())
                };
                let vid = arena.vars.alloc(N_INF);
                TruncFactor::new(vid, margin, tie)
            })
            .collect();

        // Per-team messages from neighbouring RankDiff factors (replaces TeamMessage).
        let mut lhood_lose: Vec<Gaussian> = vec![N_INF; n_teams];
        let mut lhood_win: Vec<Gaussian> = vec![N_INF; n_teams];

        // Helpers: team marginal incorporating one side of incoming RankDiff messages.
        // post_win(i) = what team i presents to the diff factor on its "winning" side.
        // post_lose(i) = what team i presents to the diff factor on its "losing" side.
        macro_rules! post_win {
            ($i:expr) => {
                team_prior[$i] * lhood_lose[$i]
            };
        }
        macro_rules! post_lose {
            ($i:expr) => {
                team_prior[$i] * lhood_win[$i]
            };
        }

        let mut step = (f64::INFINITY, f64::INFINITY);
        let mut iter = 0;

        while tuple_gt(step, 1e-6) && iter < 10 {
            step = (0.0_f64, 0.0_f64);

            // Forward sweep: diffs 0 .. n_diffs-2 (all but the last).
            for e in 0..n_diffs.saturating_sub(1) {
                let raw = post_win!(e) - post_lose!(e + 1);
                // Set diff var = raw × trunc.msg so that cavity = raw.
                arena.vars.set(trunc[e].diff, raw * trunc[e].msg);
                let d = trunc[e].propagate(&mut arena.vars);
                step = tuple_max(step, d);

                let new_ll = post_win!(e) - trunc[e].msg;
                step = tuple_max(step, lhood_lose[e + 1].delta(new_ll));
                lhood_lose[e + 1] = new_ll;
            }

            // Backward sweep: diffs n_diffs-1 .. 1 (reverse, all but the first).
            for e in (1..n_diffs).rev() {
                let raw = post_win!(e) - post_lose!(e + 1);
                arena.vars.set(trunc[e].diff, raw * trunc[e].msg);
                let d = trunc[e].propagate(&mut arena.vars);
                step = tuple_max(step, d);

                let new_lw = post_lose!(e + 1) + trunc[e].msg;
                step = tuple_max(step, lhood_win[e].delta(new_lw));
                lhood_win[e] = new_lw;
            }

            iter += 1;
        }

        // Special case: exactly 1 diff (2-team game).  The loop body is empty
        // for this case (both ranges are empty), so we run the factor once here.
        if n_diffs == 1 {
            let raw = post_win!(0) - post_lose!(1);
            arena.vars.set(trunc[0].diff, raw * trunc[0].msg);
            trunc[0].propagate(&mut arena.vars);
        }

        // Boundary updates: close the chain at both ends.
        if n_diffs > 0 {
            lhood_win[0] = post_lose!(1) + trunc[0].msg;
            lhood_lose[n_teams - 1] = post_win!(n_teams - 2) - trunc[n_diffs - 1].msg;
        }

        // Evidence = product of per-diff evidences (each cached on first propagation).
        self.evidence = trunc
            .iter()
            .map(|t| t.evidence_cached.unwrap_or(1.0))
            .product();

        // Per-team "likelihood" = product of incoming RankDiff messages.
        let m_t_ft: Vec<Gaussian> = (0..n_teams)
            .map(|si| lhood_win[si] * lhood_lose[si])
            .collect();

        // Inverse permutation: inv[orig_i] = sorted_i (O(n), avoids clone + O(n²) search).
        let mut inv = vec![0usize; n_teams];
        for (si, &orig_i) in arena.sort_buf.iter().enumerate() {
            inv[orig_i] = si;
        }

        self.likelihoods = self
            .teams
            .iter()
            .zip(self.weights.iter())
            .enumerate()
            .map(|(orig_i, (players, weights))| {
                let sorted_i = inv[orig_i];
                let m = m_t_ft[sorted_i];
                let performance = players
                    .iter()
                    .zip(weights.iter())
                    .fold(N00, |p, (player, &w)| p + (player.performance() * w));
                players
                    .iter()
                    .zip(weights.iter())
                    .map(|(player, &w)| {
                        ((m - performance.exclude(player.performance() * w)) * (1.0 / w))
                            .forget(player.beta.powi(2))
                    })
                    .collect::<Vec<_>>()
            })
            .collect::<Vec<_>>();
    }

    pub fn posteriors(&self) -> Vec<Vec<Gaussian>> {
        self.likelihoods
            .iter()
            .zip(self.teams.iter())
            .map(|(l, t)| {
                l.iter()
                    .zip(t.iter())
                    .map(|(&l, p)| l * p.prior)
                    .collect::<Vec<_>>()
            })
            .collect::<Vec<_>>()
    }
}

#[cfg(test)]
mod tests {
    use ::approx::assert_ulps_eq;

    use super::*;
    use crate::{ConstantDrift, GAMMA, Gaussian, N_INF, Player, arena::ScratchArena};

    #[test]
    fn test_1vs1() {
        let t_a = Player::new(
            Gaussian::from_ms(25.0, 25.0 / 3.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        );
        let t_b = Player::new(
            Gaussian::from_ms(25.0, 25.0 / 3.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        );

        let w = [vec![1.0], vec![1.0]];
        let g = Game::new(
            vec![vec![t_a], vec![t_b]],
            &[0.0, 1.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        let a = p[0][0];
        let b = p[1][0];

        assert_ulps_eq!(a, Gaussian::from_ms(20.794779, 7.194481), epsilon = 1e-6);
        assert_ulps_eq!(b, Gaussian::from_ms(29.205220, 7.194481), epsilon = 1e-6);

        let t_a = Player::new(
            Gaussian::from_ms(29.0, 1.0),
            25.0 / 6.0,
            ConstantDrift(GAMMA),
        );
        let t_b = Player::new(
            Gaussian::from_ms(25.0, 25.0 / 3.0),
            25.0 / 6.0,
            ConstantDrift(GAMMA),
        );

        let w = [vec![1.0], vec![1.0]];
        let g = Game::new(
            vec![vec![t_a], vec![t_b]],
            &[0.0, 1.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        let a = p[0][0];
        let b = p[1][0];

        assert_ulps_eq!(a, Gaussian::from_ms(28.896475, 0.996604), epsilon = 1e-6);
        assert_ulps_eq!(b, Gaussian::from_ms(32.189211, 6.062063), epsilon = 1e-6);

        let t_a = Player::new(Gaussian::from_ms(1.139, 0.531), 1.0, ConstantDrift(0.2125));
        let t_b = Player::new(Gaussian::from_ms(15.568, 0.51), 1.0, ConstantDrift(0.2125));

        let w = [vec![1.0], vec![1.0]];
        let g = Game::new(
            vec![vec![t_a], vec![t_b]],
            &[0.0, 1.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );

        assert_eq!(g.likelihoods[0][0], N_INF);
        assert_eq!(g.likelihoods[1][0], N_INF);
    }

    #[test]
    fn test_1vs1vs1() {
        let teams = vec![
            vec![Player::new(
                Gaussian::from_ms(25.0, 25.0 / 3.0),
                25.0 / 6.0,
                ConstantDrift(25.0 / 300.0),
            )],
            vec![Player::new(
                Gaussian::from_ms(25.0, 25.0 / 3.0),
                25.0 / 6.0,
                ConstantDrift(25.0 / 300.0),
            )],
            vec![Player::new(
                Gaussian::from_ms(25.0, 25.0 / 3.0),
                25.0 / 6.0,
                ConstantDrift(25.0 / 300.0),
            )],
        ];

        let w = [vec![1.0], vec![1.0], vec![1.0]];
        let g = Game::new(
            teams.clone(),
            &[1.0, 2.0, 0.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        let a = p[0][0];
        let b = p[1][0];

        assert_ulps_eq!(a, Gaussian::from_ms(25.000000, 6.238469), epsilon = 1e-6);
        assert_ulps_eq!(b, Gaussian::from_ms(31.311358, 6.698818), epsilon = 1e-6);

        let w = [vec![1.0], vec![1.0], vec![1.0]];
        let g = Game::new(
            teams.clone(),
            &[2.0, 1.0, 0.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        let a = p[0][0];
        let b = p[1][0];

        assert_ulps_eq!(a, Gaussian::from_ms(31.311358, 6.698818), epsilon = 1e-6);
        assert_ulps_eq!(b, Gaussian::from_ms(25.000000, 6.238469), epsilon = 1e-6);

        let w = [vec![1.0], vec![1.0], vec![1.0]];
        let g = Game::new(teams, &[1.0, 2.0, 0.0], &w, 0.5, &mut ScratchArena::new());
        let p = g.posteriors();

        let a = p[0][0];
        let b = p[1][0];
        let c = p[2][0];

        // T1 ULP shift: mu rounds to 25.0 (was 24.999999) under natural-parameter storage.
        assert_ulps_eq!(a, Gaussian::from_ms(25.0, 6.092561), epsilon = 1e-6);
        assert_ulps_eq!(b, Gaussian::from_ms(33.379314, 6.483575), epsilon = 1e-6);
        assert_ulps_eq!(c, Gaussian::from_ms(16.620685, 6.483575), epsilon = 1e-6);
    }

    #[test]
    fn test_1vs1_draw() {
        let t_a = Player::new(
            Gaussian::from_ms(25.0, 25.0 / 3.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        );
        let t_b = Player::new(
            Gaussian::from_ms(25.0, 25.0 / 3.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        );

        let w = [vec![1.0], vec![1.0]];
        let g = Game::new(
            vec![vec![t_a], vec![t_b]],
            &[0.0, 0.0],
            &w,
            0.25,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        let a = p[0][0];
        let b = p[1][0];

        assert_ulps_eq!(a, Gaussian::from_ms(24.999999, 6.469480), epsilon = 1e-6);
        assert_ulps_eq!(b, Gaussian::from_ms(24.999999, 6.469480), epsilon = 1e-6);

        let t_a = Player::new(
            Gaussian::from_ms(25.0, 3.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        );
        let t_b = Player::new(
            Gaussian::from_ms(29.0, 2.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        );

        let w = [vec![1.0], vec![1.0]];
        let g = Game::new(
            vec![vec![t_a], vec![t_b]],
            &[0.0, 0.0],
            &w,
            0.25,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        let a = p[0][0];
        let b = p[1][0];

        assert_ulps_eq!(a, Gaussian::from_ms(25.736001, 2.709956), epsilon = 1e-6);
        assert_ulps_eq!(b, Gaussian::from_ms(28.672888, 1.916471), epsilon = 1e-6);
    }

    #[test]
    fn test_1vs1vs1_draw() {
        let t_a = Player::new(
            Gaussian::from_ms(25.0, 25.0 / 3.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        );
        let t_b = Player::new(
            Gaussian::from_ms(25.0, 25.0 / 3.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        );
        let t_c = Player::new(
            Gaussian::from_ms(25.0, 25.0 / 3.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        );

        let w = [vec![1.0], vec![1.0], vec![1.0]];
        let g = Game::new(
            vec![vec![t_a], vec![t_b], vec![t_c]],
            &[0.0, 0.0, 0.0],
            &w,
            0.25,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        let a = p[0][0];
        let b = p[1][0];
        let c = p[2][0];

        // Goldens updated for natural-parameter storage: mu rounds to 25.0 (was 24.999999),
        // sigma shifts by ~3e-7 ULPs (within 1e-6 of original). Both bounded differences.
        assert_ulps_eq!(a, Gaussian::from_ms(25.0, 5.729069), epsilon = 1e-6);
        assert_ulps_eq!(b, Gaussian::from_ms(25.0, 5.707424), epsilon = 1e-6);
        assert_ulps_eq!(c, Gaussian::from_ms(25.0, 5.729069), epsilon = 1e-6);

        let t_a = Player::new(
            Gaussian::from_ms(25.0, 3.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        );
        let t_b = Player::new(
            Gaussian::from_ms(25.0, 3.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        );
        let t_c = Player::new(
            Gaussian::from_ms(29.0, 2.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        );

        let w = [vec![1.0], vec![1.0], vec![1.0]];
        let g = Game::new(
            vec![vec![t_a], vec![t_b], vec![t_c]],
            &[0.0, 0.0, 0.0],
            &w,
            0.25,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        let a = p[0][0];
        let b = p[1][0];
        let c = p[2][0];

        assert_ulps_eq!(a, Gaussian::from_ms(25.488507, 2.638208), epsilon = 1e-6);
        assert_ulps_eq!(b, Gaussian::from_ms(25.510671, 2.628751), epsilon = 1e-6);
        assert_ulps_eq!(c, Gaussian::from_ms(28.555920, 1.885689), epsilon = 1e-6);
    }

    #[test]
    fn test_2vs1vs2_mixed() {
        let t_a = vec![
            Player::new(
                Gaussian::from_ms(12.0, 3.0),
                25.0 / 6.0,
                ConstantDrift(25.0 / 300.0),
            ),
            Player::new(
                Gaussian::from_ms(18.0, 3.0),
                25.0 / 6.0,
                ConstantDrift(25.0 / 300.0),
            ),
        ];
        let t_b = vec![Player::new(
            Gaussian::from_ms(30.0, 3.0),
            25.0 / 6.0,
            ConstantDrift(25.0 / 300.0),
        )];
        let t_c = vec![
            Player::new(
                Gaussian::from_ms(14.0, 3.0),
                25.0 / 6.0,
                ConstantDrift(25.0 / 300.0),
            ),
            Player::new(
                Gaussian::from_ms(16., 3.0),
                25.0 / 6.0,
                ConstantDrift(25.0 / 300.0),
            ),
        ];

        let w = [vec![1.0, 1.0], vec![1.0], vec![1.0, 1.0]];
        let g = Game::new(
            vec![t_a, t_b, t_c],
            &[1.0, 0.0, 0.0],
            &w,
            0.25,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        assert_ulps_eq!(p[0][0], Gaussian::from_ms(13.051, 2.864), epsilon = 1e-3);
        assert_ulps_eq!(p[0][1], Gaussian::from_ms(19.051, 2.864), epsilon = 1e-3);
        assert_ulps_eq!(p[1][0], Gaussian::from_ms(29.292, 2.764), epsilon = 1e-3);
        assert_ulps_eq!(p[2][0], Gaussian::from_ms(13.658, 2.813), epsilon = 1e-3);
        assert_ulps_eq!(p[2][1], Gaussian::from_ms(15.658, 2.813), epsilon = 1e-3);
    }

    #[test]
    fn test_1vs1_weighted() {
        let w_a = vec![1.0];
        let w_b = vec![2.0];

        let t_a = vec![Player::new(
            Gaussian::from_ms(25.0, 25.0 / 3.0),
            25.0 / 6.0,
            ConstantDrift(0.0),
        )];
        let t_b = vec![Player::new(
            Gaussian::from_ms(25.0, 25.0 / 3.0),
            25.0 / 6.0,
            ConstantDrift(0.0),
        )];

        let w = [w_a, w_b];
        let g = Game::new(
            vec![t_a.clone(), t_b.clone()],
            &[1.0, 0.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        assert_ulps_eq!(
            p[0][0],
            Gaussian::from_ms(30.625173, 7.765472),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[1][0],
            Gaussian::from_ms(13.749653, 5.733840),
            epsilon = 1e-6
        );

        let w_a = vec![1.0];
        let w_b = vec![0.7];

        let w = [w_a, w_b];
        let g = Game::new(
            vec![t_a.clone(), t_b.clone()],
            &[1.0, 0.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        assert_ulps_eq!(
            p[0][0],
            Gaussian::from_ms(27.630080, 7.206676),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[1][0],
            Gaussian::from_ms(23.158943, 7.801628),
            epsilon = 1e-6
        );

        let w_a = vec![1.6];
        let w_b = vec![0.7];

        let w = [w_a, w_b];
        let g = Game::new(
            vec![t_a, t_b],
            &[1.0, 0.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        assert_ulps_eq!(
            p[0][0],
            Gaussian::from_ms(26.142438, 7.573088),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[1][0],
            Gaussian::from_ms(24.500183, 8.193278),
            epsilon = 1e-6
        );

        let w_a = vec![1.0];
        let w_b = vec![0.0];

        let t_a = vec![Player::new(
            Gaussian::from_ms(2.0, 6.0),
            1.0,
            ConstantDrift(0.0),
        )];
        let t_b = vec![Player::new(
            Gaussian::from_ms(2.0, 6.0),
            1.0,
            ConstantDrift(0.0),
        )];

        let w = [w_a, w_b];
        let g = Game::new(
            vec![t_a, t_b],
            &[1.0, 0.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        assert_ulps_eq!(
            p[0][0],
            Gaussian::from_ms(5.557067, 4.052826),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[1][0],
            Gaussian::from_ms(2.000000, 6.000000),
            epsilon = 1e-6
        );

        let w_a = vec![1.0];
        let w_b = vec![-1.0];

        let t_a = vec![Player::new(
            Gaussian::from_ms(2.0, 6.0),
            1.0,
            ConstantDrift(0.0),
        )];
        let t_b = vec![Player::new(
            Gaussian::from_ms(2.0, 6.0),
            1.0,
            ConstantDrift(0.0),
        )];

        let w = [w_a, w_b];
        let g = Game::new(
            vec![t_a, t_b],
            &[1.0, 0.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        assert_ulps_eq!(p[0][0], p[1][0], epsilon = 1e-6);
    }

    #[test]
    fn test_2vs2_weighted() {
        let t_a = vec![
            Player::new(
                Gaussian::from_ms(25.0, 25.0 / 3.0),
                25.0 / 6.0,
                ConstantDrift(0.0),
            ),
            Player::new(
                Gaussian::from_ms(25.0, 25.0 / 3.0),
                25.0 / 6.0,
                ConstantDrift(0.0),
            ),
        ];
        let w_a = vec![0.4, 0.8];

        let t_b = vec![
            Player::new(
                Gaussian::from_ms(25.0, 25.0 / 3.0),
                25.0 / 6.0,
                ConstantDrift(0.0),
            ),
            Player::new(
                Gaussian::from_ms(25.0, 25.0 / 3.0),
                25.0 / 6.0,
                ConstantDrift(0.0),
            ),
        ];
        let w_b = vec![0.9, 0.6];

        let w = [w_a, w_b];
        let g = Game::new(
            vec![t_a.clone(), t_b.clone()],
            &[1.0, 0.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        assert_ulps_eq!(
            p[0][0],
            Gaussian::from_ms(27.539023, 8.129639),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[0][1],
            Gaussian::from_ms(30.078046, 7.485372),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[1][0],
            Gaussian::from_ms(19.287197, 7.243465),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[1][1],
            Gaussian::from_ms(21.191465, 7.867608),
            epsilon = 1e-6
        );

        let w_a = vec![1.3, 1.5];
        let w_b = vec![0.7, 0.4];

        let w = [w_a, w_b];
        let g = Game::new(
            vec![t_a.clone(), t_b.clone()],
            &[1.0, 0.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        assert_ulps_eq!(
            p[0][0],
            Gaussian::from_ms(25.190190, 8.220511),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[0][1],
            Gaussian::from_ms(25.219450, 8.182783),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[1][0],
            Gaussian::from_ms(24.897589, 8.300779),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[1][1],
            Gaussian::from_ms(24.941479, 8.322717),
            epsilon = 1e-6
        );

        let w_a = vec![1.6, 0.2];
        let w_b = vec![0.7, 2.4];

        let w = [w_a, w_b];
        let g = Game::new(
            vec![t_a.clone(), t_b.clone()],
            &[1.0, 0.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        assert_ulps_eq!(
            p[0][0],
            Gaussian::from_ms(31.674697, 7.501180),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[0][1],
            Gaussian::from_ms(25.834337, 8.320970),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[1][0],
            Gaussian::from_ms(22.079819, 8.180607),
            epsilon = 1e-6
        );
        assert_ulps_eq!(
            p[1][1],
            Gaussian::from_ms(14.987953, 6.308469),
            epsilon = 1e-6
        );

        let w = [vec![1.0, 1.0], vec![1.0]];
        let g = Game::new(
            vec![
                t_a.clone(),
                vec![Player::new(
                    Gaussian::from_ms(25.0, 25.0 / 3.0),
                    25.0 / 6.0,
                    ConstantDrift(0.0),
                )],
            ],
            &[1.0, 0.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let post_2vs1 = g.posteriors();

        let w_a = vec![1.0, 1.0];
        let w_b = vec![1.0, 0.0];

        let w = [w_a, w_b];
        let g = Game::new(
            vec![t_a, t_b.clone()],
            &[1.0, 0.0],
            &w,
            0.0,
            &mut ScratchArena::new(),
        );
        let p = g.posteriors();

        assert_ulps_eq!(p[0][0], post_2vs1[0][0], epsilon = 1e-6);
        assert_ulps_eq!(p[0][1], post_2vs1[0][1], epsilon = 1e-6);
        assert_ulps_eq!(p[1][0], post_2vs1[1][0], epsilon = 1e-6);
        assert_ulps_eq!(p[1][1], t_b[1].prior, epsilon = 1e-6);
    }
}