Implement Curve448 primitives

This patch adds the necessary primitives for "curve448", defined in
RFC 7748.  Those primitives are namely: addition, doubling, scalar
multiplication of the generator or an arbitrary point, inversion, and
square root.

30 files changed, 1212 insertions, 39 deletions

diff --git a/.gitignore b/.gitignore
index b79c53f5..0afe61de 100644
--- a/.gitignore
+++ b/.gitignore
@@ -49,6 +49,7 @@ core
 /ecc-384.h
 /ecc-521.h
 /ecc-25519.h
+/ecc-448.h
 /version.h
 /nettle.aux
 /nettle.cp
diff --git a/Makefile.in b/Makefile.in
index 9f5b065a..036a3a1d 100644
--- a/Makefile.in
+++ b/Makefile.in
@@ -175,7 +175,7 @@ hogweed_SOURCES = sexp.c sexp-format.c \
 		  ecc-mod.c ecc-mod-inv.c \
 		  ecc-mod-arith.c ecc-pp1-redc.c ecc-pm1-redc.c \
 		  ecc-192.c ecc-224.c ecc-256.c ecc-384.c ecc-521.c \
-		  ecc-25519.c \
+		  ecc-25519.c ecc-448.c \
 		  ecc-size.c ecc-j-to-a.c ecc-a-to-j.c \
 		  ecc-dup-jj.c ecc-add-jja.c ecc-add-jjj.c \
 		  ecc-eh-to-a.c \
@@ -186,6 +186,7 @@ hogweed_SOURCES = sexp.c sexp-format.c \
 		  ecc-ecdsa-sign.c ecdsa-sign.c \
 		  ecc-ecdsa-verify.c ecdsa-verify.c ecdsa-keygen.c \
 		  curve25519-mul-g.c curve25519-mul.c curve25519-eh-to-x.c \
+		  curve448-mul-g.c curve448-mul.c curve448-eh-to-x.c \
 		  eddsa-compress.c eddsa-decompress.c eddsa-expand.c \
 		  eddsa-hash.c eddsa-pubkey.c eddsa-sign.c eddsa-verify.c \
 		  ed25519-sha512-pubkey.c \
@@ -196,7 +197,7 @@ OPT_SOURCES = fat-x86_64.c fat-arm.c mini-gmp.c
 HEADERS = aes.h arcfour.h arctwo.h asn1.h blowfish.h \
 	  base16.h base64.h bignum.h buffer.h camellia.h cast128.h \
 	  cbc.h ccm.h cfb.h chacha.h chacha-poly1305.h ctr.h \
-	  curve25519.h des.h dsa.h dsa-compat.h eax.h \
+	  curve25519.h curve448.h des.h dsa.h dsa-compat.h eax.h \
 	  ecc-curve.h ecc.h ecdsa.h eddsa.h \
 	  gcm.h gost28147.h gosthash94.h hmac.h \
 	  knuth-lfib.h hkdf.h \
@@ -387,6 +388,9 @@ ecc-521.h: eccdata.stamp
 ecc-25519.h: eccdata.stamp
 	./eccdata$(EXEEXT_FOR_BUILD) 255 11 6 $(NUMB_BITS) > $@T && mv $@T $@
 
+ecc-448.h: eccdata.stamp
+	./eccdata$(EXEEXT_FOR_BUILD) 448 38 6 $(NUMB_BITS) > $@T && mv $@T $@
+
 eccdata.stamp: eccdata.c
 	$(MAKE) eccdata$(EXEEXT_FOR_BUILD)
 	echo stamp > eccdata.stamp
@@ -397,6 +401,7 @@ ecc-256.$(OBJEXT): ecc-256.h
 ecc-384.$(OBJEXT): ecc-384.h
 ecc-521.$(OBJEXT): ecc-521.h
 ecc-25519.$(OBJEXT): ecc-25519.h
+ecc-448.$(OBJEXT): ecc-448.h
 
 .asm.$(OBJEXT): $(srcdir)/asm.m4 machine.m4 config.m4
 	$(M4) $(srcdir)/asm.m4 machine.m4 config.m4 $< >$*.s
@@ -650,6 +655,7 @@ distcheck: dist
 clean-here:
 	-rm -f $(TARGETS) *.$(OBJEXT) *.s *.so *.dll *.a \
 		ecc-192.h ecc-224.h ecc-256.h ecc-384.h ecc-521.h ecc-25519.h \
+		ecc-448.h \
 		aesdata$(EXEEXT_FOR_BUILD) \
 		desdata$(EXEEXT_FOR_BUILD) \
 		twofishdata$(EXEEXT_FOR_BUILD) \
diff --git a/curve448-eh-to-x.c b/curve448-eh-to-x.c
new file mode 100644
index 00000000..4bc78303
--- /dev/null
+++ b/curve448-eh-to-x.c
@@ -0,0 +1,73 @@
+/* curve448-eh-to-x.c
+
+   Copyright (C) 2017 Daiki Ueno
+   Copyright (C) 2017 Red Hat, Inc.
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <string.h>
+
+#include "curve448.h"
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+/* Transform a point on the edwards448 Edwards curve to the curve448
+   Montgomery curve, and return the x coordinate. */
+void
+curve448_eh_to_x (mp_limb_t *xp, const mp_limb_t *p, mp_limb_t *scratch)
+{
+#define vp (p + ecc->p.size)
+#define t0 scratch
+#define t1 (scratch + ecc->p.size)
+#define t2 (scratch + 2*ecc->p.size)
+
+  const struct ecc_curve *ecc = &_nettle_curve448;
+  mp_limb_t cy;
+
+  /* If u = U/W and v = V/W are the coordinates of the point on
+     edwards448 we get the curve448 x coordinate as
+
+     x = v^2 / u^2 = (V/W)^2 / (U/W)^2 = (V/U)^2
+  */
+  /* Needs a total of 9*size storage. */
+  ecc->p.invert (&ecc->p, t0, p, t1 + ecc->p.size);
+  ecc_modp_mul (ecc, t1, t0, vp);
+  ecc_modp_mul (ecc, t2, t1, t1);
+
+  cy = mpn_sub_n (xp, t2, ecc->p.m, ecc->p.size);
+  cnd_copy (cy, xp, t2, ecc->p.size);
+#undef vp
+#undef t0
+#undef t1
+#undef t2
+}
diff --git a/curve448-mul-g.c b/curve448-mul-g.c
new file mode 100644
index 00000000..a396595a
--- /dev/null
+++ b/curve448-mul-g.c
@@ -0,0 +1,74 @@
+/* curve448-mul-g.c
+
+   Copyright (C) 2017 Daiki Ueno
+   Copyright (C) 2017 Red Hat, Inc.
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <string.h>
+
+#include "curve448.h"
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+/* Intended to be compatible with NaCl's crypto_scalarmult_base. */
+void
+curve448_mul_g (uint8_t *r, const uint8_t *n)
+{
+  const struct ecc_curve *ecc = &_nettle_curve448;
+  uint8_t t[CURVE448_SIZE];
+  mp_limb_t *scratch;
+  mp_size_t itch;
+
+#define ng scratch
+#define x (scratch + 3*ecc->p.size)
+#define scratch_out (scratch + 4*ecc->p.size)
+
+  memcpy (t, n, sizeof(t));
+  t[0] &= ~3;
+  t[CURVE448_SIZE-1] = (t[CURVE448_SIZE-1] & 0x7f) | 0x80;
+
+  itch = 5*ecc->p.size + ecc->mul_g_itch;
+  scratch = gmp_alloc_limbs (itch);
+
+  mpn_set_base256_le (x, ecc->p.size, t, CURVE448_SIZE);
+
+  ecc_mul_g_eh (ecc, ng, x, scratch_out);
+  curve448_eh_to_x (x, ng, scratch_out);
+
+  mpn_get_base256_le (r, CURVE448_SIZE, x, ecc->p.size);
+  gmp_free_limbs (scratch, itch);
+#undef ng
+#undef x
+#undef scratch_out
+}
diff --git a/curve448-mul.c b/curve448-mul.c
new file mode 100644
index 00000000..afa814a4
--- /dev/null
+++ b/curve448-mul.c
@@ -0,0 +1,148 @@
+/* curve448-mul.c
+
+   Copyright (C) 2017 Daiki Ueno
+   Copyright (C) 2017 Red Hat, Inc.
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include <string.h>
+
+#include "curve448.h"
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+/* Intended to be compatible with NaCl's crypto_scalarmult. */
+void
+curve448_mul (uint8_t *q, const uint8_t *n, const uint8_t *p)
+{
+  const struct ecc_curve *ecc = &_nettle_curve448;
+  mp_size_t itch;
+  mp_limb_t *scratch;
+  int i;
+  mp_limb_t cy;
+
+  /* FIXME: Could save some more scratch space, e.g., by letting BB
+     overlap C, D, and CB overlap A, D. And possibly reusing some of
+     x2, z2, x3, z3. */
+#define x1 scratch
+#define x2 (scratch + ecc->p.size)
+#define z2 (scratch + 2*ecc->p.size)
+#define x3 (scratch + 3*ecc->p.size)
+#define z3 (scratch + 4*ecc->p.size)
+
+#define A  (scratch + 5*ecc->p.size)
+#define B  (scratch + 6*ecc->p.size)
+#define C  (scratch + 7*ecc->p.size)
+#define D  (scratch + 8*ecc->p.size)
+#define AA  (scratch + 9*ecc->p.size)
+#define BB  (scratch + 10*ecc->p.size)
+#define E  (scratch + 10*ecc->p.size) /* Overlap BB */
+#define DA  (scratch + 9*ecc->p.size) /* Overlap AA */
+#define CB  (scratch + 10*ecc->p.size) /* Overlap BB */
+
+#define a24 39081
+
+  itch = ecc->p.size * 14;
+  scratch = gmp_alloc_limbs (itch);
+
+  /* Note that 255 % GMP_NUMB_BITS == 0 isn't supported, so x1 always
+     holds at least 256 bits. */
+  mpn_set_base256_le (x1, ecc->p.size, p, CURVE448_SIZE);
+
+  /* Initialize, x2 = x1, z2 = 1 */
+  mpn_copyi (x2, x1, ecc->p.size);
+  z2[0] = 1;
+  mpn_zero (z2+1, ecc->p.size - 1);
+
+  /* Get x3, z3 from doubling. Since bit 447 is forced to 1. */
+  ecc_modp_add (ecc, A, x2, z2);
+  ecc_modp_sub (ecc, B, x2, z2);
+  ecc_modp_sqr (ecc, AA, A);
+  ecc_modp_sqr (ecc, BB, B);
+  ecc_modp_mul (ecc, x3, AA, BB);
+  ecc_modp_sub (ecc, E, AA, BB);
+  ecc_modp_addmul_1 (ecc, AA, E, a24);
+  ecc_modp_mul (ecc, z3, E, AA);
+
+  for (i = 446; i >= 2; i--)
+    {
+      int bit = (n[i/8] >> (i & 7)) & 1;
+
+      cnd_swap (bit, x2, x3, 2*ecc->p.size);
+
+      /* Formulas from RFC 7748. We compute new coordinates in
+	 memory-address order, since mul and sqr clobbers higher
+	 limbs. */
+      ecc_modp_add (ecc, A, x2, z2);
+      ecc_modp_sub (ecc, B, x2, z2);
+      ecc_modp_sqr (ecc, AA, A);
+      ecc_modp_sqr (ecc, BB, B);
+      ecc_modp_mul (ecc, x2, AA, BB);
+      ecc_modp_sub (ecc, E, AA, BB); /* Last use of BB */
+      ecc_modp_addmul_1 (ecc, AA, E, a24);
+      ecc_modp_add (ecc, C, x3, z3);
+      ecc_modp_sub (ecc, D, x3, z3);
+      ecc_modp_mul (ecc, z2, E, AA); /* Last use of E and AA */
+      ecc_modp_mul (ecc, DA, D, A);  /* Last use of D, A. FIXME: could
+					let CB overlap. */
+      ecc_modp_mul (ecc, CB, C, B);
+
+      ecc_modp_add (ecc, C, DA, CB);
+      ecc_modp_sqr (ecc, x3, C);
+      ecc_modp_sub (ecc, C, DA, CB);
+      ecc_modp_sqr (ecc, DA, C);
+      ecc_modp_mul (ecc, z3, DA, x1);
+
+      /* FIXME: Could be combined with the loop's initial cnd_swap. */
+      cnd_swap (bit, x2, x3, 2*ecc->p.size);
+    }
+  /* Do the 2 low zero bits, just duplicating x2 */
+  for ( ; i >= 0; i--)
+    {
+      ecc_modp_add (ecc, A, x2, z2);
+      ecc_modp_sub (ecc, B, x2, z2);
+      ecc_modp_sqr (ecc, AA, A);
+      ecc_modp_sqr (ecc, BB, B);
+      ecc_modp_mul (ecc, x2, AA, BB);
+      ecc_modp_sub (ecc, E, AA, BB);
+      ecc_modp_addmul_1 (ecc, AA, E, a24);
+      ecc_modp_mul (ecc, z2, E, AA);
+    }
+  ecc->p.invert (&ecc->p, x3, z2, z3 + ecc->p.size);
+  ecc_modp_mul (ecc, z3, x2, x3);
+  cy = mpn_sub_n (x2, z3, ecc->p.m, ecc->p.size);
+  cnd_copy (cy, x2, z3, ecc->p.size);
+  mpn_get_base256_le (q, CURVE448_SIZE, x2, ecc->p.size);
+
+  gmp_free_limbs (scratch, itch);
+}
diff --git a/curve448.h b/curve448.h
new file mode 100644
index 00000000..a27831e5
--- /dev/null
+++ b/curve448.h
@@ -0,0 +1,58 @@
+/* curve448.h
+
+   Copyright (C) 2017 Daiki Ueno
+   Copyright (C) 2017 Red Hat, Inc.
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#ifndef NETTLE_CURVE448_H
+#define NETTLE_CURVE448_H
+
+#include "nettle-types.h"
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+/* Name mangling */
+#define curve448_mul_g nettle_curve448_mul_g
+#define curve448_mul nettle_curve448_mul
+
+#define CURVE448_SIZE 56
+
+void
+curve448_mul_g (uint8_t *q, const uint8_t *n);
+
+void
+curve448_mul (uint8_t *q, const uint8_t *n, const uint8_t *p);
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif /* NETTLE_CURVE448_H */
diff --git a/ecc-448.c b/ecc-448.c
new file mode 100644
index 00000000..24d970e4
--- /dev/null
+++ b/ecc-448.c
@@ -0,0 +1,269 @@
+/* ecc-448.c
+
+   Arithmetic and tables for curve448,
+
+   Copyright (C) 2017 Daiki Ueno
+   Copyright (C) 2017 Red Hat, Inc.
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#if HAVE_CONFIG_H
+# include "config.h"
+#endif
+
+#include "ecc.h"
+#include "ecc-internal.h"
+
+#define USE_REDC 0
+
+#include "ecc-448.h"
+
+/* Needs 2*ecc->size limbs at rp, and 2*ecc->size additional limbs of
+   scratch space. No overlap allowed. */
+static void
+ecc_mod_pow_2k (const struct ecc_modulo *m,
+		mp_limb_t *rp, const mp_limb_t *xp,
+		unsigned k, mp_limb_t *tp)
+{
+  if (k & 1)
+    {
+      ecc_mod_sqr (m, rp, xp);
+      k--;
+    }
+  else
+    {
+      ecc_mod_sqr (m, tp, xp);
+      ecc_mod_sqr (m, rp, tp);
+      k -= 2;
+    }
+  while (k > 0)
+    {
+      ecc_mod_sqr (m, tp, rp);
+      ecc_mod_sqr (m, rp, tp);
+      k -= 2;
+    }
+}
+
+/* Computes a^{(p-3)/4} = a^{2^446-2^222-1} mod m. Needs 9 * n scratch
+   space. */
+static void
+ecc_mod_pow_446m224m1 (const struct ecc_modulo *p,
+		       mp_limb_t *rp, const mp_limb_t *ap,
+		       mp_limb_t *scratch)
+{
+#define t0 scratch
+#define t1 (scratch + 3*ECC_LIMB_SIZE)
+#define t2 (scratch + 6*ECC_LIMB_SIZE)
+
+  ecc_mod_sqr (p, rp, ap);	        /* a^2 */
+  ecc_mod_mul (p, t0, ap, rp);		/* a^3 */
+  ecc_mod_sqr (p, rp, t0);		/* a^6 */
+  ecc_mod_mul (p, t0, ap, rp);		/* a^{2^3-1} */
+  ecc_mod_pow_2k (p, rp, t0, 3, t2);	/* a^{2^6-2^3} */
+  ecc_mod_mul (p, t1, t0, rp);		/* a^{2^6-1} */
+  ecc_mod_pow_2k (p, rp, t1, 3, t2);	/* a^{2^9-2^3} */
+  ecc_mod_mul (p, t1, t0, rp);		/* a^{2^9-1} */
+  ecc_mod_pow_2k (p, t0, t1, 9, t2);	/* a^{2^18-2^9} */
+  ecc_mod_mul (p, rp, t1, t0);		/* a^{2^18-1} */
+  ecc_mod_sqr (p, t1, rp);		/* a^{2^19-2} */
+  ecc_mod_mul (p, t0, ap, t1);		/* a^{2^19-1} */
+  ecc_mod_pow_2k (p, t1, t0, 18, t2);	/* a^{2^37-2^18} */
+  ecc_mod_mul (p, t0, rp, t1);		/* a^{2^37-1} */
+  ecc_mod_pow_2k (p, t1, t0, 37, t2);	/* a^{2^74-2^37} */
+  ecc_mod_mul (p, rp, t0, t1);		/* a^{2^74-1} */
+  ecc_mod_pow_2k (p, t1, rp, 37, t2);	/* a^{2^111-2^37} */
+  ecc_mod_mul (p, rp, t0, t1);		/* a^{2^111-1} */
+  ecc_mod_pow_2k (p, t1, rp, 111, t2);	/* a^{2^222-2^111} */
+  ecc_mod_mul (p, t0, rp, t1);		/* a^{2^222-1} */
+  ecc_mod_sqr (p, t1, t0);		/* a^{2^223-2} */
+  ecc_mod_mul (p, rp, ap, t1);		/* a^{2^223-1} */
+  ecc_mod_pow_2k (p, t1, rp, 223, t2);	/* a^{2^446-2^223} */
+  ecc_mod_mul (p, rp, t0, t1);		/* a^{2^446-2^222-1} */
+#undef t0
+#undef t1
+#undef t2
+}
+
+/* Needs 9*ECC_LIMB_SIZE scratch space. */
+#define ECC_448_INV_ITCH (9*ECC_LIMB_SIZE)
+
+static void ecc_448_inv (const struct ecc_modulo *p,
+			 mp_limb_t *rp, const mp_limb_t *ap,
+			 mp_limb_t *scratch)
+{
+#define t0 scratch
+
+  ecc_mod_pow_446m224m1 (p, rp, ap, scratch); /* a^{2^446-2^222-1} */
+  ecc_mod_sqr (p, t0, rp);		      /* a^{2^447-2^223-2} */
+  ecc_mod_sqr (p, rp, t0);		      /* a^{2^448-2^224-4} */
+  ecc_mod_mul (p, t0, ap, rp);		      /* a^{2^448-2^224-3} */
+
+  mpn_copyi (rp, t0, ECC_LIMB_SIZE); /* FIXME: Eliminate copy? */
+#undef t0
+}
+
+/* First, do a canonical reduction, then check if zero */
+static int
+ecc_448_zero_p (const struct ecc_modulo *p, mp_limb_t *xp)
+{
+  mp_limb_t cy;
+  mp_limb_t w;
+  mp_size_t i;
+  cy = mpn_sub_n (xp, xp, p->m, ECC_LIMB_SIZE);
+  cnd_add_n (cy, xp, p->m, ECC_LIMB_SIZE);
+
+  for (i = 0, w = 0; i < ECC_LIMB_SIZE; i++)
+    w |= xp[i];
+  return w == 0;
+}
+
+/* Compute x such that x^2 = u/v (mod p). Returns one on success, zero
+   on failure.
+
+   To avoid a separate inversion, we use a trick of djb's, to
+   compute the candidate root as
+
+     x = (u/v)^{(p+1)/4} = u^3 v (u^5 v^3)^{(p-3)/4}.
+*/
+
+/* Needs 4*n space + scratch for ecc_mod_pow_446m224m1. */
+#define ECC_448_SQRT_ITCH (13*ECC_LIMB_SIZE)
+
+static int
+ecc_448_sqrt(const struct ecc_modulo *p, mp_limb_t *rp,
+	     const mp_limb_t *up, const mp_limb_t *vp,
+	     mp_limb_t *scratch)
+{
+#define u3v scratch
+#define u5v3 (scratch + ECC_LIMB_SIZE)
+#define u5v3p (scratch + 2*ECC_LIMB_SIZE)
+#define u2 (scratch + 2*ECC_LIMB_SIZE)
+#define u3 (scratch + 3*ECC_LIMB_SIZE)
+#define uv (scratch + 2*ECC_LIMB_SIZE)
+#define u2v2 (scratch + 3*ECC_LIMB_SIZE)
+
+#define scratch_out (scratch + 4 * ECC_LIMB_SIZE)
+
+#define x2 scratch
+#define vx2 (scratch + ECC_LIMB_SIZE)
+#define t0 (scratch + 2*ECC_LIMB_SIZE)
+
+					/* Live values */
+  ecc_mod_sqr (p, u2, up);		/* u2 */
+  ecc_mod_mul (p, u3, u2, up);		/* u3 */
+  ecc_mod_mul (p, u3v, u3, vp);		/* u3v */
+  ecc_mod_mul (p, uv, up, vp);		/* u3v, uv */
+  ecc_mod_sqr (p, u2v2, uv);		/* u3v, u2v2 */
+  ecc_mod_mul (p, u5v3, u3v, u2v2);	/* u3v, u5v3 */
+  ecc_mod_pow_446m224m1 (p, u5v3p, u5v3, scratch_out); /* u3v, u5v3p */
+  ecc_mod_mul (p, rp, u5v3p, u3v);	/* none */
+
+  /* If square root exists, have v x^2 = u */
+  ecc_mod_sqr (p, x2, rp);
+  ecc_mod_mul (p, vx2, x2, vp);
+  ecc_mod_sub (p, t0, vx2, up);
+
+  return ecc_448_zero_p (p, t0);
+
+#undef u3v
+#undef u5v3
+#undef u5v3p
+#undef u2
+#undef u3
+#undef uv
+#undef u2v2
+#undef scratch_out
+#undef x2
+#undef vx2
+#undef t0
+}
+
+const struct ecc_curve _nettle_curve448 =
+{
+  {
+    448,
+    ECC_LIMB_SIZE,
+    ECC_BMODP_SIZE,
+    0,
+    ECC_448_INV_ITCH,
+    ECC_448_SQRT_ITCH,
+
+    ecc_p,
+    ecc_Bmodp,
+    ecc_Bmodp_shifted,
+    NULL,
+    ecc_pp1h,
+
+    ecc_mod,	      /* FIXME: Implement optimized mod function */
+    ecc_mod,	      /* FIXME: Implement optimized reduce function */
+    ecc_448_inv,
+    ecc_448_sqrt,
+  },
+  {
+    446,
+    ECC_LIMB_SIZE,
+    ECC_BMODQ_SIZE,
+    0,
+    ECC_MOD_INV_ITCH (ECC_LIMB_SIZE),
+    0,
+
+    ecc_q,
+    ecc_Bmodq,
+    ecc_Bmodq_shifted,
+    NULL,
+    ecc_qp1h,
+
+    ecc_mod,	      /* FIXME: Implement optimized mod function */
+    ecc_mod,	      /* FIXME: Implement optimized reduce function */
+    ecc_mod_inv,
+    NULL,
+  },
+
+  0, /* No redc */
+  ECC_PIPPENGER_K,
+  ECC_PIPPENGER_C,
+
+  ECC_ADD_EH_ITCH (ECC_LIMB_SIZE),
+  ECC_ADD_EHH_ITCH (ECC_LIMB_SIZE),
+  ECC_DUP_EH_ITCH (ECC_LIMB_SIZE),
+  ECC_MUL_A_EH_ITCH (ECC_LIMB_SIZE),
+  ECC_MUL_G_EH_ITCH (ECC_LIMB_SIZE),
+  ECC_EH_TO_A_ITCH (ECC_LIMB_SIZE, ECC_448_INV_ITCH),
+
+  ecc_add_eh_untwisted,
+  ecc_add_ehh_untwisted,
+  ecc_dup_eh_untwisted,
+  ecc_mul_a_eh,
+  ecc_mul_g_eh,
+  ecc_eh_to_a,
+
+  ecc_b,
+  ecc_g,
+  ecc_unit,
+  ecc_table
+};
diff --git a/ecc-add-eh.c b/ecc-add-eh.c
index c07ff49a..0b0a1457 100644
--- a/ecc-add-eh.c
+++ b/ecc-add-eh.c
@@ -73,11 +73,11 @@ ecc_add_eh (const struct ecc_curve *ecc,
 #define C (scratch)
 #define D (scratch + 1*ecc->p.size)
 #define T (scratch + 2*ecc->p.size)
-#define E (scratch + 3*ecc->p.size) 
+#define E (scratch + 3*ecc->p.size)
 #define B (scratch + 4*ecc->p.size)
 #define F D
 #define G E
-  
+
   ecc_modp_mul (ecc, C, x1, x2);
   ecc_modp_mul (ecc, D, y1, y2);
   ecc_modp_add (ecc, x3, x1, y1);
@@ -91,7 +91,7 @@ ecc_add_eh (const struct ecc_curve *ecc,
   ecc_modp_add (ecc, C, D, C); /* ! */
   ecc_modp_sqr (ecc, B, z1);
   ecc_modp_sub (ecc, F, B, E);
-  ecc_modp_add (ecc, G, B, E);  
+  ecc_modp_add (ecc, G, B, E);
 
   /* x3 */
   ecc_modp_mul (ecc, B, G, T); /* ! */
@@ -105,3 +105,71 @@ ecc_add_eh (const struct ecc_curve *ecc,
   ecc_modp_mul (ecc, B, F, G);
   mpn_copyi (z3, B, ecc->p.size);
 }
+
+void
+ecc_add_eh_untwisted (const struct ecc_curve *ecc,
+		      mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
+		      mp_limb_t *scratch)
+{
+#define x1 p
+#define y1 (p + ecc->p.size)
+#define z1 (p + 2*ecc->p.size)
+
+#define x2 q
+#define y2 (q + ecc->p.size)
+
+#define x3 r
+#define y3 (r + ecc->p.size)
+#define z3 (r + 2*ecc->p.size)
+
+  /* Formulas (from djb,
+     http://www.hyperelliptic.org/EFD/g1p/auto-edwards-projective.html#doubling-dbl-2007-bl):
+
+     Computation	Operation	Live variables
+
+     C = x1*x2		mul		C
+     D = y1*y2		mul		C, D
+     T = (x1+y1)(x2+y2) - C - D		C, D, T
+     E = b*C*D		2 mul		C, E, T  (Replace C <-- D - C)
+     B = z1^2		sqr		B, C, E, T
+     F = B - E				B, C, E, F, T
+     G = B + E     			C, F, G, T
+     x3 = z1*F*T	3 mul		C, F, G, T
+     y3 = z1*G*(D-C)	2 mul		F, G
+     z3 = F*G		mul
+  */
+#define C (scratch)
+#define D (scratch + 1*ecc->p.size)
+#define T (scratch + 2*ecc->p.size)
+#define E (scratch + 3*ecc->p.size)
+#define B (scratch + 4*ecc->p.size)
+#define F D
+#define G E
+
+  ecc_modp_mul (ecc, C, x1, x2);
+  ecc_modp_mul (ecc, D, y1, y2);
+  ecc_modp_add (ecc, x3, x1, y1);
+  ecc_modp_add (ecc, y3, x2, y2);
+  ecc_modp_mul (ecc, T, x3, y3);
+  ecc_modp_sub (ecc, T, T, C);
+  ecc_modp_sub (ecc, T, T, D);
+  ecc_modp_mul (ecc, x3, C, D);
+  ecc_modp_mul (ecc, E, x3, ecc->b);
+
+  ecc_modp_sub (ecc, C, D, C);
+  ecc_modp_sqr (ecc, B, z1);
+  ecc_modp_sub (ecc, F, B, E);
+  ecc_modp_add (ecc, G, B, E);
+
+  /* x3 */
+  ecc_modp_mul (ecc, B, F, T);
+  ecc_modp_mul (ecc, x3, B, z1);
+
+  /* y3 */
+  ecc_modp_mul (ecc, B, G, z1);
+  ecc_modp_mul (ecc, y3, B, C); /* Clobbers z1 in case r == p. */
+
+  /* z3 */
+  ecc_modp_mul (ecc, B, F, G);
+  mpn_copyi (z3, B, ecc->p.size);
+}
diff --git a/ecc-add-ehh.c b/ecc-add-ehh.c
index 8fdc9ec3..027a6e77 100644
--- a/ecc-add-ehh.c
+++ b/ecc-add-ehh.c
@@ -78,7 +78,7 @@ ecc_add_ehh (const struct ecc_curve *ecc,
 #define C scratch
 #define D (scratch + ecc->p.size)
 #define T (scratch + 2*ecc->p.size)
-#define E (scratch + 3*ecc->p.size) 
+#define E (scratch + 3*ecc->p.size)
 #define A (scratch + 4*ecc->p.size)
 #define B (scratch + 5*ecc->p.size)
 #define F D
@@ -94,7 +94,7 @@ ecc_add_ehh (const struct ecc_curve *ecc,
   ecc_modp_mul (ecc, x3, C, D);
   ecc_modp_mul (ecc, E, x3, ecc->b);
   ecc_modp_add (ecc, C, D, C);	/* ! */
-  
+
   ecc_modp_mul (ecc, A, z1, z2);
   ecc_modp_sqr (ecc, B, A);
 
@@ -113,3 +113,76 @@ ecc_add_ehh (const struct ecc_curve *ecc,
   ecc_modp_mul (ecc, B, F, G);
   mpn_copyi (z3, B, ecc->p.size);
 }
+
+void
+ecc_add_ehh_untwisted (const struct ecc_curve *ecc,
+		       mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
+		       mp_limb_t *scratch)
+{
+#define x1 p
+#define y1 (p + ecc->p.size)
+#define z1 (p + 2*ecc->p.size)
+
+#define x2 q
+#define y2 (q + ecc->p.size)
+#define z2 (q + 2*ecc->p.size)
+
+#define x3 r
+#define y3 (r + ecc->p.size)
+#define z3 (r + 2*ecc->p.size)
+
+  /* Formulas (from djb,
+     http://www.hyperelliptic.org/EFD/g1p/auto-edwards-projective.html#addition-add-2007-bl):
+
+     Computation	Operation	Live variables
+
+     C = x1*x2		mul		C
+     D = y1*y2		mul		C, D
+     T = (x1+y1)(x2+y2) - C - D, mul	C, D, T
+     E = b*C*D		2 mul		C, E, T (Replace C <-- D - C)
+     A = z1*z2		mul		A, C, E, T
+     B = A^2		sqr		A, B, C, E, T
+     F = B - E				A, B, C, E, F, T
+     G = B + E     			A, C, F, G, T
+     x3 = A*F*T		2 mul		A, C, G
+     y3 = A*G*(D-C)	2 mul		F, G
+     z3 = F*G		mul
+  */
+#define C scratch
+#define D (scratch + ecc->p.size)
+#define T (scratch + 2*ecc->p.size)
+#define E (scratch + 3*ecc->p.size)
+#define A (scratch + 4*ecc->p.size)
+#define B (scratch + 5*ecc->p.size)
+#define F D
+#define G E
+
+  ecc_modp_mul (ecc, C, x1, x2);
+  ecc_modp_mul (ecc, D, y1, y2);
+  ecc_modp_add (ecc, A, x1, y1);
+  ecc_modp_add (ecc, B, x2, y2);
+  ecc_modp_mul (ecc, T, A, B);
+  ecc_modp_sub (ecc, T, T, C);
+  ecc_modp_sub (ecc, T, T, D);
+  ecc_modp_mul (ecc, x3, C, D);
+  ecc_modp_mul (ecc, E, x3, ecc->b);
+  ecc_modp_sub (ecc, C, D, C);
+
+  ecc_modp_mul (ecc, A, z1, z2);
+  ecc_modp_sqr (ecc, B, A);
+
+  ecc_modp_sub (ecc, F, B, E);
+  ecc_modp_add (ecc, G, B, E);
+
+  /* x3 */
+  ecc_modp_mul (ecc, B, F, T);
+  ecc_modp_mul (ecc, x3, B, A);
+
+  /* y3 */
+  ecc_modp_mul (ecc, B, G, C);
+  ecc_modp_mul (ecc, y3, B, A);
+
+  /* z3 */
+  ecc_modp_mul (ecc, B, F, G);
+  mpn_copyi (z3, B, ecc->p.size);
+}
diff --git a/ecc-dup-eh.c b/ecc-dup-eh.c
index 2a5c5a07..1b9a3f69 100644
--- a/ecc-dup-eh.c
+++ b/ecc-dup-eh.c
@@ -36,7 +36,7 @@
 #include "ecc.h"
 #include "ecc-internal.h"
 
-/* Double a point on an Edwards curve, in homogeneous coordinates */
+/* Double a point on a twisted Edwards curve, in homogeneous coordinates */
 void
 ecc_dup_eh (const struct ecc_curve *ecc,
 	    mp_limb_t *r, const mp_limb_t *p,
@@ -103,3 +103,56 @@ ecc_dup_eh (const struct ecc_curve *ecc,
   ecc_modp_mul (ecc, b, e, j);
   mpn_copyi (r + 2*ecc->p.size, b, ecc->p.size);
 }
+
+void
+ecc_dup_eh_untwisted (const struct ecc_curve *ecc,
+		      mp_limb_t *r, const mp_limb_t *p,
+		      mp_limb_t *scratch)
+{
+  /* Formulas (from djb,
+     http://www.hyperelliptic.org/EFD/g1p/auto-edwards-projective.html#doubling-dbl-2007-bl):
+
+     Computation	Operation	Live variables
+
+     b = (x+y)^2	sqr		b
+     c = x^2		sqr		b, c
+     d = y^2		sqr		b, c, d
+     e = c+d				b, c, d, e
+     h = z^2		sqr		b, c, d, e, h
+     j = e-2*h				b, c, d, e, j
+     x' = (b-e)*j	mul		c, d, e, j
+     y' = e*(c-d)	mul		e, j
+     z' = e*j		mul
+  */
+#define b scratch
+#define c (scratch  + ecc->p.size)
+#define d (scratch  + 2*ecc->p.size)
+#define e (scratch  + 3*ecc->p.size)
+#define j (scratch  + 4*ecc->p.size)
+
+  /* b */
+  ecc_modp_add (ecc, e, p, p + ecc->p.size);
+  ecc_modp_sqr (ecc, b, e);
+
+  /* c */
+  ecc_modp_sqr (ecc, c, p);
+  /* d */
+  ecc_modp_sqr (ecc, d, p + ecc->p.size);
+  /* h, can use r as scratch, even for in-place operation. */
+  ecc_modp_sqr (ecc, r, p + 2*ecc->p.size);
+  /* e, */
+  ecc_modp_add (ecc, e, c, d);
+  /* j */
+  ecc_modp_add (ecc, r, r, r);
+  ecc_modp_sub (ecc, j, e, r);
+
+  /* x' */
+  ecc_modp_sub (ecc, b, b, e);
+  ecc_modp_mul (ecc, r, b, j);
+  /* y' */
+  ecc_modp_sub (ecc, c, c, d); /* Redundant */
+  ecc_modp_mul (ecc, r + ecc->p.size, e, c);
+  /* z' */
+  ecc_modp_mul (ecc, b, e, j);
+  mpn_copyi (r + 2*ecc->p.size, b, ecc->p.size);
+}
diff --git a/ecc-eh-to-a.c b/ecc-eh-to-a.c
index 8173b887..adb47132 100644
--- a/ecc-eh-to-a.c
+++ b/ecc-eh-to-a.c
@@ -73,7 +73,7 @@ ecc_eh_to_a (const struct ecc_curve *ecc,
 	     is only used by ecdsa code, and ecdsa on Edwards curves
 	     makes little sense and is is only used by tests. */
 	  unsigned shift;
-	  assert (ecc->p.bit_size == 255);
+	  assert (ecc->p.bit_size == 255 || ecc->p.bit_size == 448);
 	  shift = ecc->q.bit_size - 1 - GMP_NUMB_BITS * (ecc->p.size - 1);
 	  cy = mpn_submul_1 (r, ecc->q.m, ecc->p.size,
 			     r[ecc->p.size-1] >> shift);
diff --git a/ecc-internal.h b/ecc-internal.h
index 18c1bf7d..a3116101 100644
--- a/ecc-internal.h
+++ b/ecc-internal.h
@@ -62,6 +62,9 @@
 #define ecc_dup_eh _nettle_ecc_dup_eh
 #define ecc_add_eh _nettle_ecc_add_eh
 #define ecc_add_ehh _nettle_ecc_add_ehh
+#define ecc_dup_eh_untwisted _nettle_ecc_dup_eh_untwisted
+#define ecc_add_eh_untwisted _nettle_ecc_add_eh_untwisted
+#define ecc_add_ehh_untwisted _nettle_ecc_add_ehh_untwisted
 #define ecc_mul_g _nettle_ecc_mul_g
 #define ecc_mul_a _nettle_ecc_mul_a
 #define ecc_mul_g_eh _nettle_ecc_mul_g_eh
@@ -72,6 +75,7 @@
 #define sec_tabselect _nettle_sec_tabselect
 #define sec_modinv _nettle_sec_modinv
 #define curve25519_eh_to_x _nettle_curve25519_eh_to_x
+#define curve448_eh_to_x _nettle_curve448_eh_to_x
 
 extern const struct ecc_curve _nettle_secp_192r1;
 extern const struct ecc_curve _nettle_secp_224r1;
@@ -84,6 +88,7 @@ extern const struct ecc_curve _nettle_secp_521r1;
    different coordinates). And we're not quite ready to provide
    general ecc operations over an arbitrary type of curve. */
 extern const struct ecc_curve _nettle_curve25519;
+extern const struct ecc_curve _nettle_curve448;
 
 #define ECC_MAX_SIZE ((521 + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS)
 
@@ -329,7 +334,7 @@ ecc_add_jjj (const struct ecc_curve *ecc,
 	     mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
 	     mp_limb_t *scratch);
 
-/* Point doubling on an Edwards curve, with homogeneous
+/* Point doubling on a twisted Edwards curve, with homogeneous
    cooordinates. */
 void
 ecc_dup_eh (const struct ecc_curve *ecc,
@@ -346,6 +351,21 @@ ecc_add_ehh (const struct ecc_curve *ecc,
 	     mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
 	     mp_limb_t *scratch);
 
+void
+ecc_dup_eh_untwisted (const struct ecc_curve *ecc,
+		      mp_limb_t *r, const mp_limb_t *p,
+		      mp_limb_t *scratch);
+
+void
+ecc_add_eh_untwisted (const struct ecc_curve *ecc,
+		      mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
+		      mp_limb_t *scratch);
+
+void
+ecc_add_ehh_untwisted (const struct ecc_curve *ecc,
+		       mp_limb_t *r, const mp_limb_t *p, const mp_limb_t *q,
+		       mp_limb_t *scratch);
+
 /* Computes N * the group generator. N is an array of ecc_size()
    limbs. It must be in the range 0 < N < group order, then R != 0,
    and the algorithm can work without any intermediate values getting
@@ -391,6 +411,10 @@ void
 curve25519_eh_to_x (mp_limb_t *xp, const mp_limb_t *p,
 		    mp_limb_t *scratch);
 
+void
+curve448_eh_to_x (mp_limb_t *xp, const mp_limb_t *p,
+		  mp_limb_t *scratch);
+
 /* Current scratch needs: */
 #define ECC_MOD_INV_ITCH(size) (2*(size))
 #define ECC_J_TO_A_ITCH(size) (5*(size))
@@ -415,8 +439,11 @@ curve25519_eh_to_x (mp_limb_t *xp, const mp_limb_t *p,
 #define ECC_MUL_A_EH_ITCH(size) \
   (((3 << ECC_MUL_A_EH_WBITS) + 10) * (size))
 #endif
-#define ECC_ECDSA_SIGN_ITCH(size) (12*(size))
+#define ECC_ECDSA_KEYGEN_ITCH(size) (11*(size))
+#define ECC_ECDSA_SIGN_ITCH(size) (13*(size))
 #define ECC_MOD_RANDOM_ITCH(size) (size)
 #define ECC_HASH_ITCH(size) (1+(size))
 
+#define ECC_MAX(x,y) ((x) > (y) ? (x) : (y))
+
 #endif /* NETTLE_ECC_INTERNAL_H_INCLUDED */
diff --git a/ecc-point-mul-g.c b/ecc-point-mul-g.c
index 46fceb81..02cce0d7 100644
--- a/ecc-point-mul-g.c
+++ b/ecc-point-mul-g.c
@@ -44,15 +44,14 @@
 void
 ecc_point_mul_g (struct ecc_point *r, const struct ecc_scalar *n)
 {
-  TMP_DECL(scratch, mp_limb_t, 3*ECC_MAX_SIZE + ECC_MUL_G_ITCH (ECC_MAX_SIZE));
   const struct ecc_curve *ecc = r->ecc;
   mp_limb_t size = ecc->p.size;
-  mp_size_t itch = 3*size + ecc->mul_g_itch;
+  mp_size_t itch = 3*size + ECC_MAX(ecc->mul_g_itch, ecc->h_to_a_itch);
+  mp_limb_t *scratch = gmp_alloc_limbs (itch);
 
   assert (n->ecc == ecc);
 
-  TMP_ALLOC (scratch, itch);
-
   ecc->mul_g (ecc, scratch, n->p, scratch + 3*size);
   ecc->h_to_a (ecc, 0, r->p, scratch, scratch + 3*size);
+  gmp_free_limbs (scratch, itch);
 }
diff --git a/ecc-point-mul.c b/ecc-point-mul.c
index 2be1c5c4..deb7d8ad 100644
--- a/ecc-point-mul.c
+++ b/ecc-point-mul.c
@@ -46,7 +46,7 @@ ecc_point_mul (struct ecc_point *r, const struct ecc_scalar *n,
 {
   const struct ecc_curve *ecc = r->ecc;
   mp_limb_t size = ecc->p.size;
-  mp_size_t itch = 3*size + ecc->mul_itch;
+  mp_size_t itch = 3*size + ECC_MAX(ecc->mul_itch, ecc->h_to_a_itch);
   mp_limb_t *scratch = gmp_alloc_limbs (itch);
 
   assert (n->ecc == ecc);
diff --git a/ecc-point.c b/ecc-point.c
index 31e3115a..4733b344 100644
--- a/ecc-point.c
+++ b/ecc-point.c
@@ -85,6 +85,21 @@ ecc_point_set (struct ecc_point *p, const mpz_t x, const mpz_t y)
       mpz_mul_ui (rhs, rhs, 121665);
       mpz_clear (x2);
     }
+  else if (p->ecc->p.bit_size == 448)
+    {
+      /* curve448 special case. FIXME: Do in some cleaner way? */
+      mpz_t x2, d;
+      mpz_init (x2);
+      mpz_init_set_ui (d, 39081);
+      mpz_mul (x2, x, x); /* x^2 */
+      mpz_mul (d, d, x2); /* 39081 x^2 */
+      mpz_set_ui (rhs, 1);
+      mpz_submul (rhs, d, lhs); /* 1 - 39081 x^2 y^2 */
+      /* Check that x^2 + y^2 = 1 - 39081 x^2 y^2 */
+      mpz_add (lhs, x2, lhs);	/* x^2 + y^2 */
+      mpz_clear (d);
+      mpz_clear (x2);
+    }
   else
     {
       /* Check that y^2 = x^3 - 3*x + b (mod p) */
diff --git a/eccdata.c b/eccdata.c
index 7cfc33ca..74002c1f 100644
--- a/eccdata.c
+++ b/eccdata.c
@@ -55,10 +55,8 @@ enum ecc_type
   {
     /* y^2 = x^3 - 3x + b (mod p) */
     ECC_TYPE_WEIERSTRASS,
-#if 0
     /* x^2 + y^2 = 1 - d x^2 y^2 */
     ECC_TYPE_EDWARDS,
-#endif
     /* -x^2 + y^2 = 1 - d x^2 y^2 */
     ECC_TYPE_TWISTED_EDWARDS,
   };
@@ -256,6 +254,54 @@ ecc_add (const struct ecc_curve *ecc, struct ecc_point *r,
 	  mpz_clear (y);
 	}
     }
+  else if (ecc->type == ECC_TYPE_EDWARDS)
+    {
+      mpz_t s, t, x, y;
+      mpz_init (s);
+      mpz_init (t);
+      mpz_init (x);
+      mpz_init (y);
+
+      /* t = d p_x p_y q_x q_y */
+      mpz_mul (t, ecc->b, p->x);
+      mpz_mod (t, t, ecc->p);
+      mpz_mul (t, t, p->y);
+      mpz_mod (t, t, ecc->p);
+      mpz_mul (t, t, q->x);
+      mpz_mod (t, t, ecc->p);
+      mpz_mul (t, t, q->y);
+      mpz_mod (t, t, ecc->p);
+
+      /* x' = (p_x q_y + q_x p_y) / (1 + t) */
+      mpz_mul (x, p->x, q->y);
+      mpz_mod (x, x, ecc->p);
+      mpz_addmul (x, q->x, p->y);
+      mpz_mod (x, x, ecc->p);
+      mpz_add_ui (s, t, 1);
+      mpz_invert (s, s, ecc->p);
+      mpz_mul (x, x, s);
+      mpz_mod (x, x, ecc->p);
+
+      /* y' = (p_y q_y - p_x q_x) / (1 - t) */
+      mpz_mul (y, p->y, q->y);
+      mpz_mod (y, y, ecc->p);
+      mpz_submul (y, p->x, q->x);
+      mpz_mod (y, y, ecc->p);
+      mpz_set_ui (s, 1);
+      mpz_sub (s, s, t);
+      mpz_invert (s, s, ecc->p);
+      mpz_mul (y, y, s);
+      mpz_mod (y, y, ecc->p);
+
+      mpz_swap (x, r->x);
+      mpz_swap (y, r->y);
+      r->is_zero = mpz_cmp_ui (r->x, 0) == 0 && mpz_cmp_ui (r->y, 1) == 0;
+
+      mpz_clear (s);
+      mpz_clear (t);
+      mpz_clear (x);
+      mpz_clear (y);
+    }
   else
     {
       /* Untwisted:
@@ -620,6 +666,88 @@ ecc_curve_init (struct ecc_curve *ecc, unsigned bit_size)
 		   "4cf22832ea2f0ff0df38ab61ca32112f");
       break;
 
+    case 448:
+      /* curve448, y^2 = x^3 + 156326 x^2 + x (mod p), with p = 2^{448} - 2^{224} - 1.
+
+	 According to RFC 7748, this is 4-isogenious to the Edwards
+	 curve called "edwards448"
+
+	   x^2 + y^2 = 1 - 39081 x^2 y^2 (mod p).
+
+	 And since the constant is not a square, the Edwards formulas
+	 should be "complete", with no special cases needed for
+	 doubling, neutral element, negatives, etc.
+
+	 Generator is x = 5, with y coordinate
+	 355293926785568175264127502063783334808976399387714271831880898435169088786967410002932673765864550910142774147268105838985595290606362,
+	 according to
+
+	   x = Mod(5, 2^448-2^224-1); sqrt(x^3 + 156326*x^2 + x)
+
+	 in PARI/GP. Also, in PARI notation,
+
+	   curve448 = Mod([0, 156326, 0, 1, 0], 2^448-2^224-1)
+       */
+      ecc_curve_init_str (ecc, ECC_TYPE_EDWARDS,
+			  "fffffffffffffffffffffffffffffff"
+			  "ffffffffffffffffffffffffeffffff"
+			  "fffffffffffffffffffffffffffffff"
+			  "fffffffffffffffffff",
+			  /* -39081 mod p, from PARI/GP
+			     c = Mod(-39081, p)
+			  */
+			  "fffffffffffffffffffffffffffffff"
+			  "ffffffffffffffffffffffffeffffff"
+			  "fffffffffffffffffffffffffffffff"
+			  "fffffffffffffff6756",
+			  /* Order of the subgroup is 2^446 - q_0, where
+			     q_0 = 13818066809895115352007386748515426880336692474882178609894547503885,
+			     224 bits.
+			  */
+			  "3ffffffffffffffffffffffffffffff"
+			  "fffffffffffffffffffffffff7cca23"
+			  "e9c44edb49aed63690216cc2728dc58"
+			  "f552378c292ab5844f3",
+			  "4f1970c66bed0ded221d15a622bf36d"
+			  "a9e146570470f1767ea6de324a3d3a4"
+			  "6412ae1af72ab66511433b80e18b009"
+			  "38e2626a82bc70cc05e",
+			  "693f46716eb6bc248876203756c9c76"
+			  "24bea73736ca3984087789c1e05a0c2"
+			  "d73ad3ff1ce67c39c4fdbd132c4ed7c"
+			  "8ad9808795bf230fa14");
+      ecc->ref = ecc_alloc (3);
+      ecc_set_str (&ecc->ref[0], /* 2 g */
+		   "aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa"
+		   "aaaaaaaaaaaaaaaaaaaaaaa955555555"
+		   "55555555555555555555555555555555"
+		   "5555555555555555",
+		   "ae05e9634ad7048db359d6205086c2b0"
+		   "036ed7a035884dd7b7e36d728ad8c4b8"
+		   "0d6565833a2a3098bbbcb2bed1cda06b"
+		   "daeafbcdea9386ed");
+      ecc_set_str (&ecc->ref[1], /* 3 g */
+		   "865886b9108af6455bd64316cb694333"
+		   "2241b8b8cda82c7e2ba077a4a3fcfe8d"
+		   "aa9cbf7f6271fd6e862b769465da8575"
+		   "728173286ff2f8f",
+		   "e005a8dbd5125cf706cbda7ad43aa644"
+		   "9a4a8d952356c3b9fce43c82ec4e1d58"
+		   "bb3a331bdb6767f0bffa9a68fed02daf"
+		   "b822ac13588ed6fc");
+
+      ecc_set_str (&ecc->ref[2], /* 4 g */
+		   "49dcbc5c6c0cce2c1419a17226f929ea"
+		   "255a09cf4e0891c693fda4be70c74cc3"
+		   "01b7bdf1515dd8ba21aee1798949e120"
+		   "e2ce42ac48ba7f30",
+		   "d49077e4accde527164b33a5de021b97"
+		   "9cb7c02f0457d845c90dc3227b8a5bc1"
+		   "c0d8f97ea1ca9472b5d444285d0d4f5b"
+		   "32e236f86de51839");
+
+      break;
+
     default:
       fprintf (stderr, "No known curve for size %d\n", bit_size);
       exit(EXIT_FAILURE);     
diff --git a/ecdsa-keygen.c b/ecdsa-keygen.c
index fa559a9e..aa2dfb08 100644
--- a/ecdsa-keygen.c
+++ b/ecdsa-keygen.c
@@ -47,9 +47,9 @@ ecdsa_generate_keypair (struct ecc_point *pub,
 			struct ecc_scalar *key,
 			void *random_ctx, nettle_random_func *random)
 {
-  TMP_DECL(p, mp_limb_t, 3*ECC_MAX_SIZE + ECC_MUL_G_ITCH (ECC_MAX_SIZE));
+  TMP_DECL(p, mp_limb_t, 3*ECC_MAX_SIZE + ECC_ECDSA_KEYGEN_ITCH (ECC_MAX_SIZE));
   const struct ecc_curve *ecc = pub->ecc;
-  mp_size_t itch = 3*ecc->p.size + ecc->mul_g_itch;
+  mp_size_t itch = 3*ecc->p.size + ECC_ECDSA_KEYGEN_ITCH (ecc->p.size);
 
   assert (key->ecc == ecc);
 
diff --git a/examples/ecc-benchmark.c b/examples/ecc-benchmark.c
index ea0be173..c149e24d 100644
--- a/examples/ecc-benchmark.c
+++ b/examples/ecc-benchmark.c
@@ -335,6 +335,7 @@ const struct ecc_curve * const curves[] = {
   &_nettle_curve25519,
   &_nettle_secp_256r1,
   &_nettle_secp_384r1,
+  &_nettle_curve448,
   &_nettle_secp_521r1,
 };
 
diff --git a/nettle.texinfo b/nettle.texinfo
index 9a3ca04e..5eb5752d 100644
--- a/nettle.texinfo
+++ b/nettle.texinfo
@@ -115,7 +115,7 @@ Public-key algorithms
 
 * Side-channel silence::
 * ECDSA::
-* Curve 25519::
+* Curve 25519 and Curve 448::
 
 @end detailmenu
 @end menu
@@ -4894,7 +4894,7 @@ curve'' is used as a shorthand for the bitsize of the curve's prime
 @menu
 * Side-channel silence::
 * ECDSA::
-* Curve 25519::
+* Curve 25519 and Curve 448::
 @end menu
 
 @node Side-channel silence, ECDSA, , Elliptic curves
@@ -4928,7 +4928,7 @@ accesses depend only on the size of the input data and its location in
 memory, not on the actual data bits. This implies a performance penalty
 in several of the building blocks.
 
-@node ECDSA, Curve 25519, Side-channel silence, Elliptic curves
+@node ECDSA, Curve 25519 and Curve 448, Side-channel silence, Elliptic curves
 @comment  node-name,  next,  previous,  up
 @subsubsection ECDSA
 
@@ -5032,10 +5032,11 @@ random octets and store them at @code{dst}. For advice, see
 @xref{Randomness}.
 @end deftypefun
 
-@node Curve 25519, , ECDSA, Elliptic curves
+@node Curve 25519 and Curve 448, , ECDSA, Elliptic curves
 @comment  node-name,  next,  previous,  up
-@subsubsection Curve25519
+@subsubsection Curve25519 and Curve448
 @cindex Curve 25519
+@cindex Curve 448
 
 @c FIXME: Make 2^255 pretty in all output formats. Use @sup?
 @c There are other places too (2^32, 2^130).
@@ -5110,6 +5111,40 @@ This function is intended to be compatible with the function
 @code{crypto_scalar_mult} in the NaCl library.
 @end deftypefun
 
+Similarly, Nettle also implements Curve448, an elliptic curve of
+Montgomery type, @math{y^2 = x^3 + 156326 x^2 + x @pmod{p}}, with
+@math{p = 2^448 - 2^224 - 1}.  This particular curve was proposed by
+Mike Hamburg in 2015, for fast Diffie-Hellman key exchange, and is also
+described in @cite{RFC 7748}.
+
+Nettle defines Curve 448 in @file{<nettle/curve448.h>}.
+
+@defvr Constant CURVE448_SIZE
+The octet length of the strings representing curve448 points and scalars, 56.
+@end defvr
+
+@deftypefun void curve448_mul_g (uint8_t *@var{q}, const uint8_t *@var{n})
+Computes @math{Q = N G}, where @math{G} is the group generator and
+@math{N} is an integer. The input argument @var{n} and the output
+argument @var{q} use a little-endian representation of the scalar and
+the x-coordinate, respectively. They are both of size
+@code{CURVE448_SIZE}.
+
+This function is intended to be compatible with the function
+@code{crypto_scalar_mult_base} in the NaCl library.
+@end deftypefun
+
+@deftypefun void curve448_mul (uint8_t *@var{q}, const uint8_t *@var{n}, const uint8_t *@var{p})
+Computes @math{Q = N P}, where @math{P} is an input point and @math{N}
+is an integer. The input arguments @var{n} and @var{p} and the output
+argument @var{q} use a little-endian representation of the scalar and
+the x-coordinates, respectively. They are all of size
+@code{CURVE448_SIZE}.
+
+This function is intended to be compatible with the function
+@code{crypto_scalar_mult} in the NaCl library.
+@end deftypefun
+
 @subsubsection EdDSA
 @cindex eddsa
 
diff --git a/testsuite/.gitignore b/testsuite/.gitignore
index 066bcee2..1e2a69a6 100644
--- a/testsuite/.gitignore
+++ b/testsuite/.gitignore
@@ -18,6 +18,7 @@
 /cnd-memcpy-test
 /ctr-test
 /curve25519-dh-test
+/curve448-dh-test
 /cxx-test
 /des-test
 /des3-test
diff --git a/testsuite/.test-rules.make b/testsuite/.test-rules.make
index efb7df3c..853976e1 100644
--- a/testsuite/.test-rules.make
+++ b/testsuite/.test-rules.make
@@ -232,6 +232,9 @@ dsa-keygen-test$(EXEEXT): dsa-keygen-test.$(OBJEXT)
 curve25519-dh-test$(EXEEXT): curve25519-dh-test.$(OBJEXT)
 	$(LINK) curve25519-dh-test.$(OBJEXT) $(TEST_OBJS) -o curve25519-dh-test$(EXEEXT)
 
+curve448-dh-test$(EXEEXT): curve448-dh-test.$(OBJEXT)
+	$(LINK) curve448-dh-test.$(OBJEXT) $(TEST_OBJS) -o curve448-dh-test$(EXEEXT)
+
 ecc-mod-test$(EXEEXT): ecc-mod-test.$(OBJEXT)
 	$(LINK) ecc-mod-test.$(OBJEXT) $(TEST_OBJS) -o ecc-mod-test$(EXEEXT)
 
diff --git a/testsuite/Makefile.in b/testsuite/Makefile.in
index f8f85701..3dd3c268 100644
--- a/testsuite/Makefile.in
+++ b/testsuite/Makefile.in
@@ -44,7 +44,7 @@ TS_HOGWEED_SOURCES = sexp-test.c sexp-format-test.c \
 		     rsa-sec-decrypt-test.c \
 		     rsa-compute-root-test.c \
 		     dsa-test.c dsa-keygen-test.c \
-		     curve25519-dh-test.c \
+		     curve25519-dh-test.c curve448-dh-test.c \
 		     ecc-mod-test.c ecc-modinv-test.c ecc-redc-test.c \
 		     ecc-sqrt-test.c \
 		     ecc-dup-test.c ecc-add-test.c \
diff --git a/testsuite/curve448-dh-test.c b/testsuite/curve448-dh-test.c
new file mode 100644
index 00000000..7d142d6f
--- /dev/null
+++ b/testsuite/curve448-dh-test.c
@@ -0,0 +1,100 @@
+/* curve448-dh-test.c
+
+   Copyright (C) 2017 Daiki Ueno
+   Copyright (C) 2017 Red Hat, Inc.
+
+   This file is part of GNU Nettle.
+
+   GNU Nettle is free software: you can redistribute it and/or
+   modify it under the terms of either:
+
+     * the GNU Lesser General Public License as published by the Free
+       Software Foundation; either version 3 of the License, or (at your
+       option) any later version.
+
+   or
+
+     * the GNU General Public License as published by the Free
+       Software Foundation; either version 2 of the License, or (at your
+       option) any later version.
+
+   or both in parallel, as here.
+
+   GNU Nettle is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+   General Public License for more details.
+
+   You should have received copies of the GNU General Public License and
+   the GNU Lesser General Public License along with this program.  If
+   not, see http://www.gnu.org/licenses/.
+*/
+
+#include "testutils.h"
+
+#include "curve448.h"
+
+static void
+test_g (const uint8_t *s, const uint8_t *r)
+{
+  uint8_t p[CURVE448_SIZE];
+  curve448_mul_g (p, s);
+  if (!MEMEQ (CURVE448_SIZE, p, r))
+    {
+      printf ("curve448_mul_g failure:\ns = ");
+      print_hex (CURVE448_SIZE, s);
+      printf ("\np = ");
+      print_hex (CURVE448_SIZE, p);
+      printf (" (bad)\nr = ");
+      print_hex (CURVE448_SIZE, r);
+      printf (" (expected)\n");
+      abort ();
+    }
+}
+
+static void
+test_a (const uint8_t *s, const uint8_t *b, const uint8_t *r)
+{
+  uint8_t p[CURVE448_SIZE];
+  curve448_mul (p, s, b);
+  if (!MEMEQ (CURVE448_SIZE, p, r))
+    {
+      printf ("curve448_mul failure:\ns = ");
+      print_hex (CURVE448_SIZE, s);
+      printf ("\nb = ");
+      print_hex (CURVE448_SIZE, b);
+      printf ("\np = ");
+      print_hex (CURVE448_SIZE, p);
+      printf (" (bad)\nr = ");
+      print_hex (CURVE448_SIZE, r);
+      printf (" (expected)\n");
+      abort ();
+    }
+}
+
+void
+test_main (void)
+{
+  /* From RFC 7748. */
+  test_g (H("9a8f4925d1519f5775cf46b04b5800d4ee9ee8bae8bc5565d498c28d"
+	    "d9c9baf574a9419744897391006382a6f127ab1d9ac2d8c0a598726b"),
+	  H("9b08f7cc31b7e3e67d22d5aea121074a273bd2b83de09c63faa73d2c"
+	    "22c5d9bbc836647241d953d40c5b12da88120d53177f80e532c41fa0"));
+  test_g (H("1c306a7ac2a0e2e0990b294470cba339e6453772b075811d8fad0d1d"
+	    "6927c120bb5ee8972b0d3e21374c9c921b09d1b0366f10b65173992d"),
+	  H("3eb7a829b0cd20f5bcfc0b599b6feccf6da4627107bdb0d4f345b430"
+	    "27d8b972fc3e34fb4232a13ca706dcb57aec3dae07bdc1c67bf33609"));
+
+  test_a (H("9a8f4925d1519f5775cf46b04b5800d4ee9ee8bae8bc5565d498c28d"
+	    "d9c9baf574a9419744897391006382a6f127ab1d9ac2d8c0a598726b"),
+	  H("3eb7a829b0cd20f5bcfc0b599b6feccf6da4627107bdb0d4f345b430"
+	    "27d8b972fc3e34fb4232a13ca706dcb57aec3dae07bdc1c67bf33609"),
+	  H("07fff4181ac6cc95ec1c16a94a0f74d12da232ce40a77552281d282b"
+	    "b60c0b56fd2464c335543936521c24403085d59a449a5037514a879d"));
+  test_a (H("1c306a7ac2a0e2e0990b294470cba339e6453772b075811d8fad0d1d"
+	    "6927c120bb5ee8972b0d3e21374c9c921b09d1b0366f10b65173992d"),
+	  H("9b08f7cc31b7e3e67d22d5aea121074a273bd2b83de09c63faa73d2c"
+	    "22c5d9bbc836647241d953d40c5b12da88120d53177f80e532c41fa0"),
+	  H("07fff4181ac6cc95ec1c16a94a0f74d12da232ce40a77552281d282b"
+	    "b60c0b56fd2464c335543936521c24403085d59a449a5037514a879d"));
+}
diff --git a/testsuite/ecc-add-test.c b/testsuite/ecc-add-test.c
index ad2bd292..ed4eed57 100644
--- a/testsuite/ecc-add-test.c
+++ b/testsuite/ecc-add-test.c
@@ -19,12 +19,14 @@ test_main (void)
 
       ecc_a_to_j (ecc, g, ecc->g);
 
-      if (ecc->p.bit_size == 255)
+      if (ecc->p.bit_size == 255 || ecc->p.bit_size == 448)
 	{
 	  mp_limb_t *z = xalloc_limbs (ecc_size_j (ecc));
 
-	  ASSERT (ecc->add_hh == ecc_add_eh);
-	  ASSERT (ecc->add_hhh == ecc_add_ehh);
+	  ASSERT ((ecc->p.bit_size == 255 && ecc->add_hh == ecc_add_eh)
+		  || (ecc->p.bit_size == 448 && ecc->add_hh == ecc_add_eh_untwisted));
+	  ASSERT ((ecc->p.bit_size == 255 && ecc->add_hhh == ecc_add_ehh)
+		  || (ecc->p.bit_size == 448 && ecc->add_hhh == ecc_add_ehh_untwisted));
 	  ASSERT (ecc->add_hh_itch <= ecc->add_hhh_itch);
 
 	  /* Zero point has x = 0, y = 1, z = 1 */
diff --git a/testsuite/ecc-dup-test.c b/testsuite/ecc-dup-test.c
index 0ae4444a..2499c130 100644
--- a/testsuite/ecc-dup-test.c
+++ b/testsuite/ecc-dup-test.c
@@ -14,11 +14,12 @@ test_main (void)
 
       ecc_a_to_j (ecc, g, ecc->g);
 
-      if (ecc->p.bit_size == 255)
+      if (ecc->p.bit_size == 255 || ecc->p.bit_size == 448)
 	{
 	  mp_limb_t *z = xalloc_limbs (ecc_size_j (ecc));
 
-	  ASSERT (ecc->dup == ecc_dup_eh);
+	  ASSERT ((ecc->p.bit_size == 255 && ecc->dup == ecc_dup_eh)
+		  || (ecc->p.bit_size == 448 && ecc->dup == ecc_dup_eh_untwisted));
 
 	  /* Zero point has x = 0, y = 1, z = 1 */
 	  mpn_zero (z, 3*ecc->p.size);
diff --git a/testsuite/ecc-mul-a-test.c b/testsuite/ecc-mul-a-test.c
index 245016aa..019f4d34 100644
--- a/testsuite/ecc-mul-a-test.c
+++ b/testsuite/ecc-mul-a-test.c
@@ -17,7 +17,7 @@ test_main (void)
       mp_limb_t *p = xalloc_limbs (ecc_size_j (ecc));
       mp_limb_t *q = xalloc_limbs (ecc_size_j (ecc));
       mp_limb_t *n = xalloc_limbs (size);
-      mp_limb_t *scratch = xalloc_limbs (ecc->mul_itch);
+      mp_limb_t *scratch = xalloc_limbs (ecc->mul_itch + ecc->h_to_a_itch);
       unsigned j;
       
       mpn_zero (n, size);
@@ -39,7 +39,7 @@ test_main (void)
       mpn_sub_1 (n, ecc->q.m, size, 1);
       ecc->mul (ecc, p, n, ecc->g, scratch);
       ecc->h_to_a (ecc, 0, p, p, scratch);
-      if (ecc->p.bit_size == 255)
+      if (ecc->p.bit_size == 255 || ecc->p.bit_size == 448)
 	/* For edwards curves, - (x,y ) == (-x, y). FIXME: Swap x and
 	   y, to get identical negation? */
 	mpn_sub_n (p, ecc->p.m, p, size);
diff --git a/testsuite/ecc-mul-g-test.c b/testsuite/ecc-mul-g-test.c
index 27239484..0bedfdea 100644
--- a/testsuite/ecc-mul-g-test.c
+++ b/testsuite/ecc-mul-g-test.c
@@ -17,7 +17,7 @@ test_main (void)
       mp_limb_t *p = xalloc_limbs (ecc_size_j (ecc));
       mp_limb_t *q = xalloc_limbs (ecc_size_j (ecc));
       mp_limb_t *n = xalloc_limbs (size);
-      mp_limb_t *scratch = xalloc_limbs (ecc->mul_g_itch);
+      mp_limb_t *scratch = xalloc_limbs (ecc->mul_g_itch + ecc->h_to_a_itch);
 
       mpn_zero (n, size);
 
@@ -41,7 +41,7 @@ test_main (void)
       mpn_sub_1 (n, ecc->q.m, size, 1);
       ecc->mul_g (ecc, p, n, scratch);
       ecc->h_to_a (ecc, 0, p, p, scratch);
-      if (ecc->p.bit_size == 255)
+      if (ecc->p.bit_size == 255 || ecc->p.bit_size == 448)
 	/* For edwards curves, - (x,y ) == (-x, y). FIXME: Swap x and
 	   y, to get identical negation? */
 	mpn_sub_n (p, ecc->p.m, p, size);
diff --git a/testsuite/ecdh-test.c b/testsuite/ecdh-test.c
index 0b39319d..ff4f7233 100644
--- a/testsuite/ecdh-test.c
+++ b/testsuite/ecdh-test.c
@@ -52,7 +52,8 @@ set_scalar (struct ecc_scalar *s,
 {
   mpz_t X;
   mpz_init_set_str (X, x, 0);
-  ecc_scalar_set (s, X);
+  if (!ecc_scalar_set (s, X))
+    abort ();
   mpz_clear (X);
 }
 
@@ -240,4 +241,17 @@ test_main(void)
 	   "38072138078045635808869930165213470653418146012939584392304609812494425185763",
 	   "10481077163111981870382976851703705086808805457403127024129174358161599078055",
 	   "29260211489972704256554624312266763530759418996739976957020673870747051409679");
+
+  /* NOTE: This isn't the standard way to do curve448
+     diffie-hellman, but it tests that the ecc_point interface works
+     also with curve448. */
+  test_dh ("curve448", &_nettle_curve448,
+	   "129458936807933142766404648460937163205634163580407624950524900086792185737444124895392953822100034523565454893159084960036749128566328",
+	   "23903108874160330022289088207864530114505726115081678533913226179385920277612083777349117962138808929878378666596532036566924169949084",
+	   "693683143993815499711046966874265987454661213870193324674425656110752379002105414428569086535475560314058341102862207145978150379762153",
+	   "66424594649188102315894632429895338306697492782714758296415311427244880255966850729749965592839835963032731282879151354354178946253531",
+	   "411851112596680430188999894591634506976361833537024658040418853047370769553774913299417695327870642536912872558385293694714169201128264",
+	   "337433451779159274143076131600929733721586133908369086734805607026091240174740218929467625260731556550599267570314197354864315711490353",
+	   "224725768629972498035446273711269105191383993674106563435257119903436206484342709996926420948730961128941009070083709026343858723205213",
+	   "514544926219850986487923720424370435708360925070646212523588162169142573918197583804309386017625350764529605929374479238949748203847320");
 }
diff --git a/testsuite/ecdsa-keygen-test.c b/testsuite/ecdsa-keygen-test.c
index a96c09ef..cd96782e 100644
--- a/testsuite/ecdsa-keygen-test.c
+++ b/testsuite/ecdsa-keygen-test.c
@@ -40,6 +40,22 @@ ecc_valid_p (struct ecc_point *pub)
 
       mpz_clear (x2);
     }
+  else if (pub->ecc->p.bit_size == 448)
+    {
+      /* Check that
+	 x^2 + y^2 = 1 - 39081 x^2 y^2 */
+      mpz_t x2, d;
+      mpz_init (x2);
+      mpz_init_set_ui (d, 39081);
+      mpz_mul (x2, x, x); /* x^2 */
+      mpz_mul (d, d, x2); /* 39081 x^2 */
+      mpz_set_ui (rhs, 1);
+      mpz_submul (rhs, d, lhs); /* 1 - 39081 x^2 y^2 */
+      mpz_add (lhs, x2, lhs);	/* x^2 + y^2 */
+
+      mpz_clear (d);
+      mpz_clear (x2);
+    }
   else
     {
       /* Check y^2 = x^3 - 3 x + b */
diff --git a/testsuite/testutils.c b/testsuite/testutils.c
index c9f21bab..a3ad295a 100644
--- a/testsuite/testutils.c
+++ b/testsuite/testutils.c
@@ -1673,6 +1673,7 @@ const struct ecc_curve * const ecc_curves[] = {
   &_nettle_secp_384r1,
   &_nettle_secp_521r1,
   &_nettle_curve25519,
+  &_nettle_curve448,
   NULL
 };
 
@@ -1724,7 +1725,7 @@ void
 test_ecc_mul_a (unsigned curve, unsigned n, const mp_limb_t *p)
 {
   /* For each curve, the points 2 g, 3 g and 4 g */
-  static const struct ecc_ref_point ref[6][3] = {
+  static const struct ecc_ref_point ref[7][3] = {
     { { "dafebf5828783f2ad35534631588a3f629a70fb16982a888",
 	"dd6bda0d993da0fa46b27bbc141b868f59331afa5c7e93ab" },
       { "76e32a2557599e6edcd283201fb2b9aadfd0d359cbb263da",
@@ -1785,15 +1786,22 @@ test_ecc_mul_a (unsigned curve, unsigned n, const mp_limb_t *p)
 	"1267b1d177ee69aba126a18e60269ef79f16ec176724030402c3684878f5b4d4" },
       { "203da8db56cff1468325d4b87a3520f91a739ec193ce1547493aa657c4c9f870",
 	"47d0e827cb1595e1470eb88580d5716c4cf22832ea2f0ff0df38ab61ca32112f" },
+    },
+    { { "aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa955555555555555555555555555555555555555555555555555555555",
+	"ae05e9634ad7048db359d6205086c2b0036ed7a035884dd7b7e36d728ad8c4b80d6565833a2a3098bbbcb2bed1cda06bdaeafbcdea9386ed" },
+      { "865886b9108af6455bd64316cb6943332241b8b8cda82c7e2ba077a4a3fcfe8daa9cbf7f6271fd6e862b769465da8575728173286ff2f8f",
+	"e005a8dbd5125cf706cbda7ad43aa6449a4a8d952356c3b9fce43c82ec4e1d58bb3a331bdb6767f0bffa9a68fed02dafb822ac13588ed6fc" },
+      { "49dcbc5c6c0cce2c1419a17226f929ea255a09cf4e0891c693fda4be70c74cc301b7bdf1515dd8ba21aee1798949e120e2ce42ac48ba7f30",
+	"d49077e4accde527164b33a5de021b979cb7c02f0457d845c90dc3227b8a5bc1c0d8f97ea1ca9472b5d444285d0d4f5b32e236f86de51839" },
     }
   };
-  assert (curve < 6);
+  assert (curve < 7);
   assert (n <= 4);
   if (n == 0)
     {
       /* Makes sense for curve25519 only */
       const struct ecc_curve *ecc = ecc_curves[curve];
-      assert (ecc->p.bit_size == 255);
+      assert (ecc->p.bit_size == 255 || ecc->p.bit_size == 448);
       if (!mpn_zero_p (p, ecc->p.size)
 	  || mpn_cmp (p + ecc->p.size, ecc->unit, ecc->p.size) != 0)
 	{