clutch.diffpriv.coupling_rules

Coupling rules

From stdpp Require Import namespaces.
From iris.proofmode Require Import proofmode.
From clutch.prelude Require Import stdpp_ext fin.
From clutch.prob Require Import differential_privacy.
From clutch.diffpriv Require Import lifting ectx_lifting.
From clutch.prob_lang Require Import lang notation tactics metatheory erasure.
From clutch.prob_lang.spec Require Import spec_rules.
From clutch.diffpriv Require Export primitive_laws.

Section rules.
  Context `{!diffprivGS Σ}.
  Implicit Types P Q : iProp Σ.
  Implicit Types Φ : val → iProp Σ.
  Implicit Types σ : state.
  Implicit Types e : expr.
  Implicit Types v : val.
  Implicit Types l : loc.

  #[local] Open Scope R.

helper lemma
  Lemma DPcoupl_steps_ctx_bind_r `{Countable A} (μ : distr A)
    e1' σ1' R (ε δ: nonnegreal) K :
    to_val e1' = None →
    DPcoupl μ (prim_step e1' σ1') R ε δ →
    DPcoupl μ (prim_step (fill K e1') σ1')
      (λ a '(e2', σ2'), ∃ e2'', (e2', σ2') = (fill K e2'', σ2') ∧ R a (e2'', σ2')) ε δ.
  Proof.
    intros Hcpl Hv.
    rewrite fill_dmap //= -(dret_id_right μ ) /=.
    eapply (DPcoupl_dbind' ε 0 _ δ 0); [lra|done|done|lra| |done].
    intros ? [] ?.
    apply DPcoupl_dret=>/=; [done|done|]. eauto.
  Qed.

  Lemma hoare_couple_laplace_exact (loc : Z)
    (num den : Z) K E :
    {{{ ⤇ fill K (Laplace #num #den #loc #()) }}}
      Laplace #num #den #loc #() @ E
      {{{ (z : Z), RET #z; ⤇ fill K #z }}}.
  Proof.
    iIntros (?) "Hr Hcnt".
    iApply wp_lift_prim_steps_coupl; [done|].
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1) & Hauth2 & Herr) /=".
    iDestruct (spec_auth_prog_agree with "Hauth2 Hr") as %->.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iExists _, 0%NNR, ε_now, 0%NNR, δ_now.
    repeat iSplit.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. solve_red.
    - iPureIntro. solve_red.
    - iPureIntro. apply DPcoupl_steps_ctx_bind_r => //.
      eapply DPcoupl_laplace_primstep_exact.
    - iIntros (???? (?& [=->] & (z & [=-> ->] & [=-> ->]))).
      iMod (spec_update_prog (fill K #z) with "Hauth2 Hr") as "[$ Hspec0]".
      do 2 iModIntro.
      iMod "Hclose'" as "_".
      iModIntro. iFrame.
      rewrite -wp_value.
      iDestruct ("Hcnt" with "[$Hspec0]") as "$".
  Qed.

  Lemma hoare_couple_laplace (loc loc' k k' : Z)
    (Hdist : (Z.abs (k + loc - loc') <= k')%Z)
    (num den : Z) (ε ε' : R) K E :
    IZR num / IZR den = ε →
    0 < IZR num / IZR den →
    ε' = (IZR k' * ε) →
    {{{ ⤇ fill K (Laplace #num #den #loc' #()) ∗ ↯m ε' }}}
      Laplace #num #den #loc #() @ E
      {{{ (z : Z), RET #z; ⤇ fill K #(z+k) }}}.
  Proof.
    iIntros (Hε εpos Hε').
    iIntros (?) "(Hr & Hε) Hcnt".
    iApply wp_lift_prim_steps_coupl; [done|].
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1) & Hauth2 & (Hε2 & Hδ)) /=".
    iDestruct (spec_auth_prog_agree with "Hauth2 Hr") as %->.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iDestruct (ecm_supply_ecm_inv with "Hε2 Hε") as %(ε_sp & ε_rest & -> & Hε'').
    iExists _, ε_sp, ε_rest, 0%NNR, δ_now.
    repeat iSplit.
    - iPureIntro. done.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. solve_red.
    - iPureIntro. solve_red.
    - iPureIntro. apply DPcoupl_steps_ctx_bind_r => //. rewrite Hε''.
      eapply DPcoupl_laplace_primstep => //.
    - iIntros (???? (?& [=->] & (z & [=-> ->] & [=-> ->]))).
      iMod (spec_update_prog (fill K #(z+k)) with "Hauth2 Hr") as "[$ Hspec0]".
      iMod (ecm_supply_decrease with "Hε2 Hε") as (????) "H".
      do 2 iModIntro.
      iMod "Hclose'" as "_".
      iModIntro. iFrame.
      rewrite -wp_value.
      iDestruct ("Hcnt" with "[$Hspec0]") as "$".
      have Heq : x3 = ε_rest.
      { apply nnreal_ext.
        have := f_equal nonneg H. simpl in *. subst. lra. }
      subst x3. iFrame.
      Unshelve. all: constructor.
  Qed.

  Lemma wp_couple_tapes_laplace (mean mean' k k' : Z) α α' zs zs' e Φ
    (Hdist : (Z.abs (k + mean - mean') <= k')%Z)
    (num den : Z) (ε ε' : R) E :
    IZR num / IZR den = ε →
    0 < IZR num / IZR den →
    ε' = (IZR k' * ε) →
    ▷ α ↪L (num, den, mean; zs) ∗ ▷ α' ↪Lₛ (num, den, mean'; zs') ∗ ↯m ε' -∗
    (∀ z, α ↪L (num, den, mean; zs ++ [z]) ∗ α' ↪Lₛ (num, den, mean'; zs' ++ [z + k]%Z) -∗
          WP e @ E {{ Φ }}) -∗
    WP e @ E {{ Φ }}.
  Proof.
    iIntros (Hε εpos Hε').
    iIntros "(>α & >α' & Hε) HΦ".
    iApply wp_lift_step_spec_couple.
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & Hε2 & Hδ2)".
    iDestruct "Hauth2" as "(HK&Hh2&Ht2&Htl2)/=".
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iDestruct (ecm_supply_ecm_inv with "Hε2 Hε") as %(ε'' & ε_now_rest & foo & Hε'').
    rewrite -Hε''. rewrite -Hε'' in Hε'. clear Hε''. clear ε'. rename ε'' into ε'.
    iDestruct (ghost_map_lookup with "Htl2 α'") as %?.
    iDestruct (ghost_map_lookup with "Htl1 α") as %?.

    iApply (spec_coupl_erasables_weak _ _ _ ε' ε_now_rest _ 0%NNR δ_now) => //.
    1: apply nnreal_ext ; simpl ; lra.
    1: eapply DPcoupl_laplace_statestep => //.
    { by eapply state_step_laplace_erasable. }
    { by eapply state_step_laplace_erasable. }
    simpl.
    iIntros (σ2 σ2' (z & -> & ->)).
    iApply spec_coupl_ret.

    iDestruct (ghost_map_lookup with "Htl1 α") as %?%lookup_total_correct.
    iDestruct (ghost_map_lookup with "Htl2 α'") as %?%lookup_total_correct.
    iMod (ghost_map_update ((Tape_Laplace num den mean (zs ++ [z]))) with "Htl1 α") as "[$ α]".
    iMod (ghost_map_update ((Tape_Laplace num den mean' (zs' ++ [z+k]%Z))) with "Htl2 α'") as "[$ α']".
    iMod (ecm_supply_decrease with "Hε2 Hε") as (????) "H".
    iModIntro. iMod "Hclose'" as "_". iFrame.
    simplify_eq.
    iDestruct ("HΦ" $! z with "[$α $α']") as "Hwp".
    iSplitL "H".
    2: done.
    iApply ecm_supply_eq; [|done].
      simplify_eq/=; lra.
  Qed.

  Lemma hoare_couple_laplace_choice (loc loc' T : Z)
    (dist_loc : (Z.abs (loc - loc') <= 1)%Z)
    (num den : Z) (ε ε' : R) K E :
    IZR num / IZR den = ε →
    0 < IZR num / IZR den →
    ε' = (2*ε) →
    {{{ ⤇ fill K (Laplace #num #den #loc' #()) ∗ ↯m ε' }}}
      Laplace #num #den #loc #() @ E
      {{{ (z : Z), RET #z;
          ∃ z' : Z, ⤇ fill K #z'
                 ∗
                   ( ⌜(T <= z ∧ T + 1 <= z')⌝
                     ∨
                       (⌜z < T ∧ z' < T + 1⌝ ∗ ↯m ε'))%Z
           }}}.
  Proof.
    iIntros (Hε εpos Hε').
    iIntros (?) "(Hr & Hε) Hcnt".
    iApply wp_lift_prim_steps_choice; [done|].
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1) & Hauth2 & (Hε2 & Hδ))".
    iDestruct (spec_auth_prog_agree with "Hauth2 Hr") as %->.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iDestruct (ecm_supply_ecm_inv with "Hε2 Hε") as %(? & ε_now_minus_ε' & H_ε_now & Hε'').
    set (P := (λ '(ez, _) , ∃ z : Z, ez = Val (LitV (LitInt z)) ∧ T <= z)%Z : cfg → Prop).
    set (R := (λ ρ ρ' : expr * state ,
                  ρ.2 = σ1 ∧ ρ'.2 = σ1' ∧
                  (P ρ →
                  let (ez, ez') := (ρ.1, ρ'.1) in
                  ∃ z z' : Z, ez = Val (LitV (LitInt z)) ∧
                                ez' = Val (LitV (LitInt z')) ∧
                                T <= z ∧ T + 1 <= z')%Z)).
    set (RR := (λ a '(e2', σ2'), ∃ e2'', (e2', σ2') = (fill K e2'', σ2') ∧ R a (e2'', σ2'))).
    set (R' := (λ ρ ρ' : expr * state ,
                   ρ.2 = σ1 ∧ ρ'.2 = σ1' ∧
                   (¬ P ρ →
                   let (ez, ez') := (ρ.1, ρ'.1) in
                   ∃ z z' : Z, ez = Val (LitV (LitInt z)) ∧
                                 ez' = Val (LitV (LitInt z')) ∧
                                 z < T ∧ z' < T + 1)%Z)).
    set (RR' := (λ a '(e2', σ2'), ∃ e2'', (e2', σ2') = (fill K e2'', σ2') ∧ R' a (e2'', σ2'))).
    iExists P, RR, RR', x, ε_now_minus_ε', 0%NNR, ε_now, 0%NNR, δ_now, 0%NNR.
    repeat iSplit.
    - iPureIntro. done.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. solve_red.
    - iPureIntro. solve_red.
    (* Disjointness of R and R' under P is fine *)
    - iPureIntro. intros [] [] []. intros P_ρ nP_ρ'. subst R R' P ; simpl in *. intros [h h'].
      destruct h as (e1'' & eq'' & R_ρρ').
      destruct h' as (? & eq''' & R'_ρρ').
      apply R_ρρ' in P_ρ. apply R'_ρρ' in nP_ρ'.
      destruct P_ρ as [?[?[?[?[]]]]]. destruct nP_ρ' as [?[?[?[?[]]]]].
      subst. simplify_eq. lia.

    (* If we're above (P holds) the coupling should be the shift Laplace translation. *)
    (* Shouldn't we get to know that P holds here? *)
    - iPureIntro. replace 0%R with (nonneg 0%NNR) => //.
      apply DPcoupl_steps_ctx_bind_r => //.
      subst. simpl in *.
      eapply DPcoupl_mono ; last first.
      (* for 2num/den we get a coupling that shifts the rhs by 1 *)
      1: eapply (DPcoupl_laplace_primstep loc loc' 1 (Z.abs (1 + loc - loc'))).
      all: try by intuition eauto.
      {
        rewrite Hε''. real_solver_partial. 1: lra.
        replace 2 with (IZR 2) => //. apply IZR_le.
        revert dist_loc. repeat apply Zabs_ind ; lia. }
      + simpl. intros [e σ] [e' σ'] (z & eq_ez & eq_ez'). repeat split. 1,2: simpl ; by simplify_eq.
        intros Pe. destruct Pe as (ey & eq_ey & above). simpl.
        exists z, (z + 1)%Z.
        repeat split ; simplify_eq => //.
        lia.

    (* if P is false we use the trivial coupling *)
    - iPureIntro. replace 0%R with (nonneg 0%NNR) => //. apply DPcoupl_steps_ctx_bind_r => //.
      subst. simpl in *.
      eapply DPcoupl_mono ; last first.
      1: eapply (DPcoupl_laplace_primstep loc loc'
                   (loc' - loc) (Z.abs ((loc' - loc) + loc - loc'))).
      all: try by intuition eauto.
      {
        replace (loc' - loc + loc - loc')%Z with 0%Z by lia.
        assert (Z.abs 0 = 0)%Z as ->. 2: lra. apply Zabs_ind ; lia. }
      + simpl. intros [e σ] [e' σ'] (z & eq_ez & eq_ez').
        repeat split. 1,2: simpl ; by simplify_eq.
        intros nPe. exists z, (z + (loc' - loc))%Z.
        repeat split ; simplify_eq => //.
        * subst P R R' ; simpl in *.
          destruct (decide (z < T)%Z) => //.
          exfalso. apply nPe. exists z. split => //. lia.
        * subst P R R' ; simpl in *.
          destruct (decide (z < T)%Z). 1: lia. exfalso ; apply nPe ; exists z ; split => // ; lia.

    - iIntros (e2 σ2 e2' σ2').
      destruct (decide (P (e2, σ2))) as [p | n].
      + iSplitL ; last first.
        { iIntros ([nP_ρ2 R'_ρ2]). exfalso. done. }
        iIntros (((ze2 & eqe2 & Pe2) & (e2'' & [eq_e2'' R_ρ2ρ2']))).
        unfold R in R_ρ2ρ2'.
        simpl in R_ρ2ρ2'. destruct R_ρ2ρ2' as (<- & <- & R_ρ2ρ2'). specialize (R_ρ2ρ2' p).
        destruct R_ρ2ρ2' as (z & z' & eq_e2_z & eq_e2''_z' & z_above & z'_above).
        inversion eq_e2''. simplify_eq.
        iMod (spec_update_prog (fill K #(_)) with "Hauth2 Hr") as "[$ Hspec0]".
        iMod (ecm_supply_decrease with "Hε2 Hε") as (???Herr Hε''') "H".
        do 2 iModIntro.
        iMod "Hclose'" as "_".
        iModIntro. iFrame.
        rewrite -wp_value.
        iDestruct ("Hcnt" with "[$Hspec0]") as "$".
        { iLeft. done. }
        simplify_eq. rewrite Hε'' Hε''' in Herr.
        rewrite Rplus_comm in Herr. apply Rplus_eq_reg_r in Herr. clear -Herr.
        apply nnreal_ext in Herr. subst. iFrame.

      + iSplitR.
        { iIntros ([P_ρ2 R_ρ2]). exfalso. done. }
        iIntros ((nP_ρ2 & (e2'' & [eq_e2'' R'_ρ2ρ2']))).
        unfold R' in R'_ρ2ρ2'.
        simpl in R'_ρ2ρ2'. destruct R'_ρ2ρ2' as (<- & <- & R'_ρ2ρ2'). specialize (R'_ρ2ρ2' nP_ρ2).
        destruct R'_ρ2ρ2' as (z & z' & eq_e2_z & eq_e2''_z' & z_below & z'_below).
        inversion eq_e2''. simplify_eq.
        iMod (spec_update_prog (fill K #(_)) with "Hauth2 Hr") as "[$ Hspec0]".
        do 2 iModIntro.
        iMod "Hclose'" as "_".
        iModIntro. iFrame.
        rewrite -wp_value.
        iDestruct ("Hcnt" with "[$Hspec0 Hε]") as "$".
        { iRight. iFrame. done. }
        Unshelve. all: exact 0%Z.
  Qed.

TODO: This should be generalizable to injective functions N -> M Then we can get the exact couplings with bijections as a corollary
  Lemma wp_couple_tapes (N M : nat) E e α αₛ ns nsₛ Φ (δ : R) :
    (N <= M)%nat →
    (S M - S N) / S M = δ →
    ▷ α ↪N (N; ns) ∗ ▷ αₛ ↪ₛN (M; nsₛ) ∗
    ↯ δ ∗
    (∀ (n : nat),
        ⌜ n ≤ N ⌝ -∗
        α ↪N (N; ns ++ [n]) ∗ αₛ ↪ₛN (M; nsₛ ++ [n]) -∗
        WP e @ E {{ Φ }})
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (NMpos NMδ) "(>Hα & >Hαₛ & Hδ & Hwp)".
    iMod ecm_zero as "Hε".
    iApply wp_lift_step_spec_couple.
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & Hε2 & Hδ2)".
    iDestruct "Hauth2" as "(HK&Hh2&Ht2&Htl2)/=".
    iDestruct "Hα" as (fs) "(%&Hα)".
    iDestruct "Hαₛ" as (fsₛ) "(%&Hαₛ)".
    iDestruct (ghost_map_lookup with "Ht2 Hαₛ") as %?.
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iDestruct (ec_supply_ec_inv with "Hδ2 Hδ") as %(?&?&->&<-).
    iDestruct (ecm_supply_ecm_inv with "Hε2 Hε") as %(?&?&->&hh).
    iApply spec_coupl_erasables_weak ; [done|done|..].
    { rewrite hh. simpl.
      apply ARcoupl_to_DPcoupl.
      apply (ARcoupl_state_state N M σ1 σ1' α αₛ fs fsₛ x NMpos NMδ H2 H1). }
    { by eapply state_step_erasable. }
    { by eapply state_step_erasable. }
    iIntros (σ2 σ2' (n & m & nm & -> & ->)).
    iApply spec_coupl_ret.
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?%lookup_total_correct.
    iDestruct (ghost_map_lookup with "Ht2 Hαₛ") as %?%lookup_total_correct.
    iMod (ghost_map_update ((N; fs ++ [n]) : tape) with "Ht1 Hα") as "[$ Hα]".
    iMod (ghost_map_update ((M; fsₛ ++ [m]) : tape) with "Ht2 Hαₛ") as "[$ Hαₛ]".
    iMod (ec_supply_decrease with "Hδ2 Hδ") as (????) "H".
    iMod (ecm_supply_decrease with "Hε2 Hε") as (????) "Hm".
    iModIntro. iMod "Hclose'" as "_". iFrame.
    pose proof (fin_to_nat_lt n).
    iDestruct ("Hwp" $! n with "[]") as "Hwp".
    { iPureIntro; lia. }
    rewrite -/add_ec_supply.
    iSplitL "H Hm".
    { iSplitL "Hm".
      - iApply ecm_supply_eq; [|done]. simplify_eq/=; lra.
      - iApply ec_supply_eq; [|done]. simplify_eq/=; lra.
    }
    iModIntro.
    iApply "Hwp".
    iSplitL "Hα".
    - iExists _. iFrame.
      rewrite fmap_app.
      simplify_eq. done.
    - iExists _. iFrame.
      rewrite nm.
      rewrite fmap_app.
      simplify_eq. done.
  Qed.

  Lemma wp_couple_tapes_bij N f `{Bij nat nat f} E e α αₛ ns nsₛ Φ :
    (forall n, n < S N -> f n < S N)%nat ->
    ▷ α ↪N (N; ns) ∗ ▷ αₛ ↪ₛN (N; nsₛ) ∗
      (∀ n : nat, α ↪N (N; ns ++ [n]) ∗ αₛ ↪ₛN (N; nsₛ ++ [f n]) ∗ ⌜ n ≤ N ⌝ -∗
                    WP e @ E {{ Φ }})
      ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (Hdom) "(>Hα & >Hαₛ & Hwp)".
    iApply wp_lift_step_spec_couple.
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & Hε2 & Hδ2)".
    iDestruct "Hauth2" as "(HK&Hh2&Ht2&Htl2)/=".
    iDestruct "Hα" as (fs) "(<-&Hα)".
    iDestruct "Hαₛ" as (fsₛ) "(<-&Hαₛ)".
    destruct (restr_bij_fin (S N) f) as [g [HBij Hfg]].
    { intros n Hn.
      by apply Hdom.
    }
    iDestruct (ghost_map_lookup with "Ht2 Hαₛ") as %?.
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    replace ε_now with (0 + ε_now)%NNR; last first.
    { apply nnreal_ext; simpl; lra. }
    replace δ_now with (0 + δ_now)%NNR; last first.
    { apply nnreal_ext; simpl; lra. }
    iApply spec_coupl_erasables_weak; [done|done|..].
    { eapply ARcoupl_to_DPcoupl, ARcoupl_exact.
      (* eauto unifies the wrong premise? *)
      apply Rcoupl_state_state; [apply HBij | apply H1 | apply H0 ]. }
    { by eapply state_step_erasable. }
    { by eapply state_step_erasable. }
    iIntros (σ2 σ2' (n & ? & ?)).
    iApply spec_coupl_ret.
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?%lookup_total_correct.
    iDestruct (ghost_map_lookup with "Ht2 Hαₛ") as %?%lookup_total_correct.
    simplify_map_eq.
    iMod (ghost_map_update ((N; fs ++ [n]) : tape) with "Ht1 Hα") as "[$ Hα]".
    iMod (ghost_map_update ((N; fsₛ ++ [g n]) : tape) with "Ht2 Hαₛ") as "[$ Hαₛ]".
    iModIntro. iMod "Hclose'" as "_".
    replace (0 + ε_now)%NNR with ε_now; last first.
    { apply nnreal_ext; simpl; lra. }
    replace (0 + δ_now)%NNR with δ_now; last first.
    { apply nnreal_ext; simpl; lra. }
    iFrame.
    iApply ("Hwp" $! (fin_to_nat n) with "[-]").
    iSplitL "Hα".
    { iExists _. iFrame.
      iPureIntro.
      rewrite fmap_app //. }
    iSplitL "Hαₛ".
    { iExists _. iFrame.
      iPureIntro.
      rewrite fmap_app -Hfg //. }
    iPureIntro.
    apply (fin.fin_to_nat_le n).
  Qed.

  Lemma wp_couple_tapes_rev (N M : nat) E e α αₛ ns nsₛ Φ (δ : R) :
    (M <= N)%nat →
    (S N - S M) / S N = δ →
    ▷ α ↪N (N; ns) ∗ ▷ αₛ ↪ₛN (M; nsₛ) ∗
    ↯ δ ∗
    (∀ (n m : nat),
        ⌜n = m⌝ -∗
        α ↪N (N; ns ++ [n]) ∗ αₛ ↪ₛN (M; nsₛ ++ [m]) -∗
        ⌜m ≤ M⌝ -∗
        WP e @ E {{ Φ }})
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (NMpos NMδ) "( >Hα & >Hαₛ & Hδ & Hwp)".
    iApply wp_lift_step_spec_couple.
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & Hε2 & Hδ2)".
    iDestruct "Hauth2" as "(HK&Hh2&Ht2&Htl2)".
    iDestruct "Hα" as (fs) "(<-&Hα)".
    iDestruct "Hαₛ" as (fsₛ) "(<-&Hαₛ)".
    iDestruct (ghost_map_lookup with "Ht2 Hαₛ") as %?.
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    replace ε_now with (0 + ε_now)%NNR; last first.
    { apply nnreal_ext; simpl; lra. }
    replace δ_now with (0 + δ_now)%NNR; last first.
    { apply nnreal_ext; simpl; lra. }
    iDestruct (ec_supply_ec_inv with "Hδ2 Hδ") as %(?&?&->&<-).
    iApply spec_coupl_erasables_weak ; [done|done|..].
    { by apply ARcoupl_to_DPcoupl, ARcoupl_state_state_rev. }
    { by eapply state_step_erasable. }
    { by eapply state_step_erasable. }
    iIntros (σ2 σ2' (n & m & nm & ? & ?)).
    iApply spec_coupl_ret.
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?%lookup_total_correct.
    iDestruct (ghost_map_lookup with "Ht2 Hαₛ") as %?%lookup_total_correct.
    simplify_map_eq.
    iMod (ghost_map_update ((N; fs ++ [n]) : tape) with "Ht1 Hα") as "[$ Hα]".
    iMod (ghost_map_update ((M; fsₛ ++ [m]) : tape) with "Ht2 Hαₛ") as "[$ Hαₛ]".
    iMod (ec_supply_decrease with "Hδ2 Hδ") as (????) "H".
    iModIntro. iMod "Hclose'" as "_".
    iFrame.
    iDestruct ("Hwp" with "[//] [-H Hε2] []") as "$".
    - iSplitL "Hα".
      { iExists _. iFrame.
        iPureIntro.
        rewrite fmap_app //. }
      iExists _. iFrame.
      iPureIntro.
      rewrite fmap_app //.
    - iPureIntro.
      apply (fin_to_nat_le m).
    - replace (0 + ε_now)%NNR with ε_now; last first.
      { apply nnreal_ext; simpl; lra. }
      iSplitR "H" ; iFrame ; try done.
      iApply ec_supply_eq; [|done].
      simplify_eq/=; lra.
  Qed.

  Lemma wp_rand_avoid_l {N} (m : nat) (z : Z) E (δ : R) :
    TCEq N (Z.to_nat z) →
    δ = 1 / (S N) →
    {{{ ↯ δ }}}
      rand #z @ E
    {{{ (n : nat), RET #n; ⌜n ≠ m⌝ ∗ ⌜ n ≤ N ⌝ }}}.
  Proof.
    iIntros (-> -> Φ) "Hδ Hwp".
    iApply wp_lift_prim_step_l_dret; [done|].
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & Hε2 & Hδ2)".
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iDestruct (ec_supply_ec_inv with "Hδ2 Hδ") as %(δ_cert & x1_r & -> & Hδ_cert).
    iExists _, 0%NNR, ε_now, δ_cert, x1_r.
    repeat iSplit.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. done.
    - iPureIntro. solve_red.
    - iPureIntro.
      by apply ARcoupl_to_DPcoupl, (ARcoupl_rand_no_coll_l _ (fin_force (Z.to_nat z) m)).
    - iIntros (?? (n & [= -> ->] & H1 & [=])).
      iMod (ec_supply_decrease with "Hδ2 Hδ") as (ec_del ec_rest Heq_rest Hdel_eq) "Hs".
      simplify_eq.
      do 2 iModIntro.
      iMod "Hclose'" as "_".
      iModIntro. iFrame.
      rewrite -wp_value.
      iDestruct ("Hwp" with "[]") as "$".
      + iPureIntro.
        split; eauto.
        * destruct (le_gt_dec m (Z.to_nat z)).
          -- intro H3. apply H1. apply fin_to_nat_inj.
             rewrite fin_force_to_nat_le; auto.
          -- pose proof (fin_to_nat_le n). lia.
        * apply (fin_to_nat_le).
      + iApply ec_supply_eq; [|done]. lra.
  Qed.

  Lemma wp_rand_avoid_r {N} (m : nat) (z : Z) K e E Φ (δ : R) :
    TCEq N (Z.to_nat z) →
    δ = 1 / (S N) →
    ⤇ fill K (rand #z) ∗
    ↯ δ ∗
    (∀ (n : nat),
       ⤇ fill K #n -∗ ⌜n ≠ m⌝ -∗ ⌜ n ≤ N ⌝ -∗ WP e @ E {{ Φ }})
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (-> ->) "(HK & Hδ & Hwp)".
    iApply wp_lift_step_spec_couple.
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & Hε2 & Hδ2)".
    iDestruct (spec_auth_prog_agree with "Hauth2 HK") as "->".
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iDestruct (ec_supply_ec_inv with "Hδ2 Hδ") as %(δ &?& -> & ?).
    replace ε_now with (0 + ε_now)%NNR; last first.
    { apply nnreal_ext; simpl; lra. }
    iApply (spec_coupl_erasable_steps 1 _ (dret σ1)); [done|done|..].
    { rewrite pexec_1 step_or_final_no_final //; last first.
      { apply reducible_not_final. solve_red. }
      eapply DPcoupl_steps_ctx_bind_r; [done|].
      by apply ARcoupl_to_DPcoupl, (ARcoupl_rand_no_coll_r _ (fin_force (Z.to_nat z) m)). }
    { apply dret_erasable. }
    iIntros (??? (n & [=-> ] & (y&->&[=-> ->]&?))) "!>".
    iApply spec_coupl_ret.
    iMod (spec_update_prog (fill K #_) with "Hauth2 HK") as "[$ HK]".
    iMod (ec_supply_decrease with "Hδ2 Hδ") as (????) "H".
    simplify_eq.
    iMod "Hclose'" as "_".
    iFrame.
    iDestruct ("Hwp" with "[$] []") as "Hwp".
    {
      iPureIntro.
      destruct (le_gt_dec m (Z.to_nat z)).
      - intro H3. apply H2. apply fin_to_nat_inj.
        rewrite fin_force_to_nat_le; auto.
      - pose proof (fin_to_nat_le y). lia.
    }
    iSplitL "Hε2 H".
    { iSplitL "Hε2".
      - iApply ecm_supply_eq; [|done]. simpl. lra.
      - iApply ec_supply_eq; [|done]. lra.
    }
    iApply "Hwp".
    iPureIntro.
    apply (fin_to_nat_le y).
  Qed.

rand(unit, N) ~ rand(unit, M) coupling, N <= M, under inj
  Lemma wp_couple_rand_rand_inj (N M : nat) (f: nat → nat) z w K E (δ : R) :
    (∀ n, n < S N → f n < S M)%nat →
    (∀ n1 n2, n1 < S N → n2 < S N → f n1 = f n2 → n1 = n2)%nat →
    TCEq N (Z.to_nat z) →
    TCEq M (Z.to_nat w) →
    (N <= M)%nat →
    (S M - S N) / S M = δ →
    {{{ ⤇ fill K (rand #w) ∗ ↯ δ }}}
      rand #z @ E
    {{{ (n : nat), RET #n; ⤇ fill K #(f n) ∗ ⌜ n ≤ N ⌝ }}}.
  Proof.
    iIntros (Hdom Hinj).

    set g := (λ m : fin (S N), Fin.of_nat_lt (Hdom m (fin_to_nat_lt m))).
    assert (Inj eq eq g).
    { intros m1 m2 Heq.
      assert (fin_to_nat (g m1) = f (fin_to_nat m1)) as H1.
      { rewrite /g fin_to_nat_to_fin //. }
      assert (fin_to_nat (g m2) = f (fin_to_nat m2)) as H2.
      { rewrite /g fin_to_nat_to_fin //. }
      apply fin_to_nat_inj.
      apply Hinj; [apply fin_to_nat_lt..|].
      rewrite -H1 -H2 //. by f_equal. }

    iIntros (-> -> HNM Hδ ?) "(Hr & Hδ) Hcnt".
    iApply wp_lift_prim_steps_coupl; [done|].
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & (Hε2 & Hδ2))".
    iDestruct (spec_auth_prog_agree with "Hauth2 Hr") as %->.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iDestruct (ec_supply_ec_inv with "Hδ2 Hδ") as %(δ_del & δ_rest & -> & Hδ_del).
    iExists _, 0%NNR, ε_now, δ_del, δ_rest.
    repeat iSplit.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. done.
    - iPureIntro. solve_red.
    - iPureIntro. solve_red.
    - iPureIntro. replace 0%R with (nonneg 0%NNR) => //.
      apply DPcoupl_steps_ctx_bind_r, ARcoupl_to_DPcoupl,
        (ARcoupl_rand_rand_inj _ _ g); done || lra.
    - iIntros (???? (?& [=->] & (n & [=-> ->] & [=-> ->]))).
      iMod (spec_update_prog (fill K #(g _)) with "Hauth2 Hr") as "[$ Hspec0]".
      iMod (ec_supply_decrease with "Hδ2 Hδ") as (????) "H".
      do 2 iModIntro.
      iMod "Hclose'" as "_".
      iModIntro. iFrame.
      rewrite -wp_value.
      rewrite /g fin_to_nat_to_fin.
      iDestruct ("Hcnt" with "[$Hspec0]") as "$".
      { iPureIntro. apply fin_to_nat_le. }
      have Heq : x3 = δ_rest.
      { apply nnreal_ext.
        have := f_equal nonneg H0. simpl in *. lra. }
      subst x3. iFrame.
  Qed.

rand(unit, N) ~ rand(unit, M) coupling, N <= M, under equality
  Lemma wp_couple_rand_rand_leq (N M : nat) z w K E (δ : R) :
    TCEq N (Z.to_nat z) →
    TCEq M (Z.to_nat w) →
    (N <= M)%nat →
    (S M - S N) / S M = δ →
    {{{ ⤇ fill K (rand #w) ∗ ↯ δ }}}
      rand #z @ E
    {{{ (n : nat), RET #n;
        ⌜ n ≤ N ⌝ ∗ ⤇ fill K #n }}}.
  Proof.
    iIntros (-> -> HNM <- ?) "(Hr & Hδ) Hwp".
    iApply (wp_couple_rand_rand_inj _ _ (λ x, x) with "[$]"); [lia|lia|done|done|].
    iModIntro. iIntros (?) "[? %]". iApply ("Hwp" $! n).
    iFrame.
    iPureIntro. lia.
  Qed.

rand(unit, N) ~ rand(unit, M) coupling, M <= N, along an injection
  Lemma wp_couple_rand_rand_rev_inj (N M : nat) (f : nat -> nat) z w K E (δ : R) :
    (∀ n, n < S M -> f n < S N)%nat →
    (∀ n1 n2, n1 < S M → n2 < S M → f n1 = f n2 → n1 = n2)%nat →
    TCEq N (Z.to_nat z) →
    TCEq M (Z.to_nat w) →
    (M <= N)%nat →
    (S N - S M) / S N = δ →
    {{{ ⤇ fill K (rand #w) ∗ ↯ δ }}}
      rand #z @ E
    {{{ (m : nat), RET #(f m); ⤇ fill K #m ∗ ⌜ m ≤ M ⌝ }}}.
  Proof.
    iIntros (Hdom Hinj).

    set g := (λ m : fin (S M), Fin.of_nat_lt (Hdom m (fin_to_nat_lt m))).
    assert (Inj eq eq g).
    { intros m1 m2 Heq.
      assert (fin_to_nat (g m1) = f (fin_to_nat m1)) as H1.
      { rewrite /g fin_to_nat_to_fin //. }
      assert (fin_to_nat (g m2) = f (fin_to_nat m2)) as H2.
      { rewrite /g fin_to_nat_to_fin //. }
      apply fin_to_nat_inj.
      apply Hinj; [apply fin_to_nat_lt..|].
      rewrite -H1 -H2 //. by f_equal. }

    iIntros (-> -> HNM <- ?) "(Hr & Hδ) Hwp".
    iApply wp_lift_prim_steps_coupl; [done|].
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & (Hε2 & Hδ2))".
    iDestruct (spec_auth_prog_agree with "Hauth2 Hr") as %->.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iDestruct (ec_supply_ec_inv with "Hδ2 Hδ") as %(δ_del & δ_rest & -> & Hδ_del).
    iExists _, 0%NNR, ε_now, δ_del, δ_rest.
    repeat iSplit.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. done.
    - iPureIntro. solve_red.
    - iPureIntro. solve_red.
    - iPureIntro. replace 0%R with (nonneg 0%NNR) => //.
      eapply DPcoupl_steps_ctx_bind_r; [done|].
      by eapply ARcoupl_to_DPcoupl, (ARcoupl_rand_rand_rev_inj _ _ g).
    - iIntros (???? (?& [=->] & (n & [=-> ->] & [=-> ->]))).
      iMod (spec_update_prog (fill K #_) with "Hauth2 Hr") as "[$ Hspec0]".
      iMod (ec_supply_decrease with "Hδ2 Hδ") as (????) "H".
      do 2 iModIntro.
      iMod "Hclose'" as "_".
      iModIntro. iFrame.
      rewrite /g fin_to_nat_to_fin //.
      rewrite -wp_value.
      iDestruct ("Hwp" with "[$Hspec0]") as "$".
      { iPureIntro. apply fin_to_nat_le. }
      have Heq : x3 = δ_rest.
      { apply nnreal_ext.
        have := f_equal nonneg H0. simpl in *. lra. }
      subst x3. iFrame.
  Qed.

rand(unit, N) ~ rand(unit, M) coupling, N <= M, under equality
  Lemma wp_couple_rand_rand_rev_leq (N M : nat) z w K E (ε : R) :
    TCEq N (Z.to_nat z) →
    TCEq M (Z.to_nat w) →
    (M <= N)%nat →
    (S N - S M) / S N = ε →
    {{{ ⤇ fill K (rand #w) ∗ ↯ ε }}}
      rand #z @ E
    {{{ (n : nat ), RET #n;
        ⌜ n ≤ N ⌝ ∗ ⤇ fill K #n }}}.
  Proof.
    iIntros (-> -> HNM <- ?) "(Hr & Hε) Hwp".
    iApply (wp_couple_rand_rand_rev_inj _ _ (λ x, x) with "[$]"); [lia|lia|done|done|..].
    iIntros "!>" (m) "[? %]".
    iSpecialize ("Hwp" $! m).
    iApply ("Hwp" with "[-]").
    iFrame.
    iPureIntro. lia.
  Qed.

rand(N) ~ rand(N) coupling


  Lemma hoare_couple_rand_rand_kanto (N M : nat) (S : nat -> nat -> Prop) (X Y : nat → nat -> R)
                            z w K E (ε δ : R) :
    TCEq N (Z.to_nat z) →
    TCEq M (Z.to_nat w) →
    (forall n m, 0 <= X n m) ->
    (forall n m, 0 <= Y n m <= 1) ->
    (forall h1 h2,
          (forall a, 0 <= h1 a <= 1) ->
          (forall b, 0 <= h2 b <= 1) ->
          (forall a b, S a b -> h1 a <= exp (X a b) * h2 b + Y a b) ->
          (Expval (dunifP N) h1 <=
             exp ε * Expval (dunifP M) h2 + δ)) ->
    {{{ ⤇ fill K (rand #w) ∗ ↯m ε ∗ ↯ δ }}}
      rand #z @ E
      {{{ (n m : nat), RET #n; ⤇ fill K #m ∗ ⌜ S n m ⌝ ∗ ↯m (X n m) ∗ ↯ (Y n m) }}}.
  Proof.
    iIntros (-> -> HX_nonneg HY_bounds Hcoupl Φ) "(Hr & Hε & Hδ) Hcnt".
    iApply (wp_lift_prim_steps_coupl_adv_frame_err_le_1); [done|].
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1) & Hauth2 & (Hε2 & Hδ2))".
    iDestruct (spec_auth_prog_agree with "Hauth2 Hr") as %->.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iDestruct (ecm_supply_ecm_inv with "Hε2 Hε") as %(εm2 & εm1 & -> & <-).
    iDestruct (ec_supply_ec_inv with "Hδ2 Hδ") as %(δa2 & δa1 & -> & <-).
    set (X' n m := mknonnegreal _ (HX_nonneg n m)).
    assert (forall n m, 0 <= Y n m) as HY_nonneg; [intros; apply (HY_bounds n m)|].
    set (Y' n m := mknonnegreal _ (HY_nonneg n m)).
    set (find_nat_left :=
           (λ '(e,σ), match e with
                       | Val #(LitInt n) =>
                           if bool_decide (0 <= n)%Z then
                             if bool_decide (σ = σ1) then Some (Z.to_nat n) else None
                             else None
                       | _ => None
                      end): (cfg -> option nat)).
    destruct (partial_inv_fun (fill K)) as (Kinv & HKinv).
    assert (∀ b a : expr, Kinv b = Some a → fill K a = b) as HKinvS; [intros; by apply HKinv|].
    assert (∀ b a : expr, Kinv b = None → fill K a ≠ b) as HKinvN; [intros; by apply HKinv|].
    assert (∀ e, Kinv (fill K e) = Some e) as HKinv3.
    { intro e. destruct (Kinv (fill K e)) eqn:Heq; eapply HKinv in Heq; by simplify_eq. }
    set (find_nat_right :=
           (λ '(e,σ), match Kinv e with
                      | Some (Val #(LitInt m)) =>
                          if bool_decide (0 <= m)%Z then
                             if bool_decide (σ = σ1') then Some (Z.to_nat m) else None
                             else None
                      | _ => None
                      end): (cfg -> option nat)).
    assert (forall ρ (n : nat), ρ = (Val #n, σ1) -> find_nat_left ρ = Some n) as Hfl_intro.
    { intros [] n Hρ. rewrite /find_nat_left /=. apply pair_equal_spec in Hρ as [-> ->].
      rewrite bool_decide_eq_true_2; [|lia]. rewrite bool_decide_eq_true_2; auto. by rewrite Nat2Z.id. }
    assert (forall ρ (n : nat), ρ = (fill K (Val #n), σ1') -> find_nat_right ρ = Some n) as Hfr_intro.
    { intros [] n Hρ. rewrite /find_nat_right /=. apply pair_equal_spec in Hρ as [-> ->].
      rewrite HKinv3. rewrite bool_decide_eq_true_2; [|lia]. rewrite bool_decide_eq_true_2; auto. by rewrite Nat2Z.id. }
    assert (forall e σ (n : nat), find_nat_left (e,σ) = Some n -> e = Val #n) as Hfl_e.
    { rewrite /find_nat_left. intros e σ ?. repeat case_match; try done. intros Hsome.
      apply Some_inj in Hsome. do 3 f_equal. rewrite -Hsome. rewrite Z2Nat.id //. eapply bool_decide_eq_true_1; eauto. }
    assert (forall e σ (n : nat), find_nat_left (e,σ) = Some n -> σ = σ1) as Hfl_s.
    { rewrite /find_nat_left. intros e σ ?. repeat case_match; try done. intros Hsome.
      apply Some_inj in Hsome. eapply bool_decide_eq_true_1; eauto. }
    assert (forall e σ (m : nat), find_nat_right (e,σ) = Some m -> (e,σ) = (fill K (Val #m), σ1')) as Hfr_em.
    { rewrite /find_nat_right. intros e σ ?. repeat case_match; try done. intros Hsome.
      apply Some_inj in Hsome. apply pair_equal_spec; split.
      - rewrite -(HKinvS e (#n)); eauto. do 3 f_equal. rewrite -Hsome. rewrite Z2Nat.id //. eapply bool_decide_eq_true_1; eauto.
      - eapply bool_decide_eq_true_1; eauto. }
    set (S_rel := (λ ρ1 ρ2 : cfg, ∃ n m : nat, find_nat_left ρ1 = Some n ∧ find_nat_right ρ2 = Some m)).
    set (E2_cfg := (λ ρ1 ρ2 : cfg, match find_nat_left ρ1, find_nat_right ρ2 with
                                 | Some n, Some m => if bool_decide (S n m) then X' n m else 0%NNR
                                 | _, _ => 0%NNR
                                 end)).
    set (D2_cfg := (λ ρ1 ρ2 : cfg, match find_nat_left ρ1, find_nat_right ρ2 with
                                 | Some n, Some m => if bool_decide (S n m) then Y' n m else 1%NNR
                                 | _, _ => 1%NNR
                                 end)).
    iExists S_rel, E2_cfg, D2_cfg, εm2, εm1, δa2, δa1.
    repeat iSplit.
    - iPureIntro. apply nnreal_ext. simpl. lra.
    - iPureIntro. apply nnreal_ext. simpl. lra.
    - iPureIntro. solve_red.
    - iPureIntro. solve_red.
    - iPureIntro. intros ρ1 ρ2. rewrite /D2_cfg.
      repeat case_match; simpl; [apply (HY_bounds _ _)|lra |lra |lra].
    - iPureIntro. intros h1 h2 Hh1 Hh2 Hkanto.
      simpl language.prim_step.
      rewrite fill_dmap //=.
      rewrite !head_prim_step_eq; [| solve_red ..].
      rewrite {1}/head_step {1}/head_step. simpl.
      rewrite dmap_comp.
      rewrite Expval_dmap; [rewrite Expval_dmap|..]; rewrite /compose /=.
      + specialize (Hcoupl (λ n, if bool_decide (n ≤ Z.to_nat z) then h1 (Val #n, σ1) else 0)
                      (λ m, if bool_decide (m ≤ Z.to_nat w) then h2 (fill K (Val #m), σ1') else 1)).
        simpl in Hcoupl.
        etrans; [|etrans; [apply Hcoupl|]].
        -- rewrite /Expval. apply SeriesC_le.
          ++ intros nn. rewrite bool_decide_eq_true_2; [|apply fin_to_nat_le].
             split; [|lra].
             real_solver.
          ++ apply (ex_seriesC_le _ (dunifP (Z.to_nat z))); auto.
             intros nn; split; [apply Rmult_le_pos; [auto|case_bool_decide; [apply Hh1|lra]]
| rewrite -{2}(Rmult_1_r (dunifP (Z.to_nat z) nn)); apply Rmult_le_compat_l; [auto|case_bool_decide; [apply Hh1|lra]]].
        -- real_solver.
        -- real_solver.
        -- intros a b Sab.
           assert (Hfl : find_nat_left (Val #a, σ1) = Some a) by (apply Hfl_intro; reflexivity).
           assert (Hfr : find_nat_right (fill K (Val #b), σ1') = Some b) by (apply Hfr_intro; reflexivity).
           assert (HSrel : S_rel (Val #a, σ1) (fill K (Val #b), σ1')).
           { rewrite /S_rel. eauto. }
           have Hk := Hkanto _ _ HSrel.
           rewrite /E2_cfg /D2_cfg Hfl Hfr /= in Hk.
           rewrite (bool_decide_eq_true_2 (S a b)) in Hk; [|exact Sab].
           case_bool_decide as Ha; case_bool_decide as Hb.
           ++ etrans; [apply Hk|].
              apply Rplus_le_compat.
              --- rewrite /X'. simpl. lra.
              --- rewrite /Y'. simpl. lra.
           ++ etrans; [apply Hk|].
              apply Rplus_le_compat_r.
              apply Rmult_le_compat_l; [left; apply exp_pos|].
              apply Hh2.
           ++ apply Rplus_le_le_0_compat;
                [apply Rmult_le_pos; [left; apply exp_pos|apply Hh2] | apply HY_nonneg].
           ++ apply Rplus_le_le_0_compat;
                [apply Rmult_le_pos; [left; apply exp_pos|lra] | apply HY_nonneg].
        -- apply Rplus_le_compat_r. apply Rmult_le_compat_l; [left; apply exp_pos|].
           rewrite /Expval. apply SeriesC_le.
            ++ intros nn. rewrite (bool_decide_eq_true_2 (fin_to_nat nn ≤ Z.to_nat w)); [|apply fin_to_nat_le].
               real_solver.
            ++ apply (ex_seriesC_le _ (dunifP (Z.to_nat w))); auto.
               intros nn; split; [apply Rmult_le_pos; [auto|apply Hh2]
                                | rewrite -{2}(Rmult_1_r (dunifP (Z.to_nat w) nn)); apply Rmult_le_compat_l; [auto|apply Hh2]].
      + real_solver.
      + apply (ex_seriesC_le _ (dunifP (Z.to_nat w))); auto. intros nn; split; [apply Rmult_le_pos; [auto|apply Hh2]| rewrite -{2}(Rmult_1_r (dunifP (Z.to_nat w) nn)); apply Rmult_le_compat_l; [auto|apply Hh2]].
      + real_solver.
      + apply (ex_seriesC_le _ (dunifP (Z.to_nat z))); auto. intros nn; split; [apply Rmult_le_pos; [auto|apply Hh1]| rewrite -{2}(Rmult_1_r (dunifP (Z.to_nat z) nn)); apply Rmult_le_compat_l; [auto|apply Hh1]].
    - iIntros (e2 σ2 e2' σ2') "%HS".
      iModIntro.
      iModIntro.
      destruct (Rlt_or_le (D2_cfg (e2,σ2) (e2',σ2') + δa1) 1) as [HYlt1 | ]; last first.
      {
        iLeft.
        iApply (fupd_wand_l _ _ True).
        iSplit; done.
      }
      iRight.
      rewrite /D2_cfg.
      rewrite /D2_cfg in HYlt1.
      rewrite /E2_cfg.
      assert (exists n, find_nat_left (e2, σ2) = Some n) as [n Hrwn].
      {
        destruct (find_nat_left (e2, σ2)).
        - eauto.
        - simpl in HYlt1.
          pose proof (cond_nonneg δa1).
          lra.
      }
      rewrite Hrwn.
      rewrite Hrwn in HYlt1.
      assert (exists m, find_nat_right (e2', σ2') = Some m) as [m Hrwm].
      {
        destruct (find_nat_right (e2', σ2')).
        - eauto.
        - simpl in HYlt1.
          pose proof (cond_nonneg δa1).
          lra.
      }
      rewrite Hrwm.
      rewrite Hrwm in HYlt1.
      case_bool_decide as Hnm; last first.
      {
        simpl in HYlt1.
        pose proof (cond_nonneg δa1).
        lra.
      }
      iMod (spec_update_prog (fill K (#m)) with "Hauth2 Hr") as "[? Hspec0]".
      iMod (ecm_supply_decrease with "Hε2 Hε") as (ε3 ε4 Hε3ε4 Hε3eq) "Hε4".
      iMod (ecm_supply_increase _ (X' n m) with "Hε4") as "[Hm Herrm]".
      iMod (ec_supply_decrease with "Hδ2 Hδ") as (δ3 δ4 Hδ3δ4 Hδ3eq) "Hδ4".
      iMod (ec_supply_increase _ (Y' n m) with "Hδ4") as "[Ha Herra]".
      {
        eapply (Rle_lt_trans _ (Y' n m + δa1)); auto.
        rewrite Rplus_comm.
        right.
        f_equal.
        f_equal.
        rewrite -(nonnegreal_add_cancel_l _ _ δa2).
        rewrite (nnreal_plus_comm δa2 δ4).
        apply nnreal_ext in Hδ3eq.
        rewrite -{1}Hδ3eq //.
      }
      rewrite (Hfl_e e2 σ2 n); auto.
      rewrite (Hfl_s e2 σ2 n); auto.
      rewrite (Hfr_em e2' σ2' m); auto.
      iMod "Hclose'" as "_".
      iModIntro. iFrame.
      rewrite -wp_value.
      iDestruct ("Hcnt" with "[$Hspec0 Herrm Herra]") as "$".
      {
        rewrite /X'.
        rewrite /Y'.
        by iFrame.
      }
      iSplitL "Hm".
      + iApply ecm_supply_eq; [|done].
        simplify_eq. simpl. lra.
      + iApply ec_supply_eq; [|done].
        simplify_eq. simpl. lra.
  Qed.

  Lemma hoare_couple_rand_rand_kanto_diag (N : nat) (X Y : nat -> R)
                            z K E (ε δ : R) :
    TCEq N (Z.to_nat z) →
    (forall n, 0 <= X n) ->
    (forall n, 0 <= Y n <= 1) ->
    (forall h1 h2,
          (forall a, 0 <= h1 a <= 1) ->
          (forall b, 0 <= h2 b <= 1) ->
          (forall a, h1 a <= exp (X a) * h2 a + Y a) ->
          (Expval (dunifP N) h1 <=
             exp ε * Expval (dunifP N) h2 + δ)) ->
    {{{ ⤇ fill K (rand #z) ∗ ↯m ε ∗ ↯ δ }}}
      rand #z @ E
      {{{ (n : nat), RET #n; ⤇ fill K #n ∗ ↯m (X n) ∗ ↯ (Y n) }}}.
  Proof.
    iIntros (HN HX_nonneg HY_bounds Hkanto Φ) "(Hr & Hε & Hδ) HΦ".
    iApply (hoare_couple_rand_rand_kanto N N (λ n m, n = m)
              (λ n _, X n) (λ n _, Y n)
              z z K E ε δ with "[$Hr $Hε $Hδ] [HΦ]").
    - intros n m. apply HX_nonneg.
    - intros n m. apply HY_bounds.
    - intros h1 h2 Hh1 Hh2 Hbound.
      apply Hkanto; [exact Hh1 | exact Hh2 |].
      intros a. apply (Hbound a a). reflexivity.
    - iModIntro. iIntros (n m) "(Hspec & %Hnm & Hm & Ha)".
      subst m. iApply "HΦ". iFrame.
  Qed.

  Lemma hoare_couple_rand_l_adv (N : nat) (E2 D2 : nat -> R) z E (ε δ : R) :
    TCEq N (Z.to_nat z) →
    (forall n, 0 <= E2 n) ->
    (forall n, 0 <= D2 n <= 1) ->
    Expval (dunifP N) (λ a, exp (E2 a)) <= exp ε ->
    Expval (dunifP N) (λ a, D2 a) <= δ ->
    {{{ ↯m ε ∗ ↯ δ }}}
      rand #z @ E
    {{{ (n : nat), RET #n; ⌜ (n <= N) ⌝ ∗ ↯m (E2 n) ∗ ↯ (D2 n) }}}.
  Proof.
    iIntros (-> HE2_nonneg HD2_bounds HexpE2 HD2δ1 Φ) "(Hε & Hδ) Hcnt".
    iApply (wp_lift_prim_step_l_coupl_adv_frame_err_le_1); [done|].
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1) & Hauth2 & (Hε2 & Hδ2))".
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iDestruct (ecm_supply_ecm_inv with "Hε2 Hε") as %(εm2 & εm1 & -> & <-).
    iDestruct (ec_supply_ec_inv with "Hδ2 Hδ") as %(δa2 & δa1 & -> & <-).
    assert (forall n, 0 <= D2 n) as HD2_nonneg; [intros; apply (HD2_bounds n)|].
    set (find_nat :=
           (λ '(e,σ), match e with
                       | Val #(LitInt n) =>
                           if bool_decide (0 <= n)%Z then
                             if bool_decide (σ = σ1) then Some (Z.to_nat n) else None
                             else None
                       | _ => None
                      end): (cfg -> option nat)).
    set (E2_cfg := (λ ρ : cfg, match find_nat ρ with
                                | Some n => mknonnegreal (E2 n) (HE2_nonneg n)
                                | None => 0%NNR
                                end)).
    set (D2_cfg := (λ ρ : cfg, match find_nat ρ with
                                | Some n => mknonnegreal (D2 n) (HD2_nonneg n)
                                | None => 0%NNR
                                end)).
    iExists E2_cfg, D2_cfg, εm2, εm1, δa2, δa1.
    assert (forall ρ (n : nat), ρ = (Val #n, σ1) -> find_nat ρ = Some n) as Hfl_intro.
    { intros [] n Hρ. rewrite /find_nat /=. apply pair_equal_spec in Hρ as [-> ->].
      rewrite bool_decide_eq_true_2; [|lia]. rewrite bool_decide_eq_true_2; auto. by rewrite Nat2Z.id. }
    assert (forall e σ (n : nat), find_nat (e,σ) = Some n -> e = Val #n) as Hfl_e.
    { rewrite /find_nat. intros e σ ?. repeat case_match; try done. intros Hsome.
      apply Some_inj in Hsome. do 3 f_equal. rewrite -Hsome. rewrite Z2Nat.id //. eapply bool_decide_eq_true_1; eauto. }
    assert (forall e σ (n : nat), find_nat (e,σ) = Some n -> σ = σ1) as Hfl_s.
    { rewrite /find_nat. intros e σ ?. repeat case_match; try done. intros Hsome.
      apply Some_inj in Hsome. eapply bool_decide_eq_true_1; eauto. }
    repeat iSplit.
    - iPureIntro. apply nnreal_ext. simpl. lra.
    - iPureIntro. apply nnreal_ext. simpl. lra.
    - iPureIntro. solve_red.
    - iPureIntro. intros ρ1. rewrite /D2_cfg.
      case_match; simpl; [apply (HD2_bounds _)|lra].
    - iPureIntro.
      simpl language.prim_step.
      rewrite head_prim_step_eq.
      rewrite /head_step. simpl.
      rewrite Expval_dmap; cycle 1.
      { intros. apply Rlt_le. apply exp_pos. }
      { rewrite /ex_expval. apply ex_seriesC_finite. }
      rewrite /compose /=. etrans; [|apply HexpE2]. apply Req_le. rewrite /Expval. apply SeriesC_ext.
      intros nn. rewrite /E2_cfg.
      rewrite (Hfl_intro (Val #(fin_to_nat nn), σ1) (fin_to_nat nn)); [|reflexivity].
      done.
    - iPureIntro.
      simpl language.prim_step.
      rewrite head_prim_step_eq.
      rewrite /head_step. simpl.
      rewrite Expval_dmap; cycle 1.
      { intros. apply cond_nonneg. }
      { rewrite /ex_expval. apply ex_seriesC_finite. }
      rewrite /compose /=. etrans; [|apply HD2δ1]. apply Req_le. rewrite /Expval. apply SeriesC_ext.
      intros nn. rewrite /D2_cfg.
      rewrite (Hfl_intro (Val #(fin_to_nat nn), σ1) (fin_to_nat nn)); [|reflexivity].
      done.
    - iPureIntro.
      simpl language.prim_step.
      rewrite head_prim_step_eq.
      rewrite /head_step. simpl.
      rewrite /dmap.
      apply (ex_expval_dbind _ _ (λ a, exp (E2_cfg a))).
      + intros. apply Rlt_le, exp_pos.
      + apply ex_seriesC_finite.
      + intros. apply ex_expval_dret.
    - iIntros (e2 σ2) "%Hstep".
      iModIntro. iModIntro.
      destruct (Rlt_or_le (D2_cfg (e2, σ2) + δa1 : R) 1) as [HD2lt | HD2ge]; last first.
      { iMod "Hclose'" as "_". iModIntro. iLeft. iPureIntro. simpl in HD2ge |- *. lra. }
      iRight.
      simpl language.prim_step in Hstep.
      rewrite head_prim_step_eq in Hstep.
      rewrite /head_step in Hstep. simpl in Hstep.
      apply dmap_pos in Hstep as (nn & Heq & Hnn).
      inversion Heq; subst e2 σ2.
      rewrite /E2_cfg /D2_cfg.
      assert (find_nat (Val #nn, σ1) = Some (fin_to_nat nn)) as Hfn.
      { apply Hfl_intro. reflexivity. }
      rewrite Hfn.
      set (En := mknonnegreal (E2 nn) (HE2_nonneg nn)).
      set (Dn := mknonnegreal (D2 nn) (HD2_nonneg nn)).
      rewrite /D2_cfg Hfn in HD2lt. simpl in HD2lt.
      iMod (ecm_supply_decrease with "Hε2 Hε") as (ε3 ε4 Hε3ε4 Hε3eq) "Hε4".
      iMod (ecm_supply_increase _ En with "Hε4") as "[Hm Herrm]".
      iMod (ec_supply_decrease with "Hδ2 Hδ") as (δ3 δ4 Hδ3δ4 Hδ3eq) "Hδ4".
      iMod (ec_supply_increase _ Dn with "Hδ4") as "[Ha Herra]".
      { assert (δ4 = δa1) as ->.
        { apply (f_equal nonneg) in Hδ3δ4. simpl in Hδ3δ4.
          apply nnreal_ext. simpl.
          assert (nonneg δ3 = nonneg δa2) as Hd3 by (rewrite Hδ3eq; reflexivity). lra. }
        simpl. lra. }
      iMod "Hclose'" as "_".
      iModIntro. iFrame.
      rewrite -wp_value.
      iDestruct ("Hcnt" with "[Herrm Herra]") as "$".
      { iSplit.
        - iPureIntro. apply le_INR, fin_to_nat_le.
        - rewrite /En /Dn. iFrame. }
      iSplitL "Hm".
      { iApply ecm_supply_eq; [|done].
        apply (f_equal nonneg) in Hε3ε4. simpl in Hε3ε4.
        assert (nonneg ε3 = nonneg εm2) as He3 by (rewrite Hε3eq; reflexivity).
        simpl. lra. }
      iApply ec_supply_eq; [|done].
      apply (f_equal nonneg) in Hδ3δ4. simpl in Hδ3δ4.
      assert (nonneg δ3 = nonneg δa2) as Hd3 by (rewrite Hδ3eq; reflexivity).
      simpl. lra.
  Qed.

  Lemma hoare_couple_rand_rand_subsample K E (ε δ : R) :
    (0 <= ε) → (0 <= δ <= 1) →
    {{{ ⤇ fill K (rand #1) ∗ ↯m (ln (1 + 0.5 * (exp ε - 1))) ∗ ↯ (0.5*δ) }}}
      rand #1 @ E
      {{{ (n m : nat), RET #n; ⤇ fill K #m ∗
       ⌜ n <= 1 /\ m <=1 ⌝
       ∗ ↯m ( if bool_decide (n = 0 /\ m = 0)%nat then 0 else ε)
       ∗ ↯ (if bool_decide (n = 0 /\ m = 0)%nat then 0 else δ) }}}.
  Proof.
    iIntros (Hε Hδ).
    iIntros (Ψ) "(HK & Hmerr & Haerr) HΨ".
    iApply (hoare_couple_rand_rand_kanto 1 1
      (fun n m => n <= 1 /\ m <= 1)
      (fun n m => if bool_decide (n = 0 /\ m = 0)%nat then 0 else ε)
      (fun n m => if bool_decide (n = 0 /\ m = 0)%nat then 0 else δ)
      1 1 K E (ln (1 + 0.5 * (exp ε - 1))) (0.5 * δ) with "[$HK $Hmerr $Haerr] [HΨ]" ).
    - intros n m. case_bool_decide; lra.
    - intros n m. case_bool_decide; split; lra.
    - (* Kantorovich condition *)
      intros h1 h2 Hh1 Hh2 Hh1h2.
      set (h1' := fun b : bool => h1 (fin_to_nat (bool_to_fin b))).
      set (h2' := fun b : bool => h2 (fin_to_nat (bool_to_fin b))).
      have Hr : (0 <= 1/2 <= 1) by lra.
      have Hdunif : dunifP 1 = dmap bool_to_fin fair_coin.
      { apply distr_ext => n.
        rewrite /pmf /= /dbind_pmf SeriesC_bool.
        rewrite /pmf /= /fair_coin_pmf /dret_pmf.
        inv_fin n; [simpl; lra | intros n].
        inv_fin n; [simpl; lra | intros n].
        inversion n. }
      rewrite Hdunif.
      rewrite !Expval_dmap.
      2: { intro b. apply (Hh2 _). }
      2: { apply ex_seriesC_finite. }
      2: { intro b. apply (Hh1 _). }
      2: { apply ex_seriesC_finite. }
      rewrite !Expval_fair_coin.
      rewrite /compose /h1' /h2'.
      fold h1' h2'.
      (* Prove pointwise bound from Hh1h2 for all bool pairs *)
      have Hpointwise : ∀ b1 b2 : bool,
          h1' b1 <= exp (match b1, b2 with false, false => 0 | _, _ => ε end) *
                    h2' b2 + match b1, b2 with false, false => 0 | _, _ => δ end.
      { intros [] []; rewrite /h1' /h2' /=.
        - apply (Hh1h2 1%nat 1%nat). split; simpl; lra.
        - apply (Hh1h2 1%nat 0%nat). split; simpl; lra.
        - apply (Hh1h2 0%nat 1%nat). split; simpl; lra.
        - apply (Hh1h2 0%nat 0%nat). split; simpl; lra. }
      (* Use biased_coin_kanto_exp_bound with h1', h2' *)
      have key := biased_coin_kanto_exp_bound (1/2) Hr ε δ ltac:(lra) Hε (proj1 Hδ) h1' h2'
                    (fun b => Hh1 (fin_to_nat (bool_to_fin b)))
                    (fun b => Hh2 (fin_to_nat (bool_to_fin b)))
                    Hpointwise.
      (* Connect biased_coin (1/2) to fair_coin *)
      rewrite /h1' /h2' in key.
      etrans; [| etrans; [apply key | ]].
      { rewrite SeriesC_bool /pmf /= /biased_coin_pmf. lra. }
      { rewrite SeriesC_bool /pmf /= /biased_coin_pmf.
        replace (1 / 2 : R) with (0.5 : R) by lra. ring_simplify. lra. }
    - iModIntro.
       iIntros (n m) "H".
       iApply "HΨ". iFrame.
  Qed.

Converting an (ε,0) multiplicative-DP requirement into an (ε',δ) budget with ε' = ½ε² + √(2 ln(1/δ)) · ε (the pure-DP → zCDP → approximate-DP conversion, with ρ = ε²/2).
No program step is taken: we use wp_lift_step_spec_couple and instantiate spec_coupl with two trivial dirac distributions, so the only content is the credit reframe. The reframe's soundness obligation is exactly conversion_obligation_framed (framed by the ambient credit x3/y3 left in the supply). If δ ≥ 1, the ↯ δ resource is already contradictory, so no upper bound on δ is needed.
  Lemma wp_eps_to_delta_conv (ε δ : R) E e Φ :
    0 <= ε -> 0 < δ ->
    ↯m (/2 * ε^2 + sqrt (2 * ln (/δ)) * ε) ∗ ↯ δ ∗
    (↯m ε -∗ WP e @ E {{ Φ }})
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    intros Heps Hd.
    iIntros "(Hεm & Hδ & Hwp)".
    (* If δ ≥ 1 then ↯ δ is already contradictory. *)
    destruct (Rlt_le_dec δ 1) as [Hδ1 | Hδ1].
    2:{ iExFalso. iApply (ec_contradict with "Hδ"); lra. }
    assert (Hd1 : 0 < δ <= 1) by lra.
    iApply wp_lift_step_spec_couple.
    iIntros (σ1 e1' σ1' ε_now δ_now) "(Hσ & Hspec & Hεm2 & Hδ2)".
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose".
    iDestruct (ecm_supply_ecm_inv with "Hεm2 Hεm") as %(x1 & x3 & Hnow & Hx1).
    iDestruct (ec_supply_ec_inv with "Hδ2 Hδ") as %(y1 & y3 & Hδnow & Hy1).
    iApply spec_coupl_rec.
    iExists (fun σ2 (ρ : cfg) => σ2 = σ1 /\ ρ = (e1', σ1')),
            (dret σ1), (dret (e1', σ1')),
            (fun _ _ => (x3 + mknonnegreal ε Heps)%NNR),
            (fun _ _ => y3).
    iSplit; [iPureIntro; apply dret_erasable|].
    iSplit; [iPureIntro; apply dret_rewritable|].
    iSplit; [iPureIntro; exists y3; intros; simpl; lra|].
    iSplit.
    { iPureIntro. intros h1 h2 Hh1 Hh2 Hpre.
      rewrite !Expval_dret.
      specialize (Hpre σ1 (e1', σ1') (conj eq_refl eq_refl)).
      assert (Hεeq : nonneg ε_now = nonneg x3 + (/2 * ε^2 + sqrt (2 * ln (/δ)) * ε)).
      { rewrite Hnow /=. rewrite Hx1. lra. }
      assert (Hδeq : nonneg δ_now = nonneg y3 + δ).
      { rewrite Hδnow /=. rewrite Hy1. lra. }
      rewrite Hεeq Hδeq.
      apply (conversion_obligation_framed (nonneg x3) (nonneg y3) δ ε
               (cond_nonneg x3) (cond_nonneg y3) Hd1 Heps
               (h1 σ1) (h2 (e1', σ1')) (Hh1 σ1) (Hh2 (e1', σ1'))).
      exact Hpre. }
    iIntros (σ2 e2' σ2') "%HS". destruct HS as [-> HS2]. injection HS2 as -> ->.
    iApply spec_coupl_ret.
    iMod (ecm_supply_decrease with "Hεm2 Hεm") as (x1' x3' Hnow' Hx1') "Hεm2".
    iMod (ec_supply_decrease with "Hδ2 Hδ") as (y1' y3' Hδnow' Hy1') "Hδ2".
    iMod (ecm_supply_increase x3' (mknonnegreal ε Heps) with "Hεm2") as "[Hεm2 Hεfrag]".
    iModIntro. iMod "Hclose" as "_".
    iFrame "Hσ Hspec".
    assert (Hx3 : nonneg x3 = nonneg x3').
    { pose proof (f_equal nonneg Hnow) as A. pose proof (f_equal nonneg Hnow') as B.
      simpl in A, B. rewrite Hx1 in A. rewrite Hx1' in B. lra. }
    assert (Hy3 : nonneg y3' = nonneg y3).
    { pose proof (f_equal nonneg Hδnow) as A. pose proof (f_equal nonneg Hδnow') as B.
      simpl in A, B. rewrite Hy1 in A. rewrite Hy1' in B. lra. }
    iModIntro.
    iSplitR "Hwp Hεfrag".
    - rewrite /err_interp /mult_ec_supply /add_ec_supply.
      iSplitL "Hεm2".
      + iApply (ecm_supply_eq with "Hεm2"). simpl. rewrite Hx3. lra.
      + iApply (ec_supply_eq with "Hδ2"). simpl. exact Hy3.
    - iApply "Hwp". iApply (ecm_eq with "Hεfrag"). done.
  Qed.

  (*
    There should be an easier proof of this using wp_couple_rand_rand_inj,
    but that uses an injective function nat -> nat as opposed to fin (S N) -> fin (S N)
  *)
  Lemma wp_couple_rand_rand N f `{Bij nat nat f} z K E :
    TCEq N (Z.to_nat z) →
    (forall n:nat, (n < S N)%nat -> (f n < S N)%nat) ->
    {{{ ⤇ fill K (rand #z) }}}
      rand #z @ E
    {{{ (n : nat), RET #n; ⌜ n ≤ N ⌝ ∗ ⤇ fill K #(f n) }}}.
  Proof.
    iIntros (H0 Hdom Ψ) "Hr HΨ".
    destruct (restr_bij_fin (S N) f Hdom) as [ff [Hbij Hff]].
    iApply wp_lift_prim_steps_coupl; [done|].
    iIntros (σ1 e1' σ1' ε δ) "[Hσ [Hs Hε]]".
    iDestruct (spec_auth_prog_agree with "Hs Hr") as %->.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose".
    iExists _, 0%NNR, ε, 0%NNR, δ.
    repeat iSplit.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. solve_red.
    - iPureIntro. solve_red.
    - iPureIntro. replace 0%R with (nonneg 0%NNR) => //.
      apply (DPcoupl_steps_ctx_bind_r _ _ _
        (λ a b, ∃ n : fin (S N), a = (Val #n, σ1) ∧ b = (Val #(f n), σ1'))) => //.
      apply ARcoupl_to_DPcoupl, ARcoupl_exact.
      eapply Rcoupl_mono.
      { apply (Rcoupl_rand_rand _ ff). destruct H0; done. }
      intros [ea σa] [eb σb] [n [Hna Hnb]].
      simplify_eq. exists n. rewrite Hff. split; reflexivity.
    - iIntros (e2 σ2 e2' σ2' (b & [= ->] & (n & [= -> ->] & [= -> ->]))).
      iMod (spec_update_prog with "Hs Hr") as "[$ Hr]".
      do 2 iModIntro.
      iMod "Hclose" as "_".
      iModIntro. iFrame.
      iApply wp_value.
      iApply "HΨ".
      iFrame.
      iPureIntro. apply fin_to_nat_le.
  Qed.

coupling rand and rand but avoid certain values
  Lemma wp_couple_rand_rand_avoid N (l:list _) z K E :
    TCEq N (Z.to_nat z) →
    NoDup l ->
    {{{ ↯ (length l/(N+1)) ∗
          ⤇ fill K (rand #z) }}}
      rand #z @ E
      {{{ (n : fin (S N)), RET #n; ⌜n∉l⌝ ∗ ⤇ fill K #n }}}.
  Proof.
    iIntros (H0 Hl Ψ) "[Hδ Hr] HΨ".
    iApply wp_lift_prim_steps_coupl; [done|].
    iIntros (σ1 e1' σ1' ε δ) "[Hσ [Hs [Hε2 Hδ2]]]".
    iDestruct (spec_auth_prog_agree with "Hs Hr") as %->.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose".
    iDestruct (ec_supply_ec_inv with "Hδ2 Hδ") as %(x & x1 & -> & H).
    iExists _, 0%NNR, ε, x, x1.
    repeat iSplit.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. solve_red.
    - iPureIntro. solve_red.
    - iPureIntro.
      eapply DPcoupl_steps_ctx_bind_r; first done.
      apply ARcoupl_to_DPcoupl, ARcoupl_rand_rand_avoid_list; first done.
      + rewrite S_INR. lra.
      + by apply TCEq_eq.
    - iIntros (e2 σ2 e2' σ2' (b & [= ->] & (?&?&[= -> ->] & [= -> ->]))) "!> !>".
      iMod (spec_update_prog with "Hs Hr") as "[$ Hr]".
      iMod (ec_supply_decrease with "Hδ2 Hδ") as (x2 x3 H1 ?) "Hnew".
      replace (x3) with (x1); last first.
      { apply nnreal_ext. inversion H1. lra. }
      iMod "Hclose" as "_".
      iModIntro. iFrame.
      iApply wp_value.
      iApply "HΨ".
      iFrame.
      by iPureIntro.
  Qed.

  Local Lemma length_remove_dups `{EqDecision A} (l:list A):
    length (remove_dups l) <= length l.
  Proof.
    induction l; first done.
    simpl.
    case_match; simpl.
    all: rewrite -!/(INR (S _)); rewrite !S_INR; lra.
  Qed.

  Lemma wp_couple_rand_rand_avoid' N (l:list _) z K E :
    TCEq N (Z.to_nat z) →
    {{{ ↯ (length l/(N+1)) ∗
          ⤇ fill K (rand #z) }}}
      rand #z @ E
      {{{ (n : fin (S N)), RET #n; ⌜n∉l⌝ ∗ ⤇ fill K #n }}}.
  Proof.
    set (l':=remove_dups l).
    iIntros (H0 Ψ) "[Hε Hr] HΨ".
    iApply (wp_couple_rand_rand_avoid with "[-HΨ]").
    - apply NoDup_remove_dups.
    - iFrame.
      iApply (ec_weaken with "[$]").
      split.
      + apply Rcomplements.Rdiv_le_0_compat.
        * apply pos_INR.
        * pose proof pos_INR N. lra.
      + rewrite !Rdiv_def.
        apply Rmult_le_compat_r.
        * apply Rlt_le. apply Rinv_0_lt_compat. pose proof pos_INR N. lra.
        * apply length_remove_dups.
    - iModIntro. iIntros (?) "[%?]".
      iApply "HΨ". iFrame.
      iPureIntro.
      by rewrite -elem_of_remove_dups.
  Qed.

fragmented state rand N ~ state rand M, N>=M, under injective function from M to N
  Lemma wp_couple_fragmented_rand_rand_inj {M N} (f: nat → nat) {_ : Inj (=) (=) f}
    ns nsₛ α αₛ e E Φ:
    (M <= N)%nat →
    (forall n : nat, (n < S M)%nat -> (f n < S N)%nat) ->
    ▷ α ↪N (N; ns) ∗ ▷ αₛ ↪ₛN (M; nsₛ) ∗
    (∀ (n : nat),
       ⌜ n ≤ N ⌝ -∗
       if bool_decide (∃ m:nat, m ≤ M /\ f m = n) then
         ∀ m : nat, α ↪N (N; ns ++ [f m]) ∗ αₛ ↪ₛN (M; nsₛ ++ [m]) ∗ ⌜ f m ≤ N ⌝ ∗ ⌜ m ≤ M ⌝ -∗
              WP e @ E {{ Φ }}
       else
         α ↪N (N; ns ++ [n]) ∗ αₛ ↪ₛN (M; nsₛ) ∗ ⌜ n ≤ N ⌝ -∗
         WP e @ E {{ Φ }})
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (Hineq Hdom) "(>Hα & >Hαₛ & Hwp)".
    edestruct (restr_inj_fin (S M) (S N) f (le_n_S M N Hineq) Hdom) as [g [HgInj HgEq]].
    iDestruct "Hα" as (fs) "(<-&Hα)".
    iDestruct "Hαₛ" as (fsₛ) "(<-&Hαₛ)".
    iApply wp_lift_step_spec_couple.
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & Hε2 & Hδ2)".
    iDestruct "Hauth2" as "(HK&Hh2&Ht2&Htl2)/=".
    iDestruct (ghost_map_lookup with "Ht2 Hαₛ") as %?.
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?.
    simplify_map_eq.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    replace (ε_now) with (0 + ε_now)%NNR; last (apply nnreal_ext; simpl; lra).
    replace (δ_now) with (0 + δ_now)%NNR; last (apply nnreal_ext; simpl; lra).
    iApply spec_coupl_erasables_weak; [done|done|..].
    { by apply ARcoupl_to_DPcoupl, ARcoupl_exact, Rcoupl_fragmented_rand_rand_inj. }
    { by eapply state_step_erasable. }
    { eapply erasable_dbind_predicate.
      - solve_distr_mass.
      - by eapply state_step_erasable.
      - apply dret_erasable. }
    iIntros (?? [n H']).
    case_bool_decide in H'.
    - destruct Hf as [m' <-].
      destruct H' as (m & ? & ? & Hfm).
      simplify_eq.
      iMod (ghost_map_update ((N; fs ++ [g _]) : tape) with "Ht1 Hα") as "[Ht1 Hα]".
      iMod (ghost_map_update ((M; fsₛ ++ [_]) : tape) with "Ht2 Hαₛ") as "[Ht2 Hαₛ]".
      iModIntro.
      iApply spec_coupl_ret.
      iMod "Hclose'" as "_".
      iSpecialize ("Hwp" $! (f m')).
      rewrite bool_decide_eq_true_2.
      2: { exists m'.
           split; auto.
           apply fin_to_nat_le. }
      iSpecialize ("Hwp" $! _ m').
      iDestruct ("Hwp" with "[$Hα $Hαₛ]") as "Hwp".
      { iPureIntro.
        split; [rewrite fmap_app /= HgEq // |].
        split; [rewrite fmap_app /= // |].
        split; auto.
        - apply Nat.lt_succ_r, Hdom, fin_to_nat_lt.
        - apply fin_to_nat_le.
      }
      assert (0 + ε_now = ε_now)%NNR as ->.
      { apply nnreal_ext; simpl; lra. }
      assert (0 + δ_now = δ_now)%NNR as ->.
      { apply nnreal_ext; simpl; lra. }
      by iFrame.
    - destruct H' as [??]. simplify_eq.
      iMod (ghost_map_update ((N; fs ++ [n]) : tape) with "Ht1 Hα") as "[Ht1 Hα]".
      iModIntro.
      iApply spec_coupl_ret.
      iMod "Hclose'" as "_".
      iSpecialize ("Hwp" $! (fin_to_nat n)).
      rewrite bool_decide_eq_false_2 //.
      2: {
        intros [m [Hm1 Hm2]].
        apply Hf.
        assert (m < S M)%nat as Hm3.
        { lia. }
        exists (nat_to_fin Hm3).
        apply fin_to_nat_inj.
        rewrite HgEq -Hm2.
        rewrite fin_to_nat_to_fin //.
      }
      iDestruct ("Hwp" with "[]") as "Hwp".
      { iPureIntro. apply fin_to_nat_le. }
      assert (0 + ε_now = ε_now)%NNR as ->.
      { apply nnreal_ext; simpl; lra. }
      assert (0 + δ_now = δ_now)%NNR as ->.
      { apply nnreal_ext; simpl; lra. }
      iFrame.
      iApply "Hwp".
      iModIntro.
      iSplitL "Hα".
      { iFrame. rewrite fmap_app //. }
      iSplitL "Hαₛ".
      { iFrame. auto. }
      iPureIntro. apply fin_to_nat_le.
      Unshelve.
      apply Nat.lt_succ_r, Hdom, fin_to_nat_lt.
  Qed.

fragmented state rand N ~ state rand M, N>=M, under equality
  Lemma wp_couple_fragmented_rand_rand_leq (M N : nat) ns nsₛ α αₛ e E Φ:
    (M <= N)%nat →
    ▷ α ↪N (N; ns) ∗ ▷ αₛ ↪ₛN (M; nsₛ) ∗
    (∀ (n : nat),
       ⌜ n ≤ N ⌝ -∗
       if (bool_decide (n < S M)%nat)
         then
           α ↪N (N; ns ++ [n]) ∗ αₛ ↪ₛN (M; nsₛ ++ [n]) ∗ ⌜ n ≤ M ⌝ -∗
           WP e @ E {{ Φ }}
         else
           α ↪N (N; ns ++ [n]) ∗ αₛ ↪ₛN (M; nsₛ) ∗ ⌜ n ≤ N ⌝ -∗
           WP e @ E {{ Φ }}
    )
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (Hineq) "(>Hα & >Hαₛ & Hwp)".

    (*
    assert (∀ x : fin (S M), fin_to_nat x < S N)*)

    iApply (wp_couple_fragmented_rand_rand_inj (λ x, x) with "[$Hα $Hαₛ Hwp]"); [done| intros; lia |].
    iIntros (n) "%Hn".
    case_bool_decide as H1.
    - destruct H1 as [n' [Hn' <-]].
      iIntros (?) "(?&?&%&%)".
      iSpecialize ("Hwp" $! m).
      rewrite bool_decide_eq_true_2; last by lia.
      iApply "Hwp"; auto.
      by iFrame.
    - iSpecialize ("Hwp" $! n).
      case_match.
      2: { iIntros "(?&?&%)".
           iApply "Hwp"; auto.
           by iFrame. }
      exfalso. apply H1.
      exists n. split; auto.
      apply bool_decide_eq_true_1 in H.
      lia.
  Qed.

fragmented state rand N ~ state rand M, M>=N, under injective function from N to M

  Lemma wp_couple_fragmented_rand_rand_inj_rev {M N} (f: nat -> nat) {_: Inj (=) (=) f}
    ns nsₛ α αₛ e E Φ:
    (N <= M)%nat →
    (forall n, n < S N -> f n < S M)%nat ->
    ▷ α ↪N (N; ns) ∗ ▷ αₛ ↪ₛN (M; nsₛ) ∗
    (∀ (m : nat),
        ⌜ m ≤ M ⌝ -∗
        if bool_decide (∃ n:nat, n ≤ N /\ f n = m) then
          ∀ n, α ↪N (N; ns ++ [n]) ∗ αₛ ↪ₛN (M; nsₛ ++ [f n]) ∗ ⌜n ≤ N⌝ ∗ ⌜f n ≤ M⌝ -∗
               WP e @ E {{ Φ }}
        else
          α ↪N (N; ns) ∗ αₛ ↪ₛN (M; nsₛ++[m]) ∗ ⌜m ≤ M⌝ -∗ WP e @ E {{ Φ }})
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (Hineq Hdom) "(>Hα & >Hαₛ & Hwp)".
    edestruct (restr_inj_fin (S N) (S M) f (le_n_S N M Hineq) Hdom) as [g [HgInj HgEq]].
    iDestruct "Hα" as (fs) "(<-&Hα)".
    iDestruct "Hαₛ" as (fsₛ) "(<-&Hαₛ)".
    iApply wp_lift_step_spec_couple.
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & Hε2 & Hδ2)".
    iDestruct "Hauth2" as "(HK&Hh2&Ht2&Htl2)/=".
    iDestruct (ghost_map_lookup with "Ht2 Hαₛ") as %?.
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    replace (ε_now) with (0+ε_now)%NNR; last (apply nnreal_ext; simpl; lra).
    replace (δ_now) with (0+δ_now)%NNR; last (apply nnreal_ext; simpl; lra).
    iApply spec_coupl_erasables_weak; [done|done|..].
    { by apply ARcoupl_to_DPcoupl, ARcoupl_exact, Rcoupl_swap, Rcoupl_fragmented_rand_rand_inj. }
    { eapply erasable_dbind_predicate.
      - solve_distr_mass.
      - by eapply state_step_erasable.
      - apply dret_erasable. }
    { by eapply state_step_erasable. }
    iIntros (?? [m H']).
    case_bool_decide in H'.
    - destruct Hf as [m' <-].
      destruct H' as (n & ? & ? & Hfn).
      simplify_eq.
      iMod (ghost_map_update ((N; fs ++ [_]) : tape) with "Ht1 Hα") as "[Ht1 Hα]".
      iMod (ghost_map_update ((M; fsₛ ++ [_]) : tape) with "Ht2 Hαₛ") as "[Ht2 Hαₛ]".
      iModIntro. iApply spec_coupl_ret. iMod "Hclose'" as "_".
      iSpecialize ("Hwp" $! (f m')).
      rewrite bool_decide_eq_true_2.
      2: { exists m'.
           split; auto.
           apply fin_to_nat_le. }
      iSpecialize ("Hwp" $! _ m').
      iDestruct ("Hwp" with "[$Hα $Hαₛ]") as "Hwp".
      { iPureIntro.
        split; [rewrite fmap_app /= // |].
        split; [rewrite fmap_app /= -HgEq // |].
        split; auto.
        - apply fin_to_nat_le.
        - apply Nat.lt_succ_r, Hdom, fin_to_nat_lt.
      }
      assert (0 + ε_now = ε_now)%NNR as ->.
      { apply nnreal_ext; simpl; lra. }
      assert (0 + δ_now = δ_now)%NNR as ->.
      { apply nnreal_ext; simpl; lra. }
      by iFrame.
    - destruct H' as [??]. simplify_eq.
      iMod (ghost_map_update ((M; fsₛ ++ [m]) : tape) with "Ht2 Hαₛ") as "[Ht2 Hαₛ]".
      iModIntro. iApply spec_coupl_ret. iMod "Hclose'" as "_".
      iSpecialize ("Hwp" $! m).
      rewrite bool_decide_eq_false_2 //.
      2: {
        intros [n [Hn1 Hn2]].
        apply Hf.
        assert (n < S N)%nat as Hn3 by lia.
        exists (nat_to_fin Hn3).
        apply fin_to_nat_inj.
        rewrite HgEq -Hn2.
        rewrite fin_to_nat_to_fin //.
      }
      iDestruct ("Hwp" with "[]") as "Hwp"; [iPureIntro; apply fin_to_nat_le | ].
      assert (0 + ε_now = ε_now)%NNR as ->.
      { apply nnreal_ext; simpl; lra. }
      assert (0 + δ_now = δ_now)%NNR as ->.
      { apply nnreal_ext; simpl; lra. }
      iFrame.
      iApply "Hwp".
      iSplitL "Hα"; [by iFrame |].
      iSplitL; [|iPureIntro; apply fin_to_nat_le ].
      iFrame.
      rewrite fmap_app /= //.
      Unshelve.
      apply Nat.lt_succ_r, Hdom, fin_to_nat_lt.
  Qed.

fragmented state rand N ~ state rand M, M>=N, under injective function from N to M, but with errors for rejection sampling!
  Lemma wp_couple_fragmented_rand_rand_inj_rev' {M N} (f : nat -> nat) {_: Inj (=) (=) f}
    ns nsₛ α αₛ e E Φ (δ : R) :
    0 <= δ →
    (N < M)%nat →
    (forall n, n < S N -> f n < S M)%nat ->
    ▷ α ↪N (N; ns) ∗ ▷ αₛ ↪ₛN (M; nsₛ) ∗ ↯ δ ∗
    (∀ (m : nat),
       ⌜ m ≤ M ⌝ -∗
       if bool_decide (∃ n:nat, n ≤ N /\ f n = m) then
         ∀ n, α ↪N (N; ns ++ [n]) ∗ αₛ ↪ₛN (M; nsₛ ++ [f n]) ∗ ⌜ n ≤ N ⌝ ∗ ⌜ f n ≤ M ⌝ -∗
              WP e @ E {{ Φ }}
     else
       ∀ (δ' : R),
         ⌜δ' = (S M / (S M - S N) * δ)%R⌝ ∗
         α ↪N (N; ns) ∗ αₛ ↪ₛN (M; nsₛ++[m]) ∗ ↯ δ' ∗ ⌜ m ≤ M ⌝ -∗
         WP e @ E {{ Φ }})
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (Hδ Hineq Hdom) "(>Hα & >Hαₛ & Hδ & Hwp)".
    edestruct (restr_inj_fin (S N) (S M) f (le_n_S N M (Nat.lt_le_incl _ _ Hineq)) Hdom) as [g [HgInj HgEq]].
    iDestruct "Hα" as (fs) "(<-&Hα)".
    iDestruct "Hαₛ" as (fsₛ) "(<-&Hαₛ)".
    iApply wp_lift_step_spec_couple.
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & Hε2 & Hδ2)".
    iDestruct "Hauth2" as "(HK&Hh2&Ht2&Htl2)".
    iDestruct (ghost_map_lookup with "Ht2 Hαₛ") as %?.
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iDestruct (ec_supply_bound with "[$][$]") as %Hle.

    set δ' := mknonnegreal _ Hδ.

    set δ_now1 := nnreal_minus δ_now δ' Hle.
    set δ_now2 := (δ_now + δ' * nnreal_nat (S N) / nnreal_nat (S M - S N))%NNR.
    set (E2 σ := if bool_decide (∃ x, σ = state_upd_tapes <[αₛ:=(M; fsₛ ++ [g x])]> σ1')
                 then δ_now1 else δ_now2).
    assert (∀ σ, E2 σ <= Rmax δ_now1 δ_now2)%R.
    { intros ?. rewrite /E2. apply Rmax_Rle. case_bool_decide; by [left| right]. }

    iApply (spec_coupl_erasables_exp E2 (Rmax δ_now1 δ_now2) 0%NNR ε_now δ_now).
    { eapply ARcoupl_exact, Rcoupl_swap, Rcoupl_fragmented_rand_rand_inj; done || lia. }
    { apply erasable_dbind_predicate.
      - solve_distr_mass.
      - by eapply state_step_erasable.
      - apply dret_erasable. }
    { by eapply state_step_erasable. }
    { done. }
    { simpl. erewrite state_step_unfold; [|done].
      (* TODO: cleanup *)
      rewrite /Expval.
      erewrite (SeriesC_ext _
                  (λ b : state,
                      if bool_decide (b ∈ (λ x, state_upd_tapes <[αₛ:=(M; fsₛ ++ [x])]> σ1')
                                        <$> (fin_enum (S M)))
                     then /(S M) * E2 b else 0)%R); last first.
      { intros n.
        case_bool_decide as Hin; last first.
        { apply Rmult_eq_0_compat_r. rewrite /dmap/dbind/dbind_pmf/pmf/=.
          apply SeriesC_0. intros. apply Rmult_eq_0_compat_l.
          rewrite /dret_pmf. case_bool_decide; last lra.
          subst. exfalso. apply Hin. erewrite list_elem_of_fmap.
          exists x; split; first done. replace (fin_enum (S M)) with (enum (fin (S M))) by done.
          apply elem_of_enum. }
        rewrite list_elem_of_fmap in Hin. destruct Hin as [y [-> ?]].
        apply Rmult_eq_compat_r. rewrite /dmap/dbind/dbind_pmf/pmf/=.
        rewrite SeriesC_scal_l.
        replace (SeriesC _) with 1%R; first lra.
        symmetry.
        rewrite /dret_pmf.
        erewrite (SeriesC_ext _ (λ x, if bool_decide (x = y) then 1 else 0))%R;
          first apply SeriesC_singleton.
        intros.
        case_bool_decide as H2i.
        - apply state_upd_tapes_same in H2i. simplify_eq.
          case_bool_decide; done.
        - case_bool_decide; last done.
          subst. done. }
      { trans (SeriesC (λ x, if bool_decide (∃ y, g y = x) then / S M * δ_now1 else / S M * δ_now2))%R.
        - rewrite Rplus_0_l.
          set (h σ := match decide (∃ x, σ = state_upd_tapes <[αₛ:=(M; fsₛ ++ [x])]> σ1') with
                    | left Hproof => Some (epsilon Hproof)
                    | _ => None
                    end).
          etrans; last eapply (SeriesC_le_inj _ h).
          + apply Req_le_sym. apply SeriesC_ext.
should be ok
            intros s. rewrite /h. case_match eqn:Heqn; last first.
            { rewrite bool_decide_eq_false_2; first (simpl;lra).
              erewrite list_elem_of_fmap.
              intros [? [->?]]. apply n.
              naive_solver. }
            pose proof epsilon_correct _ e0 as H'.
            rewrite bool_decide_eq_true_2; last first.
            { destruct e0 as [x ?]. subst. rewrite list_elem_of_fmap.
              eexists _. split; first done.
              replace (fin_enum _) with (enum (fin (S M))) by done.
              apply elem_of_enum. }
            rewrite !S_INR.
            rewrite /E2.
            simpl in *. subst.
            case_bool_decide as H1'.
            * rewrite bool_decide_eq_true_2.
              { rewrite /δ_now1. simpl; lra. }
              destruct H1' as [y ?]. exists y. rewrite H3. done.
            * rewrite bool_decide_eq_false_2.
              { rewrite /δ_now2; simpl; lra. }
              intros [x H2'].
              apply H1'. rewrite H' in H2'. apply state_upd_tapes_same in H2'. simplify_eq.
              naive_solver.
          + intros. case_bool_decide; apply Rmult_le_pos; try done.
            all: rewrite <-Rdiv_1_l; apply Rcomplements.Rdiv_le_0_compat; try lra.
            all: apply pos_INR_S.
          + intros n1 n2 m. rewrite /h. do 2 case_match; try done.
            intros.
            pose proof epsilon_correct _ e0.
            pose proof epsilon_correct _ e1. simpl in *. simplify_eq.
            rewrite H7 H8. by repeat f_equal.
          + apply ex_seriesC_finite.
        - eset (diff:=elements (((list_to_set (enum (fin(S M)))):gset _ )∖ ((list_to_set(g<$>enum (fin(S N)))):gset _))).
          erewrite (SeriesC_ext _
                      (λ x : fin (S M), (if bool_decide (x ∈ g<$> enum (fin(S N))) then / S M * δ_now1 else 0%R) +
                                         if bool_decide (x ∈ diff ) then / S M * δ_now2 else 0%R
                   ))%R; last first.
          {
annoying lemma again
            intros n. rewrite /diff.
            case_bool_decide as H1'.
            - destruct H1' as [? H1']. rewrite bool_decide_eq_true_2; last first.
              + subst. apply list_elem_of_fmap_2. apply elem_of_enum.
              + subst. rewrite bool_decide_eq_false_2; first lra.
                rewrite elem_of_elements.
                rewrite not_elem_of_difference; right.
                rewrite elem_of_list_to_set. apply list_elem_of_fmap_2; apply elem_of_enum.
            - rewrite bool_decide_eq_false_2; last first.
              { rewrite list_elem_of_fmap. intros [?[??]].
                subst. apply H1'. naive_solver. }
              rewrite bool_decide_eq_true_2; first lra.
              rewrite elem_of_elements. rewrite elem_of_difference.
              split; rewrite elem_of_list_to_set; first apply elem_of_enum.
              rewrite list_elem_of_fmap. intros [?[??]].
              subst. apply H1'. naive_solver.
          }
        rewrite SeriesC_plus; try apply ex_seriesC_finite.
        rewrite !SeriesC_list_2; last first.
        { apply NoDup_fmap_2; [done|apply NoDup_enum]. }
        { rewrite /diff. eapply NoDup_elements. }
        rewrite length_fmap. rewrite fin.length_enum_fin.
        rewrite /diff.
        replace (length _) with (S M - S N)%nat; last first.
        { erewrite <-size_list_to_set; last apply NoDup_elements.
          erewrite list_to_set_elements.
          rewrite size_difference.
          - rewrite !size_list_to_set; [|apply NoDup_fmap; [auto|apply NoDup_enum]|apply NoDup_enum]; auto.
            rewrite length_fmap.
            rewrite !fin.length_enum_fin. done.
          - intros ??. apply elem_of_list_to_set. apply elem_of_enum.
        }
        rewrite /δ_now1 /δ_now2. simpl. rewrite -/(INR (S N)) -/(INR (S M)). rewrite !S_INR.
        rewrite !Rmult_assoc.
        rewrite minus_INR; last lia.
        cut ((N+1)/ (M + 1) * δ_now - (N+1)/(M+1) *δ+
               (M-N)/ (M + 1) * δ_now + ((N + 1)/(M+1) * ((M-N)/ (M - N))) * δ <= δ_now)%R; first lra.
        rewrite Rdiv_diag; last first.
        { assert (N < M)%R; real_solver. }
        cut ((N + 1) / (M + 1) * δ_now+ (M - N) / (M + 1) * δ_now <= δ_now)%R; first lra.
        cut ((M + 1) / (M + 1) * δ_now <= δ_now)%R; first lra.
        rewrite Rdiv_diag; first lra.
        pose proof pos_INR M. lra. }
      Unshelve. all : eapply gset_fin_set. }

    iIntros (?? [m H']).
    case_bool_decide in H' as H1'.
    - destruct H' as (n&?&?&?).
      destruct H1' as [n' <-].
      assert (n' = n) as -> by (by apply (inj _)).
      simplify_eq.
      iApply spec_coupl_ret.
      iMod (ghost_map_update ((N; fs ++ [n]) : tape) with "Ht1 Hα") as "[$ Hα]".
      iMod (ghost_map_update ((M; fsₛ ++ [g n]) : tape) with "Ht2 Hαₛ") as "[$ Hαₛ]".
      iModIntro. iMod "Hclose'" as "_".
      iSpecialize ("Hwp" $! (f n)).
      rewrite bool_decide_eq_true_2.
      2: { exists n.
           split; auto.
           apply fin_to_nat_le. }

      iSpecialize ("Hwp" $! _ (n)). iFrame.
      iDestruct ("Hwp" with "[$Hα $Hαₛ]") as "Hwp".
      { iPureIntro.
        split; [rewrite fmap_app /= // |].
        split; [rewrite fmap_app /= HgEq // |].
        split; [apply fin_to_nat_le | ].
        apply Nat.lt_succ_r, Hdom, fin_to_nat_lt.
      }
      replace (δ_now) with (δ' + δ_now1)%NNR; last first.
      { apply nnreal_ext. simpl. lra. }
      iMod (ec_supply_decrease with "[$] [$]") as (????) "H".
      iFrame.
      rewrite /E2 bool_decide_eq_true_2; [|eauto].
      iApply ec_supply_eq; [|done].
      simplify_eq /=. lra.

    - destruct H' as [??]. simplify_eq.
      replace (E2 _) with (δ_now2); last first.
      { rewrite /E2. rewrite bool_decide_eq_false_2 //.
        intros [? H2']. apply state_upd_tapes_same in H2'. simplify_eq. naive_solver. }
      destruct (Rle_or_lt 1 δ_now2).
      { iModIntro. by iApply spec_coupl_ret_err_ge_1. }
      iModIntro.
      iApply spec_coupl_ret.
      iMod (ghost_map_update ((M; fsₛ ++ [m]) : tape) with "Ht2 Hαₛ") as "[? Hαₛ]".
      iMod "Hclose'" as "_".
      iSpecialize ("Hwp" $! m).
      rewrite bool_decide_eq_false_2 //.
      2: {
        intros [n [Hn1 Hn2]].
        apply H1'.
        assert (n < S N)%nat as Hn3 by lia.
        exists (nat_to_fin Hn3).
        apply fin_to_nat_inj.
        rewrite HgEq -Hn2.
        rewrite fin_to_nat_to_fin //.
      }
      rewrite !S_INR /=.
      iFrame.
      iMod (ec_supply_increase with "[$Hδ2]") as "[$ Hδ']".
      { by eapply Rle_lt_trans. }
      iCombine "Hδ Hδ'" as "Hδ".
      iApply ("Hwp" $! _ with "[$Hα $Hαₛ $Hδ]").
      iPureIntro.
      split.
      1:{
        simpl. rewrite -/(INR (S N)). rewrite S_INR.
        replace (INR M + 1 - (INR N + 1))%R with (INR M - INR N)%R by lra.
        rewrite -{1}(Rmult_1_l δ).
        rewrite Rmult_assoc (Rmult_comm δ).
        rewrite -Rmult_plus_distr_r. apply Rmult_eq_compat_r.
        rewrite Rdiv_def.
        replace (1)%R with ((INR M - INR N)*/(INR M - INR N))%R at 1; last first.
        { apply Rinv_r. apply lt_INR in Hineq. lra. }
        rewrite minus_INR; [|real_solver].
        rewrite -Rmult_plus_distr_r. lra. }
      split; auto.
      split; [ | apply fin_to_nat_le ].
      rewrite fmap_app //.
      Unshelve.
      + apply Nat.lt_succ_r, Hdom, fin_to_nat_lt.
      + apply fin_to_nat_le.
  Qed.

 Lemma wp_couple_fragmented_rand_rand_leq_rev' {M N : nat} ns nsₛ α αₛ e E Φ (ε : R) :
       0 <= ε →
       (N < M)%nat →
       ▷ α ↪N (N; ns) ∗ ▷ αₛ ↪ₛN (M; nsₛ) ∗ ↯ ε ∗
       (∀ (m : nat),
          ⌜ m ≤ M ⌝ -∗
          if bool_decide (m < S N)%nat
            then
              α ↪N (N; ns ++ [m]) ∗ αₛ ↪ₛN (M; nsₛ ++ [m]) ∗ ⌜ m ≤ N ⌝ -∗
              WP e @ E {{ Φ }}
            else
              ∀ (ε' : R),
                ⌜ε' = (S M / (S M - S N) * ε)%R⌝ ∗
                α ↪N (N; ns) ∗ αₛ ↪ₛN (M; nsₛ++[m]) ∗ ↯ ε' ∗ ⌜ m ≤ M ⌝ -∗
                WP e @ E {{ Φ }}
         )
       ⊢ WP e @ E {{ Φ }}.
     Proof.
       iIntros (Hε Hineq) "(>Hα & >Hαₛ & Hε & Hwp)".
       iApply (wp_couple_fragmented_rand_rand_inj_rev' (λ x, x) with "[$Hα $Hαₛ $Hε Hwp]"); [done|done| .. ].
       { intros; lia. }
       iIntros (n) "%Hn".
       case_bool_decide as H1.
       - destruct H1 as [n' [Hn' <-]].
         iIntros (?) "(?&?&%&%)".
         iSpecialize ("Hwp" $! n with "[//]").
         rewrite bool_decide_eq_true_2; last by lia.
         iSpecialize ("Hwp" with "[-]"); iFrame.
         done.
       - iSpecialize ("Hwp" $! n with "[//]").
         rewrite bool_decide_eq_false_2; [iFrame |].
         intro. apply H1.
         exists n. lia.
     Qed.

wp_couple_exp
  Lemma wp_couple_exp (M N p:nat)
    (f : list nat → nat)
    (Hdom: forall (l : list nat), Forall (λ x, (x < S N)%nat) l -> (f l < S M)%nat)
    (Hinj: forall (l1 l2:list nat),
        Forall (λ x, (x < S N)%nat) l1 ->
        Forall (λ x, (x < S N)%nat) l2 ->
        length l1 = p -> length l2 = p -> f l1 = f l2 -> l1 = l2) ns nsₛ α αₛ e E Φ:
    (S N ^ p = S M)%nat->
    ▷ α ↪N (N; ns) ∗ ▷ αₛ ↪ₛN (M; nsₛ) ∗
    (∀ (l : list nat) (m: nat),
       ⌜Forall (λ x, (x < S N)%nat) l⌝ -∗
       ⌜(m < S M)%nat ⌝ -∗
       ⌜length l = p /\ f l = m⌝ -∗
       α ↪N (N; ns ++ l) -∗ αₛ ↪ₛN (M; nsₛ ++ [m]) -∗
       WP e @ E {{ Φ }}
    )
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (Heq) "(>Hα & >Hαₛ & Hwp)".
    iDestruct "Hα" as (fs) "(<-&Hα)".
    iDestruct "Hαₛ" as (fsₛ) "(<-&Hαₛ)".
    destruct (restr_list_inj_fixed_length (S N) (S M) p f) as [g [Hg1 Hg2]]; auto.
    iApply wp_lift_step_spec_couple.
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & Hε2 & Hδ2)".
    iDestruct "Hauth2" as "(HK&Hh2&Ht2&Htl2)/=".
    iDestruct (ghost_map_lookup with "Ht2 Hαₛ") as %?.
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?.
    simplify_map_eq.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    replace (ε_now) with (0+ε_now)%NNR; last (apply nnreal_ext; simpl; lra).
    replace (δ_now) with (0+δ_now)%NNR; last (apply nnreal_ext; simpl; lra).
    iApply spec_coupl_erasables_weak; [done|done|..].
    - apply ARcoupl_to_DPcoupl, ARcoupl_exact. apply Rcoupl_state_state_exp.
      all: exact.
    - eapply iterM_state_step_erasable; done.
    - by eapply state_step_erasable.
    - iIntros (σ2 σ2') "%K".
      destruct K as (xs' & z & Hlen & -> & -> & <-).
      iMod (ghost_map_update ((N; fs ++ xs') : tape) with "Ht1 Hα") as "[Ht1 Hα]".
      iMod (ghost_map_update ((M; fsₛ ++ [g xs']) : tape) with "Ht2 Hαₛ") as "[Ht2 Hαₛ]".
      iModIntro.
      iApply spec_coupl_ret.
      iMod "Hclose'".
      iSpecialize ("Hwp" $! (fin_to_nat <$> xs') (g xs') with "[][][]").
      + iPureIntro.
        apply fin_forall_leq.
      + iPureIntro. apply fin_to_nat_lt.
      + iPureIntro; split; auto.
        rewrite length_fmap //.
      + iFrame.
        replace (0+_)%NNR with (ε_now).
        2:{ apply nnreal_ext. simpl; lra. }
        replace (0+_)%NNR with (δ_now).
        2:{ apply nnreal_ext. simpl; lra. }
        iFrame.
        iApply ("Hwp" with "[Hα][-]").
        * rewrite -fmap_app. iFrame. done.
        * iFrame. rewrite fmap_app. done.
  Qed.

  Lemma wp_couple_exp_decoder (M N p : nat ) ns nsₛ α αₛ e E Φ:
    (S N ^ p = S M)%nat ->
    ▷ α ↪N (N; ns) ∗ ▷ αₛ ↪ₛN (M; nsₛ) ∗
    (∀ (l : list nat) (m: nat),
       ⌜Forall (λ x, (x < S N)%nat) l⌝ -∗
       ⌜(m < S M)%nat ⌝ -∗
       ⌜length l = p⌝ -∗
       α ↪N (N; ns ++ l) -∗ αₛ ↪ₛN (M; nsₛ ++ [@decoder_nat N l]) -∗
       WP e @ E {{ Φ }}
    )
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (Heq) "(>Hα & >Hαₛ & Hwp)".
    set (f := (λ l : list nat, if bool_decide (length l = p) then @decoder_nat N l else 0%nat )).
    iApply (wp_couple_exp M N p f); auto.
    {
      rewrite -Heq /f.
      intros l Hl.
      case_bool_decide; simplify_eq.
      - by apply fin.decoder_aux_ineq.
      - lia.
    }
    {
      rewrite /f.
      intros ?????<-.
      case_bool_decide; [|done].
      rewrite bool_decide_eq_true_2; auto.
      intro.
      by apply (@fin.decoder_aux_inj N).
    }
    iFrame.
    rewrite /f.
    iIntros (??) "% % (%&%) Hα Hαs".
    case_bool_decide; [|done].
    iApply ("Hwp" with "[//] [//] [//] Hα [Hαs]").
    simplify_eq.
    iFrame.
  Qed.

  Lemma wp_couple_exp_decoder_lr (M N p : nat ) ns nsₛ α αₛ e E Φ:
    (S N ^ p = S M)%nat ->
    ▷ α ↪N (N; ns) ∗ ▷ αₛ ↪ₛN (M; nsₛ) ∗
    (∀ (l : list nat) (m: nat),
       ⌜Forall (λ x, (x < S N)%nat) l⌝ -∗
       ⌜(m < S M)%nat ⌝ -∗
       ⌜length l = p⌝ -∗
       α ↪N (N; ns ++ l) -∗ αₛ ↪ₛN (M; nsₛ ++ [@decoder_nat_lr N l]) -∗
       WP e @ E {{ Φ }}
    )
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (Heq) "(>Hα & >Hαₛ & Hwp)".
    set (f := (λ l : list nat, if bool_decide (length l = p) then @decoder_nat_lr N l else 0%nat )).
    iApply (wp_couple_exp M N p f); auto.
    {
      rewrite -Heq /f.
      intros l Hl.
      case_bool_decide; simplify_eq.
      - by apply fin.decoder_lr_aux_ineq.
      - lia.
    }
    {
      rewrite /f.
      intros ?????<-.
      case_bool_decide; [|done].
      rewrite bool_decide_eq_true_2; auto.
      intro.
      by apply (@fin.decoder_lr_aux_inj N).
    }
    iFrame.
    rewrite /f.
    iIntros (??) "% % (%&%) Hα Hαs".
    case_bool_decide; [|done].
    iApply ("Hwp" with "[//] [//] [//] Hα [Hαs]").
    simplify_eq.
    iFrame.
  Qed.

  Lemma wp_couple_exp_rev (M N p:nat)
    (f:(list nat) -> nat)
    (Hdom: forall (l : list nat), Forall (λ x, (x < S N)%nat) l -> (f l < S M)%nat)
    (Hinj: forall (l1 l2:list nat),
        Forall (λ x, (x < S N)%nat) l1 ->
        Forall (λ x, (x < S N)%nat) l2 ->
        length l1 = p -> length l2 = p -> f l1 = f l2 -> l1 = l2) ns nsₛ α αₛ e E Φ:
    (S N ^ p = S M)%nat->
    ▷ α ↪N (M; ns) ∗ ▷ αₛ ↪ₛN (N; nsₛ) ∗
    (∀ (l : list nat) (m: nat),
       ⌜Forall (λ x, (x < S N)%nat) l⌝ -∗
       ⌜(m < S M)%nat ⌝ -∗
       ⌜length l = p /\ f l = m⌝ -∗
       α ↪N (M; ns ++ [m]) -∗ αₛ ↪ₛN (N; nsₛ ++ l) -∗
       WP e @ E {{ Φ }}
    )
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (Heq) "(>Hα & >Hαₛ & Hwp)".
    iApply wp_lift_step_spec_couple.
    iDestruct "Hα" as (fs) "(<-&Hα)".
    iDestruct "Hαₛ" as (fsₛ) "(<-&Hαₛ)".
    destruct (restr_list_inj_fixed_length (S N) (S M) p f) as [g [Hg1 Hg2]]; auto.
    iIntros (σ1 e1' σ1' ε_now δ_now) "((Hh1 & Ht1 & Htl1) & Hauth2 & Hε2 & Hδ2)".
    iDestruct "Hauth2" as "(HK&Hh2&Ht2&Htl2)/=".
    iDestruct (ghost_map_lookup with "Ht2 Hαₛ") as %?.
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?.
    simplify_map_eq.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    replace (ε_now) with (0+ε_now)%NNR; last (apply nnreal_ext; simpl; lra).
    replace δ_now with (0 + δ_now)%NNR at 2 ; [|apply nnreal_ext; simpl; lra].
    iApply spec_coupl_erasables_weak; [done|done|..].
    - apply ARcoupl_to_DPcoupl, ARcoupl_exact. apply Rcoupl_swap. apply Rcoupl_state_state_exp.
      all: exact.
    - by eapply state_step_erasable.
    - eapply iterM_state_step_erasable; done.
    - iIntros (σ2 σ2') "%K".
      destruct K as (xs' & z & Hlen & -> & -> & <-).
      iMod (ghost_map_update ((M; fs ++ [g xs']) : tape) with "Ht1 Hα") as "[Ht1 Hα]".
      iMod (ghost_map_update ((N; fsₛ ++ xs') : tape) with "Ht2 Hαₛ") as "[Ht2 Hαₛ]".
      iModIntro.
      iApply spec_coupl_ret.
      iSpecialize ("Hwp" $! (fin_to_nat <$> xs') (g xs') with "[][][]").
      + iPureIntro.
        apply fin_forall_leq.
      + iPureIntro. apply fin_to_nat_lt.
      + iPureIntro; split; auto.
        rewrite length_fmap //.
      + iFrame.
        replace (0+_)%NNR with (ε_now).
        2:{ apply nnreal_ext. simpl; lra. }
        iFrame.
        iMod "Hclose'".
        iApply ("Hwp" with "[Hα][-]").
        * iFrame. rewrite fmap_app //.
        * iFrame. rewrite -fmap_app //.
  Qed.

  Lemma wp_couple_exp_decoder_rev (M N p:nat) ns nsₛ α αₛ e E Φ:
    (S N ^ p = S M)%nat->
    ▷ α ↪N (M; ns) ∗ ▷ αₛ ↪ₛN (N; nsₛ) ∗
    (∀ (l : list nat) (m: nat),
       ⌜Forall (λ x, (x < S N)%nat) l⌝ -∗
       ⌜(m < S M)%nat ⌝ -∗
       ⌜length l = p⌝ -∗
       α ↪N (M; ns ++ [@decoder_nat N l]) -∗ αₛ ↪ₛN (N; nsₛ ++ l) -∗
       WP e @ E {{ Φ }}
    )
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (Heq) "(>Hα & >Hαₛ & Hwp)".
    set (f := (λ l : list nat, if bool_decide (length l = p) then @decoder_nat N l else 0%nat )).
    iApply (wp_couple_exp_rev M N p f); auto.
    {
      rewrite -Heq /f.
      intros l Hl.
      case_bool_decide; simplify_eq.
      - by apply fin.decoder_aux_ineq.
      - lia.
    }
    {
      rewrite /f.
      intros ?????<-.
      case_bool_decide; [|done].
      rewrite bool_decide_eq_true_2; auto.
      intro.
      by apply (@fin.decoder_aux_inj N).
    }
    iFrame.
    rewrite /f.
    iIntros (??) "% % (%&%) Hα Hαs".
    case_bool_decide; [|done].
    simplify_eq.
    iApply ("Hwp" with "[//] [//] [//] Hα [Hαs]").
    iFrame.
  Qed.

  Lemma wp_couple_exp_decoder_lr_rev (M N p:nat) ns nsₛ α αₛ e E Φ:
    (S N ^ p = S M)%nat->
    ▷ α ↪N (M; ns) ∗ ▷ αₛ ↪ₛN (N; nsₛ) ∗
    (∀ (l : list nat) (m: nat),
       ⌜Forall (λ x, (x < S N)%nat) l⌝ -∗
       ⌜(m < S M)%nat ⌝ -∗
       ⌜length l = p⌝ -∗
       α ↪N (M; ns ++ [@decoder_nat_lr N l]) -∗ αₛ ↪ₛN (N; nsₛ ++ l) -∗
       WP e @ E {{ Φ }}
    )
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (Heq) "(>Hα & >Hαₛ & Hwp)".
    set (f := (λ l : list nat, if bool_decide (length l = p) then @decoder_nat_lr N l else 0%nat )).
    iApply (wp_couple_exp_rev M N p f); auto.
    {
      rewrite -Heq /f.
      intros l Hl.
      case_bool_decide; simplify_eq.
      - by apply fin.decoder_lr_aux_ineq.
      - lia.
    }
    {
      rewrite /f.
      intros ?????<-.
      case_bool_decide; [|done].
      rewrite bool_decide_eq_true_2; auto.
      intro.
      by apply (@fin.decoder_lr_aux_inj N).
    }
    iFrame.
    rewrite /f.
    iIntros (??) "% % (%&%) Hα Hαs".
    case_bool_decide; [|done].
    simplify_eq.
    iApply ("Hwp" with "[//] [//] [//] Hα [Hαs]").
    iFrame.
  Qed.

Exact couplings

  Lemma wp_couple_tape_rand N f `{Bij nat nat f} K E α z ns Φ e :
    TCEq N (Z.to_nat z) →
    (∀ n, n < S N -> f n < S N)%nat →
    ▷ α ↪N (N; ns) ∗ ⤇ fill K (rand #z) ∗
    (∀ n : nat, α ↪N (N; ns ++ [n]) ∗ ⤇ fill K #(f n) ∗ ⌜ n ≤ N ⌝ -∗ WP e @ E {{ Φ }})
    ⊢ WP e @ E {{ Φ }}.
  Proof.
    iIntros (H0 Hdom) "(>Hα & Hj & Hwp)".
    iDestruct "Hα" as (fs) "(<-&Hα)".
    destruct (restr_bij_fin (S N) f Hdom) as [ff [Hbij Hff]].
    iApply wp_lift_step_spec_couple.
    iIntros (σ1 e1' σ1' ε δ) "[[Hh1 [Ht1 Htl1]] [Hauth2 Herr]]".
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?.
    iDestruct (spec_auth_prog_agree with "Hauth2 Hj") as %-> .
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    replace ε with (0 + ε)%NNR at 2 ; [|apply nnreal_ext; simpl; lra].
    replace δ with (0 + δ)%NNR at 2 ; [|apply nnreal_ext; simpl; lra].
    iApply (spec_coupl_erasable_steps 1 _ (state_step σ1 α)); [done|done|..].
    { rewrite pexec_1 step_or_final_no_final; last first.
      { apply reducible_not_final. solve_red. }
      apply DPcoupl_steps_ctx_bind_r => /=; [done|].
      apply ARcoupl_to_DPcoupl, ARcoupl_exact, Rcoupl_pos_R, (Rcoupl_state_rand N ff); eauto.
      rewrite -H0 //.
    }
    { by eapply state_step_erasable. }
    iIntros(σ2 e2' σ2' (? & [= ->] & (? & -> & [= -> ->]) & ? & ?)).
    iApply spec_coupl_ret.
    iMod (spec_update_prog (fill K #(ff _)) with "Hauth2 Hj") as "[$ Hspec0]".
    iDestruct (ghost_map_lookup with "Ht1 Hα") as %?%lookup_total_correct.
    iMod (ghost_map_update ((_; fs ++ [_]) : tape) with "Ht1 Hα") as "[$ Hα]".
    iModIntro. iMod "Hclose'" as "_".
    iSpecialize ("Hwp" with "[Hα Hspec0]").
    {
      iFrame.
      iSplitR.
      - iPureIntro.
        rewrite fmap_app //.
      - rewrite Hff. iFrame.
        iPureIntro. apply fin_to_nat_le.
    }
    iFrame.
    done.
  Qed.

rand(unit, N) ~ state_step(α', N) coupling

  Lemma wp_couple_rand_tape N f `{Bij nat nat f} z E α ns :
    TCEq N (Z.to_nat z) →
    (∀ n, n < S N -> f n < S N)%nat →
    {{{ ▷ α ↪ₛN (N; ns) }}}
      rand #z @ E
    {{{ (n : nat), RET #n; α ↪ₛN (N; ns ++ [f n]) ∗ ⌜ n ≤ N ⌝ }}}.
  Proof.
    iIntros (H0 Hdom ?) ">Hαs Hwp".
    iDestruct "Hαs" as (fs) "(<-&Hαs)".
    destruct (restr_bij_fin (S N) f Hdom) as [ff [Hbij Hff]].
    iApply wp_lift_prim_step_l_erasable; [done|].
    iIntros (σ1 e1' σ1' ε δ) "[(Hh1 & Ht1 & ?) [Hauth2 [Hε2 Hδ2]]]".
    iDestruct (spec_auth_lookup_tape with "Hauth2 Hαs") as %?.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose'".
    iExists _, (state_step σ1' α), 0%NNR, ε, 0%NNR, δ.
    repeat iSplit.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. solve_red.
    - iPureIntro. apply ARcoupl_to_DPcoupl, ARcoupl_exact.
      eapply (Rcoupl_rand_state _ ff); eauto.
      rewrite -H0//.
    - iPureIntro. by eapply state_step_erasable.
    - iIntros (e2 σ2 σ2' (n & [= -> ->] & ->)).
      iMod (spec_auth_update_tape (_; fs ++ [ff _]) with "Hauth2 Hαs") as "[Htapes Hαs]".
      do 2 iModIntro.
      iMod "Hclose'" as "_".
      iFrame.
      iApply wp_value.
      iApply ("Hwp" $! _ with "[$Hαs]").
      iPureIntro.
      rewrite fmap_app -Hff.
      split; auto.
      apply fin_to_nat_le.
  Qed.

  Lemma wp_couple_rand_rand_lbl N f `{Bij nat nat f} z K E α :
    TCEq N (Z.to_nat z) →
    (∀ n, n < S N -> f n < S N)%nat →
    {{{ α ↪ₛN (N; []) ∗ ⤇ fill K (rand(#lbl:α) #z) }}}
      rand #z @ E
      {{{ (n : nat), RET #n; α ↪ₛN (N; []) ∗ ⤇ fill K #(f n) ∗ ⌜ n ≤ N ⌝ }}}.
  Proof.
    iIntros (-> Hdom ?) "(Hα & Hspec) Hwp".
    iApply wp_spec_update.
    iApply (wp_couple_rand_tape with "[$Hα]"); auto.
    iIntros "!>" (n) "[Hα %Hn]".
    simpl.
    iDestruct (read_spec_tape_head with "Hα") as (x xs) "[Hα [%Hrw Hret]]" .
    iMod (step_rand with "[$]") as "[? ?]".
    iModIntro.
    iApply ("Hwp" with "[-]").
    iSpecialize ("Hret" with "[$]").
    rewrite Hrw.
    by iFrame.
  Qed.

  Lemma wp_couple_rand_lbl_rand_lbl N f `{Bij nat nat f} z K E α α' :
    TCEq N (Z.to_nat z) →
    (∀ n, n < S N -> f n < S N)%nat →
    {{{ ▷ α ↪N (N; []) ∗ ▷ α' ↪ₛN (N; []) ∗ ⤇ fill K (rand(#lbl:α') #z) }}}
      rand(#lbl:α) #z @ E
    {{{ (n : nat), RET #n; α ↪N (N; []) ∗ α' ↪ₛN (N; []) ∗ ⤇ fill K #(f n) ∗ ⌜ n ≤ N ⌝ }}}.
  Proof.
    iIntros (???) "(>Hα & >Hαs & Hr) HΨ".
    iMod ec_zero.
    iApply (wp_couple_tapes_bij N f); auto.
    iFrame.
    iIntros (n) "(Hα & Hαs & %) /=".
    iDestruct (read_spec_tape_head with "Hαs") as (x xs) "[Hαs [%Hrw Hret]]" .
    iMod (step_rand with "[$Hr $Hαs]") as "[? ?]".
    iApply (wp_rand_tape with "Hα").
    iIntros "!> (Hα & %)".
    iApply ("HΨ").
    iSpecialize ("Hret" with "[$]").
    rewrite Hrw.
    by iFrame.
  Qed.

rand(α, N) ~ rand(α, N) wrong bound coupling

  Lemma wp_couple_rand_lbl_rand_lbl_wrong N M f `{Bij nat nat f} z K E α α' xs ys :
    TCEq N (Z.to_nat z) →
    N ≠ M →
    (∀ n, n < S N -> f n < S N)%nat →
    {{{ ▷ α ↪N (M; xs) ∗ ▷ α' ↪ₛN (M; ys) ∗ ⤇ fill K (rand(#lbl:α') #z) }}}
      rand(#lbl:α) #z @ E
    {{{ (n : nat), RET #n; α ↪N (M; xs) ∗ α' ↪ₛN (M; ys) ∗ ⤇ fill K #(f n) ∗ ⌜ n ≤ N ⌝ }}}.
  Proof.
    iIntros (-> ? Hdom Ψ) "(>Hα & >Hαs & Hr) Hwp".
    iApply wp_lift_prim_steps_coupl; [done|].
    iIntros (σ1 e1' σ1' ε δ) "[[Hh [Ht Htl]] [Hs Herr]]".
    iDestruct "Hα" as (fs) "(<-&Hα)".
    iDestruct "Hαs" as (fsₛ) "(<-&Hαs)".
    destruct (restr_bij_fin _ f Hdom) as [g [HBij Hfg]]; auto.
    iDestruct (ghost_map_lookup with "Ht Hα") as %?.
    iDestruct (spec_auth_lookup_tape with "Hs Hαs") as %?.
    iDestruct (spec_auth_prog_agree with "Hs Hr") as %->.
    iApply fupd_mask_intro; [set_solver|]; iIntros "Hclose".
    iExists (λ ρ2 ρ2',
               ∃ (n : fin _), ρ2 = (Val #n, σ1) ∧ ρ2' = (fill K #(f n), σ1')),
             0%NNR, ε, 0%NNR, δ.
    repeat iSplit.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. apply nnreal_ext; simpl; lra.
    - iPureIntro. solve_red.
    - iPureIntro. solve_red.
    - iPureIntro.
      rewrite /= fill_dmap //.
      rewrite -(dret_id_right (prim_step _ _)) /=.
      apply ARcoupl_to_DPcoupl, ARcoupl_exact.
      apply Rcoupl_dmap.
      eapply Rcoupl_mono; [by eapply (Rcoupl_rand_lbl_rand_lbl_wrong _ _ g)|].
      intros [] [] (b & [=] & [=])=>/=.
      simplify_eq. rewrite Hfg. eauto.
    - iIntros (e2 σ2 e2' σ2' (b & [= -> ->] & [= -> ->])) "!>".
      iModIntro.
      iMod (spec_update_prog with "Hs Hr") as "[$ Hr]".
      iFrame.
      iMod "Hclose" as "_".
      iApply wp_value.
      iApply ("Hwp" with "[-]").
      iFrame.
      iPureIntro.
      split; auto.
      split; auto.
      apply fin_to_nat_le.
  Qed.

End rules.