diff --git a/src/tea3/cli.py b/src/tea3/cli.py
index 5ac7749..31876bd 100644
--- a/src/tea3/cli.py
+++ b/src/tea3/cli.py
@@ -3,6 +3,7 @@ from tea3.tea3model import Tea3Model
 from tea3.variable_search import run_exhaustive
 from tea3.sbox import run_sbox
 from tea3.variable_xor import run_variable_xor
+from tea3.f31f32 import run_f31f32
 
 def prompt_int(message: str, lo: int, hi: int) -> int:
     while True:
@@ -155,8 +156,9 @@ def main():
     print("  2) Exhaustive variable-change search")
     print("  3) S box analysis")
     print("  4) variable XOR")
+    print("  5) F31, F32 analysis")
 
-    mode = prompt_choice("Your choice (1, 2, 3 or 4): ", {1, 2, 3, 4})
+    mode = prompt_choice("Your choice (1, 2, 3, 4 or 5): ", {1, 2, 3, 4, 5})
 
     if mode == 1:
         run_classic_cli()
@@ -164,8 +166,10 @@ def main():
         run_exhaustive_cli()
     elif mode == 3:
         run_sbox()
-    else:
+    elif mode == 4:
         run_variable_xor_cli()
+    elif mode == 5:
+        run_f31f32()
 
 
 main()
diff --git a/src/tea3/f31f32.py b/src/tea3/f31f32.py
new file mode 100644
index 0000000..a0dcbca
--- /dev/null
+++ b/src/tea3/f31f32.py
@@ -0,0 +1,63 @@
+from sage.all import BooleanPolynomialRing
+from sage.crypto.sbox import SBox
+
+from tea3.tea3model import F31, F32
+
+
+def _vars():
+    R = BooleanPolynomialRing(16, [f"u{i}" for i in range(16)])
+    u = R.gens()
+    return list(u[:8]), list(u[8:])
+
+
+def _poly_masks(p, vars_):
+    idx = {v: i for i, v in enumerate(vars_)}
+    masks = []
+    for monom in p:
+        if monom:
+            mask = 0
+            for v in monom.variables():
+                mask |= 1 << idx[v]
+            masks.append(mask)
+    return int(p.constant_coefficient()), masks
+
+
+def _truth_table(vecfun):
+    x, y = _vars()
+    polys = vecfun(x, y)
+    vars_ = x + y
+    coords = [_poly_masks(p, vars_) for p in polys]
+
+    tt = [0] * (1 << 16)
+    for n in range(1 << 16):
+        out = 0
+        for j, (c, masks) in enumerate(coords):
+            bit = c
+            for mask in masks:
+                if (n & mask) == mask:
+                    bit ^= 1
+            out |= bit << j
+        tt[n] = out
+
+    return tt
+
+
+def _sbox(vecfun):
+    return SBox(_truth_table(vecfun), big_endian=False)
+
+def run_f31f32():
+    f31 = _sbox(F31)
+    f32 = _sbox(F32)
+    print("\nF31")
+    print("differential uniformity:", f31.differential_uniformity())
+    print("max DDT coefficient:", f31.maximal_difference_probability_absolute())
+    print("max difference probability:", f31.maximal_difference_probability())
+    print("max LAT coefficient:", f31.maximal_linear_bias_absolute())
+    print("relative bias:", f31.maximal_linear_bias_relative())
+
+    print("\nF32")
+    print("differential uniformity:", f32.differential_uniformity())
+    print("max DDT coefficient:", f32.maximal_difference_probability_absolute())
+    print("max difference probability:", f32.maximal_difference_probability())
+    print("max LAT coefficient:", f32.maximal_linear_bias_absolute())
+    print("relative bias:", f32.maximal_linear_bias_relative())