diff -Nru /n/sources/plan9/sys/src/libauthsrv/decaf.mp /sys/src/libauthsrv/decaf.mp
--- /n/sources/plan9/sys/src/libauthsrv/decaf.mp	Wed Dec 31 19:00:00 1969
+++ /sys/src/libauthsrv/decaf.mp	Wed Jun 26 15:59:52 2024
@@ -0,0 +1,49 @@
+# negate r when n > (p-1)/2
+decaf_neg(p, n, r) {
+	mod(p) m = -r;
+	r = n > (p-1)>>1 ? m : r;
+}
+
+#	field F_p
+#	curve a*x**2+y**2==1+d*x**2*y**2
+#	input X,Y,Z,T (extended coordinates)
+decaf_encode(p,a,d, X,Y,Z,T, s) mod(p) {
+	r = misqrt((a-d)*(Z+Y)*(Z-Y), p);
+	u = (a-d)*r;
+	decaf_neg(p, -2*u*Z, r);
+	s = u*(r*(a*Z*X-d*Y*T)+Y)/a;
+	decaf_neg(p, s, s);
+}
+
+#	field F_p
+#	curve a*x**2+y**2==1+d*x**2*y**2
+#	input s
+#	output in extended coordinates
+decaf_decode(p,a,d, s, ok,X,Y,Z,T) {
+	if(s > (p-1)>>1){
+		ok = 0;
+	} else mod(p) {
+		ss = s^2;
+		Z = 1+a*ss;
+		u = Z^2 - 4*d*ss;
+		v = u*ss;
+		if(v == 0)
+			ok = 1;
+		else {
+			ok = msqrt(v, p);
+			if(ok != 0){
+				v = 1/ok;
+				ok = 1;
+			}
+		}
+		if(ok != 0) {
+			decaf_neg(p, u*v, v);
+			w = v * s * (2-Z);
+			if(s == 0)
+				w = w + 1;
+			X = 2*s;
+			Y = w * Z;
+			T = w * X;
+		}
+	}
+}
diff -Nru /n/sources/plan9/sys/src/libauthsrv/ed448.mp /sys/src/libauthsrv/ed448.mp
--- /n/sources/plan9/sys/src/libauthsrv/ed448.mp	Wed Dec 31 19:00:00 1969
+++ /sys/src/libauthsrv/ed448.mp	Wed Jun 26 16:01:31 2024
@@ -0,0 +1,10 @@
+# Edwards Curve Ed448-Goldilocks
+
+# x^2+y^2 = 1-39081x^2y^2
+ed448_curve(p,a,d,x,y) {
+	p = 2^448 - 2^224 - 1;
+	a = 1;
+	d = -39081;
+	x = 117812161263436946737282484343310064665180535357016373416879082147939404277809514858788439644911793978499419995990477371552926308078495;
+	y = 19;
+}
diff -Nru /n/sources/plan9/sys/src/libauthsrv/edwards.mp /sys/src/libauthsrv/edwards.mp
--- /n/sources/plan9/sys/src/libauthsrv/edwards.mp	Wed Dec 31 19:00:00 1969
+++ /sys/src/libauthsrv/edwards.mp	Wed Jun 26 16:00:48 2024
@@ -0,0 +1,40 @@
+# Edwards curve arithmetic
+edwards_add(p,a,d, X1,Y1,Z1,T1, X2,Y2,Z2,T2, X3,Y3,Z3,T3) mod(p) {
+	A = X1*X2;
+	B = Y1*Y2;
+	C = d*T1*T2;
+	D = Z1*Z2;
+	E = (X1+Y1)*(X2+Y2);
+	E = E - A - B;
+	F = D - C;
+	G = D + C;
+	H = B - a*A;
+	X3 = E*F;
+	Y3 = G*H;
+	Z3 = F*G;
+	T3 = E*H;
+}
+edwards_sel(s, X1,Y1,Z1,T1, X2,Y2,Z2,T2, X3,Y3,Z3,T3){
+	X3 = s != 0 ? X1 : X2;
+	Y3 = s != 0 ? Y1 : Y2;
+	Z3 = s != 0 ? Z1 : Z2;
+	T3 = s != 0 ? T1 : T2;
+}
+edwards_new(x,y,z,t, X,Y,Z,T) {
+	X = x;
+	Y = y;
+	Z = z;
+	T = t;
+}
+edwards_scale(p,a,d, s, X1,Y1,Z1,T1, X3,Y3,Z3,T3) {
+	X2,Y2,Z2,T2 = edwards_new(X1,Y1,Z1,T1);
+	X4,Y4,Z4,T4 = edwards_new( 0, 1, 1, 0);
+	X3,Y3,Z3,T3 = edwards_sel(s % 2, X2,Y2,Z2,T2, X4,Y4,Z4,T4);
+	k = s >> 1; j = p >> 1;
+	while(j != 0){
+		X2,Y2,Z2,T2 = edwards_add(p,a,d, X2,Y2,Z2,T2, X2,Y2,Z2,T2);
+		X4,Y4,Z4,T4 = edwards_add(p,a,d, X2,Y2,Z2,T2, X3,Y3,Z3,T3);
+		X3,Y3,Z3,T3 = edwards_sel(k % 2, X4,Y4,Z4,T4, X3,Y3,Z3,T3);
+		k = k >> 1; j = j >> 1;
+	}
+}
diff -Nru /n/sources/plan9/sys/src/libauthsrv/elligator2.mp /sys/src/libauthsrv/elligator2.mp
--- /n/sources/plan9/sys/src/libauthsrv/elligator2.mp	Wed Dec 31 19:00:00 1969
+++ /sys/src/libauthsrv/elligator2.mp	Wed Jun 26 16:01:03 2024
@@ -0,0 +1,30 @@
+#elligator2:
+#	curve a*x^2+y^2==1+d*x^2*y^2
+#	input r0
+#	n is any non-square
+#
+elligator2(p,a,d, n, r0, X,Y,Z,T) mod(p) {
+	r = n*r0*r0;
+	D = (d*r+a-d)*(d*r-a*r-d);
+	N = (r+1)*(a-2*d);
+	ND = N*D;
+	if(ND == 0) {
+		c = 1;
+		e = 0;
+	} else {
+		e = msqrt(ND, p);
+		if(e != 0) {
+			c = 1;
+			e = 1/e;
+		} else {
+			c = -1;
+			e = n*r0*misqrt(n*ND, p);
+		}
+	}
+	s = c*N*e;
+	t = -c*N*(r-1)*((a-2*d)*e)^2-1;
+	X = 2*s*t;
+	Y = (1-a*s*s)*(1+a*s*s);
+	Z = (1+a*s*s)*t;
+	T = (2*s)*(1-a*s*s);
+}
diff -Nru /n/sources/plan9/sys/src/libauthsrv/mkfile /sys/src/libauthsrv/mkfile
--- /n/sources/plan9/sys/src/libauthsrv/mkfile	Wed Oct 23 12:16:16 2002
+++ /sys/src/libauthsrv/mkfile	Wed Jun 26 15:56:34 2024
@@ -21,10 +21,24 @@
 HFILES=\
 	/sys/include/authsrv.h\
 
+MPCFILES=\
+	msqrt.mpc\
+	decaf.mpc\
+	edwards.mpc\
+	elligator2.mpc\
+	spake2ee.mpc\
+	ed448.mpc\
+
 UPDATE=\
 	mkfile\
 	$HFILES\
 	${OFILES:%.$O=%.c}\
+	${MPCFILES:.mpc=%.mp}\
 	${LIB:/$objtype/%=/386/%}\
 
+CLEANFILES=$MPCFILES
+
 </sys/src/cmd/mksyslib
+
+%.mpc:D:	%.mp
+	mpc $stem.mp > $target
diff -Nru /n/sources/plan9/sys/src/libauthsrv/msqrt.mp /sys/src/libauthsrv/msqrt.mp
--- /n/sources/plan9/sys/src/libauthsrv/msqrt.mp	Wed Dec 31 19:00:00 1969
+++ /sys/src/libauthsrv/msqrt.mp	Wed Jun 26 16:05:24 2024
@@ -0,0 +1,100 @@
+# derived from: http://eli.thegreenplace.net/2009/03/07/computing-square-roots-in-python
+
+# Compute the Legendre symbol a|p using Euler's criterion.
+# p is a prime, a is relatively prime to p (if p divides a,
+# then a|p = 0)
+legendresymbol(a, p, r) {
+	pm1 = p-1;
+	mod(p) r = a^(pm1>>1);
+	if(r == pm1)
+		r = -1;
+}
+
+# Find a quadratic residue (mod p) of 'a'. p must be an
+# odd prime.
+#
+# Solve the congruence of the form:
+#	x^2 = a (mod p)
+# And returns x. Node that p - x is also a root.
+#
+# 0 is returned if no square root exists for these
+# a and p.
+#
+# The Tonelli-Shanks algorithm is used (except
+# for some simple cases in which the solution is known
+# from an identity).
+msqrt(a, p, r) {
+	if(legendresymbol(a, p) != 1)
+		r = 0;
+	else if(a == 0)
+		r = 0;
+	else if(p == 2)
+		r = a;
+	else if(p%4 == 3){
+		e = p+1 >> 2;
+		mod(p) r = a^e;
+	} else {
+		# Partition p-1 to s * 2^e for an odd s (i.e.
+		# reduce all the powers of 2 from p-1)
+		s = p-1;
+		e = 0;
+		while(s%2 == 0){
+			s = s >> 1;
+			e = e + 1;
+		}
+
+		# Find some 'n' with a legendre symbol n|p = -1.
+		# Shouldn't take long.
+		n = 2;
+		while(legendresymbol(n, p) != -1)
+			n = n + 1;
+
+		# x is a guess of the square root that gets better
+		# with each iteration.
+		# b is the "fudge factor" - by now much we're off
+		# with the guess. The invariant x^2 == a*b (mod p)
+		# is maintained throughout the loop.
+		# g is used for successive powers of n to update
+		# both a and b
+		# e is the exponent - decreases with each update
+		mod(p){
+			x = a^(s+1 >> 1);
+			b = a^s;
+			g = n^s;
+		}
+		while(1==1){
+			t = b;
+			m = 0;
+			while(m < e){
+				if(t == 1)
+					break;
+				t = t*t % p;
+				m = m + 1;
+			}
+			if(m == 0){
+				r = x;
+				break;
+			}
+			t = 2^(e-m-1);
+			mod(p){
+				gs = g^t;
+				g = gs*gs;
+				x = x*gs;
+				b = b*g;
+			}
+			e = m;
+		}
+	}
+}
+
+# modular inverse square-root
+misqrt(a, p, r) {
+	if((p % 4) == 3){
+		e = ((p-3)>>2);
+		mod(p) r = a^e;
+	} else {
+		r = msqrt(a, p);
+		if(r != 0)
+			mod(p) r = 1/r;
+	}
+}
diff -Nru /n/sources/plan9/sys/src/libauthsrv/spake2ee.mp /sys/src/libauthsrv/spake2ee.mp
--- /n/sources/plan9/sys/src/libauthsrv/spake2ee.mp	Wed Dec 31 19:00:00 1969
+++ /sys/src/libauthsrv/spake2ee.mp	Wed Jun 26 16:01:17 2024
@@ -0,0 +1,35 @@
+#
+# this implements a variant of SPAKE2 Elligator edition described in:
+# https://www.mail-archive.com/curves@moderncrypto.org/msg00412.html
+#
+
+# derive points PM or PN from a (password) hash
+spake2ee_h2P(p,a,d, h, PX,PY,PZ,PT){
+	# find a small non-square for elligator
+	n = 2;
+	while(legendresymbol(n, p) != -1)
+		n = n + 1;
+	PX,PY,PZ,PT = elligator2(p,a,d, n, h%p);
+}
+
+# Ya = xa*G+PM, Yb = xb*G+PN
+spake2ee_1(p,a,d, x, GX,GY, PX,PY,PZ,PT, y){
+	mod(p) X,Y,Z,T = edwards_scale(p,a,d, x, GX,GY,1,GX*GY);
+	X,Y,Z,T = edwards_add(p,a,d, X,Y,Z,T, PX,PY,PZ,PT);
+	y = decaf_encode(p,a,d, X,Y,Z,T);
+}
+
+# Z = xa*(Yb-PN)
+#   = xa*(xb*G+PN-PN)
+#   = xa*xb*G
+#   = xb*xa*G
+#   = xb*(xa*G+PM-PM)
+#   = xb*(Ya-PM)
+spake2ee_2(p,a,d, PX,PY,PZ,PT, x, y, ok, z){
+	ok, X,Y,Z,T = decaf_decode(p,a,d, y);
+	if(ok != 0){
+		mod(p) X,Y,Z,T = edwards_add(p,a,d, X,Y,Z,T, -PX,PY,PZ,-PT);
+		X,Y,Z,T = edwards_scale(p,a,d, x, X,Y,Z,T);
+		z = decaf_encode(p,a,d, X,Y,Z,T);
+	}
+}