ec.c

00001 /**
00002  * @file ec.c
00003  * @brief ECC (Elliptic Curve Cryptography)
00004  *
00005  * @section License
00006  *
00007  * Copyright (C) 2010-2017 Oryx Embedded SARL. All rights reserved.
00008  *
00009  * This file is part of CycloneCrypto Open.
00010  *
00011  * This program is free software; you can redistribute it and/or
00012  * modify it under the terms of the GNU General Public License
00013  * as published by the Free Software Foundation; either version 2
00014  * of the License, or (at your option) any later version.
00015  *
00016  * This program is distributed in the hope that it will be useful,
00017  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00018  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00019  * GNU General Public License for more details.
00020  *
00021  * You should have received a copy of the GNU General Public License
00022  * along with this program; if not, write to the Free Software Foundation,
00023  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
00024  *
00025  * @author Oryx Embedded SARL (www.oryx-embedded.com)
00026  * @version 1.7.6
00027  **/
00028 
00029 //Switch to the appropriate trace level
00030 #define TRACE_LEVEL CRYPTO_TRACE_LEVEL
00031 
00032 //Dependencies
00033 #include <stdlib.h>
00034 #include "crypto.h"
00035 #include "ec.h"
00036 #include "debug.h"
00037 
00038 //Check crypto library configuration
00039 #if (EC_SUPPORT == ENABLED)
00040 
00041 //EC Public Key OID (1.2.840.10045.2.1)
00042 const uint8_t EC_PUBLIC_KEY_OID[7] = {0x2A, 0x86, 0x48, 0xCE, 0x3D, 0x02, 0x01};
00043 
00044 
00045 /**
00046  * @brief Initialize EC domain parameters
00047  * @param[in] params Pointer to the EC domain parameters to be initialized
00048  **/
00049 
00050 void ecInitDomainParameters(EcDomainParameters *params)
00051 {
00052    //Initialize structure
00053    params->type = EC_CURVE_TYPE_NONE;
00054    params->mod = NULL;
00055 
00056    //Initialize EC domain parameters
00057    mpiInit(&params->p);
00058    mpiInit(&params->a);
00059    mpiInit(&params->b);
00060    ecInit(&params->g);
00061    mpiInit(&params->q);
00062 }
00063 
00064 
00065 /**
00066  * @brief Release EC domain parameters
00067  * @param[in] params Pointer to the EC domain parameters to free
00068  **/
00069 
00070 void ecFreeDomainParameters(EcDomainParameters *params)
00071 {
00072    //Release previously allocated resources
00073    mpiFree(&params->p);
00074    mpiFree(&params->a);
00075    mpiFree(&params->b);
00076    ecFree(&params->g);
00077    mpiFree(&params->q);
00078 }
00079 
00080 
00081 /**
00082  * @brief Load EC domain parameters
00083  * @param[out] params Pointer to the structure to be initialized
00084  * @param[in] curveInfo Elliptic curve parameters
00085  * @return Error code
00086  **/
00087 
00088 error_t ecLoadDomainParameters(EcDomainParameters *params, const EcCurveInfo *curveInfo)
00089 {
00090    error_t error;
00091 
00092    //Debug message
00093    TRACE_DEBUG("Loading %s EC domain parameters...\r\n", curveInfo->name);
00094 
00095    //Curve type
00096    params->type = curveInfo->type;
00097 
00098    //Import prime modulus
00099    MPI_CHECK(mpiReadRaw(&params->p, curveInfo->p, curveInfo->pLen));
00100    //Import parameter a
00101    MPI_CHECK(mpiReadRaw(&params->a, curveInfo->a, curveInfo->aLen));
00102    //Import parameter b
00103    MPI_CHECK(mpiReadRaw(&params->b, curveInfo->b, curveInfo->bLen));
00104    //Import the x-coordinate of the base point G
00105    MPI_CHECK(mpiReadRaw(&params->g.x, curveInfo->gx, curveInfo->gxLen));
00106    //Import the y-coordinate of the base point G
00107    MPI_CHECK(mpiReadRaw(&params->g.y, curveInfo->gy, curveInfo->gyLen));
00108    //Import base point order q
00109    MPI_CHECK(mpiReadRaw(&params->q, curveInfo->q, curveInfo->qLen));
00110 
00111    //Normalize base point G
00112    MPI_CHECK(mpiSetValue(&params->g.z, 1));
00113 
00114    //Fast modular reduction
00115    params->mod = curveInfo->mod;
00116 
00117    //Debug message
00118    TRACE_DEBUG("  p:\r\n");
00119    TRACE_DEBUG_MPI("    ", &params->p);
00120    TRACE_DEBUG("  a:\r\n");
00121    TRACE_DEBUG_MPI("    ", &params->a);
00122    TRACE_DEBUG("  b:\r\n");
00123    TRACE_DEBUG_MPI("    ", &params->b);
00124    TRACE_DEBUG("  Gx:\r\n");
00125    TRACE_DEBUG_MPI("    ", &params->g.x);
00126    TRACE_DEBUG("  Gy:\r\n");
00127    TRACE_DEBUG_MPI("    ", &params->g.y);
00128    TRACE_DEBUG("  q:\r\n");
00129    TRACE_DEBUG_MPI("    ", &params->q);
00130 
00131 end:
00132    //Return status code
00133    return error;
00134 }
00135 
00136 
00137 /**
00138  * @brief Initialize elliptic curve point
00139  * @param[in,out] r Pointer to the EC point to be initialized
00140  **/
00141 
00142 void ecInit(EcPoint *r)
00143 {
00144    //Initialize structure
00145    mpiInit(&r->x);
00146    mpiInit(&r->y);
00147    mpiInit(&r->z);
00148 }
00149 
00150 
00151 /**
00152  * @brief Release an elliptic curve point
00153  * @param[in,out] r Pointer to the EC point to initialize to free
00154  **/
00155 
00156 void ecFree(EcPoint *r)
00157 {
00158    //Release previously allocated resources
00159    mpiFree(&r->x);
00160    mpiFree(&r->y);
00161    mpiFree(&r->z);
00162 }
00163 
00164 
00165 /**
00166  * @brief Copy EC point
00167  * @param[out] r Destination EC point
00168  * @param[in] s Source EC point
00169  * @return Error code
00170  **/
00171 
00172 error_t ecCopy(EcPoint *r, const EcPoint *s)
00173 {
00174    error_t error;
00175 
00176    //R and S are the same instance?
00177    if(r == s)
00178       return NO_ERROR;
00179 
00180    //Copy coordinates
00181    MPI_CHECK(mpiCopy(&r->x, &s->x));
00182    MPI_CHECK(mpiCopy(&r->y, &s->y));
00183    MPI_CHECK(mpiCopy(&r->z, &s->z));
00184 
00185 end:
00186    //Return status code
00187    return error;
00188 }
00189 
00190 
00191 /**
00192  * @brief Convert an octet string to an EC point
00193  * @param[in] params EC domain parameters
00194  * @param[out] r EC point resulting from the conversion
00195  * @param[in] data Pointer to the octet string
00196  * @param[in] length Length of the octet string
00197  * @return Error code
00198  **/
00199 
00200 error_t ecImport(const EcDomainParameters *params,
00201    EcPoint *r, const uint8_t *data, size_t length)
00202 {
00203    error_t error;
00204    size_t k;
00205 
00206    //Get the length in octets of the prime
00207    k = mpiGetByteLength(&params->p);
00208 
00209    //Check the length of the octet string
00210    if(length != (k * 2 + 1))
00211       return ERROR_DECODING_FAILED;
00212 
00213    //Compressed point representation is not supported
00214    if(data[0] != 0x04)
00215       return ERROR_ILLEGAL_PARAMETER;
00216 
00217    //Convert the x-coordinate to a multiple precision integer
00218    error = mpiReadRaw(&r->x, data + 1, k);
00219    //Any error to report?
00220    if(error)
00221       return error;
00222 
00223    //Convert the y-coordinate to a multiple precision integer
00224    error = mpiReadRaw(&r->y, data + k + 1, k);
00225    //Any error to report?
00226    if(error)
00227       return error;
00228 
00229    //Successful processing
00230    return NO_ERROR;
00231 }
00232 
00233 
00234 /**
00235  * @brief Convert an EC point to an octet string
00236  * @param[in] params EC domain parameters
00237  * @param[in] a EC point to be converted
00238  * @param[out] data Pointer to the octet string
00239  * @param[out] length Length of the resulting octet string
00240  * @return Error code
00241  **/
00242 
00243 error_t ecExport(const EcDomainParameters *params,
00244    const EcPoint *a, uint8_t *data, size_t *length)
00245 {
00246    error_t error;
00247    size_t k;
00248 
00249    //Get the length in octets of the prime
00250    k = mpiGetByteLength(&params->p);
00251 
00252    //Point compression is not used
00253    data[0] = 0x04;
00254 
00255    //Convert the x-coordinate to an octet string
00256    error = mpiWriteRaw(&a->x, data + 1, k);
00257    //Conversion failed?
00258    if(error)
00259       return error;
00260 
00261    //Convert the y-coordinate to an octet string
00262    error = mpiWriteRaw(&a->y, data + k + 1, k);
00263    //Conversion failed?
00264    if(error)
00265       return error;
00266 
00267    //Return the total number of bytes that have been written
00268    *length = k * 2 + 1;
00269    //Successful processing
00270    return NO_ERROR;
00271 }
00272 
00273 
00274 /**
00275  * @brief Compute projective representation
00276  * @param[in] params EC domain parameters
00277  * @param[out] r Projective representation of the point
00278  * @param[in] s Affine representation of the point
00279  * @return Error code
00280  **/
00281 
00282 error_t ecProjectify(const EcDomainParameters *params, EcPoint *r, const EcPoint *s)
00283 {
00284    error_t error;
00285 
00286    //Copy point
00287    EC_CHECK(ecCopy(r, s));
00288    //Map the point to projective space
00289    MPI_CHECK(mpiSetValue(&r->z, 1));
00290 
00291 end:
00292    //Return status code
00293    return error;
00294 }
00295 
00296 
00297 /**
00298  * @brief Recover affine representation
00299  * @param[in] params EC domain parameters
00300  * @param[out] r Affine representation of the point
00301  * @param[in] s Projective representation of the point
00302  * @return Error code
00303  **/
00304 
00305 error_t ecAffinify(const EcDomainParameters *params, EcPoint *r, const EcPoint *s)
00306 {
00307    error_t error;
00308    Mpi a;
00309    Mpi b;
00310 
00311    //Point at the infinity?
00312    if(!mpiCompInt(&s->z, 0))
00313       return ERROR_INVALID_PARAMETER;
00314 
00315    //Initialize multiple precision integers
00316    mpiInit(&a);
00317    mpiInit(&b);
00318 
00319    //Compute a = 1/Sz mod p
00320    MPI_CHECK(mpiInvMod(&a, &s->z, &params->p));
00321 
00322    //Set Rx = a^2 * Sx mod p
00323    EC_CHECK(ecSqrMod(params, &b, &a));
00324    EC_CHECK(ecMulMod(params, &r->x, &b, &s->x));
00325 
00326    //Set Ry = a^3 * Sy mod p
00327    EC_CHECK(ecMulMod(params, &b, &b, &a));
00328    EC_CHECK(ecMulMod(params, &r->y, &b, &s->y));
00329 
00330    //Set Rz = 1
00331    MPI_CHECK(mpiSetValue(&r->z, 1));
00332 
00333 end:
00334    //Release multiple precision integers
00335    mpiFree(&a);
00336    mpiFree(&b);
00337 
00338    //Return status code
00339    return error;
00340 }
00341 
00342 
00343 /**
00344  * @brief Check whether the affine point S is on the curve
00345  * @param[in] params EC domain parameters
00346  * @param[in] s Affine representation of the point
00347  * @return TRUE if the affine point S is on the curve, else FALSE
00348  **/
00349 
00350 bool_t ecIsPointAffine(const EcDomainParameters *params, const EcPoint *s)
00351 {
00352    error_t error;
00353    Mpi t1;
00354    Mpi t2;
00355 
00356    //Initialize multiple precision integers
00357    mpiInit(&t1);
00358    mpiInit(&t2);
00359 
00360    //Compute t1 = (Sx^3 + a * Sx + b) mod p
00361    EC_CHECK(ecSqrMod(params, &t1, &s->x));
00362    EC_CHECK(ecMulMod(params, &t1, &t1, &s->x));
00363    EC_CHECK(ecMulMod(params, &t2, &params->a, &s->x));
00364    EC_CHECK(ecAddMod(params, &t1, &t1, &t2));
00365    EC_CHECK(ecAddMod(params, &t1, &t1, &params->b));
00366 
00367    //Compute t2 = Sy^2
00368    EC_CHECK(ecSqrMod(params, &t2, &s->y));
00369 
00370    //Check whether the point is on the elliptic curve
00371    if(mpiComp(&t1, &t2))
00372       error = ERROR_FAILURE;
00373 
00374 end:
00375    //Release multiple precision integers
00376    mpiFree(&t1);
00377    mpiFree(&t2);
00378 
00379    //Return TRUE if the affine point S is on the curve, else FALSE
00380    return error ? FALSE : TRUE;
00381 }
00382 
00383 
00384 /**
00385  * @brief Point doubling
00386  * @param[in] params EC domain parameters
00387  * @param[out] r Resulting point R = 2S
00388  * @param[in] s Point S
00389  * @return Error code
00390  **/
00391 
00392 error_t ecDouble(const EcDomainParameters *params, EcPoint *r, const EcPoint *s)
00393 {
00394    error_t error;
00395    Mpi t1;
00396    Mpi t2;
00397    Mpi t3;
00398    Mpi t4;
00399    Mpi t5;
00400 
00401    //Initialize multiple precision integers
00402    mpiInit(&t1);
00403    mpiInit(&t2);
00404    mpiInit(&t3);
00405    mpiInit(&t4);
00406    mpiInit(&t5);
00407 
00408    //Set t1 = Sx
00409    MPI_CHECK(mpiCopy(&t1, &s->x));
00410    //Set t2 = Sy
00411    MPI_CHECK(mpiCopy(&t2, &s->y));
00412    //Set t3 = Sz
00413    MPI_CHECK(mpiCopy(&t3, &s->z));
00414 
00415    //Point at the infinity?
00416    if(!mpiCompInt(&t3, 0))
00417    {
00418       //Set R = (1, 1, 0)
00419       MPI_CHECK(mpiSetValue(&r->x, 1));
00420       MPI_CHECK(mpiSetValue(&r->y, 1));
00421       MPI_CHECK(mpiSetValue(&r->z, 0));
00422    }
00423    else
00424    {
00425       //SECP K1 elliptic curve?
00426       if(params->type == EC_CURVE_TYPE_SECP_K1)
00427       {
00428          //Compute t5 = t1^2
00429          EC_CHECK(ecSqrMod(params, &t5, &t1));
00430          //Compute t4 = 3 * t5
00431          EC_CHECK(ecAddMod(params, &t4, &t5, &t5));
00432          EC_CHECK(ecAddMod(params, &t4, &t4, &t5));
00433       }
00434       //SECP R1 elliptic curve?
00435       else if(params->type == EC_CURVE_TYPE_SECP_R1)
00436       {
00437          //Compute t4 = t3^2
00438          EC_CHECK(ecSqrMod(params, &t4, &t3));
00439          //Compute t5 = t1 - t4
00440          EC_CHECK(ecSubMod(params, &t5, &t1, &t4));
00441          //Compute t4 = t1 + t4
00442          EC_CHECK(ecAddMod(params, &t4, &t1, &t4));
00443          //Compute t5 = t4 * t5
00444          EC_CHECK(ecMulMod(params, &t5, &t4, &t5));
00445          //Compute t4 = 3 * t5
00446          EC_CHECK(ecAddMod(params, &t4, &t5, &t5));
00447          EC_CHECK(ecAddMod(params, &t4, &t4, &t5));
00448       }
00449       else
00450       {
00451          //Compute t4 = t3^4
00452          EC_CHECK(ecSqrMod(params, &t4, &t3));
00453          EC_CHECK(ecSqrMod(params, &t4, &t4));
00454          //Compute t4 = a * t4
00455          EC_CHECK(ecMulMod(params, &t4, &t4, &params->a));
00456          //Compute t5 = t1^2
00457          EC_CHECK(ecSqrMod(params, &t5, &t1));
00458          //Compute t4 = t4 + 3 * t5
00459          EC_CHECK(ecAddMod(params, &t4, &t4, &t5));
00460          EC_CHECK(ecAddMod(params, &t4, &t4, &t5));
00461          EC_CHECK(ecAddMod(params, &t4, &t4, &t5));
00462       }
00463 
00464       //Compute t3 = t3 * t2
00465       EC_CHECK(ecMulMod(params, &t3, &t3, &t2));
00466       //Compute t3 = 2 * t3
00467       EC_CHECK(ecAddMod(params, &t3, &t3, &t3));
00468       //Compute t2 = t2^2
00469       EC_CHECK(ecSqrMod(params, &t2, &t2));
00470       //Compute t5 = t1 * t2
00471       EC_CHECK(ecMulMod(params, &t5, &t1, &t2));
00472       //Compute t5 = 4 * t5
00473       EC_CHECK(ecAddMod(params, &t5, &t5, &t5));
00474       EC_CHECK(ecAddMod(params, &t5, &t5, &t5));
00475       //Compute t1 = t4^2
00476       EC_CHECK(ecSqrMod(params, &t1, &t4));
00477       //Compute t1 = t1 - 2 * t5
00478       EC_CHECK(ecSubMod(params, &t1, &t1, &t5));
00479       EC_CHECK(ecSubMod(params, &t1, &t1, &t5));
00480       //Compute t2 = t2^2
00481       EC_CHECK(ecSqrMod(params, &t2, &t2));
00482       //Compute t2 = 8 * t2
00483       EC_CHECK(ecAddMod(params, &t2, &t2, &t2));
00484       EC_CHECK(ecAddMod(params, &t2, &t2, &t2));
00485       EC_CHECK(ecAddMod(params, &t2, &t2, &t2));
00486       //Compute t5 = t5 - t1
00487       EC_CHECK(ecSubMod(params, &t5, &t5, &t1));
00488       //Compute t5 = t4 * t5
00489       EC_CHECK(ecMulMod(params, &t5, &t4, &t5));
00490       //Compute t2 = t5 - t2
00491       EC_CHECK(ecSubMod(params, &t2, &t5, &t2));
00492 
00493       //Set Rx = t1
00494       MPI_CHECK(mpiCopy(&r->x, &t1));
00495       //Set Ry = t2
00496       MPI_CHECK(mpiCopy(&r->y, &t2));
00497       //Set Rz = t3
00498       MPI_CHECK(mpiCopy(&r->z, &t3));
00499    }
00500 
00501 end:
00502    //Release multiple precision integers
00503    mpiFree(&t1);
00504    mpiFree(&t2);
00505    mpiFree(&t3);
00506    mpiFree(&t4);
00507    mpiFree(&t5);
00508 
00509    //Return status code
00510    return error;
00511 }
00512 
00513 
00514 /**
00515  * @brief Point addition (helper routine)
00516  * @param[in] params EC domain parameters
00517  * @param[out] r Resulting point R = S + T
00518  * @param[in] s First operand
00519  * @param[in] t Second operand
00520  * @return Error code
00521  **/
00522 
00523 error_t ecAdd(const EcDomainParameters *params, EcPoint *r, const EcPoint *s, const EcPoint *t)
00524 {
00525    error_t error;
00526    Mpi t1;
00527    Mpi t2;
00528    Mpi t3;
00529    Mpi t4;
00530    Mpi t5;
00531    Mpi t6;
00532    Mpi t7;
00533 
00534    //Initialize multiple precision integers
00535    mpiInit(&t1);
00536    mpiInit(&t2);
00537    mpiInit(&t3);
00538    mpiInit(&t4);
00539    mpiInit(&t5);
00540    mpiInit(&t6);
00541    mpiInit(&t7);
00542 
00543    //Set t1 = Sx
00544    MPI_CHECK(mpiCopy(&t1, &s->x));
00545    //Set t2 = Sy
00546    MPI_CHECK(mpiCopy(&t2, &s->y));
00547    //Set t3 = Sz
00548    MPI_CHECK(mpiCopy(&t3, &s->z));
00549    //Set t4 = Tx
00550    MPI_CHECK(mpiCopy(&t4, &t->x));
00551    //Set t5 = Ty
00552    MPI_CHECK(mpiCopy(&t5, &t->y));
00553 
00554    //Check whether Tz != 1
00555    if(mpiCompInt(&t->z, 1))
00556    {
00557       //Compute t6 = Tz
00558       MPI_CHECK(mpiCopy(&t6, &t->z));
00559       //Compute t7 = t6^2
00560       EC_CHECK(ecSqrMod(params, &t7, &t6));
00561       //Compute t1 = t1 * t7
00562       EC_CHECK(ecMulMod(params, &t1, &t1, &t7));
00563       //Compute t7 = t6 * t7
00564       EC_CHECK(ecMulMod(params, &t7, &t6, &t7));
00565       //Compute t2 = t2 * t7
00566       EC_CHECK(ecMulMod(params, &t2, &t2, &t7));
00567    }
00568 
00569    //Compute t7 = t3^2
00570    EC_CHECK(ecSqrMod(params, &t7, &t3));
00571    //Compute t4 = t4 * t7
00572    EC_CHECK(ecMulMod(params, &t4, &t4, &t7));
00573    //Compute t7 = t3 * t7
00574    EC_CHECK(ecMulMod(params, &t7, &t3, &t7));
00575    //Compute t5 = t5 * t7
00576    EC_CHECK(ecMulMod(params, &t5, &t5, &t7));
00577    //Compute t4 = t1 - t4
00578    EC_CHECK(ecSubMod(params, &t4, &t1, &t4));
00579    //Compute t5 = t2 - t5
00580    EC_CHECK(ecSubMod(params, &t5, &t2, &t5));
00581 
00582    //Check whether t4 == 0
00583    if(!mpiCompInt(&t4, 0))
00584    {
00585       //Check whether t5 == 0
00586       if(!mpiCompInt(&t5, 0))
00587       {
00588          //Set R = (0, 0, 0)
00589          MPI_CHECK(mpiSetValue(&r->x, 0));
00590          MPI_CHECK(mpiSetValue(&r->y, 0));
00591          MPI_CHECK(mpiSetValue(&r->z, 0));
00592       }
00593       else
00594       {
00595          //Set R = (1, 1, 0)
00596          MPI_CHECK(mpiSetValue(&r->x, 1));
00597          MPI_CHECK(mpiSetValue(&r->y, 1));
00598          MPI_CHECK(mpiSetValue(&r->z, 0));
00599       }
00600    }
00601    else
00602    {
00603       //Compute t1 = 2 * t1 - t4
00604       EC_CHECK(ecAddMod(params, &t1, &t1, &t1));
00605       EC_CHECK(ecSubMod(params, &t1, &t1, &t4));
00606       //Compute t2 = 2 * t2 - t5
00607       EC_CHECK(ecAddMod(params, &t2, &t2, &t2));
00608       EC_CHECK(ecSubMod(params, &t2, &t2, &t5));
00609 
00610       //Check whether Tz != 1
00611       if(mpiCompInt(&t->z, 1))
00612       {
00613          //Compute t3 = t3 * t6
00614          EC_CHECK(ecMulMod(params, &t3, &t3, &t6));
00615       }
00616 
00617       //Compute t3 = t3 * t4
00618       EC_CHECK(ecMulMod(params, &t3, &t3, &t4));
00619       //Compute t7 = t4^2
00620       EC_CHECK(ecSqrMod(params, &t7, &t4));
00621       //Compute t4 = t4 * t7
00622       EC_CHECK(ecMulMod(params, &t4, &t4, &t7));
00623       //Compute t7 = t1 * t7
00624       EC_CHECK(ecMulMod(params, &t7, &t1, &t7));
00625       //Compute t1 = t5^2
00626       EC_CHECK(ecSqrMod(params, &t1, &t5));
00627       //Compute t1 = t1 - t7
00628       EC_CHECK(ecSubMod(params, &t1, &t1, &t7));
00629       //Compute t7 = t7 - 2 * t1
00630       EC_CHECK(ecAddMod(params, &t6, &t1, &t1));
00631       EC_CHECK(ecSubMod(params, &t7, &t7, &t6));
00632       //Compute t5 = t5 * t7
00633       EC_CHECK(ecMulMod(params, &t5, &t5, &t7));
00634       //Compute t4 = t2 * t4
00635       EC_CHECK(ecMulMod(params, &t4, &t2, &t4));
00636       //Compute t2 = t5 - t4
00637       EC_CHECK(ecSubMod(params, &t2, &t5, &t4));
00638 
00639       //Compute t2 = t2 / 2
00640       if(mpiIsEven(&t2))
00641       {
00642          MPI_CHECK(mpiShiftRight(&t2, 1));
00643       }
00644       else
00645       {
00646          MPI_CHECK(mpiAdd(&t2, &t2, &params->p));
00647          MPI_CHECK(mpiShiftRight(&t2, 1));
00648       }
00649 
00650       //Set Rx = t1
00651       MPI_CHECK(mpiCopy(&r->x, &t1));
00652       //Set Ry = t2
00653       MPI_CHECK(mpiCopy(&r->y, &t2));
00654       //Set Rz = t3
00655       MPI_CHECK(mpiCopy(&r->z, &t3));
00656    }
00657 
00658 end:
00659    //Release multiple precision integers
00660    mpiFree(&t1);
00661    mpiFree(&t2);
00662    mpiFree(&t3);
00663    mpiFree(&t4);
00664    mpiFree(&t5);
00665    mpiFree(&t6);
00666    mpiFree(&t7);
00667 
00668    //Return status code
00669    return error;
00670 }
00671 
00672 
00673 /**
00674  * @brief Point addition
00675  * @param[in] params EC domain parameters
00676  * @param[out] r Resulting point R = S + T
00677  * @param[in] s First operand
00678  * @param[in] t Second operand
00679  * @return Error code
00680  **/
00681 
00682 error_t ecFullAdd(const EcDomainParameters *params, EcPoint *r, const EcPoint *s, const EcPoint *t)
00683 {
00684    error_t error;
00685 
00686    //Check whether Sz == 0
00687    if(!mpiCompInt(&s->z, 0))
00688    {
00689       //Set R = T
00690       MPI_CHECK(mpiCopy(&r->x, &t->x));
00691       MPI_CHECK(mpiCopy(&r->y, &t->y));
00692       MPI_CHECK(mpiCopy(&r->z, &t->z));
00693    }
00694    //Check whether Tz == 0
00695    else if(!mpiCompInt(&t->z, 0))
00696    {
00697       //Set R = S
00698       MPI_CHECK(mpiCopy(&r->x, &s->x));
00699       MPI_CHECK(mpiCopy(&r->y, &s->y));
00700       MPI_CHECK(mpiCopy(&r->z, &s->z));
00701    }
00702    else
00703    {
00704       //Compute R = S + T
00705       EC_CHECK(ecAdd(params, r, s, t));
00706 
00707       //Check whether R == (0, 0, 0)
00708       if(!mpiCompInt(&r->x, 0) && !mpiCompInt(&r->y, 0) && !mpiCompInt(&r->z, 0))
00709       {
00710          //Compute R = 2 * S
00711          EC_CHECK(ecDouble(params, r, s));
00712       }
00713    }
00714 
00715 end:
00716    //Return status code
00717    return error;
00718 }
00719 
00720 
00721 /**
00722  * @brief Point subtraction
00723  * @param[in] params EC domain parameters
00724  * @param[out] r Resulting point R = S - T
00725  * @param[in] s First operand
00726  * @param[in] t Second operand
00727  * @return Error code
00728  **/
00729 
00730 error_t ecFullSub(const EcDomainParameters *params, EcPoint *r, const EcPoint *s, const EcPoint *t)
00731 {
00732    error_t error;
00733    EcPoint u;
00734 
00735    //Initialize EC point
00736    ecInit(&u);
00737 
00738    //Set Ux = Tx and Uz = Tz
00739    MPI_CHECK(mpiCopy(&u.x, &t->x));
00740    MPI_CHECK(mpiCopy(&u.z, &t->z));
00741    //Set Uy = p - Ty
00742    MPI_CHECK(mpiSub(&u.y, &params->p, &t->y));
00743 
00744    //Compute R = S + U
00745    EC_CHECK(ecFullAdd(params, r, s, &u));
00746 
00747 end:
00748    //Release EC point
00749    ecFree(&u);
00750 
00751    //Return status code
00752    return error;
00753 }
00754 
00755 
00756 /**
00757  * @brief Scalar multiplication
00758  * @param[in] params EC domain parameters
00759  * @param[out] r Resulting point R = d.S
00760  * @param[in] d An integer d such as 0 <= d < p
00761  * @param[in] s EC point
00762  * @return Error code
00763  **/
00764 
00765 error_t ecMult(const EcDomainParameters *params, EcPoint *r, const Mpi *d, const EcPoint *s)
00766 {
00767    error_t error;
00768    uint_t i;
00769    Mpi h;
00770 
00771    //Initialize multiple precision integer
00772    mpiInit(&h);
00773 
00774    //Check whether d == 0
00775    if(!mpiCompInt(d, 0))
00776    {
00777       //Set R = (1, 1, 0)
00778       MPI_CHECK(mpiSetValue(&r->x, 1));
00779       MPI_CHECK(mpiSetValue(&r->y, 1));
00780       MPI_CHECK(mpiSetValue(&r->z, 0));
00781    }
00782    //Check whether d == 1
00783    else if(!mpiCompInt(d, 1))
00784    {
00785       //Set R = S
00786       MPI_CHECK(mpiCopy(&r->x, &s->x));
00787       MPI_CHECK(mpiCopy(&r->y, &s->y));
00788       MPI_CHECK(mpiCopy(&r->z, &s->z));
00789    }
00790    //Check whether Sz == 0
00791    else if(!mpiCompInt(&s->z, 0))
00792    {
00793       //Set R = (1, 1, 0)
00794       MPI_CHECK(mpiSetValue(&r->x, 1));
00795       MPI_CHECK(mpiSetValue(&r->y, 1));
00796       MPI_CHECK(mpiSetValue(&r->z, 0));
00797    }
00798    else
00799    {
00800       //Check whether Sz != 1
00801       if(mpiCompInt(&s->z, 1))
00802       {
00803          //Normalize S
00804          EC_CHECK(ecAffinify(params, r, s));
00805          EC_CHECK(ecProjectify(params, r, r));
00806       }
00807       else
00808       {
00809          //Set R = S
00810          MPI_CHECK(mpiCopy(&r->x, &s->x));
00811          MPI_CHECK(mpiCopy(&r->y, &s->y));
00812          MPI_CHECK(mpiCopy(&r->z, &s->z));
00813       }
00814 
00815 //Left-to-right binary method
00816 #if 0
00817       for(i = mpiGetBitLength(d) - 1; i >= 1; i--)
00818       {
00819          //Point doubling
00820          EC_CHECK(ecDouble(params, r, r));
00821 
00822          if(mpiGetBitValue(d, i - 1))
00823          {
00824             //Compute R = R + S
00825             EC_CHECK(ecFullAdd(params, r, r, s));
00826          }
00827       }
00828 //Fast left-to-right binary method
00829 #else
00830       //Precompute h = 3 * d
00831       MPI_CHECK(mpiAdd(&h, d, d));
00832       MPI_CHECK(mpiAdd(&h, &h, d));
00833 
00834       //Scalar multiplication
00835       for(i = mpiGetBitLength(&h) - 2; i >= 1; i--)
00836       {
00837          //Point doubling
00838          EC_CHECK(ecDouble(params, r, r));
00839 
00840          //Check whether h(i) == 1 and k(i) == 0
00841          if(mpiGetBitValue(&h, i) && !mpiGetBitValue(d, i))
00842          {
00843             //Compute R = R + S
00844             EC_CHECK(ecFullAdd(params, r, r, s));
00845          }
00846          //Check whether h(i) == 0 and k(i) == 1
00847          else if(!mpiGetBitValue(&h, i) && mpiGetBitValue(d, i))
00848          {
00849             //Compute R = R - S
00850             EC_CHECK(ecFullSub(params, r, r, s));
00851          }
00852       }
00853 #endif
00854    }
00855 
00856 end:
00857    //Release multiple precision integer
00858    mpiFree(&h);
00859 
00860    //Return status code
00861    return error;
00862 }
00863 
00864 
00865 /**
00866  * @brief An auxiliary function for the twin multiplication
00867  * @param[in] t An integer T such as 0 <= T <= 31
00868  * @return Output value
00869  **/
00870 
00871 uint_t ecTwinMultF(uint_t t)
00872 {
00873    if(18 <= t && t < 22)
00874       return 9;
00875    else if(14 <= t && t < 18)
00876       return 10;
00877    else if(22 <= t && t < 24)
00878       return 11;
00879    else if(4 <= t && t < 12)
00880       return 14;
00881    else
00882       return 12;
00883 }
00884 
00885 
00886 /**
00887  * @brief Twin multiplication
00888  * @param[in] params EC domain parameters
00889  * @param[out] r Resulting point R = d0.S + d1.T
00890  * @param[in] d0 An integer d such as 0 <= d0 < p
00891  * @param[in] s EC point
00892  * @param[in] d1 An integer d such as 0 <= d1 < p
00893  * @param[in] t EC point
00894  * @return Error code
00895  **/
00896 
00897 error_t ecTwinMult(const EcDomainParameters *params, EcPoint *r,
00898    const Mpi *d0, const EcPoint *s, const Mpi *d1, const EcPoint *t)
00899 {
00900    error_t error;
00901    int_t k;
00902    uint_t m;
00903    uint_t m0;
00904    uint_t m1;
00905    uint_t c0;
00906    uint_t c1;
00907    uint_t h0;
00908    uint_t h1;
00909    int_t u0;
00910    int_t u1;
00911    EcPoint spt;
00912    EcPoint smt;
00913 
00914    //Initialize EC points
00915    ecInit(&spt);
00916    ecInit(&smt);
00917 
00918    //Precompute SpT = S + T
00919    EC_CHECK(ecFullAdd(params, &spt, s, t));
00920    //Precompute SmT = S - T
00921    EC_CHECK(ecFullSub(params, &smt, s, t));
00922 
00923    //Let m0 be the bit length of d0
00924    m0 = mpiGetBitLength(d0);
00925    //Let m1 be the bit length of d1
00926    m1 = mpiGetBitLength(d1);
00927    //Let m = MAX(m0, m1)
00928    m = MAX(m0, m1);
00929 
00930    //Let c be a 2 x 6 binary matrix
00931    c0 = mpiGetBitValue(d0, m - 4);
00932    c0 |= mpiGetBitValue(d0, m - 3) << 1;
00933    c0 |= mpiGetBitValue(d0, m - 2) << 2;
00934    c0 |= mpiGetBitValue(d0, m - 1) << 3;
00935    c1 = mpiGetBitValue(d1, m - 4);
00936    c1 |= mpiGetBitValue(d1, m - 3) << 1;
00937    c1 |= mpiGetBitValue(d1, m - 2) << 2;
00938    c1 |= mpiGetBitValue(d1, m - 1) << 3;
00939 
00940    //Set R = (1, 1, 0)
00941    MPI_CHECK(mpiSetValue(&r->x, 1));
00942    MPI_CHECK(mpiSetValue(&r->y, 1));
00943    MPI_CHECK(mpiSetValue(&r->z, 0));
00944 
00945    //Calculate both multiplications at the same time
00946    for(k = m; k >= 0; k--)
00947    {
00948       //Compute h(0) = 16 * c(0,1) + 8 * c(0,2) + 4 * c(0,3) + 2 * c(0,4) + c(0,5)
00949       h0 = c0 & 0x1F;
00950       //Check whether c(0,0) == 1
00951       if(c0 & 0x20)
00952          h0 = 31 - h0;
00953 
00954       //Compute h(1) = 16 * c(1,1) + 8 * c(1,2) + 4 * c(1,3) + 2 * c(1,4) + c(1,5)
00955       h1 = c1 & 0x1F;
00956       //Check whether c(1,0) == 1
00957       if(c1 & 0x20)
00958          h1 = 31 - h1;
00959 
00960       //Compute u(0)
00961       if(h0 < ecTwinMultF(h1))
00962          u0 = 0;
00963       else if(c0 & 0x20)
00964          u0 = -1;
00965       else
00966          u0 = 1;
00967 
00968       //Compute u(1)
00969       if(h1 < ecTwinMultF(h0))
00970          u1 = 0;
00971       else if(c1 & 0x20)
00972          u1 = -1;
00973       else
00974          u1 = 1;
00975 
00976       //Update c matrix
00977       c0 <<= 1;
00978       c0 |= mpiGetBitValue(d0, k - 5);
00979       c0 ^= u0 ? 0x20 : 0x00;
00980       c1 <<= 1;
00981       c1 |= mpiGetBitValue(d1, k - 5);
00982       c1 ^= u1 ? 0x20 : 0x00;
00983 
00984       //Point doubling
00985       EC_CHECK(ecDouble(params, r, r));
00986 
00987       //Check u(0) and u(1)
00988       if(u0 == -1 && u1 == -1)
00989       {
00990          //Compute R = R - SpT
00991          EC_CHECK(ecFullSub(params, r, r, &spt));
00992       }
00993       else if(u0 == -1 && u1 == 0)
00994       {
00995          //Compute R = R - S
00996          EC_CHECK(ecFullSub(params, r, r, s));
00997       }
00998       else if(u0 == -1 && u1 == 1)
00999       {
01000          //Compute R = R - SmT
01001          EC_CHECK(ecFullSub(params, r, r, &smt));
01002       }
01003       else if(u0 == 0 && u1 == -1)
01004       {
01005          //Compute R = R - T
01006          EC_CHECK(ecFullSub(params, r, r, t));
01007       }
01008       else if(u0 == 0 && u1 == 1)
01009       {
01010          //Compute R = R + T
01011          EC_CHECK(ecFullAdd(params, r, r, t));
01012       }
01013       else if(u0 == 1 && u1 == -1)
01014       {
01015          //Compute R = R + SmT
01016          EC_CHECK(ecFullAdd(params, r, r, &smt));
01017       }
01018       else if(u0 == 1 && u1 == 0)
01019       {
01020          //Compute R = R + S
01021          EC_CHECK(ecFullAdd(params, r, r, s));
01022       }
01023       else if(u0 == 1 && u1 == 1)
01024       {
01025          //Compute R = R + SpT
01026          EC_CHECK(ecFullAdd(params, r, r, &spt));
01027       }
01028    }
01029 
01030 end:
01031    //Release EC points
01032    ecFree(&spt);
01033    ecFree(&smt);
01034 
01035    //Return status code
01036    return error;
01037 }
01038 
01039 
01040 /**
01041  * @brief Fast modular addition
01042  * @param[in] params EC domain parameters
01043  * @param[out] r Resulting integer R = (A + B) mod p
01044  * @param[in] a An integer such as 0 <= A < p
01045  * @param[in] b An integer such as 0 <= B < p
01046  * @return Error code
01047  **/
01048 
01049 error_t ecAddMod(const EcDomainParameters *params, Mpi *r, const Mpi *a, const Mpi *b)
01050 {
01051    error_t error;
01052 
01053    //Compute R = A + B
01054    MPI_CHECK(mpiAdd(r, a, b));
01055 
01056    //Compute R = (A + B) mod p
01057    if(mpiComp(r, &params->p) >= 0)
01058    {
01059        MPI_CHECK(mpiSub(r, r, &params->p));
01060    }
01061 
01062 end:
01063    //Return status code
01064    return error;
01065 }
01066 
01067 
01068 /**
01069  * @brief Fast modular subtraction
01070  * @param[in] params EC domain parameters
01071  * @param[out] r Resulting integer R = (A - B) mod p
01072  * @param[in] a An integer such as 0 <= A < p
01073  * @param[in] b An integer such as 0 <= B < p
01074  * @return Error code
01075  **/
01076 
01077 error_t ecSubMod(const EcDomainParameters *params, Mpi *r, const Mpi *a, const Mpi *b)
01078 {
01079    error_t error;
01080 
01081    //Compute R = A - B
01082    MPI_CHECK(mpiSub(r, a, b));
01083 
01084    //Compute R = (A - B) mod p
01085    if(mpiCompInt(r, 0) < 0)
01086    {
01087        MPI_CHECK(mpiAdd(r, r, &params->p));
01088    }
01089 
01090 end:
01091    //Return status code
01092    return error;
01093 }
01094 
01095 
01096 /**
01097  * @brief Fast modular multiplication
01098  * @param[in] params EC domain parameters
01099  * @param[out] r Resulting integer R = (A * B) mod p
01100  * @param[in] a An integer such as 0 <= A < p
01101  * @param[in] b An integer such as 0 <= B < p
01102  * @return Error code
01103  **/
01104 
01105 error_t ecMulMod(const EcDomainParameters *params, Mpi *r, const Mpi *a, const Mpi *b)
01106 {
01107    error_t error;
01108 
01109    //Compute R = A * B
01110    MPI_CHECK(mpiMul(r, a, b));
01111 
01112    //Compute R = (A * B) mod p
01113    if(params->mod != NULL)
01114    {
01115       MPI_CHECK(params->mod(r, &params->p));
01116    }
01117    else
01118    {
01119       MPI_CHECK(mpiMod(r, r, &params->p));
01120    }
01121 
01122 end:
01123    //Return status code
01124    return error;
01125 }
01126 
01127 
01128 /**
01129  * @brief Fast modular squaring
01130  * @param[in] params EC domain parameters
01131  * @param[out] r Resulting integer R = (A ^ 2) mod p
01132  * @param[in] a An integer such as 0 <= A < p
01133  * @return Error code
01134  **/
01135 
01136 error_t ecSqrMod(const EcDomainParameters *params, Mpi *r, const Mpi *a)
01137 {
01138    error_t error;
01139 
01140    //Compute R = A ^ 2
01141    MPI_CHECK(mpiMul(r, a, a));
01142 
01143    //Compute R = (A ^ 2) mod p
01144    if(params->mod != NULL)
01145    {
01146       MPI_CHECK(params->mod(r, &params->p));
01147    }
01148    else
01149    {
01150       MPI_CHECK(mpiMod(r, r, &params->p));
01151    }
01152 
01153 end:
01154    //Return status code
01155    return error;
01156 }
01157 
01158 #endif
01159