Make bn_sqra_low() closer to bn_mula_low(). It should also fix #172.

include/low/relic_bn_low.h

@@ -202,7 +202,7 @@ void bn_rshd_low(dig_t *c, const dig_t *a, int size, int digits);
  * @param[in] a				- the digit vector to multiply.
  * @param[in] digit			- the digit to multiply.
  * @param[in] size			- the number of digits to multiply.
- * @return the carry of the addition.
+ * @return the carry of the last addition.
  */
 dig_t bn_mula_low(dig_t *c, const dig_t *a, dig_t digit, int size);
 
@@ -251,8 +251,9 @@ void bn_muld_low(dig_t *c, const dig_t *a, int sa, const dig_t *b, int sb,
  * @param[out] c			- the result.
  * @param[in] a				- the digit vector to square.
  * @param[in] size			- the number of digitss to square.
+ * @return the carry of the last addition.
  */
-void bn_sqra_low(dig_t *c, const dig_t *a, int size);
+dig_t bn_sqra_low(dig_t *c, const dig_t *a, int size);
 
 /**
  * Squares a digit vector. Computes c = a * a.

src/bn/relic_bn_sqr.c

@@ -170,9 +170,11 @@ void bn_sqr_basic(bn_t c, const bn_t a) {
 		bn_zero(t);
 		t->used = digits;
 
-		for (i = 0; i < a->used; i++) {
-			bn_sqra_low(t->dp + (2 * i), a->dp + i, a->used - i);
+		for (i = 0; i < a->used - 1; i++) {
+			t->dp[a->used + i + 1] =
+					bn_sqra_low(t->dp + 2 * i, a->dp + i, a->used - i);
 		}
+		bn_sqra_low(t->dp + 2 * i, a->dp + i, 1);
 
 		t->sign = RLC_POS;
 		bn_trim(t);

src/fp/relic_fp_mul.c

@@ -74,6 +74,8 @@ static void fp_mul_karat_imp(dv_t c, const fp_t a, const fp_t b, int size,
 		/* a0b0 = a0 * b0 and a1b1 = a1 * b1 */
 		if (level <= 1) {
 #if FP_MUL == BASIC
+			dv_zero(a0b0, h);
+			dv_zero(a1b1, h);
 			for (int i = 0; i < h; i++) {
 				carry = bn_mula_low(a0b0 + i, a, *(b + i), h);
 				*(a0b0 + i + h) = carry;
@@ -111,6 +113,7 @@ static void fp_mul_karat_imp(dv_t c, const fp_t a, const fp_t b, int size,
 		if (level <= 1) {
 			/* t = (a1 + a0)*(b1 + b0) */
 #if FP_MUL == BASIC
+			dv_zero(t, h1 + 1);
 			for (int i = 0; i < h1 + 1; i++) {
 				carry = bn_mula_low(t + i, a1, *(b1 + i), h1 + 1);
 				*(t + i + h1 + 1) = carry;

src/fp/relic_fp_sqr.c

@@ -77,12 +77,14 @@ static void fp_sqr_karat_imp(dv_t c, const fp_t a, int size, int level) {
 		if (level <= 1) {
 			/* a0a0 = a0 * a0 and a1a1 = a1 * a1 */
 #if FP_SQR == BASIC
-			for (i = 0; i < h; i++) {
-				bn_sqra_low(a0a0 + (2 * i), a + i, h - i);
+			for (i = 0; i < h - 1; i++) {
+				a0a0[h + i + 1] = bn_sqra_low(a0a0 + 2 * i, a + i, h - i);
 			}
-			for (i = 0; i < h1; i++) {
-				bn_sqra_low(a1a1 + (2 * i), a + h + i, h1 - i);
+			bn_sqra_low(a0a0 + 2 * i, a + i, 1);
+			for (i = 0; i < h1 - 1; i++) {
+				a1a1[h1 + i + 1] = bn_sqra_low(a1a1 + 2 * i, a + h + i, h1 - i);
 			}
+			bn_sqra_low(a1a1 + 2 * i, a + h + i, 1);
 #elif FP_SQR == COMBA || FP_SQR == INTEG
 			bn_sqrn_low(a0a0, a, h);
 			bn_sqrn_low(a1a1, a + h, h1);
@@ -113,9 +115,10 @@ static void fp_sqr_karat_imp(dv_t c, const fp_t a, int size, int level) {
 		if (level <= 1) {
 			/* a1a1 = (a1 + a0)*(a1 + a0) */
 #if FP_SQR == BASIC
-			for (i = 0; i < h1 + 1; i++) {
-				bn_sqra_low(t1 + (2 * i), t0 + i, h1 + 1 - i);
+			for (i = 0; i < h1; i++) {
+				t1[h1 + i + 2] = bn_sqra_low(t1 + (2 * i), t0 + i, h1 + 1 - i);
 			}
+			bn_sqra_low(t1 + (2 * i), t0 + i, 1);
 #elif FP_SQR == COMBA || FP_SQR == INTEG
 			bn_sqrn_low(t1, t0, h1 + 1);
 #elif FP_SQR == MULTP
@@ -167,9 +170,11 @@ void fp_sqr_basic(fp_t c, const fp_t a) {
 		dv_new(t);
 		dv_zero(t, 2 * RLC_FP_DIGS);
 
-		for (i = 0; i < RLC_FP_DIGS; i++) {
-			bn_sqra_low(t + (2 * i), a + i, RLC_FP_DIGS - i);
+		for (i = 0; i < RLC_FP_DIGS - 1; i++) {
+			t[RLC_FP_DIGS + i + 1] =
+					bn_sqra_low(t + 2 * i, a + i, RLC_FP_DIGS - i);
 		}
+		bn_sqra_low(t + 2 * i, a + i, 1);
 
 		fp_rdc(c, t);
 	}

src/low/easy/relic_bn_sqr_low.c

@@ -78,39 +78,33 @@
 /* Public definitions                                                         */
 /*============================================================================*/
 
-void bn_sqra_low(dig_t *c, const dig_t *a, int size) {
+dig_t bn_sqra_low(dig_t *c, const dig_t *a, int size) {
 	int i;
 	dig_t c0, c1;
-	dig_t digit;
 	dbl_t r, r0, r1;
 
-	digit = *a;
-
 	/* Accumulate this column with the square of a->dp[i]. */
-	r = (dbl_t)(*c) + (dbl_t)(digit) * (dbl_t)(digit);
-
-	*c = (dig_t)r;
+	r = (dbl_t)(*c) + (dbl_t)(a[0]) * (dbl_t)(a[0]);
+	c[0] = (dig_t)r;
 
 	/* Update the carry. */
 	c0 = (dig_t)(r >> (dbl_t)RLC_DIG);
 	c1 = 0;
 
-	c++;
-	a++;
-	for (i = 0; i < size - 1; i++, a++, c++) {
-		r = (dbl_t)(digit) * (dbl_t)(*a);
+	for (i = 1; i < size; i++) {
+		r = (dbl_t)(a[0]) * (dbl_t)(a[i]);
 		r0 = r + r;
-		r1 = r0 + (dbl_t)(*c) + (dbl_t)(c0);
-		*c = (dig_t)r1;
+		r1 = r0 + (dbl_t)(c[i]) + (dbl_t)(c0);
+		c[i] = (dig_t)r1;
 
 		/* Accumulate the old delayed carry. */
 		c0 = (dig_t)((r1 >> (dbl_t)RLC_DIG) + c1);
 		/* Compute the new delayed carry. */
 		c1 = (r0 < r) || (r1 < r0) || (c0 < c1);
 	}
-	*c += c0;
-	c1 += (*c++ < c0);
-	*c += c1;
+	c[size] += c0;
+	c1 += (c[size] < c0);
+	return c1;
 }
 
 void bn_sqrn_low(dig_t *c, const dig_t *a, int size) {

src/low/gmp/relic_bn_sqr_low.c

@@ -40,15 +40,16 @@
 /* Public definitions                                                         */
 /*============================================================================*/
 
-void bn_sqra_low(dig_t *c, const dig_t *a, int size) {
-	dig_t carry, digit = *a;
+dig_t bn_sqra_low(dig_t *c, const dig_t *a, int size) {
+	dig_t carry, digit = a[0];
 
-	carry = mpn_addmul_1(c, a, size, digit);
+	carry = mpn_addmul_1(c, a, size, a[0]);
 	mpn_add_1(c + size, c + size, size, carry);
 	if (size > 1) {
 		carry = mpn_addmul_1(c + 1, a + 1, size - 1, digit);
-		mpn_add_1(c + size, c + size, size, carry);
+		return mpn_add_1(c + size, c + size, size, carry);
 	}
+	return 0;
 }
 
 void bn_sqrn_low(dig_t *c, const dig_t *a, int size) {

src/low/gmp_sec/relic_bn_sqr_low.c

@@ -41,15 +41,16 @@
 /* Public definitions                                                         */
 /*============================================================================*/
 
-void bn_sqra_low(dig_t *c, const dig_t *a, int size) {
-	dig_t carry, digit = *a;
+dig_t bn_sqra_low(dig_t *c, const dig_t *a, int size) {
+	dig_t carry, digit = a[0];
 
-	carry = bn_mula_low(c, a, digit, size);
-	bn_add1_low(c + size, c + size, carry, size);
+	carry = mpn_addmul_1(c, a, size, a[0]);
+	mpn_add_1(c + size, c + size, size, carry);
 	if (size > 1) {
-		carry = bn_mula_low(c + 1, a + 1, digit, size - 1);
-		bn_add1_low(c + size, c + size, carry, size);
+		carry = mpn_addmul_1(c + 1, a + 1, size - 1, digit);
+		return mpn_add_1(c + size, c + size, size, carry);
 	}
+	return 0;
 }
 
 void bn_sqrn_low(dig_t *c, const dig_t *a, int size) {