diff --git a/ecdsa.c b/ecdsa.c
index c757c31ff3..bf44463394 100644
--- a/ecdsa.c
+++ b/ecdsa.c
@@ -290,7 +290,11 @@ void point_jacobian_add(const curve_point *p1, jacobian_curve_point *p2, const e
 	bn_add(&xz, &p2->x);
 	// xz = x1' + x2
 
-	is_doubling = bn_is_zero(&h) | bn_is_equal(&h, prime);
+	// check for h == 0 % prime.  Note that h never normalizes to
+	// zero, since h = x1' + 2*prime - x2 > 0 and a positive
+	// multiple of prime is always normalized to prime by
+	// bn_fast_mod.
+	is_doubling = bn_is_equal(&h, prime);
 
 	bn_multiply(&p1->y, &yz, prime);        // yz = y1' = y1*z2^3;
 	bn_subtractmod(&yz, &p2->y, &r, prime);
diff --git a/test_curves.py b/test_curves.py
index 25e07dfb17..a05e2638fe 100755
--- a/test_curves.py
+++ b/test_curves.py
@@ -129,7 +129,7 @@ def test_mult_half(curve, r):
     y = int2bn(x)
     lib.bn_mult_half(y, int2bn(curve.p))
     y = bn2int(y)
-    if y > curve.p:
+    if y >= curve.p:
         y -= curve.p
     half = ecdsa.numbertheory.inverse_mod(2, curve.p)
     assert y == (x * half) % curve.p
@@ -156,6 +156,17 @@ def test_subtract2(r):
     assert z == z_
 
 
+def test_add(curve, r):
+    x = r.randrange(0, 2 ** 256)
+    y = r.randrange(0, 2 ** 256)
+    z_ = x + y
+    z = int2bn(x)
+    lib.bn_add(z, int2bn(y))
+    z = bn2int(z)
+
+    assert z == z_
+
+
 def test_addmod(curve, r):
     x = r.randrange(0, 2 ** 256)
     y = r.randrange(0, 2 ** 256)
@@ -163,7 +174,8 @@ def test_addmod(curve, r):
     z = int2bn(x)
     lib.bn_addmod(z, int2bn(y), int2bn(curve.p))
     z = bn2int(z)
-
+    if z >= curve.p:
+        z = z - curve.p
     assert z == z_