diff --git a/include/relic_ep.h b/include/relic_ep.h
index d4a1416d6..e223d5f79 100644
--- a/include/relic_ep.h
+++ b/include/relic_ep.h
@@ -65,6 +65,8 @@ enum {
     EP_B12,
     /* New family from Fotiadis-Martindale family with embedding degree 16. */
     EP_N16,
+	/* Fotiadis-Martindale famile with embedding degree 16. */
+	EP_FM16,
     /* Kachisa-Schaefer-Scott family with embedding degree 16. */
     EP_K16,
     /* Fotiadis-Martindale family with embedding degree 18. */
@@ -175,6 +177,8 @@ enum {
     SG18_P638,
 	/** New family with embeeding degree 16. */
 	N16_P765,
+	/* Fotiadis-Moartindale with embedding degree 16. */
+	FM16_P765,
 	/** Kachisa-Schaefer-Scott with embedding degree 16. */
 	K16_P766,
 	/** 1536-bit supersingular curve. */
diff --git a/include/relic_fp.h b/include/relic_fp.h
index e5bf72632..9c3b7ce4a 100644
--- a/include/relic_fp.h
+++ b/include/relic_fp.h
@@ -158,6 +158,8 @@ enum {
     SG18_638,
 	/** 765-bit prime for new family with embedding degree 16. */
 	N16_765,
+	/** 765-bit prime for FM curve with embeddind degree 16. */
+	FM16_765,
 	/** 766-bit prime for KSS curve with embedding degree 16. */
 	K16_766,
 	/** 1024-bit prime for CTIDH. */
diff --git a/src/ep/relic_ep_map.c b/src/ep/relic_ep_map.c
index 55d2b5f78..f02c080c7 100644
--- a/src/ep/relic_ep_map.c
+++ b/src/ep/relic_ep_map.c
@@ -240,11 +240,12 @@ void ep_map_swift(ep_t p, const uint8_t *msg, size_t len) {
 	/* enough space for two field elements plus extra bytes for uniformity */
 	const size_t len_per_elm = (FP_PRIME + ep_param_level() + 7) / 8;
 	uint8_t s, *pseudo_random_bytes = RLC_ALLOCA(uint8_t, 2 * len_per_elm + 1);
-	fp_t c, t, u, v, w, y, x1, y1, z1;
+	fp_t a, c, t, u, v, w, y, x1, y1, z1;
 	ctx_t *ctx = core_get();
 	bn_t k;
 
 	bn_null(k);
+	fp_null(a);
 	fp_null(c);
 	fp_null(t);
 	fp_null(u);
@@ -276,77 +277,153 @@ void ep_map_swift(ep_t p, const uint8_t *msg, size_t len) {
 		fp_prime_conv(t, k);
 		s = pseudo_random_bytes[2 * len_per_elm] & 1;
 
+		fp_copy(a, ep_curve_get_a());
+
 		if ((ep_curve_opt_b() == RLC_ZERO) && (ctx->mod8 == 1)) {
 			/* This is the approach due to Koshelev introduced in
 			 * https://eprint.iacr.org/2021/1034.pdf */
 			
-			/* Compute t^2 = 3c*sqrt(a)*(2c^3*x^6 - 3*c^2*x^4 - 3*c*x^2 + 2).*/
-			/* Compute w = 3*c. */
-			fp_set_dig(c, -fp_prime_get_qnr());
-			fp_neg(c, c);
-			fp_dbl(w, c);
-			fp_add(w, w, c);
-
-			/* Compute x^2, x^4 and x^6 in sequence. */
-			fp_sqr(z1, u);
-			fp_sqr(y1, z1);
-			fp_mul(t, z1, y1);
-
-			fp_dbl(t, t);
-			fp_mul(t, t, c);
-			fp_mul(t, t, c);
-			fp_mul(t, t, c);
-
-			fp_mul(v, y1, c);
-			fp_mul(v, v, w);
-			fp_sub(t, t, v);
-
-			/* v = -3*c*x^2. */
-			fp_mul(v, w, z1);
-			fp_neg(v, v);
-			fp_add(t, t, v);
-			fp_add_dig(t, t, 2);
-
-			/* Assume a = 1 for simplicitly. */
-			fp_mul(t, t, w);
-			fp_mul(t, t, ctx->ep_map_c[6]);
-			dig_t c1 = fp_is_sqr(t);
-			/* If t is not square, compute u = 1/(uc), t = sqrt(t/c)/(c*u^3)*/
-			fp_inv(v, c);
-			fp_inv(x1, u);
-			fp_mul(y1, t, v);
-			/* If t is a square, extract its square root. */
-			dv_copy_cond(t, y1, RLC_FP_DIGS, !c1);
-			fp_srt(t, t);
-			fp_mul(y1, t, v);
-			fp_sqr(y, x1);
-			fp_mul(y, y, x1);
-			fp_mul(y1, y1, y);
-			fp_mul(x1, x1, v);
-			dv_copy_cond(u, x1, RLC_FP_DIGS, !c1);
-			dv_copy_cond(t, y1, RLC_FP_DIGS, !c1);
-
-			/* Compute x = sqrt(a)*(c*x^2 - 2)/(-3*c*x^2). */
-			fp_sqr(z1, u);
-			fp_mul(v, w, z1);
-			fp_neg(v, v);
-			fp_inv(v, v);
-			fp_mul(p->x, z1, c);
-			fp_sub_dig(p->x, p->x, 2);
-			fp_mul(p->x, p->x, v);
-			fp_mul(p->x, p->x, ctx->ep_map_c[6]);
-
-			/* Compute y = y*2*sqrt(a)/(3^2*c^2*x^3). */
-			fp_mul(z1, z1, u);
-			fp_sqr(w, w);
-			fp_mul(w, w, z1);
-			fp_inv(w, w);
-			fp_dbl(p->y, ctx->ep_map_c[6]);
-			fp_mul(p->y, p->y, t);
-			fp_mul(p->y, p->y, w);
-			fp_set_dig(p->z, 1);
-			p->coord = BASIC;
-			ep_mul_cof(p, p);
+			if (fp_is_sqr(a)) {
+				/* Compute t^2 = 3c*sqrt(a)*(2c^3*x^6 - 3*c^2*x^4 - 3*c*x^2 + 2).*/
+				/* Compute w = 3*c. */
+				fp_set_dig(c, -fp_prime_get_qnr());
+				fp_neg(c, c);
+				fp_dbl(w, c);
+				fp_add(w, w, c);
+
+				/* Compute x^2, x^4 and x^6 in sequence. */
+				fp_sqr(z1, u);
+				fp_sqr(y1, z1);
+				fp_mul(t, z1, y1);
+
+				fp_dbl(t, t);
+				fp_mul(t, t, c);
+				fp_mul(t, t, c);
+				fp_mul(t, t, c);
+
+				fp_mul(v, y1, c);
+				fp_mul(v, v, w);
+				fp_sub(t, t, v);
+
+				/* v = -3*c*x^2. */
+				fp_mul(v, w, z1);
+				fp_neg(v, v);
+				fp_add(t, t, v);
+				fp_add_dig(t, t, 2);
+
+				/* Assume a = 1 for simplicitly. */
+				fp_mul(t, t, w);
+				fp_mul(t, t, ctx->ep_map_c[6]);
+				dig_t c1 = fp_is_sqr(t);
+				/* If t is not square, compute u = 1/(uc), t = sqrt(t/c)/(c*u^3)*/
+				fp_inv(v, c);
+				fp_inv(x1, u);
+				fp_mul(y1, t, v);
+				/* If t is a square, extract its square root. */
+				dv_copy_cond(t, y1, RLC_FP_DIGS, !c1);
+				fp_srt(t, t);
+				fp_mul(y1, t, v);
+				fp_sqr(y, x1);
+				fp_mul(y, y, x1);
+				fp_mul(y1, y1, y);
+				fp_mul(x1, x1, v);
+				dv_copy_cond(u, x1, RLC_FP_DIGS, !c1);
+				dv_copy_cond(t, y1, RLC_FP_DIGS, !c1);
+
+				/* Compute x = sqrt(a)*(c*x^2 - 2)/(-3*c*x^2). */
+				fp_sqr(z1, u);
+				fp_mul(v, w, z1);
+				fp_neg(v, v);
+				fp_inv(v, v);
+				fp_mul(p->x, z1, c);
+				fp_sub_dig(p->x, p->x, 2);
+				fp_mul(p->x, p->x, v);
+				fp_mul(p->x, p->x, ctx->ep_map_c[6]);
+
+				/* Compute y = y*2*sqrt(a)/(3^2*c^2*x^3). */
+				fp_mul(z1, z1, u);
+				fp_sqr(w, w);
+				fp_mul(w, w, z1);
+				fp_inv(w, w);
+				fp_dbl(p->y, ctx->ep_map_c[6]);
+				fp_mul(p->y, p->y, t);
+				fp_mul(p->y, p->y, w);
+				fp_set_dig(p->z, 1);
+				p->coord = BASIC;
+			} else {
+				/* Compute c = 3*a^2, t^2 = 6a(9u^5 − 14au^3 + 3cu).*/
+				fp_neg(a, a);
+				fp_sqr(c, a);
+				fp_dbl(t, c);
+				fp_add(c, c, t);
+				fp_dbl(t, c);
+				fp_add(t, t, c);
+				fp_mul(t, t, u);
+
+				fp_sqr(v, u);
+				fp_mul(w, v, u);
+				fp_mul(x1, w, a);
+				fp_mul_dig(x1, x1, 14);
+				fp_sub(t, t, x1);
+
+				fp_mul(w, w, v);
+				fp_dbl(x1, w);
+				fp_add(w, w, x1);
+				fp_dbl(x1, w);
+				fp_add(w, w, x1);
+				fp_add(t, t, w);
+				fp_mul(t, t, a);
+				fp_dbl(t, t);
+				fp_dbl(x1, t);
+				fp_add(t, t, x1);
+				dig_t c1 = fp_is_sqr(t);
+				/* If t is not square, compute u = a/u, t = a*sqrt(a*t)/u^3*/
+				fp_inv(x1, u);
+				fp_mul(y1, t, a);
+				/* If t is a square, extract its square root. */
+				dv_copy_cond(t, y1, RLC_FP_DIGS, !c1);
+				fp_srt(t, t);
+				fp_mul(y1, t, a);
+				fp_sqr(y, x1);
+				fp_mul(y, y, x1);
+				fp_mul(y1, y1, y);
+				fp_mul(x1, x1, a);
+				dv_copy_cond(u, x1, RLC_FP_DIGS, !c1);
+				dv_copy_cond(t, y1, RLC_FP_DIGS, !c1);
+
+				/* Compute x = 2^4*i*3*a^2*u / (3*(3*u^2 - a))^2. */
+				fp_copy(y, ctx->ep_map_c[6]);
+				fp_mul(c, c, u);
+				fp_mul(x1, c, y);
+				fp_dbl(x1, x1);
+				fp_dbl(x1, x1);
+				fp_dbl(x1, x1);
+				fp_dbl(p->x, x1);
+				fp_sqr(v, u);
+				fp_dbl(z1, v);
+				fp_add(z1, z1, v);
+				fp_sub(z1, z1, a);
+				fp_dbl(p->z, z1);
+				fp_add(p->z, p->z, z1);
+
+				/* Compute y = 3*2*(i-1)*a*(3^2*u^2 + a)*t / (3*(3*u^2 - a))^3. */
+				fp_sub_dig(y, y, 1);
+				fp_mul(y1, y, a);
+				fp_dbl(y1, y1);
+				fp_dbl(p->y, y1);
+				fp_add(p->y, p->y, y1);
+				fp_mul(p->y, p->y, t);
+				fp_dbl(y1, v);
+				fp_add(y1, y1, v);
+				fp_dbl(v, y1);
+				fp_add(y1, y1, v);
+				fp_add(y1, y1, a);
+				fp_mul(p->y, p->y, y1);
+
+				/* Multiply by cofactor. */
+				p->coord = JACOB;
+				ep_norm(p, p);
+			}
 		} else if ((ep_curve_opt_b() == RLC_ZERO) && (ctx->mod8 != 1)) {
 			/* This is the approach due to Koshelev introduced in
 			 * https://eprint.iacr.org/2021/1604.pdf */
@@ -435,7 +512,6 @@ void ep_map_swift(ep_t p, const uint8_t *msg, size_t len) {
 			fp_mul(c, ctx->ep_map_c[5], x1);
 			fp_copy(p->x, v);
 			fp_sqr(p->y, y);
-			p->coord = BASIC;
 			/* We use zp as temporary, but there is no problem with \psi. */
 			int index = 0;
 			fp_copy(y1, u);
@@ -486,9 +562,8 @@ void ep_map_swift(ep_t p, const uint8_t *msg, size_t len) {
 			fp_mul(z1, z1, y);
 			dv_copy_cond(p->x, t, RLC_FP_DIGS, c3);
 			dv_copy_cond(p->y, z1, RLC_FP_DIGS, c3);
-			/* Multiply by cofactor. */
+			p->coord = BASIC;
 			fp_set_dig(p->z, 1);
-			ep_mul_cof(p, p);
 		} else {
 			/* This is the SwiftEC case per se. */
 			if (ep_curve_opt_a() != RLC_ZERO) {
@@ -555,16 +630,18 @@ void ep_map_swift(ep_t p, const uint8_t *msg, size_t len) {
 					fp_copy(p->y, t);
 					fp_set_dig(p->z, 1);
 					p->coord = BASIC;
-					ep_mul_cof(p, p);
 				}
 			}
 		}
+		/* Multiply by cofactor. */
+		ep_mul_cof(p, p);
 	}
 	RLC_CATCH_ANY {
 		RLC_THROW(ERR_CAUGHT);
 	}
 	RLC_FINALLY {
 		bn_free(k);
+		fp_free(a);
 		fp_free(c);
 		fp_free(t);
 		fp_free(u);
diff --git a/src/ep/relic_ep_mul_cof.c b/src/ep/relic_ep_mul_cof.c
index 1da814637..a9a9c1d9f 100644
--- a/src/ep/relic_ep_mul_cof.c
+++ b/src/ep/relic_ep_mul_cof.c
@@ -71,7 +71,7 @@ void ep_mul_cof(ep_t r, const ep_t p) {
 				}
 				break;
 			case EP_N16:
-				/* if (u % 2) == 0, compute = (u * (u**3+1)//2)*P
+				/* if (u % 2) == 0, compute = (u * (u**3+1)/2)*P
     			 * else Compute (u * (u**3+1))*P */
 				fp_prime_get_par(k);
 				bn_sqr(l, k);
@@ -87,6 +87,24 @@ void ep_mul_cof(ep_t r, const ep_t p) {
 					ep_mul_basic(r, p, k);
 				}
 				break;
+			case EP_FM16:
+				/* Compute (u/2)*P + [u^3]*phi([u/2]P) */
+				fp_prime_get_par(k);
+				bn_sqr(l, k);
+				bn_mul(l, l, k);
+				if (bn_is_even(k)) {
+					bn_hlv(k, k);
+				}
+				if (bn_bits(k) < RLC_DIG) {
+					ep_mul_dig(r, p, k->dp[0]);
+				} else {
+					ep_mul_basic(r, p, k);
+				}
+				ep_psi(v, r);
+				ep_mul_basic(v, v, l);
+				ep_add(r, r, v);
+				ep_norm(r, r);
+				break;
 			case EP_K16:
 				/* Compute 1250*(P + [(u+1)/2]phi(P)) */
 				fp_prime_get_par(k);
diff --git a/src/ep/relic_ep_param.c b/src/ep/relic_ep_param.c
index eec143808..b5ef9dcdb 100644
--- a/src/ep/relic_ep_param.c
+++ b/src/ep/relic_ep_param.c
@@ -659,6 +659,7 @@
 /** @} */
 #endif
 
+#if defined(EP_ENDOM) && FP_PRIME == 765
 /**
  * Parameters for a 765-bit pairing-friendly prime curve.
  */
@@ -671,6 +672,20 @@
 #define N16_P765_H		"26597655A836F22BC9CF003FFFFFFFFFFFFFA30FAB330D5A7F0000000000000000000000384F01000000000000000000"
 /** @} */
 
+/**
+ * Parameters for a 765-bit pairing-friendly prime curve.
+ */
+/** @{ */
+#define FM16_P765_A		"1"
+#define FM16_P765_B		"0"
+#define FM16_P765_X		"82CD0FCADA08D0804B28E07F7DB2A82B87CDB5A711F9D1242A5F5C7AC4EE055EE696D96973E038E253D643226B111E7842EE367B341E4AF3DD871AD8EFA11EE08BC5507CBEC91F4871A040549486E5C192EB0BD5B3396D0502042D1A4983AF0"
+#define FM16_P765_Y		"12DD98F37AD410B2F5F5D43B8F6CCAE4ECB2F13308530AE334B36E0C27792B99C9ABD810DF8F824A1F3503445B65C16B76B769E9258BFA1782C2D8FF2C1FDE165D8F6B1B5F0E6DCD934F1D93C646524DE73F46A6CA43CB3BDA9EAC7F7AFBD16A"
+#define FM16_P765_R		"97806EA97185E5A71C9C881E055495FC0100000000000000000000000000000000000000000000000000000000000001"
+#define FM16_P765_H		"25E01BAA5C617969C72722078155257F004000000000000000000000000000000000000038223FF01000000000000000"
+/** @} */
+#endif
+
+#if defined(EP_ENDOM) && FP_PRIME == 766
 /**
  * Parameters for a 766-bit pairing-friendly prime curve.
  */
@@ -682,6 +697,7 @@
 #define K16_P766_R		"1B6C1BFC8E56CCE359E1D8A9B94553D096A506CE2ECF4A33C5D526AC5F3B61CB0A6D76FCD8487EDEE0B0F9BA2DFA29D5AB0B164B8792C233ED1E6EB350BA9F4D37112A98DE816BEB1EA8DDB1"
 #define K16_P766_H		"2327FFFFFFFFE8905E7E6E0003E7E080C57EE9EF4"
 /** @} */
+#endif
 
 #if defined(EP_SUPER) && FP_PRIME == 1536
 /**
@@ -699,7 +715,7 @@
 
 #if defined(EP_PLAIN) && FP_PRIME == 3072
 /**
- * Parameters for a 1536-bit supersingular elliptic curve.
+ * Parameters for a 3072-bit supersingular elliptic curve.
  */
 /** @{ */
 #define K1_P3072_A		"-4"
@@ -1138,6 +1154,11 @@ void ep_param_set(int param) {
 				endom = 1;
 				pairf = EP_N16;
 				break;
+			case FM16_P765:
+				ASSIGN(FM16_P765, FM16_765);
+				endom = 1;
+				pairf = EP_FM16;
+				break;
 #endif
 #if defined(EP_ENDOM) && FP_PRIME == 766
 			case K16_P766:
@@ -1206,6 +1227,8 @@ void ep_param_set(int param) {
 					bn_sub_dig(lamb, lamb, 1);
 					break;
 				case EP_N16:
+				case EP_FM16:
+					/* lambda = -z^4. */
 					bn_sqr(lamb, lamb);
 					bn_sqr(lamb, lamb);
 					bn_neg(lamb, lamb);
@@ -1415,6 +1438,7 @@ int ep_param_set_any_endom(void) {
 #endif
 #elif FP_PRIME == 765
 	ep_param_set(N16_P765);
+	//ep_param_set(FM16_P765);
 #elif FP_PRIME == 766
 	ep_param_set(K16_P766);
 #else
@@ -1531,7 +1555,8 @@ int ep_param_set_any_pairf(void) {
 	extension = 3;
 #endif
 #elif FP_PRIME == 765
-	ep_param_set(N16_P765);
+	//ep_param_set(N16_P765);
+	ep_param_set(FM16_P765);
 	type = RLC_EP_MTYPE;
 	extension = 4;
 #elif FP_PRIME == 766
@@ -1690,6 +1715,9 @@ void ep_param_print(void) {
 		case N16_P765:
 			util_banner("Curve N16-P765:", 0);
 			break;
+		case FM16_P765:
+			util_banner("Curve FM16-P765:", 0);
+			break;
 		case K16_P766:
 			util_banner("Curve K16-P766:", 0);
 			break;
@@ -1752,6 +1780,8 @@ int ep_param_level(void) {
 		case B12_P455:
 			return 140;
 		case NIST_P384:
+		case N16_P765:
+		case FM16_P765:
 		case K16_P766:
 		case K18_P638:
 		case B24_P509:
@@ -1778,6 +1808,7 @@ int ep_param_embed(void) {
 		case EP_B12:
 			return 12;
 		case EP_N16:
+		case EP_FM16:
 		case EP_K16:
 			return 16;
 		case EP_K18:
diff --git a/src/epx/relic_ep4_curve.c b/src/epx/relic_ep4_curve.c
index 5b228c439..a7f7cf6cc 100644
--- a/src/epx/relic_ep4_curve.c
+++ b/src/epx/relic_ep4_curve.c
@@ -153,6 +153,28 @@
 /** @} */
 #endif
 
+#if defined(EP_ENDOM) && FP_PRIME == 765
+/** @{ */
+#define FM16_P765_A0	"0"
+#define FM16_P765_A1	"0"
+#define FM16_P765_A2	"1"
+#define FM16_P765_A3	"0"
+#define FM16_P765_B0	"0"
+#define FM16_P765_B1	"0"
+#define FM16_P765_B2	"0"
+#define FM16_P765_B3	"0"
+#define FM16_P765_X0	"EA52AB67B443C778FDAD4A69F92C71AF91DFA14F1E9311F9931580427EA8AE3FB8FDBEA439156418AECEAF3779118AA3063FD5AE91CCAA71CB2A80142809E4329A79F14BD9EE4B544AE5DD4E165B7859E9A9E53FE20CC118A377FAC7EE03383"
+#define FM16_P765_X1	"FF84E6B156239F35EE8CE26A50005A921F08ABE1188932030B8D967D4D4BB3141AA2740E0D47BE45DCC0A999FFAB378143100C8FB5EEA72123D9B6F4369AD068936ECE2819BF1CC377A76A26AD024FE8A5BA7D814E2C1B10194317FBE167A0"
+#define FM16_P765_X2	"9D76B11D4F9C55C1FE9F26B63238C3FAF0910E728D3ADC4D7257F79125D273B85124D6BD0EAA25FBB4B7524F3D17F3879242A2CBECA4F46C90028E439CB95C506EEC3077612132F9168A865B0D8DF11B1F14BBF3DC8BF97F136A724B635F6AC"
+#define FM16_P765_X3	"BC47D8993B5921DD0C1159BD25EF396A2FFCBDB7AE2F3CD9BA1C3F08BABEB4080F308F393B38EC9DAA97FFB076BE5A11CE4EAB250A803CC14E5E8C4240B3396E84461D4D18E97A7F958B1B62FC3F9B660993C8D26BAFB491BD654DBFEB9D0AD"
+#define FM16_P765_Y0	"11AD3FC23DBE9ECF9683D753985679D485E3ADFF0EB31727F30C51034ED9CA51D1BBE4680B4A29802FA21867D4AFFA3183BD5A9F4537C9B54C08270D3FEC78A16F1B204344167D955B9B9D21620FEECD08A164B5E83A06F8364B9FE668DD6D16"
+#define FM16_P765_Y1	"8679BA060D85D894CFF27F6535CC817058AC5F7B33D7A9EB131D5DD4C05B6E8205A1ED28928D930F5B4B020D5317873E257AFAAD153FE8CADDC751DA135F3CA0333F4A4CB43445A7B997E3352BB9F199B1FFFA9435993A166F31A3D81FD1E3C"
+#define FM16_P765_Y2	"8E6C8D9A79BD0D7187D0D2C1285FB705138053984D865E1E5FDAB315D100AFEE4C403F592849ECD2146E8FED537E92BD606F3B856F4530D3F61E6F064F34580C19E6E9798C45C2FEC328C565649DDDC4E9E7EF6DC323CC3CC401AA27F817264"
+#define FM16_P765_Y3	"E90F6A9653575AA7A158E8414F47C7CB777889FA45C550A76455C999C307B2F92DEFF4F129625FA4E3C2D0403EEB9D06F238862875E97B5DABDC51EAF5095DF38089EA8D165F767CF8E11C9DD34CE98B44A96BA1125ED81ABB137A822894FF6"
+#define FM16_P765_R		"97806EA97185E5A71C9C881E055495FC0100000000000000000000000000000000000000000000000000000000000001"
+#define FM16_P765_H		"6822EA3AA20627C198952D3A862F7922329CB37E53F890BF58E136C8578F4018ADB15BC58B2906276B272ED80976F28344FE7D36BE178B8A5589BE16EB483451CED64B6A44730FAA989051E1B87AC870364BCD14BE0893D560147CFFA3D313EE4A40AD71F8DA71980AD969768E3208C77240C6765CC9FEC7E67B23C3DEE0FC9B3FA681E45258CF2399008513635184F09CEB78E2D98FF766926C33D0CD33C59FDB1177ADFDB42CE45C73F65053289168A22A122C0DF279C7A10CC55B3CA63256BC71D953DC548E52B03790D30877280233B3D5FEC4E128CAAFDD23395D77E82F519F3D794D6BBD50A96BD483C9BA59537F39A1F892FAFD7C87426883477201546AA2D3774A597A2EB03DDB76B17C8AC8DBD425932DBC7FBAD3E701C7E8BC6488EC12F0A3756D5D3A3EFB06BD20ABDE389CAC64A7E7CCB9CB1EB996E6F75F1367D92B703BEFE41C003BF0200000002"
+/** @} */
+#endif
 
 #if defined(EP_ENDOM) && FP_PRIME == 766
 /** @{ */
@@ -403,6 +425,9 @@ void ep4_curve_set_twist(int type) {
 			case N16_P765:
 				ASSIGN(N16_P765);
 				break;
+			case FM16_P765:
+				ASSIGN(FM16_P765);
+				break;
 #elif FP_PRIME == 766
 			case K16_P766:
 				ASSIGN(K16_P766);
diff --git a/src/epx/relic_ep4_mul.c b/src/epx/relic_ep4_mul.c
index b1653d038..1181bbbfd 100644
--- a/src/epx/relic_ep4_mul.c
+++ b/src/epx/relic_ep4_mul.c
@@ -66,6 +66,8 @@ static void ep4_psi(ep4_t r, const ep4_t p) {
 				ep4_frb(r, p, 5);
 				ep4_neg(r, r);
 				break;
+			case EP_FM16:
+				/* u = p mod r */
 			default:
 				ep4_frb(r, p, 1);
 				break;
diff --git a/src/epx/relic_ep4_mul_cof.c b/src/epx/relic_ep4_mul_cof.c
index ed3285e24..228232624 100644
--- a/src/epx/relic_ep4_mul_cof.c
+++ b/src/epx/relic_ep4_mul_cof.c
@@ -124,7 +124,6 @@ static void ep4_mul_cof_k16(ep4_t r, const ep4_t p) {
 		ep4_free(t4);
 		ep4_free(t5);
 		bn_free(x);
-
 	}
 }
 
@@ -198,7 +197,89 @@ static void ep4_mul_cof_n16(ep4_t r, const ep4_t p) {
 		ep4_free(t4);
 		ep4_free(t5);
 		bn_free(x);
+	}
+}
+
+/**
+ * Multiplies a point by the cofactor in a FM16 curve.
+ *
+ * @param[out] r			- the result.
+ * @param[in] p				- the point to multiply.
+ */
+static void ep4_mul_cof_fm16(ep4_t r, const ep4_t p) {
+	bn_t x;
+	ep4_t t0, t1, t2, t3;
+
+	ep4_null(t0);
+	ep4_null(t1);
+	ep4_null(t2);
+	ep4_null(t3);
+	bn_null(x);
+
+	RLC_TRY {
+		bn_new(x);
+		ep4_new(t0);
+		ep4_new(t1);
+		ep4_new(t2);
+		ep4_new(t3);
+
+		fp_prime_get_par(x);
+		if (bn_is_even(x)) {
+			bn_hlv(x, x);
+		} else {
+			RLC_THROW(ERR_NO_VALID);
+		}
+
+		/* [1 + u^7/2, u^6/2, u^5/2, u^4/2, u^3(2-u)/2, u^2(2-u)/2, u(2-u)/2, (2-u)/2] */
 
+		/* Compute [1 - u/2]\psi^7(P). */
+		ep4_mul_basic(t1, p, x);
+		ep4_sub(t0, p, t1);
+		ep4_frb(t0, t0, 7);
+		bn_dbl(x, x);
+		/* Compute [u - u^2/2]\psi^6(P). */
+		ep4_mul_basic(t2, t1, x);
+		ep4_dbl(t1, t1);
+		ep4_sub(t3, t1, t2);
+		ep4_frb(t3, t3, 6);
+		ep4_add(t0, t0, t3);
+		/* Compute [u^2 - u^3/2]\psi^5(P). */
+		ep4_mul_basic(t1, t2, x);
+		ep4_dbl(t2, t2);
+		ep4_sub(t3, t2, t1);
+		ep4_frb(t3, t3, 5);
+		ep4_add(t0, t0, t3);
+		/* Compute [u^3 - u^4/2]\psi^4(P). */
+		ep4_mul_basic(t2, t1, x);
+		ep4_dbl(t1, t1);
+		ep4_sub(t3, t1, t2);
+		ep4_frb(t3, t3, 4);
+		ep4_add(t0, t0, t3);
+		/* Compute [u^4/2]\psi^3(P). */
+		ep4_frb(t3, t2, 3);
+		ep4_add(t0, t0, t3);
+		/* Compute [u^5/2]\psi^2(P). */
+		ep4_mul_basic(t2, t2, x);
+		ep4_frb(t3, t2, 2);
+		ep4_add(t0, t0, t3);
+		/* Compute [u^6/2]\psi(P). */
+		ep4_mul_basic(t2, t2, x);
+		ep4_frb(t3, t2, 1);
+		ep4_add(t0, t0, t3);
+		/* Compute [1 + u^7/2]P. */
+		ep4_mul_basic(t2, t2, x);
+		ep4_add(t0, t0, t2);
+		ep4_add(t0, t0, p);
+
+		ep4_norm(r, t0);
+	} RLC_CATCH_ANY {
+		RLC_THROW(ERR_CAUGHT);
+	} RLC_FINALLY {
+		ep4_free(t0);
+		ep4_free(t1);
+		ep4_free(t2);
+		ep4_free(t3);
+		bn_free(x);
 	}
 }
 
@@ -280,6 +361,9 @@ void ep4_mul_cof(ep4_t r, const ep4_t p) {
 			case EP_N16:
 				ep4_mul_cof_n16(r, p);
 				break;
+			case EP_FM16:
+				ep4_mul_cof_fm16(r, p);
+				break;
 			case EP_B24:
 				ep4_mul_cof_b24(r, p);
 				break;
diff --git a/src/fp/relic_fp_param.c b/src/fp/relic_fp_param.c
index 4e9fd7a61..3f5108cda 100644
--- a/src/fp/relic_fp_param.c
+++ b/src/fp/relic_fp_param.c
@@ -579,6 +579,17 @@ void fp_param_set(int param) {
 				bn_neg(t0, t0);
 				fp_prime_set_pairf(t0, EP_N16);
 				break;
+			case FM16_765:
+				/* u = 2^48-2^44-2^38+2^31 */
+				bn_set_2b(t0, 48);
+				bn_set_2b(t1, 44);
+				bn_sub(t0, t0, t1);
+				bn_set_2b(t1, 38);
+				bn_sub(t0, t0, t1);
+				bn_set_2b(t1, 31);
+				bn_add(t0, t0, t1);
+				fp_prime_set_pairf(t0, EP_FM16);
+				break;
 #elif FP_PRIME == 766
 			case K16_766:
 				/* u = 2^78-2^76-2^28+2^14+2^7+1 */
@@ -787,6 +798,7 @@ int fp_param_set_any_tower(void) {
 #endif
 #elif FP_PRIME == 765
 	fp_param_set(N16_765);
+	//fp_param_set(FM16_765);
 #elif FP_PRIME == 766
 	fp_param_set(K16_766);
 #elif FP_PRIME == 1536
diff --git a/src/fp/relic_fp_prime.c b/src/fp/relic_fp_prime.c
index c18319e98..0d37b8e92 100644
--- a/src/fp/relic_fp_prime.c
+++ b/src/fp/relic_fp_prime.c
@@ -437,6 +437,25 @@ void fp_prime_set_pairf(const bn_t x, int pairf) {
 				bn_div_dig(p, p, 4);
 				fp_prime_set_dense(p);
 				break;
+			case EP_FM16:
+				/* p = (x^16 + x^10 + 5*x^8 + x^2 + 4*x + 4)/4 */
+				bn_sqr(t1, t0);
+				bn_mul(p, t1, t0);
+				bn_sqr(p, p);
+				bn_add_dig(p, p, 1);
+				bn_mul(p, p, t1);
+				bn_add_dig(p, p, 5);
+				bn_mul(p, p, t1);
+				bn_mul(p, p, t1);
+				bn_mul(p, p, t1);
+				bn_add_dig(p, p, 1);
+				bn_mul(p, p, t0);
+				bn_add_dig(p, p, 4);
+				bn_mul(p, p, t0);
+				bn_add_dig(p, p, 4);
+				bn_div_dig(p, p, 4);
+				fp_prime_set_dense(p);
+				break;
 			case EP_K16:
 				/* p = (u^10 + 2*u^9 + 5*u^8 + 48*u^6 + 152*u^5 + 240*u^4 +
 						625*u^2 + 2398*u + 3125) div 980 */
diff --git a/src/fpx/relic_fp4_mul.c b/src/fpx/relic_fp4_mul.c
index 16904cc59..1f9eaefbe 100644
--- a/src/fpx/relic_fp4_mul.c
+++ b/src/fpx/relic_fp4_mul.c
@@ -174,6 +174,7 @@ void fp4_mul_frb(fp4_t c, const fp4_t a, int i, int j) {
 			for (int k = 0; k < j; k++) {
 	        	fp2_mul(c[0], c[0], t);
 				fp2_mul(c[1], c[1], t);
+				fp4_mul_art(c, c);
 				/* If constant in base field, then second component is zero. */
 				if (core_get()->frb4 == 1) {
 					fp4_mul_art(c, c);
diff --git a/src/fpx/relic_fpx_cyc.c b/src/fpx/relic_fpx_cyc.c
index cc23f9415..f6afe9ae7 100644
--- a/src/fpx/relic_fpx_cyc.c
+++ b/src/fpx/relic_fpx_cyc.c
@@ -601,7 +601,7 @@ void fp12_exp_cyc(fp12_t c, const fp12_t a, const bn_t b) {
 		return fp12_set_dig(c, 1);
 	}
 
-	if (1) {//(bn_bits(b) > RLC_DIG) && ((bn_ham(b) << 3) > bn_bits(b))) {
+	if ((bn_bits(b) > RLC_DIG) && ((bn_ham(b) << 3) > bn_bits(b))) {
 		fp12_t r, s, t[1 << (RLC_WIDTH - 2)];
 		int8_t naf[RLC_FP_BITS + 1], *k, w = RLC_WIDTH;
 
@@ -1057,6 +1057,9 @@ static void fp16_gls(fp16_t c, const fp16_t a) {
 				fp16_frb(c, a, 5);
 				fp16_inv_cyc(c, c);
 				break;
+			case EP_FM16:
+				fp16_frb(c, a, 1);
+				break;
 		}
 	}
 	RLC_CATCH_ANY {
diff --git a/src/fpx/relic_fpx_field.c b/src/fpx/relic_fpx_field.c
index cb05e31ce..7203b179b 100644
--- a/src/fpx/relic_fpx_field.c
+++ b/src/fpx/relic_fpx_field.c
@@ -82,7 +82,7 @@ void fp2_field_init(void) {
 		fp_zero(t0[0]);
 		fp_set_dig(t0[1], 1);
 		/* If it does not work, attempt (u + 2), otherwise double. */
-		/* We cannot used QR test here because Frobenius constants below. */
+		/* We cannot used QR test here due to Frobenius constants below. */
 		if (fp2_srt(t1, t0)) {
 			ctx->qnr2 = 2;
 			fp_set_dig(t0[0], ctx->qnr2);
diff --git a/src/low/x64-asm-12l/macro.s b/src/low/x64-asm-12l/macro.s
index a33c14129..dfd99f58b 100644
--- a/src/low/x64-asm-12l/macro.s
+++ b/src/low/x64-asm-12l/macro.s
@@ -47,6 +47,7 @@
 #define P11 0x3C410B7E6EC19106
 #define U0	0xC18CA908C52344BB
 #elif FP_PRIME == 765
+/* AFG16-765 */
 #define P0	0x0000000000000001
 #define P1	0x00000000384F0100
 #define P2	0x7D00000000000000
@@ -57,9 +58,24 @@
 #define P7	0xEE9C1E7F21BD9E92
 #define P8	0x249F514A2A836FBF
 #define P9	0x8866F5670199231B
-#define P10 0xB2847B1232833CC3
-#define P11 0x16FAB993B0C96754
+#define P10	0xB2847B1232833CC3
+#define P11	0x16FAB993B0C96754
 #define U0	0xFFFFFFFFFFFFFFFF
+/* FM16-765
+#define P0	0x1000EFC080000001
+#define P1	0x0000000038223FF0
+#define P2	0x0000000000000000
+#define P3	0x0140000000000000
+#define P4	0x93D3AA2586A9BB7B
+#define P5	0x2C9088558A226AF0
+#define P6	0x7071D6BA0697D5A1
+#define P7	0xFE00400021385D1A
+#define P8	0x1629227BB6527E4E
+#define P9	0xD4A66E04AA631EEA
+#define P10 0xBC5664C6F237BCB4
+#define P11 0x166A30BEAF4CE221
+#define U0	0xD000EFC07FFFFFFF
+*/
 #endif
 
 #if defined(__APPLE__)
@@ -193,13 +209,31 @@
 
 .macro _RDCN0 i, j, k, R0, R1, R2 A, P
 	movq	8*\i(\A), %rax
+#if U0 == 0xFFFFFFFFFFFFFF
+	.if \j != 2
+		mulq	8*\j(\P)
+		addq	%rax, \R0
+		adcq	%rdx, \R1
+		adcq	$0, \R2
+	.endif
+#else
 	mulq	8*\j(\P)
 	addq	%rax, \R0
 	adcq	%rdx, \R1
 	adcq	$0, \R2
+#endif
 	.if \j > 1
 		_RDCN0 "(\i + 1)", "(\j - 1)", \k, \R0, \R1, \R2, \A, \P
 	.else
+#if U0 == 0xFFFFFFFFFFFFFFFF
+		addq	8*\k(\A), \R0
+		adcq	$0, \R1
+		adcq	$0, \R2
+		negq	\R0
+		movq	\R0, 8*\k(\A)
+		adcq	$0, \R1
+		adcq	$0, \R2
+#else
 		addq	8*\k(\A), \R0
 		adcq	$0, \R1
 		adcq	$0, \R2
@@ -210,6 +244,7 @@
 		addq	%rax , \R0
 		adcq	%rdx , \R1
 		adcq	$0   , \R2
+#endif
 		xorq	\R0, \R0
 	.endif
 .endm
@@ -219,11 +254,21 @@
 .endm
 
 .macro _RDCN1 i, j, k, l, R0, R1, R2 A, P
+#if U0 == 0xFFFFFFFFFFFFFF
+	.if \j != 2
+		movq	8*\i(\A), %rax
+		mulq	8*\j(\P)
+		addq	%rax, \R0
+		adcq	%rdx, \R1
+		adcq	$0, \R2
+	.endif
+#else
 	movq	8*\i(\A), %rax
 	mulq	8*\j(\P)
 	addq	%rax, \R0
 	adcq	%rdx, \R1
 	adcq	$0, \R2
+#endif
 	.if \j > \l
 		_RDCN1 "(\i + 1)", "(\j - 1)", \k, \l, \R0, \R1, \R2, \A, \P
 	.else
@@ -242,6 +287,12 @@
 // r8, r9, r10, r11, r12, r13, r14, r15, rbp, rbx, rsp, //rsi, rdi, //rax, rcx, rdx
 .macro FP_RDCN_LOW C, R0, R1, R2, A, P
 	xorq	\R1, \R1
+#if U0 == 0xFFFFFFFFFFFFFFFF
+	movq	0(\A), \R0
+	negq	\R0
+	movq	\R0 , 0(\A)
+	adcq	$0   , \R1
+#else
 	movq	$U0, %rcx
 
 	movq	0(\A), \R0
@@ -251,9 +302,9 @@
 	mulq	0(\P)
 	addq	%rax , \R0
 	adcq	%rdx , \R1
+#endif
 	xorq    \R2  , \R2
 	xorq    \R0  , \R0
-
 	RDCN0	0, 1, \R1, \R2, \R0, \A, \P
 	RDCN0	0, 2, \R2, \R0, \R1, \A, \P
 	RDCN0	0, 3, \R0, \R1, \R2, \A, \P
diff --git a/src/pc/relic_pc_util.c b/src/pc/relic_pc_util.c
index fddd35d44..515d36d3b 100644
--- a/src/pc/relic_pc_util.c
+++ b/src/pc/relic_pc_util.c
@@ -173,6 +173,16 @@ int g1_is_valid(const g1_t a) {
 					g1_neg(u, u);
 					r = g1_on_curve(a) && (g1_cmp(u, a) == RLC_EQ);
 					break;
+				case EP_FM16:
+					/* Check that P == (u**4)*\phi(P). */
+					fp_prime_get_par(n);
+					g1_mul_any(u, a, n);
+					g1_mul_any(u, u, n);
+					g1_mul_any(u, u, n);
+					g1_mul_any(u, u, n);
+					ep_psi(u, u);
+					r = g1_on_curve(a) && (g1_cmp(u, a) == RLC_EQ);
+					break;
 				case EP_K18:
 					/* Check that [a_0]P + [a_1]\psi(P)) == O, for
 					 * a_0 = 19a_1 + 1, a_1 = (x/7)^3 */
@@ -384,6 +394,13 @@ int g2_is_valid(const g2_t a) {
 				g2_frb(t, t, 4);
 				r = g2_on_curve(a) && (g2_cmp(w, t) == RLC_EQ);
 				break;
+			case EP_FM16:
+				/* Check that u*Q == psi(Q). */
+				fp_prime_get_par(n);
+				g2_mul_any(u, a, n);
+				g2_frb(v, a, 1);
+				r = g2_on_curve(a) && (g2_cmp(u, v) == RLC_EQ);
+				break;
 			case EP_K18:
 				/* Check that P + u*psi2P + 2*psi3P == \mathcal{O}. */
 				fp_prime_get_par(n);
@@ -571,6 +588,14 @@ int gt_is_valid(const gt_t a) {
 				r = (gt_cmp(w, t) == RLC_EQ);
 				r &= fp16_test_cyc((void *)a);
 				break;
+			case EP_FM16:
+				/* Check that u*Q == psi(Q). */
+				fp_prime_get_par(n);
+				gt_exp(u, a, n);
+				gt_frb(v, a, 1);
+				r = (gt_cmp(u, v) == RLC_EQ);
+				r &= fp16_test_cyc((void *)a);
+				break;
 			case EP_K18:
 				/* Check that P + u*psi2P + 2*psi3P == \mathcal{O}. */
 				gt_frb(u, a, 2);
diff --git a/src/pp/relic_pp_exp_k16.c b/src/pp/relic_pp_exp_k16.c
index 48c71360c..6f64ef199 100644
--- a/src/pp/relic_pp_exp_k16.c
+++ b/src/pp/relic_pp_exp_k16.c
@@ -334,6 +334,100 @@ static void pp_exp_new(fp16_t c, fp16_t a) {
 	}
 }
 
+/**
+ * Computes the final exponentiation of a pairing defined over a FM16 curve.
+ *
+ * @param[out] c			- the result.
+ * @param[in] a				- the extension field element to exponentiate.
+ */
+static void pp_exp_fm16(fp16_t c, fp16_t a) {
+	fp16_t t0, t1, t2, t3, t4, t5, t6, t7;
+	bn_t x, x_;
+
+	bn_null(x);
+	bn_null(x_);
+	fp16_null(t0);
+	fp16_null(t1);
+	fp16_null(t2);
+	fp16_null(t3);
+	fp16_null(t4);
+	fp16_null(t5);
+	fp16_null(t6);
+	fp16_null(t7);
+
+	RLC_TRY {
+		bn_new(x);
+		bn_new(x_);
+		fp16_new(t0);
+		fp16_new(t1);
+		fp16_new(t2);
+		fp16_new(t3);
+		fp16_new(t4);
+		fp16_new(t5);
+		fp16_new(t6);
+		fp16_new(t7);
+
+		fp_prime_get_par(x);
+
+		/* First, compute m = f^(p^8 - 1). */
+		fp16_conv_cyc(c, a);
+
+		/* Now compute m^((p^8 + 1) / r). */
+		bn_hlv(x_, x);
+
+		fp16_exp_cyc(t1, c, x_);
+		fp16_exp_cyc(t1, t1, x_);
+		fp16_exp_cyc(t2, t1, x);
+		fp16_exp_cyc(t3, t2, x);
+		fp16_exp_cyc(t4, t3, x);
+		fp16_exp_cyc(t5, t4, x);
+		fp16_exp_cyc(t6, t5, x);
+		fp16_exp_cyc(t7, t6, x);
+
+		fp16_mul(t0, t1, c);
+		fp16_mul(t0, t0, t7);
+		fp16_frb(t7, t0, 7);
+		fp16_exp_cyc(t0, t0, x);
+		fp16_frb(t1, t0, 6);
+		fp16_mul(t7, t7, t1);
+		fp16_exp_cyc(t0, t0, x);
+		fp16_frb(t1, t0, 5);
+		fp16_mul(t7, t7, t1);
+		fp16_exp_cyc(t0, t0, x);
+		fp16_frb(t1, t0, 4);
+		fp16_mul(t7, t7, t1);
+		fp16_exp_cyc(t0, t0, x);
+		fp16_frb(t1, t0, 3);
+		fp16_mul(t7, t7, t1);
+		fp16_exp_cyc(t0, t0, x);
+		fp16_frb(t1, t0, 2);
+		fp16_mul(t7, t7, t1);
+		fp16_exp_cyc(t0, t0, x);
+		fp16_frb(t1, t0, 1);
+		fp16_mul(t7, t7, t1);
+		fp16_exp_cyc(t0, t0, x);
+		fp16_mul(t7, t7, t0);
+		fp16_mul(c, c, t7);
+	}
+	RLC_CATCH_ANY {
+		RLC_THROW(ERR_CAUGHT);
+	}
+	RLC_FINALLY {
+		bn_free(x);
+		bn_free(x_);
+		fp16_free(t0);
+		fp16_free(t1);
+		fp16_free(t2);
+		fp16_free(t3);
+		fp16_free(t4);
+		fp16_free(t5);
+		fp16_free(t6);
+		fp16_free(t7);
+	}
+}
+ 
+
+
 
 /*============================================================================*/
 /* Public definitions                                                         */
@@ -347,5 +441,8 @@ void pp_exp_k16(fp16_t c, fp16_t a) {
 		case EP_N16:
 			pp_exp_new(c, a);
 			break;
+		case EP_FM16:
+			pp_exp_fm16(c, a);
+			break;
 	}
 }
diff --git a/src/pp/relic_pp_map_k16.c b/src/pp/relic_pp_map_k16.c
index cce7b6aeb..837076802 100644
--- a/src/pp/relic_pp_map_k16.c
+++ b/src/pp/relic_pp_map_k16.c
@@ -489,6 +489,7 @@ void pp_map_oatep_k16(fp16_t r, const ep_t p, const ep4_t q) {
 
 		if (!ep_is_infty(_p[0]) && !ep4_is_infty(_q[0])) {
 			switch (ep_curve_is_pairf()) {
+				case EP_FM16:
 				case EP_N16:
 					/* r = f_{|a|,Q}(P). */
 					pp_mil_k16(r, t, _q, _p, 1, a);
@@ -559,6 +560,7 @@ void pp_map_sim_oatep_k16(fp16_t r, const ep_t *p, const ep4_t *q, int m) {
 
 		if (j > 0) {
 			switch (ep_curve_is_pairf()) {
+				case EP_FM16:
 				case EP_N16:
 					/* r = f_{|a|,Q}(P). */
 					pp_mil_k16(r, t, _q, _p, j, a);
diff --git a/test/test_fpx.c b/test/test_fpx.c
index 7a84680a1..26dbe8259 100644
--- a/test/test_fpx.c
+++ b/test/test_fpx.c
@@ -9074,7 +9074,8 @@ int main(void) {
 		}
 	}
 
-	if (fp_prime_get_qnr() && fp_prime_get_cnr() && (ep_param_embed() >= 12)) {
+	if (fp_prime_get_qnr() && fp_prime_get_cnr() &&
+			(ep_param_embed() >= 12) && (ep_param_embed() != 16)) {
 		util_banner("Dodecic extension:", 0);
 		util_banner("Utilities:", 1);