source: git/libpolys/coeffs/modulop.h @ 2a4a23

spielwiese
Last change on this file since 2a4a23 was 2a4a23, checked in by Hans Schoenemann <hannes@…>, 6 years ago
fix: Z/p-arithmetic for all architectures
  • Property mode set to 100644
File size: 4.8 KB
Line 
1#ifndef MODULOP_H
2#define MODULOP_H
3/****************************************
4*  Computer Algebra System SINGULAR     *
5****************************************/
6/*
7* ABSTRACT: numbers modulo p (<=32749)
8*/
9#include "misc/auxiliary.h"
10
11
12// define if a*b is with mod instead of tables
13//#define HAVE_GENERIC_MULT
14// define if 1/b is from  tables
15//#define HAVE_INVTABLE
16// define if an if should be used
17//#define HAVE_GENERIC_ADD
18
19//#undef HAVE_GENERIC_ADD
20//#undef HAVE_GENERIC_MULT
21//#undef HAVE_INVTABLE
22
23//#define HAVE_GENERIC_ADD 1
24//#define HAVE_GENERIC_MULT 1
25//#define HAVE_INVTABLE 1
26
27// enable large primes (32749 < p < 2^31-)
28#define NV_OPS
29#define NV_MAX_PRIME 32749
30#define FACTORY_MAX_PRIME 536870909
31
32struct n_Procs_s; typedef struct  n_Procs_s  *coeffs;
33struct snumber; typedef struct snumber *   number;
34
35BOOLEAN npInitChar(coeffs r, void* p);
36
37// inline number npMultM(number a, number b, int npPrimeM)
38// // return (a*b)%n
39// {
40//    double ab;
41//    long q, res;
42//
43//    ab = ((double) ((int)a)) * ((double) ((int)b));
44//    q  = (long) (ab/((double) npPrimeM));  // q could be off by (+/-) 1
45//    res = (long) (ab - ((double) q)*((double) npPrimeM));
46//    res += (res >> 31) & npPrimeM;
47//    res -= npPrimeM;
48//    res += (res >> 31) & npPrimeM;
49//    return (number)res;
50// }
51#ifdef HAVE_GENERIC_MULT
52static inline number npMultM(number a, number b, const coeffs r)
53{
54  return (number)
55    ((((unsigned long) a)*((unsigned long) b)) % ((unsigned long) r->ch));
56}
57static inline void npInpMultM(number &a, number b, const coeffs r)
58{
59  a=(number)
60    ((((unsigned long) a)*((unsigned long) b)) % ((unsigned long) r->ch));
61}
62#else
63static inline number npMultM(number a, number b, const coeffs r)
64{
65  long x = (long)r->npLogTable[(long)a]+ r->npLogTable[(long)b];
66  #ifdef HAVE_GENERIC_ADD
67  if (x>r->npPminus1M) x-=r->npPminus1M;
68  #else
69    x-=r->npPminus1M;
70    #if SIZEOF_LONG == 8
71     x += (x >> 63) & r->npPminus1M;
72    #else
73     x += (x >> 31) & r->npPminus1M;
74    #endif
75  #endif
76  return (number)(long)r->npExpTable[x];
77}
78static inline void npInpMultM(number &a, number b, const coeffs r)
79{
80  long x = (long)r->npLogTable[(long)a]+ r->npLogTable[(long)b];
81  #ifdef HAVE_GENERIC_ADD
82  if (x>r->npPminus1M) x-=r->npPminus1M;
83  #else
84    x-=r->npPminus1M;
85    #if SIZEOF_LONG == 8
86     x += (x >> 63) & r->npPminus1M;
87    #else
88     x += (x >> 31) & r->npPminus1M;
89    #endif
90  #endif
91  a=(number)(long)r->npExpTable[x];
92}
93#endif
94
95#if 0
96inline number npAddAsm(number a, number b, int m)
97{
98  number r;
99    asm ("addl %2, %1; cmpl %3, %1; jb 0f; subl %3, %1; 0:"
100         : "=&r" (r)
101         : "%0" (a), "g" (b), "g" (m)
102         : "cc");
103  return r;
104}
105inline number npSubAsm(number a, number b, int m)
106{
107  number r;
108  asm ("subl %2, %1; jnc 0f; addl %3, %1; 0:"
109        : "=&r" (r)
110        : "%0" (a), "g" (b), "g" (m)
111        : "cc");
112  return r;
113}
114#endif
115#ifdef HAVE_GENERIC_ADD
116static inline number npAddM(number a, number b, const coeffs r)
117{
118  unsigned long R = (unsigned long)a + (unsigned long)b;
119  return (number)(R >= r->ch ? R - r->ch : R);
120}
121static inline void npInpAddM(number &a, number b, const coeffs r)
122{
123  unsigned long R = (unsigned long)a + (unsigned long)b;
124  a=(number)(R >= r->ch ? R - r->ch : R);
125}
126static inline number npSubM(number a, number b, const coeffs r)
127{
128  return (number)((long)a<(long)b ?
129                       r->ch-(long)b+(long)a : (long)a-(long)b);
130}
131#else
132static inline number npAddM(number a, number b, const coeffs r)
133{
134   unsigned long res = (long)((unsigned long)a + (unsigned long)b);
135   res -= r->ch;
136#if SIZEOF_LONG == 8
137   res += ((long)res >> 63) & r->ch;
138#else
139   res += ((long)res >> 31) & r->ch;
140#endif
141   return (number)res;
142}
143static inline void npInpAddM(number &a, number b, const coeffs r)
144{
145   unsigned long res = (long)((unsigned long)a + (unsigned long)b);
146   res -= r->ch;
147#if SIZEOF_LONG == 8
148   res += ((long)res >> 63) & r->ch;
149#else
150   res += ((long)res >> 31) & r->ch;
151#endif
152   a=(number)res;
153}
154static inline number npSubM(number a, number b, const coeffs r)
155{
156   long res = ((long)a - (long)b);
157#if SIZEOF_LONG == 8
158   res += (res >> 63) & r->ch;
159#else
160   res += (res >> 31) & r->ch;
161#endif
162   return (number)res;
163}
164#endif
165
166static inline number npNegM(number a, const coeffs r)
167{
168  return (number)((long)(r->ch)-(long)(a));
169}
170
171static inline BOOLEAN npIsZeroM (number  a, const coeffs)
172{
173  return 0 == (long)a;
174}
175
176// inline number npMultM(number a, number b, int npPrimeM)
177// {
178//   return (number)(((long)a*(long)b) % npPrimeM);
179// }
180
181// The folloing is reused inside gnumpc.cc, gnumpfl.cc and longrat.cc
182long    npInt         (number &n, const coeffs r);
183
184// The following is currently used in OPAE.cc, OPAEQ.cc and OPAEp.cc for setting their SetMap...
185nMapFunc npSetMap(const coeffs src, const coeffs dst); // FIXME! BUG?
186
187#define npEqualM(A,B,r)  ((A)==(B))
188
189
190#endif
Note: See TracBrowser for help on using the repository browser.