1 files changed, 50 insertions, 0 deletions
diff --git a/src/zmodpow.c b/src/zmodpow.c
new file mode 100644
index 0000000..42bed2f
--- /dev/null
+++ b/src/zmodpow.c
@@ -0,0 +1,50 @@
+/* See LICENSE file for copyright and license details. */
+#include "internals"
+
+#include <errno.h>
+
+#define tb  libzahl_tmp_pow_b
+#define tc  libzahl_tmp_pow_c
+#define td  libzahl_tmp_pow_d
+
+
+void
+zmodpow(z_t a, z_t b, z_t c, z_t d)
+{
+	size_t i, n;
+
+	if (zsignum(c) <= 0) {
+		if (zzero(c)) {
+			if (zzero(b)) {
+				errno = EDOM; /* Indeterminate form: 0:th power of 0 */
+				FAILURE_JUMP();
+			}
+			zsetu(a, 1);
+		} else if (zzero(b) || zzero(d)) {
+			errno = EDOM; /* Undefined form: Division by 0 */
+			FAILURE_JUMP();
+		} else {
+			SET_SIGNUM(a, 0);
+		}
+		return;
+	} else if (zzero(d)) {
+		errno = EDOM; /* Undefined form: Division by 0 */
+		FAILURE_JUMP();
+	} else if (zzero(b)) {
+		SET_SIGNUM(a, 0);
+		return;
+	}
+
+	n = zbits(c);
+
+	zmod(tb, b, d);
+	zset(tc, c);
+	zset(td, d);
+	zsetu(a, 1);
+
+	for (i = 0; i < n; i++) {
+		if (zbtest(tc, i))
+			zmodmul(a, a, tb, td);
+		zmodsqr(tb, tb, td);
+	}
+}