# 公式
本文未采用规范文档中的算法，而是采用《Fast Prime Field Elliptic Curve Cryptography with 256 Bit Primes》这篇文章中的算法：
Input:
if (U1 == U2) then
if (S1 != S2) return POINT_AT_INFINITY
else return POINT_DOUBLE (X1, Y1, Z1)
abort
end
1
2
3
4
5
2
3
4
5
Output:
根据公式实现代码就简单了：
# 算法实现
void gm_point_add(gm_point_t * r, const gm_point_t * a, const gm_point_t * b) {
// U1 = X1 * (Z2)^2
// U2 = X2 * (Z1)^2
// S1 = Y1 * (Z2)^3
// S2 = Y2 * (Z1)^3
// if U1==U2 && S1==S2 return double(a)
// H = U2 - U1
// R = S2 - S1
// X3 = R^2 - H^3 - 2U1 * H^2
// Y3 = R * ( U1 * H^2 - X3) - S1 * H^3
// Z3 = H * Z1 * Z2
gm_bn_t U1, U2;
gm_bn_t S1, S2;
gm_bn_t H;
gm_bn_t R;
gm_bn_t X3, Y3, Z3;
if (gm_is_at_infinity(a)) {
gm_point_copy(r, b);
return;
}
if (gm_is_at_infinity(b)) {
gm_point_copy(r, a);
return;
}
// Z1 ^ 2
gm_bn_sqr(H, a->Z, GM_BN_P);
// Z2 ^ 2
gm_bn_sqr(R, b->Z, GM_BN_P);
// U1 = X1 * (Z2)^2
gm_bn_mont_mul(U1, a->X, R, GM_BN_P);
// U2 = X2 * (Z1)^2
gm_bn_mont_mul(U2, b->X, H, GM_BN_P);
// Z1 ^ 3
gm_bn_mont_mul(H, H, a->Z, GM_BN_P);
// Z2 ^ 3
gm_bn_mont_mul(R, R, b->Z, GM_BN_P);
// S1 = Y1 * (Z2)^3
gm_bn_mont_mul(S1, a->Y, R, GM_BN_P);
// S2 = Y2 * (Z1)^3
gm_bn_mont_mul(S2, b->Y, H, GM_BN_P);
if(gm_bn_cmp(U1, U2) == 0 && gm_bn_cmp(S1, S2) == 0) {
gm_point_double(r, a);
return;
}
// H = U2 - U1
gm_bn_sub(H, U2, U1, GM_BN_P);
// R = S2 - S1
gm_bn_sub(R, S2, S1, GM_BN_P);
// Z3 = H * Z1 * Z2
gm_bn_mont_mul(Z3, H, a->Z, GM_BN_P);
gm_bn_mont_mul(Z3, Z3, b->Z, GM_BN_P);
// H ^ 2
gm_bn_sqr(U2, H, GM_BN_P);
// H ^ 3
gm_bn_mont_mul(H, U2, H, GM_BN_P);
// U1 * H ^2
gm_bn_mont_mul(Y3, U1, U2, GM_BN_P);
// 2 * U1 * H ^2
gm_bn_add(U2, Y3, Y3, GM_BN_P);
// R^2
gm_bn_sqr(X3, R, GM_BN_P);
// X3 = R^2 - H^3 - 2U1 * H^2
gm_bn_sub(X3, X3, H, GM_BN_P);
gm_bn_sub(X3, X3, U2, GM_BN_P);
// Y3 = R * ( U1 * H^2 - X3) - S1 * H^3
// ( U1 * H^2 - X3)
gm_bn_sub(Y3, Y3, X3, GM_BN_P);
// R * ( U1 * H^2 - X3)
gm_bn_mont_mul(Y3, Y3, R, GM_BN_P);
// S1 * H^3
gm_bn_mont_mul(S1, S1, H, GM_BN_P);
// Y3 = R * ( U1 * H^2 - X3) - S1 * H^3
gm_bn_sub(Y3, Y3, S1, GM_BN_P);
// output
gm_bn_copy(r->X, X3);
gm_bn_copy(r->Y, Y3);
gm_bn_copy(r->Z, Z3);
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
# 单元测试
test_ec函数增加：
else if(strcmp(alg, "point_add") == 0){ // 10万
int i = 0;
for (i = 0; i < 100000; i++) {
gm_point_add(&r, &r, &p2);
}
}
1
2
3
4
5
6
2
3
4
5
6
main函数增加：
TEST_EC_ALG("point_add",
"EB04AAE0D53FBA1E3611D5B9ED6EFA3EE5BA57C41AA7A09DDC5816AF09057757CE6FA0678392F4716E45F58E7322C76D5997B1FE44C36D8A5A59B146EE162B93",
"AA307A6575D8348037CEC6F860F6312317B34C81838834EDB008F54A2E590FDC593293D89FE9C933E6CE7E91CD4ABF81EC3C26395622B65754A8C0EE8FB354E9",
"3B94D813424CC514F0C05A6D0D7A84AD321AFD86DEC60BA7A4CEDEA4F82C67E5A9ABC4ACE4B15ADE7B5175F8E9E0E33C2B2A88107C67BE15596CD83CBBF4E244");
1
2
3
4
2
3
4
算法效率：
saintdeMacBook-Pro:bn saint$ time ./a.out point_add
r = 3B94D813424CC514F0C05A6D0D7A84AD321AFD86DEC60BA7A4CEDEA4F82C67E5A9ABC4ACE4B15ADE7B5175F8E9E0E33C2B2A88107C67BE15596CD83CBBF4E244
test result: ok
real 0m2.081s
user 0m1.986s
sys 0m0.020s
1
2
3
4
5
6
7
2
3
4
5
6
7
未经本人同意，禁止转载！