mirror of
https://github.com/lcn2/calc.git
synced 2025-08-16 01:03:29 +03:00
a31078bbec
Some folks might think: “you still use RCS”?!? And we will say, hey, at least we switched from SCCS to RCS back in … I think it was around 1994 ... at least we are keeping up! :-) :-) :-) Logs say that SCCS version 18 became RCS version 19 on 1994 March 18. RCS served us well. But now it is time to move on. And so we are switching to git. Calc releases produce a lot of file changes. In the 125 releases of calc since 1996, when I started managing calc releases, there have been 15473 file mods!
273 lines
6.4 KiB
C
273 lines
6.4 KiB
C
/*
|
|
* c_pmodm127 - calculate q mod 2^(2^127-1)
|
|
*
|
|
* Copyright (C) 2004-2007 Landon Curt Noll
|
|
*
|
|
* Calc is open software; you can redistribute it and/or modify it under
|
|
* the terms of the version 2.1 of the GNU Lesser General Public License
|
|
* as published by the Free Software Foundation.
|
|
*
|
|
* Calc is distributed in the hope that it will be useful, but WITHOUT
|
|
* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
|
|
* or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General
|
|
* Public License for more details.
|
|
*
|
|
* A copy of version 2.1 of the GNU Lesser General Public License is
|
|
* distributed with calc under the filename COPYING-LGPL. You should have
|
|
* received a copy with calc; if not, write to Free Software Foundation, Inc.
|
|
* 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
|
|
*
|
|
* Under source code control: 2004/07/28 22:12:25
|
|
* File existed as early as: 2004
|
|
*
|
|
* Share and enjoy! :-) http://www.isthe.com/chongo/tech/comp/calc/
|
|
*/
|
|
|
|
|
|
#if defined(CUSTOM)
|
|
|
|
#include <stdio.h>
|
|
|
|
#include "have_const.h"
|
|
#include "value.h"
|
|
#include "custom.h"
|
|
#include "zmath.h"
|
|
|
|
#include "have_unused.h"
|
|
|
|
/* 2^255 */
|
|
STATIC HALF h255[] = {
|
|
#if BASEB == 32
|
|
(HALF)0x00000000, (HALF)0x00000000, (HALF)0x00000000, (HALF)0x00000000,
|
|
(HALF)0x00000000, (HALF)0x00000000, (HALF)0x00000000, (HALF)0x80000000
|
|
#else /* BASEB == 32 */
|
|
(HALF)0x0000, (HALF)0x0000, (HALF)0x0000, (HALF)0x0000,
|
|
(HALF)0x0000, (HALF)0x0000, (HALF)0x0000, (HALF)0x0000,
|
|
(HALF)0x0000, (HALF)0x0000, (HALF)0x0000, (HALF)0x0000,
|
|
(HALF)0x0000, (HALF)0x0000, (HALF)0x0000, (HALF)0x8000
|
|
#endif /* BASEB == 32 */
|
|
};
|
|
ZVALUE p255 = {
|
|
h255, 8, 0
|
|
};
|
|
|
|
|
|
/* static declarations */
|
|
S_FUNC void zmod5_or_zmod(ZVALUE *zp);
|
|
STATIC BOOL havelastmod = FALSE;
|
|
STATIC ZVALUE lastmod[1];
|
|
STATIC ZVALUE lastmodinv[1];
|
|
|
|
|
|
/*
|
|
* c_pmodm127 - calculate q mod 2^(2^127-1)
|
|
*
|
|
* given:
|
|
* count = 1;
|
|
* vals[0] real number; (q - potential factor)
|
|
*
|
|
* returns:
|
|
* result real number; (q mod 2^(2^127-1))
|
|
*/
|
|
/*ARGSUSED*/
|
|
VALUE
|
|
c_pmodm127(char UNUSED *name, int UNUSED count, VALUE **vals)
|
|
{
|
|
VALUE result; /* what we will return */
|
|
ZVALUE q; /* test factor */
|
|
ZVALUE temp; /* temp calculation value */
|
|
int i; /* exponent value */
|
|
|
|
/*
|
|
* arg check
|
|
*/
|
|
result.v_type = V_NULL;
|
|
if (vals[0]->v_type != V_NUM) {
|
|
math_error("Non-numeric argument for pmodm127");
|
|
/*NOTREACHED*/
|
|
}
|
|
if (qisfrac(vals[0]->v_num)) {
|
|
math_error("Non-integer argument for pmodm127");
|
|
/*NOTREACHED*/
|
|
}
|
|
if (qisneg(vals[0]->v_num) || qiszero(vals[0]->v_num)) {
|
|
math_error("argument for pmodm127 <= 0");
|
|
/*NOTREACHED*/
|
|
}
|
|
|
|
/*
|
|
* look at the numerator
|
|
*/
|
|
q = vals[0]->v_num->num;
|
|
|
|
/*
|
|
* setup lastmod with q
|
|
*/
|
|
if (havelastmod && zcmp(q, *lastmod)) {
|
|
zfree(*lastmod);
|
|
zfree(*lastmodinv);
|
|
havelastmod = FALSE;
|
|
}
|
|
if (!havelastmod) {
|
|
zcopy(q, lastmod);
|
|
zbitvalue(2 * q.len * BASEB, &temp);
|
|
zquo(temp, q, lastmodinv, 0);
|
|
zfree(temp);
|
|
havelastmod = TRUE;
|
|
}
|
|
|
|
/*
|
|
* start with 2^255
|
|
*/
|
|
result.v_num = qalloc();
|
|
zcopy(p255, &result.v_num->num);
|
|
result.v_type = V_NUM;
|
|
|
|
/*
|
|
* compute 2^(2^127-1) mod q by modular exponentation
|
|
*
|
|
* We implement the following calc code in C:
|
|
*
|
|
* (* given q, our test factor, as the arg to this function *)
|
|
* result = 2^255;
|
|
* for (i=8; i < 127; ++i) {
|
|
* result %= q; (* mod *)
|
|
* result *= result; (* square *)
|
|
* result <<= 1; (* times 2 *)
|
|
* }
|
|
* result %= q; (* result is now 2^(2^127-1) % q *)
|
|
*/
|
|
for (i=8; i<127; ++i) {
|
|
#if 0
|
|
/* use of zmod is a bit slower than zmod5_or_zmod */
|
|
(void) zmod(result.v_num->num, *lastmod, &temp, 0);
|
|
zfree(result.v_num->num);
|
|
result.v_num->num = temp;
|
|
#else
|
|
zmod5_or_zmod(&result.v_num->num); /* mod */
|
|
#endif
|
|
#if 0
|
|
/* use of zmul is slightly slower than zsquare */
|
|
zmul(result.v_num->num, result.v_num->num, &temp); /* square */
|
|
#else
|
|
zsquare(result.v_num->num, &temp); /* square */
|
|
#endif
|
|
/*
|
|
* We could manually shift here, but this would o speed
|
|
* up the operation only a very tiny bit at the expense
|
|
* of a bunch of special code.
|
|
*/
|
|
zfree(result.v_num->num);
|
|
zshift(temp, 1, &result.v_num->num); /* times 2 */
|
|
zfree(temp);
|
|
}
|
|
zmod5_or_zmod(&result.v_num->num); /* result = 2^(2^127-1) % q */
|
|
|
|
/*
|
|
* cleanup and return result
|
|
*/
|
|
return result;
|
|
}
|
|
|
|
|
|
/*
|
|
* zmod5_or_zmod - fast integer modulo the modulus computation
|
|
*
|
|
* "borrowed" from ../zmod.c
|
|
*
|
|
* Given the address of a positive integer whose word count does not
|
|
* exceed twice that of the modulus stored at lastmod, to evaluate and store
|
|
* at that address the value of the integer modulo the modulus.
|
|
*
|
|
* Unlike the static routine in ../zmod.c, we will call the zmod and return
|
|
* the result of the zmod5_or_zmod conditions do not apply to the argument
|
|
* and saved mod.
|
|
*/
|
|
S_FUNC void
|
|
zmod5_or_zmod(ZVALUE *zp)
|
|
{
|
|
LEN len, modlen, j;
|
|
ZVALUE tmp1, tmp2;
|
|
ZVALUE z1, z2, z3;
|
|
HALF *a, *b;
|
|
FULL f;
|
|
HALF u;
|
|
|
|
int subcount = 0;
|
|
|
|
if (zrel(*zp, *lastmod) < 0)
|
|
return;
|
|
modlen = lastmod->len;
|
|
len = zp->len;
|
|
z1.v = zp->v + modlen - 1;
|
|
z1.len = len - modlen + 1;
|
|
z1.sign = z2.sign = z3.sign = 0;
|
|
if (z1.len > modlen + 1) {
|
|
/* in ../zmod.c we did a math_error("Bad call to zmod5!!!"); */
|
|
/* here we just call the slower zmod and return the result */
|
|
(void) zmod(*zp, *lastmod, &tmp1, 0);
|
|
/* replace zp with tmp1 mod result */
|
|
zfree(*zp);
|
|
*zp = tmp1;
|
|
return;
|
|
}
|
|
z2.v = lastmodinv->v + modlen + 1 - z1.len;
|
|
z2.len = lastmodinv->len - modlen - 1 + z1.len;
|
|
zmul(z1, z2, &tmp1);
|
|
z3.v = tmp1.v + z1.len;
|
|
z3.len = tmp1.len - z1.len;
|
|
if (z3.len > 0) {
|
|
zmul(z3, *lastmod, &tmp2);
|
|
j = modlen;
|
|
a = zp->v;
|
|
b = tmp2.v;
|
|
u = 0;
|
|
len = modlen;
|
|
while (j-- > 0) {
|
|
f = (FULL) *a - (FULL) *b++ - (FULL) u;
|
|
*a++ = (HALF) f;
|
|
u = - (HALF) (f >> BASEB);
|
|
}
|
|
if (z1.len > 1) {
|
|
len++;
|
|
if (tmp2.len > modlen)
|
|
f = (FULL) *a - (FULL) *b - (FULL) u;
|
|
else
|
|
f = (FULL) *a - (FULL) u;
|
|
*a++ = (HALF) f;
|
|
}
|
|
while (len > 0 && *--a == 0)
|
|
len--;
|
|
zp->len = len;
|
|
zfree(tmp2);
|
|
}
|
|
zfree(tmp1);
|
|
while (len > 0 && zrel(*zp, *lastmod) >= 0) {
|
|
subcount++;
|
|
if (subcount > 2) {
|
|
math_error("Too many subtractions in zmod5_or_zmod");
|
|
/*NOTREACHED*/
|
|
}
|
|
j = modlen;
|
|
a = zp->v;
|
|
b = lastmod->v;
|
|
u = 0;
|
|
while (j-- > 0) {
|
|
f = (FULL) *a - (FULL) *b++ - (FULL) u;
|
|
*a++ = (HALF) f;
|
|
u = - (HALF) (f >> BASEB);
|
|
}
|
|
if (len > modlen) {
|
|
f = (FULL) *a - (FULL) u;
|
|
*a++ = (HALF) f;
|
|
}
|
|
while (len > 0 && *--a == 0)
|
|
len--;
|
|
zp->len = len;
|
|
}
|
|
if (len == 0)
|
|
zp->len = 1;
|
|
}
|
|
|
|
#endif /* CUSTOM */
|