# GMP library usage

Mehari Zeraberuk mehzera75 at gmail.com
Wed Jun 27 10:36:28 CEST 2012

```I have implemented modified Montgomery algorithm for scalar multiplication
on elliptic curves over prime field. That means, I have generated an
elliptic curve and then I want to perform scalar multiplication to obtain a
different point taking a base point.

It works fine without any error. The problem that I have noticed is when I
compare my results with sage notebook() output with the same input. That
is, both my program routine written using GMP library and the sage notebook
have taken the same elliptic curve and the same base point. I have computed
15 outputs but these all outputs don't match exactly.

Let say, I have a base point P in given elliptic curve E. And what I have
computed is 2P, 3P, 4P, ..., 16P.

When I compare my results with sage's notebook result:
2P and 3P are the same in both. When I go forward they are not the same but
the values exist in each other. Example,value of 4P in my program is the
same as the value of 6P in sage. And the value of 6P in my program is the
value 4P in sage. This relation exists
in between [5P,7p], [9P,13P], [11P,15P] of my program and sage's output.

The only value which doesn't exist in sage's output is the value of 16P of
my program.

Here is the output:
Sages output first then my program's output:
[2]P (394338732385753539704329114170058740226583411304 : 0 : 1)
[2]P (394338732385753539704329114170058740226583411304:0:1)

[3]P (61173517351975174025225734916610194804006447529 : 0 : 1)
[3]P (61173517351975174025225734916610194804006447529:0:1)

[4]P (77045660826527804546356207752200062408330735107 : 0 : 1)
[4]P (445461892534884558282398070150296628962757233304:0:1)-----------|
|
[5]P (212407178542199801410776275004692921476663966405 : 0 : 1)       |
[5]P (1102091121897397878843851737345661014474664901937:0:1)          |
|
[6]P (445461892534884558282398070150296628962757233304 : 0 : 1)-------|
[6]P (77045660826527804546356207752200062408330735107:0:1)

[7]P (1102091121897397878843851737345661014474664901937 : 0 : 1)
[7]P (212407178542199801410776275004692921476663966405:0:1)

[8]P (1063735419006570896603586319780523366738061424268 : 0 : 1)
[8]P (1157112365715925853006705955346605399358986947254:0:1)

[9]P (942403039234610257670757092808579821116134722128 : 0 : 1)
[9]P (455741917462293985843834200468460237371525116124:0:1)

[10]P (262694984029940552988483558698403039057677978822 : 0 : 1)
[10]P (540417159514270199435927221390544500908169034214:0:1)

[11]P (727583468425696204778422848236071880936637020034 : 0 : 1)
[11]P (704236662606484567083005298243909411654391423965:0:1)

[12]P (413415304392460041897472137129378692656644760201 : 0 : 1)
[12]P (262694984029940552988483558698403039057677978822:0:1)

[13]P (455741917462293985843834200468460237371525116124 : 0 : 1)
[13]P (942403039234610257670757092808579821116134722128:0:1)

[14]P (1157112365715925853006705955346605399358986947254 : 0 : 1)
[14]P (413415304392460041897472137129378692656644760201:0:1)

[15]P (704236662606484567083005298243909411654391423965 : 0 : 1)
[15]P (727583468425696204778422848236071880936637020034:0:1)

[16]P (540417159514270199435927221390544500908169034214 : 0 : 1)
[16]P (160022303035305699800324895708510435712240396897:0:1)

It is obvious that modified Montgomery algorithm is working properly form
the outputs I have got above.
Thus, I guess there exists some sort of error in my code though it computes
properly without any error.
I guess the problem is around passing structure pointers but I couldn't
figure out where the problem is.
what I have used from gmp librayr to compute scalar multiplication.
-------------- next part --------------
void EC_GFp_xECDBL(EC_GROUP *group, EC_POINT *dbl_point, const EC_POINT *point){

mpz_t R1, R2, R3, R4, R5, R6;

mpz_init(R1);
mpz_init(R2);
mpz_init(R3);
mpz_init(R4);
mpz_init(R5);
mpz_init(R6);

mpz_set(R1, point->X);
mpz_set(R2, point->Z);

mpz_powm_ui(R3, R1, 2, group->field);
mpz_powm_ui(R4, R2, 2, group->field);

EC_mulmod(&R5, (const mpz_t *)&R4, (const mpz_t *)&group->a, (const mpz_t *)&group->field);
EC_mulmod(&R6, (const mpz_t *)&R4, (const mpz_t *)&group->b, (const mpz_t *)&group->field);
EC_mulmod(&R4, (const mpz_t *)&R4, (const mpz_t *)&R6, (const mpz_t *)&group->field);
EC_mulmod(&R2, (const mpz_t *)&R1, (const mpz_t *)&R2, (const mpz_t *)&group->field);
EC_submod(&R1, (const mpz_t *)&R3, (const mpz_t *)&R5,(const mpz_t *)&group->field);
EC_addmod(&R5, (const mpz_t *)&R3, (const mpz_t *)&R5, (const mpz_t *)&group->field);
mpz_powm_ui(R1, R1, 2, group->field);

EC_mulmod(&R3, (const mpz_t *)&R6, (const mpz_t *)&R2, (const mpz_t *)&group->field);
EC_addmod(&R3, (const mpz_t *)&R3, (const mpz_t *)&R3, (const mpz_t *)&group->field);
EC_addmod(&R3, (const mpz_t *)&R3, (const mpz_t *)&R3, (const mpz_t *)&group->field);
EC_addmod(&R3, (const mpz_t *)&R3, (const mpz_t *)&R3, (const mpz_t *)&group->field);
EC_submod(&R1, (const mpz_t *)&R1, (const mpz_t *)&R3, (const mpz_t *)&group->field);
EC_mulmod(&R5, (const mpz_t *)&R5, (const mpz_t *)&R2, (const mpz_t *)&group->field);
EC_addmod(&R4, (const mpz_t *)&R5, (const mpz_t *)&R4, (const mpz_t *)&group->field);
EC_addmod(&R4, (const mpz_t *)&R4, (const mpz_t *)&R4, (const mpz_t *)&group->field);
EC_addmod(&R4, (const mpz_t *)&R4, (const mpz_t *)&R4, (const mpz_t *)&group->field);

mpz_set(dbl_point->X, R1);
mpz_set(dbl_point->Z, R4);

mpz_clear(R1);
mpz_clear(R2);
mpz_clear(R3);
mpz_clear(R4);
mpz_clear(R5);
mpz_clear(R6);
}

void EC_GFp_xECADDDBL(EC_GROUP *group, EC_POINT *point1, EC_POINT *point2, const mpz_t *x){

mpz_t R1, R2, R3, R4, R5, R6, R7;

mpz_init(R1);
mpz_init(R2);
mpz_init(R3);
mpz_init(R4);
mpz_init(R5);
mpz_init(R6);
mpz_init(R7);

mpz_set(R1, point1->X);
mpz_set(R2, point1->Z);
mpz_set(R3, point2->X);
mpz_set(R4, point2->Z);

EC_mulmod(&R6, (const mpz_t *)&R1, (const mpz_t *)&R4, (const mpz_t *)&group->field);//mpz_mul(R6, R1, R4); mpz_mod(R6, R6, P);
EC_mulmod(&R1, (const mpz_t *)&R1, (const mpz_t *)&R3, (const mpz_t *)&group->field);//mpz_mul(R1, R1, R3); mpz_mod(R1, R1, P);
EC_mulmod(&R4, (const mpz_t *)&R2, (const mpz_t *)&R4, (const mpz_t *)&group->field);//mpz_mul(R4, R2, R4); mpz_mod(R4, R4, P);
EC_mulmod(&R2, (const mpz_t *)&R3, (const mpz_t *)&R2, (const mpz_t *)&group->field);//mpz_mul(R2, R3, R2); mpz_mod(R2, R2, P);
EC_submod(&R3, (const mpz_t *)&R6, (const mpz_t *)&R2, (const mpz_t *)&group->field);//mpz_sub(R3, R6, R2); mpz_mod(R3, R3, p);

mpz_powm_ui(R3, R3, 2, group->field);

EC_mulmod(&R5, x, (const mpz_t *)&R3, (const mpz_t *)&group->field);//mpz_mul(R5, x, R5);mpz_mod(R5, R5, P);
EC_mulmod(&R7, (const mpz_t *)&group->a, (const mpz_t *)&R4, (const mpz_t *)&group->field);//mpz_mul(R7, a, R7); mpz_mod(R7, R7, P);
EC_addmod(&R1, (const mpz_t *)&R1, (const mpz_t *)&R7, (const mpz_t *)&group->field);// mpz_add(R1, R1, R1); mpz_mod(R1, R1, P);
EC_addmod(&R2, (const mpz_t *)&R2, (const mpz_t *)&R6, (const mpz_t *)&group->field);//mpz_add(R1, R1, R1); mpz_mod(R1, R1, P);

EC_mulmod(&R1, (const mpz_t *)&R1, (const mpz_t *)&R2, (const mpz_t *)&group->field);
mpz_powm_ui(R2, R4, 2, group->field);
EC_mulmod(&R7, (const mpz_t *)&group->b, (const mpz_t *)&R2, (const mpz_t *)&group->field);
EC_addmod(&R1, (const mpz_t *)&R1, (const mpz_t *)&R7, (const mpz_t *)&group->field);
EC_addmod(&R1, (const mpz_t *)&R1, (const mpz_t *)&R1, (const mpz_t *)&group->field);
EC_submod(&R5, (const mpz_t *)&R1, (const mpz_t *)&R5, (const mpz_t *)&group->field);//mpz_sub(R5, R1, R5); mpz_mod(R5, R5, p);
EC_addmod(&R5, (const mpz_t *)&R7, (const mpz_t *)&R5, (const mpz_t *)&group->field);//mpz_add(R5, R7, R5); mpz_mod(R5, R5, P);
EC_addmod(&R5, (const mpz_t *)&R7, (const mpz_t *)&R5, (const mpz_t *)&group->field);//mpz_add(R5, R7, R5); mpz_mod(R5, R5, P);
EC_mulmod(&R2, (const mpz_t *)&group->a, (const mpz_t *)&R2, (const mpz_t *)&group->field);//mpz_mul(R2, a, R2); mpz_mod(R2, R2, P);

mpz_powm_ui(R1, R6, 2, group->field);

EC_addmod(&R1, (const mpz_t *)&R1, (const mpz_t *)&R2, (const mpz_t *)&group->field);////mpz_add(R1, R1, R2); mpz_mod(R1, R1, P);
EC_addmod(&R2, (const mpz_t *)&R2, (const mpz_t *)&R2, (const mpz_t *)&group->field);//mpz_add(R2, R2, R2); mpz_mod(R2, R2, P);
EC_submod(&R2, (const mpz_t *)&R1, (const mpz_t *)&R2, (const mpz_t *)&group->field);//mpz_sub(R2, R1, R2); mpz_mod(R2, R2, p);

mpz_powm_ui(R2, R2, 2, group->field);

EC_mulmod(&R1, (const mpz_t *)&R6, (const mpz_t *)&R1, (const mpz_t *)&group->field);//mpz_mul(R1, R6, R1); mpz_mod(R1, R1, P);
EC_mulmod(&R7, (const mpz_t *)&R4, (const mpz_t *)&R7, (const mpz_t *)&group->field);//mpz_mul(R7, R4, R7); mpz_mod(R7, R7, P);
EC_addmod(&R1, (const mpz_t *)&R1, (const mpz_t *)&R7, (const mpz_t *)&group->field);//mpz_add(R1, R1, R7); mpz_mod(R1, R1, P);
EC_mulmod(&R7, (const mpz_t *)&R6, (const mpz_t *)&R7, (const mpz_t *)&group->field);//mpz_mul(R7, R6, R7); mpz_mod(R7, R7, P);
EC_addmod(&R7, (const mpz_t *)&R7, (const mpz_t *)&R7, (const mpz_t *)&group->field);//mpz_add(R7, R7, R7); mpz_mod(R7, R7, P);
EC_addmod(&R7, (const mpz_t *)&R7, (const mpz_t *)&R7, (const mpz_t *)&group->field);//mpz_add(R7, R7, R7); mpz_mod(R7, R7, P);
EC_addmod(&R7, (const mpz_t *)&R7, (const mpz_t *)&R7, (const mpz_t *)&group->field);//mpz_add(R7, R7, R7); mpz_mod(R7, R7, P);
EC_submod(&R7, (const mpz_t *)&R2, (const mpz_t *)&R7, (const mpz_t *)&group->field);//mpz_sub(R7, R2, R7); mpz_mod(R7, R7, p);
EC_mulmod(&R6, (const mpz_t *)&R4, (const mpz_t *)&R1, (const mpz_t *)&group->field);//mpz_mul(R6, R4, R1); mpz_mod(R6, R6, P);
EC_addmod(&R6, (const mpz_t *)&R6, (const mpz_t *)&R6, (const mpz_t *)&group->field);//mpz_add(R6, R6, R6); mpz_mod(R6, R6, P);
EC_addmod(&R6, (const mpz_t *)&R6, (const mpz_t *)&R6, (const mpz_t *)&group->field);//mpz_add(R6, R6, R6); mpz_mod(R6, R6, P);

mpz_set(point1->X, R5);
mpz_set(point1->Z, R3);
mpz_set(point2->X, R7);
mpz_set(point2->Z, R6);

mpz_clear(R1);
mpz_clear(R2);
mpz_clear(R3);
mpz_clear(R4);
mpz_clear(R5);
mpz_clear(R6);
mpz_clear(R7);
}

void ec_GFp_montgomery_point_multiply(EC_GROUP *group, EC_POINT *r, const mpz_t scalar, EC_POINT *point){

mpz_t x;
int inf = 0;
long i;

EC_POINT dbl_point, scalar_point;

ec_GFp_point_init(&dbl_point);//initiailize point
ec_GFp_point_init(&scalar_point);

mpz_init(x);

if ((scalar == NULL) || (mpz_sgn(scalar) == 0)  || (point == NULL)){// || EC_POINT_is_at_infinity((const EC_GROUP*)group,(const EC_POINT*) point)){
EC_POINT_set_to_infinity(group, r);
inf = 1;
goto endp_mont;
}

i = mpz_size(scalar)-1;
word = mpz_getlimbn(scalar, i);

mpz_set(x, point->X);

EC_GFp_xECDBL(group, &dbl_point, point);//double base point

if(mpz_cmp_ui(scalar, 2) == 0){//it is double and no need of further process
EC_POINT_Invert(group, &dbl_point);
mpz_set(r->X, dbl_point.X);
mpz_set(r->Z, dbl_point.Z);

goto endp_mont;
}

/*If top most bit was at word break, go to the next word*/
i--;
}
ec_GFp_point_copy(&scalar_point, point);
while(i >= 0){
word = mpz_getlimbn(scalar, i);//assign one word  or limb from scalar to the variable word
EC_GFp_xECADDDBL(group, &scalar_point, &dbl_point, (const mpz_t *)&x);
}else{
EC_GFp_xECADDDBL(group, &dbl_point, &scalar_point, (const mpz_t *)&x);
}
}
i--;
}
//Value of Y could be recovered at this stage when a need arises
EC_POINT_Invert(group, &scalar_point);
mpz_set(r->X, scalar_point.X);
mpz_set(r->Z, scalar_point.Z);

endp_mont:if(inf){
mpz_set_ui(r->X,0);
mpz_set_ui(r->Z,0);
}
ec_GFp_point_clear(&dbl_point);
ec_GFp_point_clear(&scalar_point);
mpz_clear(x);
}

void EC_GFp_scalar_PointMultiplication(EC_GROUP group, EC_POINT base_point){

int i = 0;
mpz_t scalar;
mpz_init(scalar);
EC_POINT result;
EC_GROUP *gr;
EC_POINT *point, *res;
ec_GFp_point_init(&result);

point = &base_point;
res = &result;
gr = &group;
for(i = 0; i < 15; i++){
mpz_set_ui(scalar, i + 2);
gmp_printf("\n[%Zd]P ", scalar);
ec_GFp_montgomery_point_multiply(gr, res, scalar, point);
print_point(&result);
}
ec_GFp_point_clear(&result);
mpz_clear(scalar);
}
void EC_GFp_Test(){
EC_GROUP group;
EC_POINT base_point;
ec_GFp_point_init(&base_point);
ec_GFp_group_init(&group);
mpz_set_str(group.field, "1200265344304484529820484878267404952334866443931", 10);
mpz_set_str(group.a, "33438906192665209100362084073678594186987856027", 10);
mpz_set_str(group.b, "1200265344304484529820484878267404952334866443538", 10);

mpz_set_str(base_point.X, "625634078674629860449874299851267027831111903985", 10);
mpz_set_str(base_point.Y, "165579142211495461916769326936888490114303428494", 10);
mpz_set_ui(base_point.Z, 1);

EC_GFp_scalar_PointMultiplication(group, base_point);
ec_GFp_point_clear(&base_point);
ec_GFp_group_clear(&group);
}
```