diff --git a/src/pc/relic_pc_exp.c b/src/pc/relic_pc_exp.c
index daa5d517d..a844f5951 100644
--- a/src/pc/relic_pc_exp.c
+++ b/src/pc/relic_pc_exp.c
@@ -107,14 +107,14 @@ static void gt_psi(gt_t c, const gt_t a) {
 #define RLC_GT_TABLE		(1 << (RLC_WIDTH - 2))
 
 /**
- * Exponentiates an element from G_T in constant time.
+ * Exponentiates an element from G_T using the w-NAF method.
  *
  * @param[out] c			- the result.
  * @param[in] a				- the element to exponentiate.
  * @param[in] b				- the exponent.
  * @param[in] f				- the maximum Frobenius power.
  */
-void gt_exp_gls_imp(gt_t c, const gt_t a, const bn_t b, size_t f) {
+void gt_exp_gls_naf(gt_t c, const gt_t a, const bn_t b, size_t f) {
 	int8_t *naf  = RLC_ALLOCA(int8_t, f * (RLC_FP_BITS + 1));
 	int8_t n0, *s = RLC_ALLOCA(int8_t, f);
 	gt_t q, *t = RLC_ALLOCA(gt_t, f * RLC_GT_TABLE);
@@ -255,17 +255,141 @@ void gt_exp_gls_imp(gt_t c, const gt_t a, const bn_t b, size_t f) {
 	}
 }
 
-#if FP_PRIME < 1536
+/**
+ * Exponentiates an element from G_T using the variable-time GLS-SAC method.
+ *
+ * @param[out] c			- the result.
+ * @param[in] a				- the element to exponentiate.
+ * @param[in] b				- the exponent.
+ * @param[in] f				- the maximum Frobenius power.
+ */
+static void gt_exp_gls_sac(gt_t c, const gt_t a, const bn_t b, size_t d,
+		size_t f) {
+	size_t l, s = (1 << (f / d - 1));
+	bn_t n, *_b = RLC_ALLOCA(bn_t, f), u;
+	int8_t col, *e = RLC_ALLOCA(int8_t, d);
+	int8_t *sac = RLC_ALLOCA(int8_t, d * f * RLC_FP_BITS);
+	gt_t *q = RLC_ALLOCA(gt_t, f), *t = RLC_ALLOCA(gt_t, d * s);
+
+	if (sac == NULL || e == NULL || t == NULL || _b == NULL || q == NULL) {
+		RLC_THROW(ERR_NO_MEMORY);
+		return;
+	}
+
+	bn_null(n);
+	bn_null(u);
+
+	RLC_TRY {
+		bn_new(n);
+		bn_new(u);
+		for (int i = 0; i < f; i++) {
+			bn_null(_b[i]);
+			gt_null(q[i]);
+			bn_new(_b[i]);
+			gt_new(q[i]);
+		}
+		for (size_t i = 0; i < d; i++) {
+			for (int j = 0; j < s; j++) {
+				gt_null(t[i * s + j]);
+				gt_new(t[i * s + j]);
+			}
+		}
+
+		gt_get_ord(n);
+		fp_prime_get_par(u);
+		if (ep_curve_is_pairf() == EP_SG18) {
+			/* Compute base -3*u for the recoding below. */
+			bn_dbl(n, u);
+			bn_add(u, u, n);
+			bn_neg(u, u);
+		}
+		bn_mod(_b[0], b, n);
+		bn_rec_frb(_b, f, _b[0], u, n, ep_curve_is_pairf() == EP_BN);
+
+		gt_copy(q[0], a);
+		for (size_t i = 1; i < f; i++) {
+			gt_psi(q[i], q[i - 1]);
+		}
+		for (size_t i = 0; i < f; i++) {
+			gt_copy(c, q[i]);
+			if (bn_sign(_b[i]) == RLC_NEG) {
+				gt_inv(c, q[i]);
+			}
+			bn_abs(_b[i], _b[i]);
+		}
+		for (size_t i = 0; i < d; i++) {
+			e[i] = bn_is_even(_b[i * f / d]);
+			bn_add_dig(_b[i * f / d], _b[i * f / d], e[i]);
+		}
+		
+		for (size_t i = 0; i < d; i++) {
+			gt_copy(t[i * s], q[i * f / d]);
+			for (size_t j = 1; j < s; j++) {
+				l = util_bits_dig(j);
+				gt_mul(t[i * s + j], t[i * s + (j ^ (1 << (l - 1)))], q[l + i * f / d]);
+			}
+			l = RLC_FP_BITS;
+			bn_rec_sac(sac + i * f * RLC_FP_BITS, &l, _b + i * f / d, f / d, bn_bits(n));
+		}
+
+		gt_set_unity(c);
+		for (int j = l - 1; j >= 0; j--) {
+			gt_sqr(c, c);
+			for (size_t i = 0; i < d; i++) {
+				col = 0;
+				for (int k = f / d - 1; k > 0; k--) {
+					col <<= 1;
+					col += sac[i * f * RLC_FP_BITS + k * l + j];
+				}
+			
+				if (sac[i * f * RLC_FP_BITS + j]) {
+					gt_inv(q[1], t[i * s + col]);
+					gt_mul(c, c, q[1]);
+				} else {
+					gt_mul(c, c, t[i * s + col]);
+				}	
+			}
+		}
+
+		for (size_t i = 0; i < d; i++) {
+			if (e[i]) {
+				gt_inv(q[1], q[i * f / d]);
+				gt_mul(c, c, q[1]);
+			}
+		}
+	}
+	RLC_CATCH_ANY {
+		RLC_THROW(ERR_CAUGHT);
+	}
+	RLC_FINALLY {
+		bn_free(n);
+		bn_free(u);
+		for (int i = 0; i < f; i++) {
+			bn_free(_b[i]);
+			gt_free(q[i]);
+		}
+		for (size_t i = 0; i < d; i++) {
+			for (int j = 0; j < s; j++) {
+				gt_free(t[i * d + j]);
+			}
+		}
+		RLC_FREE(e);
+		RLC_FREE(_b);
+		RLC_FREE(q);
+		RLC_FREE(t);
+		RLC_FREE(sac);
+	}
+}
 
 /**
- * Exponentiates an element from G_T in constant time.
+ * Exponentiates an element from G_T in constant time using the GLS-SAC method.
  *
  * @param[out] c			- the result.
  * @param[in] a				- the element to exponentiate.
  * @param[in] b				- the exponent.
  * @param[in] f				- the maximum Frobenius power.
  */
-static void gt_exp_reg_gls(gt_t c, const gt_t a, const bn_t b, size_t d,
+static void gt_exp_reg_sac(gt_t c, const gt_t a, const bn_t b, size_t d,
 		size_t f) {
 	size_t l, s = (1 << (f / d - 1));
 	bn_t n, *_b = RLC_ALLOCA(bn_t, f), u;
@@ -380,8 +504,6 @@ static void gt_exp_reg_gls(gt_t c, const gt_t a, const bn_t b, size_t d,
 	}
 }
 
-#else
-
 /**
  * Exponentiates an element from G_T in constant time.
  *
@@ -530,8 +652,6 @@ void gt_exp_reg_gls(gt_t c, const gt_t a, const bn_t b, size_t f) {
 	}
 }
 
-#endif
-
 /*============================================================================*/
 /* Public definitions                                                         */
 /*============================================================================*/
@@ -652,7 +772,11 @@ void gt_exp(gt_t c, const gt_t a, const bn_t b) {
 #if FP_PRIME == 1536 || FP_PRIME == 544
 	RLC_CAT(RLC_GT_LOWER, exp_cyc)(c, a, b);
 #elif FP_PRIME < 1536
-	gt_exp_gls_imp(c, a, b, ep_curve_frdim());
+	if (ep_curve_embed() == 18) {
+		gt_exp_gls_sac(c, a, b, 1, ep_curve_frdim());
+	} else {
+		gt_exp_gls_naf(c, a, b, ep_curve_frdim());
+	}
 #else
 	RLC_CAT(RLC_GT_LOWER, exp)(c, a, b);
 #endif
@@ -670,7 +794,7 @@ void gt_exp_sec(gt_t c, const gt_t a, const bn_t b) {
 #if FP_PRIME < 1536
 	size_t d = ep_curve_frdim();
 	d = (d > 4 ? d / 4 : 1);
-	gt_exp_reg_gls(c, a, b, d, ep_curve_frdim());
+	gt_exp_reg_sac(c, a, b, d, ep_curve_frdim());
 #elif FP_PRIME == 1536
 	gt_exp_reg_gls(c, a, b, 1);
 #else