diff --git a/src/Init/Data/Fin/Basic.lean b/src/Init/Data/Fin/Basic.lean
index 611bb8ed08..7a757b6cc7 100644
--- a/src/Init/Data/Fin/Basic.lean
+++ b/src/Init/Data/Fin/Basic.lean
@@ -21,16 +21,16 @@ def elim0.{u} {α : Sort u} : Fin 0 → α
 variable {n : Nat}
 
 protected def ofNat {n : Nat} (a : Nat) : Fin (succ n) :=
-  ⟨a % succ n, Nat.modLt _ (Nat.zeroLtSucc _)⟩
+  ⟨a % succ n, Nat.mod_lt _ (Nat.zeroLtSucc _)⟩
 
 protected def ofNat' {n : Nat} (a : Nat) (h : n > 0) : Fin n :=
-  ⟨a % n, Nat.modLt _ h⟩
+  ⟨a % n, Nat.mod_lt _ h⟩
 
 private theorem mlt {b : Nat} : {a : Nat} → a < n → b % n < n
-  | 0,   h => Nat.modLt _ h
+  | 0,   h => Nat.mod_lt _ h
   | a+1, h =>
     have n > 0 from Nat.ltTrans (Nat.zeroLtSucc _) h;
-    Nat.modLt _ this
+    Nat.mod_lt _ this
 
 protected def add : Fin n → Fin n → Fin n
   | ⟨a, h⟩, ⟨b, _⟩ => ⟨(a + b) % n, mlt h⟩
@@ -106,8 +106,8 @@ instance : OfNat (Fin (noindex! (n+1))) i where
 theorem vneOfNe {i j : Fin n} (h : i ≠ j) : val i ≠ val j :=
   fun h' => absurd (eqOfVeq h') h
 
-theorem modnLt : ∀ {m : Nat} (i : Fin n), m > 0 → (i % m).val < m
-  | m, ⟨a, h⟩, hp =>  Nat.ltOfLeOfLt (modLe _ _) (modLt _ hp)
+theorem modn_lt : ∀ {m : Nat} (i : Fin n), m > 0 → (i % m).val < m
+  | m, ⟨a, h⟩, hp =>  Nat.ltOfLeOfLt (mod_le _ _) (mod_lt _ hp)
 
 end Fin
 
diff --git a/src/Init/Data/Nat.lean b/src/Init/Data/Nat.lean
index cf7be93c6e..7afb1ce9e5 100644
--- a/src/Init/Data/Nat.lean
+++ b/src/Init/Data/Nat.lean
@@ -6,5 +6,6 @@ Authors: Leonardo de Moura
 prelude
 import Init.Data.Nat.Basic
 import Init.Data.Nat.Div
+import Init.Data.Nat.Gcd
 import Init.Data.Nat.Bitwise
 import Init.Data.Nat.Control
diff --git a/src/Init/Data/Nat/Basic.lean b/src/Init/Data/Nat/Basic.lean
index a95f30ee47..56bb2c9199 100644
--- a/src/Init/Data/Nat/Basic.lean
+++ b/src/Init/Data/Nat/Basic.lean
@@ -185,6 +185,16 @@ theorem subLt : ∀ {n m : Nat}, 0 < n → 0 < m → n - m < n
       show n - m < succ n from
       lt_succ_of_le (subLe n m)
 
+theorem sub_succ (n m : Nat) : n - succ m = pred (n - m) :=
+  rfl
+
+theorem succ_sub_succ (n m : Nat) : succ n - succ m = n - m :=
+  succ_sub_succ_eq_sub n m
+
+protected theorem sub_self : ∀ (n : Nat), n - n = 0
+  | 0        => by rw Nat.sub_zero
+  | (succ n) => by rw [succ_sub_succ, Nat.sub_self n]
+
 protected theorem ltOfLtOfLe {n m k : Nat} : n < m → m ≤ k → n < k :=
   Nat.leTrans
 
@@ -231,7 +241,7 @@ theorem eqZeroOfLeZero {n : Nat} (h : n ≤ 0) : n = 0 :=
 theorem ltOfSuccLt {n m : Nat} : succ n < m → n < m :=
   leOfSuccLe
 
-theorem ltOfSuccLtSucc {n m : Nat} : succ n < succ m → n < m :=
+theorem lt_of_succ_lt_succ {n m : Nat} : succ n < succ m → n < m :=
   leOfSuccLeSucc
 
 theorem ltOfSuccLe {n m : Nat} (h : succ n ≤ m) : n < m :=
diff --git a/src/Init/Data/Nat/Div.lean b/src/Init/Data/Nat/Div.lean
index 9955e8dc96..da983c0054 100644
--- a/src/Init/Data/Nat/Div.lean
+++ b/src/Init/Data/Nat/Div.lean
@@ -8,11 +8,11 @@ import Init.WF
 import Init.Data.Nat.Basic
 namespace Nat
 
-private def divRecLemma {x y : Nat} : 0 < y ∧ y ≤ x → x - y < x :=
+private def div_rec_lemma {x y : Nat} : 0 < y ∧ y ≤ x → x - y < x :=
   fun ⟨ypos, ylex⟩ => subLt (Nat.ltOfLtOfLe ypos ylex) ypos
 
 private def div.F (x : Nat) (f : ∀ x₁, x₁ < x → Nat → Nat) (y : Nat) : Nat :=
-  if h : 0 < y ∧ y ≤ x then f (x - y) (divRecLemma h) y + 1 else zero
+  if h : 0 < y ∧ y ≤ x then f (x - y) (div_rec_lemma h) y + 1 else zero
 
 @[extern "lean_nat_div"]
 protected def div (a b : @& Nat) : Nat :=
@@ -20,18 +20,18 @@ protected def div (a b : @& Nat) : Nat :=
 
 instance : Div Nat := ⟨Nat.div⟩
 
-private theorem divDefAux (x y : Nat) : x / y = if h : 0 < y ∧ y ≤ x then (x - y) / y + 1 else 0 :=
+private theorem div_eq_aux (x y : Nat) : x / y = if h : 0 < y ∧ y ≤ x then (x - y) / y + 1 else 0 :=
   congrFun (WellFounded.fixEq ltWf div.F x) y
 
-theorem divDef (x y : Nat) : x / y = if 0 < y ∧ y ≤ x then (x - y) / y + 1 else 0 :=
-  difEqIf (0 < y ∧ y ≤ x) ((x - y) / y + 1) 0 ▸ divDefAux x y
+theorem div_eq (x y : Nat) : x / y = if 0 < y ∧ y ≤ x then (x - y) / y + 1 else 0 :=
+  difEqIf (0 < y ∧ y ≤ x) ((x - y) / y + 1) 0 ▸ div_eq_aux x y
 
 private theorem div.induction.F.{u}
         (C : Nat → Nat → Sort u)
         (ind  : ∀ x y, 0 < y ∧ y ≤ x → C (x - y) y → C x y)
         (base : ∀ x y, ¬(0 < y ∧ y ≤ x) → C x y)
         (x : Nat) (f : ∀ (x₁ : Nat), x₁ < x → ∀ (y : Nat), C x₁ y) (y : Nat) : C x y :=
-  if h : 0 < y ∧ y ≤ x then ind x y h (f (x - y) (divRecLemma h) y) else base x y h
+  if h : 0 < y ∧ y ≤ x then ind x y h (f (x - y) (div_rec_lemma h) y) else base x y h
 
 theorem div.inductionOn.{u}
       {motive : Nat → Nat → Sort u}
@@ -42,7 +42,7 @@ theorem div.inductionOn.{u}
   WellFounded.fix Nat.ltWf (div.induction.F motive ind base) x y
 
 private def mod.F (x : Nat) (f : ∀ x₁, x₁ < x → Nat → Nat) (y : Nat) : Nat :=
-  if h : 0 < y ∧ y ≤ x then f (x - y) (divRecLemma h) y else x
+  if h : 0 < y ∧ y ≤ x then f (x - y) (div_rec_lemma h) y else x
 
 @[extern "lean_nat_mod"]
 protected def mod (a b : @& Nat) : Nat :=
@@ -50,11 +50,11 @@ protected def mod (a b : @& Nat) : Nat :=
 
 instance : Mod Nat := ⟨Nat.mod⟩
 
-private theorem modDefAux (x y : Nat) : x % y = if h : 0 < y ∧ y ≤ x then (x - y) % y else x :=
+private theorem mod_eq_aux (x y : Nat) : x % y = if h : 0 < y ∧ y ≤ x then (x - y) % y else x :=
   congrFun (WellFounded.fixEq ltWf mod.F x) y
 
-theorem modDef (x y : Nat) : x % y = if 0 < y ∧ y ≤ x then (x - y) % y else x :=
-  difEqIf (0 < y ∧ y ≤ x) ((x - y) % y) x ▸ modDefAux x y
+theorem mod_eq (x y : Nat) : x % y = if 0 < y ∧ y ≤ x then (x - y) % y else x :=
+  difEqIf (0 < y ∧ y ≤ x) ((x - y) % y) x ▸ mod_eq_aux x y
 
 theorem mod.inductionOn.{u}
       {motive : Nat → Nat → Sort u}
@@ -64,28 +64,28 @@ theorem mod.inductionOn.{u}
       : motive x y :=
   div.inductionOn x y ind base
 
-theorem modZero (a : Nat) : a % 0 = a :=
+theorem mod_zero (a : Nat) : a % 0 = a :=
   have (if 0 < 0 ∧ 0 ≤ a then (a - 0) % 0 else a) = a from
     have h : ¬ (0 < 0 ∧ 0 ≤ a) from fun ⟨h₁, _⟩ => absurd h₁ (Nat.ltIrrefl _)
     ifNeg h
-  (modDef a 0).symm ▸ this
+  (mod_eq a 0).symm ▸ this
 
-theorem modEqOfLt {a b : Nat} (h : a < b) : a % b = a :=
+theorem mod_eq_of_lt {a b : Nat} (h : a < b) : a % b = a :=
   have (if 0 < b ∧ b ≤ a then (a - b) % b else a) = a from
     have h' : ¬(0 < b ∧ b ≤ a) from fun ⟨_, h₁⟩ => absurd h₁ (Nat.notLeOfGt h)
     ifNeg h'
-  (modDef a b).symm ▸ this
+  (mod_eq a b).symm ▸ this
 
-theorem modEqSubMod {a b : Nat} (h : a ≥ b) : a % b = (a - b) % b :=
+theorem mod_eq_sub_mod {a b : Nat} (h : a ≥ b) : a % b = (a - b) % b :=
   match eqZeroOrPos b with
   | Or.inl h₁ => h₁.symm ▸ (Nat.sub_zero a).symm ▸ rfl
-  | Or.inr h₁ => (modDef a b).symm ▸ ifPos ⟨h₁, h⟩
+  | Or.inr h₁ => (mod_eq a b).symm ▸ ifPos ⟨h₁, h⟩
 
-theorem modLt (x : Nat) {y : Nat} : y > 0 → x % y < y := by
+theorem mod_lt (x : Nat) {y : Nat} : y > 0 → x % y < y := by
    refine mod.inductionOn (motive := fun x y => y > 0 → x % y < y) x y ?k1 ?k2
    case k1 =>
      intro x y ⟨_, h₁⟩ h₂ h₃
-     rw [modEqSubMod h₁]
+     rw [mod_eq_sub_mod h₁]
      exact h₂ h₃
    case k2 =>
      intro x y h₁ h₂
@@ -94,15 +94,34 @@ theorem modLt (x : Nat) {y : Nat} : y > 0 → x % y < y := by
      | Or.inl h₁ => exact absurd h₂ h₁
      | Or.inr h₁ =>
        have hgt : y > x from gtOfNotLe h₁
-       have heq : x % y = x from modEqOfLt hgt
+       have heq : x % y = x from mod_eq_of_lt hgt
        rw [← heq] at hgt
        exact hgt
 
-theorem modLe (x y : Nat) : x % y ≤ x := by
+theorem mod_le (x y : Nat) : x % y ≤ x := by
   match Nat.ltOrGe x y with
-  | Or.inl h₁ => rw [modEqOfLt h₁]; apply Nat.leRefl
+  | Or.inl h₁ => rw [mod_eq_of_lt h₁]; apply Nat.leRefl
   | Or.inr h₁ => match eqZeroOrPos y with
-    | Or.inl h₂ => rw [h₂, Nat.modZero x]; apply Nat.leRefl
-    | Or.inr h₂ => exact Nat.leTrans (Nat.leOfLt (Nat.modLt _ h₂)) h₁
+    | Or.inl h₂ => rw [h₂, Nat.mod_zero x]; apply Nat.leRefl
+    | Or.inr h₂ => exact Nat.leTrans (Nat.leOfLt (mod_lt _ h₂)) h₁
+
+@[simp] theorem zero_mod (b : Nat) : 0 % b = 0 := by
+  rw [mod_eq]
+  have ¬ (0 < b ∧ b ≤ 0) by
+    intro ⟨h₁, h₂⟩
+    exact absurd (Nat.ltOfLtOfLe h₁ h₂) (Nat.ltIrrefl 0)
+  simp [this]
+
+@[simp] theorem mod_self (n : Nat) : n % n = 0 := by
+  rw [mod_eq_sub_mod (Nat.leRefl _), Nat.sub_self, zero_mod]
+
+theorem mod_one (x : Nat) : x % 1 = 0 := by
+  have h : x % 1 < 1 from mod_lt x decide!
+  have (y : Nat) → y < 1 → y = 0 by
+    intro y
+    cases y with
+    | zero   => intro h; rfl
+    | succ y => intro h; apply absurd (Nat.lt_of_succ_lt_succ h) (Nat.notLtZero y)
+  exact this _ h
 
 end Nat
diff --git a/src/Init/Data/Nat/Gcd.lean b/src/Init/Data/Nat/Gcd.lean
new file mode 100644
index 0000000000..d9f3386125
--- /dev/null
+++ b/src/Init/Data/Nat/Gcd.lean
@@ -0,0 +1,36 @@
+/-
+Copyright (c) 2021 Microsoft Corporation. All rights reserved.
+Released under Apache 2.0 license as described in the file LICENSE.
+Authors: Leonardo de Moura
+-/
+prelude
+import Init.Data.Nat.Div
+
+namespace Nat
+
+private def gcdF (x : Nat) : (∀ x₁, x₁ < x → Nat → Nat) → Nat → Nat :=
+  match x with
+  | 0      => fun _ y => y
+  | succ x => fun f y => f (y % succ x) (mod_lt _ (zeroLtSucc  _)) (succ x)
+
+@[extern "lean_nat_gcd"]
+def gcd (a b : @& Nat) : Nat :=
+  WellFounded.fix ltWf gcdF a b
+
+@[simp] theorem gcd_zero_left (y : Nat) : gcd 0 y = y :=
+  rfl
+
+theorem gcd_succ (x y : Nat) : gcd (succ x) y = gcd (y % succ x) (succ x) :=
+  rfl
+
+@[simp] theorem gcd_one_left (n : Nat) : gcd 1 n = 1 := by
+  rw [gcd_succ, mod_one]
+  rfl
+
+@[simp] theorem gcd_zero_right (n : Nat) : gcd n 0 = n := by
+  cases n <;> simp [gcd_succ]
+
+@[simp] theorem gcd_self (n : Nat) : gcd n n = n := by
+  cases n <;> simp [gcd_succ]
+
+end Nat
\ No newline at end of file
diff --git a/src/Init/Data/UInt.lean b/src/Init/Data/UInt.lean
index ccc9e5056d..9a6273b2cb 100644
--- a/src/Init/Data/UInt.lean
+++ b/src/Init/Data/UInt.lean
@@ -343,5 +343,5 @@ def USize.decLe (a b : USize) : Decidable (a ≤ b) :=
 instance (a b : USize) : Decidable (a < b) := USize.decLt a b
 instance (a b : USize) : Decidable (a ≤ b) := USize.decLe a b
 
-theorem USize.modnLt {m : Nat} : ∀ (u : USize), m > 0 → USize.toNat (u % m) < m
-  | ⟨u⟩, h => Fin.modnLt u h
+theorem USize.modn_lt {m : Nat} : ∀ (u : USize), m > 0 → USize.toNat (u % m) < m
+  | ⟨u⟩, h => Fin.modn_lt u h
diff --git a/src/Std/Data/HashMap.lean b/src/Std/Data/HashMap.lean
index 98a5147365..ef0c50eb48 100644
--- a/src/Std/Data/HashMap.lean
+++ b/src/Std/Data/HashMap.lean
@@ -29,7 +29,7 @@ namespace HashMapImp
 variable {α : Type u} {β : Type v}
 
 def mkIdx {n : Nat} (h : n > 0) (u : USize) : { u : USize // u.toNat < n } :=
-  ⟨u % n, USize.modnLt _ h⟩
+  ⟨u % n, USize.modn_lt _ h⟩
 
 @[inline] def reinsertAux (hashFn : α → USize) (data : HashMapBucket α β) (a : α) (b : β) : HashMapBucket α β :=
   let ⟨i, h⟩ := mkIdx data.property (hashFn a)
diff --git a/src/Std/Data/HashSet.lean b/src/Std/Data/HashSet.lean
index a57adbea91..e15f92db37 100644
--- a/src/Std/Data/HashSet.lean
+++ b/src/Std/Data/HashSet.lean
@@ -28,7 +28,7 @@ namespace HashSetImp
 variable {α : Type u}
 
 def mkIdx {n : Nat} (h : n > 0) (u : USize) : { u : USize // u.toNat < n } :=
-  ⟨u % n, USize.modnLt _ h⟩
+  ⟨u % n, USize.modn_lt _ h⟩
 
 @[inline] def reinsertAux (hashFn : α → USize) (data : HashSetBucket α) (a : α) : HashSetBucket α :=
   let ⟨i, h⟩ := mkIdx data.property (hashFn a)
diff --git a/src/include/lean/lean.h b/src/include/lean/lean.h
index 9f013f9c3e..d715e0d2da 100644
--- a/src/include/lean/lean.h
+++ b/src/include/lean/lean.h
@@ -1433,6 +1433,7 @@ static inline lean_obj_res lean_nat_lxor(b_lean_obj_arg a1, b_lean_obj_arg a2) {
 lean_obj_res lean_nat_shiftl(b_lean_obj_arg a1, b_lean_obj_arg a2);
 lean_obj_res lean_nat_shiftr(b_lean_obj_arg a1, b_lean_obj_arg a2);
 lean_obj_res lean_nat_pow(b_lean_obj_arg a1, b_lean_obj_arg a2);
+lean_obj_res lean_nat_gcd(b_lean_obj_arg a1, b_lean_obj_arg a2);
 
 /* Integers */
 
diff --git a/src/runtime/object.cpp b/src/runtime/object.cpp
index 277d5e5c9c..f94fb96897 100644
--- a/src/runtime/object.cpp
+++ b/src/runtime/object.cpp
@@ -1315,6 +1315,20 @@ extern "C" lean_obj_res lean_nat_pow(b_lean_obj_arg a1, b_lean_obj_arg a2) {
         return mpz_to_nat(mpz_value(a1).pow(lean_unbox(a2)));
 }
 
+extern "C" lean_obj_res lean_nat_gcd(b_lean_obj_arg a1, b_lean_obj_arg a2) {
+    if (lean_is_scalar(a1)) {
+      if (lean_is_scalar(a2))
+        return mpz_to_nat(gcd(mpz::of_size_t(lean_unbox(a1)), mpz::of_size_t(lean_unbox(a2))));
+      else
+        return mpz_to_nat(gcd(mpz::of_size_t(lean_unbox(a1)), mpz_value(a2)));
+    } else {
+      if (lean_is_scalar(a2))
+        return mpz_to_nat(gcd(mpz_value(a1), mpz::of_size_t(lean_unbox(a2))));
+      else
+        return mpz_to_nat(gcd(mpz_value(a1), mpz_value(a2)));
+    }
+}
+
 // =======================================
 // Integers
 
diff --git a/tests/lean/gcd.lean b/tests/lean/gcd.lean
new file mode 100644
index 0000000000..365b97b8b6
--- /dev/null
+++ b/tests/lean/gcd.lean
@@ -0,0 +1,16 @@
+open Nat
+
+#eval   gcd 8 12
+#reduce gcd 8 12
+
+#eval   gcd 0 9
+#reduce gcd 0 9
+
+#eval   gcd 4 0
+#reduce gcd 4 0
+
+#eval   gcd 1 10
+#reduce gcd 1 10
+
+#eval gcd (2*3*7*11*17) (2*7*23)
+#eval gcd (2*7*23) (2*3*7*11*17)
diff --git a/tests/lean/gcd.lean.expected.out b/tests/lean/gcd.lean.expected.out
new file mode 100644
index 0000000000..3098b7af35
--- /dev/null
+++ b/tests/lean/gcd.lean.expected.out
@@ -0,0 +1,10 @@
+4
+4
+9
+9
+4
+4
+1
+1
+14
+14