source: git/factory/ffops.cc @ fbb0173

spielwiese
Last change on this file since fbb0173 was 650f2d8, checked in by Mohamed Barakat <mohamed.barakat@…>, 13 years ago
renamed assert.h -> cf_assert.h in factory
  • Property mode set to 100644
File size: 1.8 KB
RevLine 
[493c477]1/* emacs edit mode for this file is -*- C++ -*- */
[341696]2/* $Id$ */
[2dd068]3
[84250a6]4#include <config.h>
5
[071a57]6#include <string.h>
7
[650f2d8]8#include "cf_assert.h"
[84250a6]9
[2dd068]10#include "cf_defs.h"
11#include "ffops.h"
12
[a43bd1]13int ff_prime = 0;
14int ff_halfprime = 0;
[02af91]15bool ff_big = false;
[cfa45e]16short * ff_invtab = new short [32767];
[2dd068]17
18void ff_setprime ( const int p )
19{
[0eb89c]20    if ( p != ff_prime ) {
[806c18]21        ff_prime = p;
22        ff_halfprime = ff_prime / 2;
23        if ( ! ff_big )
24            memset( ff_invtab, 0, ff_prime*sizeof(short) );
[0eb89c]25    }
[2dd068]26}
27
28int ff_newinv ( const int a )
29{
[b7c8bd]30    int p, q, r1, r2, y1, y2;
31    if (a < 2)
[806c18]32          return (ff_invtab[a] = a);
[b7c8bd]33    r1 = p = ff_prime;
34    q = r1 / a;
35    y1 = -q;
36    r1 -= a * q;
37    if (r1 == 1)
38    {
39      y1 += p;
40      ff_invtab[y1] = a;
41      return (ff_invtab[a] = y1);
42    }
43    r2 = a;
44    y2 = 1;
45    for (;;)
46    {
47      q = r2 / r1;
48      y2 -= y1 * q;
49      r2 -= r1 * q;
50      if (r2 == 1)
51      {
52        if (y2 < 0)
53          y2 += p;
54        ff_invtab[y2] = a;
55        return (ff_invtab[a] = y2);
56      }
57      q = r1 / r2;
58      y1 -= y2 * q;
59      r1 -= r2 * q;
60      if (r1 == 1)
61      {
62        if (y1 < 0)
63          y1 += p;
64        ff_invtab[y1] = a;
65        return (ff_invtab[a] = y1);
66      }
[02af91]67    }
68}
69
70int ff_biginv ( const int a )
71{
[b7c8bd]72    int p, q, r1, r2, y1, y2;
73    if (a < 2)
74      return a;
75    r1 = p = ff_prime;
76    q = r1 / a;
77    y1 = -q;
78    r1 -= a * q;
79    if (r1 == 1)
80      return p + y1;
81    r2 = a;
82    y2 = 1;
83    for (;;)
84    {
85      q = r2 / r1;
86      y2 -= y1 * q;
87      r2 -= r1 * q;
88      if (r2 == 1)
89      {
90        if (y2 > 0)
91          return y2;
92        else
93          return p + y2;
94      }
95      q = r1 / r2;
96      y1 -= y2 * q;
97      r1 -= r2 * q;
98      if (r1 == 1)
99      {
100        if (y1 > 0)
101          return y1;
102        else
103          return p + y1;
104      }
[2dd068]105    }
106}
Note: See TracBrowser for help on using the repository browser.