diff --git a/theories/analysis/RealSeries.ec b/theories/analysis/RealSeries.ec
index ac61c29b7d..15c0cd39a8 100644
--- a/theories/analysis/RealSeries.ec
+++ b/theories/analysis/RealSeries.ec
@@ -120,7 +120,7 @@ pose R := big predT (fun x => `|s2 x|) (undup I).
 exists (M+R) => J uniq_J.
 rewrite (@bigEM (mem I)) addrC &(ler_add).
 - rewrite (@eq_bigr _ _ (fun x => `|s1 x|)) //= 1:/#.
-  by rewrite -big_filter; smt(filter_uniq).
+  rewrite -big_filter. by apply /sum_s1 /filter_uniq.
 - have P : perm_eq (filter (mem I) J) (filter (mem J) (undup I)).
     by apply: uniq_perm_eq; smt(filter_uniq undup_uniq mem_filter mem_undup).
   rewrite -big_filter (eq_big_perm P) big_filter big_mkcond.
diff --git a/theories/crypto/SplitRO.ec b/theories/crypto/SplitRO.ec
index 9182e2e078..bbd552d338 100644
--- a/theories/crypto/SplitRO.ec
+++ b/theories/crypto/SplitRO.ec
@@ -150,18 +150,12 @@ section PROOFS.
   proof.
     proc; inline *.
     swap{2} 5 -3; swap{2} 6 -2; sp 0 2.
-    seq 1 2 : (#pre /\ r{1} = ofpair (r{2}, r0{2})).
+    alias{2} 3 rr = ofpair (r, r0).
+    seq 1 3 : (#pre /\ r{1} = ofpair (r{2}, r0{2})).
     + conseq />.
-      alias{2} 2 one = 1. swap{2} 2 1.
-      alias{2} 3 rr = ofpair (r, r0).
-      kill{2} 4 ! 1; first by auto.
       transitivity{2}
-        {
-          r <$ sampleto1 x;
-          r0 <$ sampleto2 x;
-          rr <- ofpair (r, r0);
-        }
-        (={x} /\ (x0 = x /\ x1 = x){2} ==> r{1} = rr{2})
+        { r <$ sampleto1 x; r0 <$ sampleto2 x; rr <- ofpair (r, r0); }
+        (={x} ==> r{1} = rr{2})
         (={x} /\ (x0 = x /\ x1 = x){2} ==> rr{1} = ofpair(r, r0){2}); 1,2,4: by auto => /#.
       rndsem*{2} 0.
       by auto => *; rewrite -dmap_dprodE sample_spec /#.
diff --git a/theories/datatypes/FSet.ec b/theories/datatypes/FSet.ec
index 335924ccd1..f2a8bff4f2 100644
--- a/theories/datatypes/FSet.ec
+++ b/theories/datatypes/FSet.ec
@@ -466,6 +466,9 @@ proof. by rewrite fsetDUl fsetDv fsetU0. qed.
 lemma fsetDKv (A B : 'a fset) : (A `&` B) `\` B = fset0.
 proof. by rewrite fsetDIl fsetDv fsetI0. qed.
 
+lemma fsetDID (A B : 'a fset) : A `\` B = A `\` (B `&` A).
+proof. by rewrite fsetDIr fsetDv fsetU0. qed.
+
 (* -------------------------------------------------------------------- *)
 lemma subsetIl (A B : 'a fset) : (A `&` B) \subset A.
 proof. by apply/subsetP=> x; rewrite inE; case. qed.
@@ -747,6 +750,17 @@ proof.
   by apply/filter_uniq/uniq_elems.
 qed.
 
+lemma foldU (a : 'a) (f : 'a -> 'b -> 'b) (z : 'b) (A : 'a fset):
+  (forall a a' b, f a (f a' b) = f a' (f a b)) =>
+  ! mem A a =>
+  fold f z (A `|` fset1 a) = f a (fold f z A).
+proof.
+move=> assoc_f. rewrite (foldC a _ _ (A `|` fset1 a) assoc_f). 
+- by apply in_fsetU; right; exact in_fset1. 
+rewrite fsetDK fsetDID fset1I.
+move=> />. by rewrite fsetD0.
+qed.
+
 (* -------------------------------------------------------------------- *)
 
 op rangeset (m n : int) = oflist (range m n).
diff --git a/theories/datatypes/List.ec b/theories/datatypes/List.ec
index 7d8d0c2aa0..3a2e99e959 100644
--- a/theories/datatypes/List.ec
+++ b/theories/datatypes/List.ec
@@ -1605,7 +1605,7 @@ by rewrite revK size_rev opprD opprK addrA.
 qed.
 
 (* -------------------------------------------------------------------- *)
-(*                        Duplicate-freenes                             *)
+(*                        Duplicate-freeness                            *)
 (* -------------------------------------------------------------------- *)
 op uniq (s : 'a list) =
   with s = []      => true
@@ -2164,6 +2164,17 @@ lemma mapiP x0 (f : int -> 'a -> 'b) (s : 'a list) y :
     exists n, (0 <= n && n < size s) /\ y = f n (nth x0 s n).
 proof. exact: mapi_recP. qed.
 
+lemma mapi_cat (f : int -> 'a -> 'b) l l' :
+  mapi f (l ++ l') = mapi f l ++ mapi_rec f l' (size l).
+proof.
+elim/last_ind: l l' => //.
+smt(size_rcons cats1 cat_rcons).
+qed.
+
+lemma mapi_rcons (f : int -> 'a -> 'b) l x :
+  mapi f (rcons l x) = rcons (mapi f l) (f (size l) x).
+proof. smt(mapi_cat cats1). qed.
+
 (* -------------------------------------------------------------------- *)
 (*                          Element Replacement                         *)
 (* -------------------------------------------------------------------- *)
diff --git a/theories/datatypes/Xreal.ec b/theories/datatypes/Xreal.ec
index a395b40172..f36ce48080 100644
--- a/theories/datatypes/Xreal.ec
+++ b/theories/datatypes/Xreal.ec
@@ -956,6 +956,15 @@ proof.
 qed.
 
 (* -------------------------------------------------------------------- *)
+lemma Ep_bool_gen (d : bool distr) (f : bool -> xreal) :
+  Ep d f = mu1 d false ** f false + mu1 d true ** f true.
+proof.
+have -> := Ep_fin [false; true] d f //; 1: by move => [] //.
+rewrite /big.
+by have -> : map (d ** f) (filter predT [false; true]) =
+  [(d ** f) false; (d ** f) true].
+qed.
+
 lemma Ep_dbool (f : bool -> xreal) :
   Ep {0,1} f = of_reald 0.5 ** f true + of_reald 0.5 ** f false.
 proof.
diff --git a/theories/distributions/Distr.ec b/theories/distributions/Distr.ec
index f4c8ddf94e..0d3811a497 100644
--- a/theories/distributions/Distr.ec
+++ b/theories/distributions/Distr.ec
@@ -64,6 +64,20 @@ inductive isdistr (m : 'a -> real) =
        (forall x, 0%r <= m x)
      & (forall s, uniq s => big predT m s <= 1%r).
 
+lemma isdistr_summable_equiv (m : 'a -> real) :
+  isdistr m <=> (forall x, 0%r <= m x) /\ summable m /\ sum m <= 1%r.
+proof.
+rewrite/bounded. split => [ [ * ] | [ ? [ * ] ] ]; last first.
+- split; 1: assumption.
+  move => ? ?. apply (ler_trans (sum m)) => //.
+  by apply ler_big_sum => //. 
+have sumM : summable m.
+- exists 1%r.
+  have <- : (fun i => `|m i|) = m; smt().
+do 2 ! (split => //).
+by apply (lerfin_sum sumM).
+qed.
+
 lemma eq_isdistr (d1 d2 : 'a -> real) :
   d1 == d2 => isdistr d1 = isdistr d2.
 proof. by move=> /fun_ext=> ->. qed.
diff --git a/theories/distributions/Geometric.ec b/theories/distributions/Geometric.ec
new file mode 100644
index 0000000000..82477e1617
--- /dev/null
+++ b/theories/distributions/Geometric.ec
@@ -0,0 +1,139 @@
+require import AllCore Real RealSeq RealSeries Xreal StdOrder StdBigop List.
+import RField RealOrder Xreal Distr DBool Bigreal.BRA.
+
+section Parametrized.
+(* The main parameter, generally the probability of success *)
+declare op p : { real | 0%r < p <= 1%r } as in01_p.
+
+(* Probability mass function. Note that this is the definition
+   which starts at [n = 0]. *)
+op g_mass i = if 0 <= i then p * (1%r - p) ^ i else 0%r.
+op geometric = mk g_mass.
+
+local lemma summable_mass : summable g_mass.
+proof.
+apply (summable_from_bounded _ (fun i => Some i)); 1: by split => /#. 
+exists 1%r => n.
+rewrite pmap_some map_id.
+rewrite /partial (@eq_big_seq _ (fun i => p * (1%r - p) ^ i)).
+- smt(mem_range in01_p expr_ge0).
+case (0 <= n) => nat_n; 2: smt(big_geq).
+rewrite -mulr_sumr.
+apply ler_pdivl_mull; 1: smt(in01_p).
+have -> : inv p * 1%r = (1%r / (1%r - (1%r - p))) by smt().
+by apply Bigreal.sum_pow_le; smt(in01_p).
+qed.
+
+local lemma sum_mass : sum g_mass = 1%r. 
+proof.
+rewrite (sumEw _ (fun i => Some i) (fun i => 0 <= i) _ _ summable_mass); 2: smt(). 
+- by split => /#.
+have -> : (fun n => big predT g_mass (pmap (fun i => Some i) (range 0 n))) =
+  (fun n => if 0 <= n then (1%r - (1%r - p) ^ n) else 0%r).
+- apply/fun_ext => n. rewrite pmap_some map_id.
+  rewrite (@eq_big_seq _ (fun i => p * (1%r - p) ^ i)).
+  + smt(mem_range in01_p expr_ge0).
+  case (0 <= n) => nat_n; 2: smt(big_geq).
+  rewrite -mulr_sumr Bigreal.sum_pow; 1,2: smt(in01_p). 
+  by field; smt(in01_p).
+rewrite (lim_eq 0 _ (fun n => (1%r - (1%r - p) ^ n))); 1: by auto => />.
+rewrite limD; 1: exact cnvC.
+- by apply/cnvN/cnv_pow; smt(in01_p).
+by rewrite limC limN lim_pow; 1: smt(in01_p).
+qed.
+
+lemma isdistr_geometric : isdistr g_mass.
+proof.
+apply isdistr_summable_equiv.
+split; 1: smt(in01_p expr_ge0).
+split; 1: exact summable_mass.
+by apply/lerr_eq/sum_mass.
+qed.
+
+lemma geometric_ll : is_lossless geometric.
+proof.
+rewrite /is_lossless mu_mass (massK _ isdistr_geometric).
+exact sum_mass.
+qed.
+end section Parametrized.
+
+abstract theory ModuleBased.
+op p : { real | 0%r < p < 1%r } as in01_p.
+
+clone FixedBiased as Bernoulli
+  with op p <- p
+proof*.
+realize in01_p by exact in01_p. 
+
+op bernoulli = Bernoulli.dbiased.
+
+module M = {
+  (* Direct sampling *)
+  proc sample() = {
+      var i;
+      i <$ geometric p;
+      return i;
+  }
+  
+  (* Rejection sampling from a bernoulli *)
+  proc rej() = {
+    var i <- -1;
+    var b <- false;
+    while (!b) {
+      i <- i + 1;
+      b <$ bernoulli;
+    }
+    return i;
+  }
+}.
+
+section Proofs.
+lemma rej_ll : islossless M.rej.
+proof.
+proc.
+sp. conseq (: true ==> _) => //.
+while (true) (b2i (!b)) 1 p; 1..3: smt().
+- move=> H. seq 2 : true 1%r 1%r 0%r 0%r (true) => //.
+  + auto; smt(Bernoulli.dbiased_ll).
+- auto; smt(Bernoulli.dbiased_ll).
+split; 1: smt(in01_p).
+move => z. 
+rnd. wp. skip.
+rewrite Bernoulli.dbiasedE.
+by auto => />. 
+qed.
+
+
+local lemma rej_bound &m k :
+  Pr[ M.rej() @ &m : res = k ] <= g_mass p k.
+proof.
+byehoare => //.
+proc.
+while ((-1 <= i /\ (b => 0 <= i)) `|` (if b then b2r (i = k) else
+  if i < k then p * (1%r - p) ^ (k - i - 1) else 0%r)%xr); 1,3: by auto => /#. 
+wp. skip => /> &h. rewrite Ep_bool_gen Bernoulli.dbiased1E Bernoulli.dbiased1E.
+apply xle_cxr_r => /> *.
+have -> /= : -1 <= i{h} + 1 by smt().
+have -> /= : 0 <= i{h} + 1 by smt().
+case (i{h} + 1 < k) => * />; first by smt(Domain.exprS).
+smt(Domain.expr0).
+qed.
+
+lemma rej_distr &m k :
+  Pr[ M.rej() @ &m : res = k ] = mu1 (geometric p) k.
+proof.
+rewrite Pr[mu1_le_eq_mu1 (geometric p)] => //.
+- exact rej_ll.
+move => ?. rewrite (muK _ _ (isdistr_geometric _ _)). 
+- smt(in01_p).
+exact rej_bound. 
+qed.
+
+equiv rej_equiv : M.rej ~ M.sample : true ==> ={res}.
+bypr (res{1}) (res{2}) => // *.
+rewrite rej_distr.
+byphoare; 2,3: by auto.
+proc. rnd. skip. by auto => />.
+qed.
+end section Proofs.
+end ModuleBased.