My Project
kChinese.cc
Go to the documentation of this file.
1 #include "misc/auxiliary.h"
2 
3 #include "misc/intvec.h"
4 #include "misc/options.h"
6 #include "polys/matpol.h"
7 #include "polys/simpleideals.h"
8 #include "coeffs/longrat.h"
9 #include "Singular/feOpt.h"
10 #include "kernel/polys.h"
11 #include "kernel/mod2.h"
12 
13 #ifdef HAVE_VSPACE
14 
15 #define mpz_isNeg(A) ((A)->_mp_size<0)
16 number nlRInit (long i);
17 
19 #include "kernel/ideals.h"
20 #include <sys/types.h>
21 #include <sys/wait.h>
22 
23 // send a number via a string s
24 static char *send_number(char *s, number n)
25 {
26  long *d=(long*)s;
27  if (SR_HDL(n) & SR_INT)
28  {
29  *d=(long)n;
30  s+=SIZEOF_LONG;
31  }
32  else
33  {
34  *d=n->s*2;/* n->s in 0..3: 0..6, use +8 for neg. numbers */
35  s+=SIZEOF_LONG;
36  if (mpz_isNeg(n->z)) { *d+=8; mpz_neg(n->z,n->z); }
37  size_t l;
38  d=(long*)s;
39  s+=SIZEOF_LONG;
40  mpz_export(s,&l,-1,sizeof(mp_limb_t),0,0,n->z);
41  *d=l;
42  s+=l*sizeof(mp_limb_t);
43  if (n->s!=3)
44  {
45  d=(long*)s;
46  s+=SIZEOF_LONG;
47  mpz_export(s,&l,-1,sizeof(mp_limb_t),0,0,n->n);
48  *d=l;
49  s+=l*sizeof(mp_limb_t);
50  }
51  }
52  return s;
53 }
54 
55 static char * get_number(char *s, number *n)
56 {
57  // format: last bit 1: imm. number (long)
58  // otherwise: 0,2: size(long) mpz, size(long) mpz
59  // 6: size(long) mpz
60  // 8,10: size(long) -mpz, size(long) mpz
61  // 14: size(long) -mpz
62  long *d=(long*)s;
63  s+=SIZEOF_LONG;
64  if (((*d)&1)==1) // immidiate number
65  {
66  *n=(number)(*d);
67  }
68  else
69  {
70  *n=nlRInit(0);
71  BOOLEAN neg=(*d>=8);
72  if (neg) *d-=8;
73  (*n)->s=(*d)/2;
74  d=(long*)s;
75  s+=SIZEOF_LONG;
76  size_t l=*d;
77  mpz_realloc2((*n)->z,l*sizeof(mp_limb_t)*8);
78  mpz_import((*n)->z,l,-1,sizeof(mp_limb_t),0,0,s);
79  if (neg) mpz_neg((*n)->z,(*n)->z);
80  s+=l*sizeof(mp_limb_t);
81  if ((*n)->s!=3)
82  {
83  d=(long*)s;
84  s+=SIZEOF_LONG;
85  l=*d;
86  mpz_init2((*n)->n,l*sizeof(mp_limb_t)*8);
87  mpz_import((*n)->n,l,-1,sizeof(mp_limb_t),0,0,s);
88  s+=l*sizeof(mp_limb_t);
89  }
90  }
91  return s;
92 }
93 
94 static long size_number(number n)
95 {
96  long ll=SIZEOF_LONG;
97  if (SR_HDL(n) & SR_INT)
98  {
99  return SIZEOF_LONG;
100  }
101  else
102  {
103  if (n->s==3)
104  {
105  ll+=SIZEOF_LONG*2; /* n->s, mpz size */
106  long l=mpz_size1(n->z);
107  ll+=l*sizeof(mp_limb_t);
108  }
109  else
110  {
111  ll+=SIZEOF_LONG*3; /* n->s, mpz size(n->z) mpz size(n->n)*/
112  size_t l=mpz_size1(n->z);
113  ll+=l*sizeof(mp_limb_t);
114  l=mpz_size1(n->n);
115  ll+=l*sizeof(mp_limb_t);
116  }
117  }
118  return ll;
119 }
120 
121 static char* send_mon(char *s, poly m, const ring r)
122 {
123  // format: number exp[0..r->ExpL_Size]
124  s=send_number(s,p_GetCoeff(m,r));
125  memcpy(s,m->exp,r->ExpL_Size*sizeof(long));
126  s+=r->ExpL_Size*sizeof(long);
127  return s;
128 }
129 
130 static char* get_mon(char *s, poly *m, const ring r)
131 {
132  (*m)=p_Init(r);
133  s=get_number(s,&p_GetCoeff(*m,r));
134  memcpy((*m)->exp,s,r->ExpL_Size*sizeof(long));
135  s+=r->ExpL_Size*sizeof(long);
136  return s;
137 }
138 
139 static long size_mon(poly m, const ring r)
140 {
141  long ll=size_number(p_GetCoeff(m,r));
142  ll+=r->ExpL_Size*sizeof(long);
143  return ll;
144 }
145 
146 static char* send_poly(char *s, int ind, poly p, const ring r)
147 {
148  // format: index(long) length(long) mon...
149  //p_Write(p,r);PrintLn();
150  long *d=(long*)s;
151  *d=ind;
152  s+=SIZEOF_LONG;
153  long l=pLength(p);
154  d=(long*)s;
155  *d=l;
156  s+=SIZEOF_LONG;
157  while(p!=NULL)
158  {
159  s=send_mon(s,p,r);
160  pIter(p);
161  }
162  return s;
163 }
164 
165 static char* get_poly(char *s,int &ind, poly *p,const ring r)
166 {
167  long *d=(long*)s;
168  ind=*d;
169  s+=SIZEOF_LONG;
170  d=(long*)s;
171  long l=*d;
172  s+=SIZEOF_LONG;
173  for(long i=0;i<l;i++)
174  {
175  poly m;
176  s=get_mon(s,&m,r);
177  pNext(m)=*p;
178  *p=m;
179  }
180  *p=pReverse(*p);
181  return s;
182 }
183 
184 static long size_poly(poly p, const ring r)
185 {
186  long l=SIZEOF_LONG*2;
187  while(p!=NULL)
188  {
189  l+=size_mon(p,r);
190  pIter(p);
191  }
192  return l;
193 }
194 
195 ideal id_ChineseRemainder_0(ideal *xx, number *q, int rl, const ring r)
196 {
197  int cnt=0;int rw=0; int cl=0;
198  // find max. size of xx[.]:
199  for(int j=rl-1;j>=0;j--)
200  {
201  int i=IDELEMS(xx[j])*xx[j]->nrows;
202  if (i>cnt) cnt=i;
203  if (xx[j]->nrows >rw) rw=xx[j]->nrows; // for lifting matrices
204  if (xx[j]->ncols >cl) cl=xx[j]->ncols; // for lifting matrices
205  }
206  if (rw*cl !=cnt)
207  {
208  WerrorS("format mismatch in CRT");
209  return NULL;
210  }
211  int cpus=(int)(long)feOptValue(FE_OPT_CPUS);
212  if ((cpus==1) || (2*cpus>=cnt))
213  /* at least 2 polys for each process, or switch to seriell version */
214  return id_ChineseRemainder(xx,q,rl,r);
215  ideal result=idInit(cnt,xx[0]->rank);
216  result->nrows=rw; // for lifting matrices
217  result->ncols=cl; // for lifting matrices
218  int parent_pid=getpid();
219  using namespace vspace;
220  vmem_init();
221  // Create a queue of int
222  VRef<Queue<int> > queue = vnew<Queue<int> >();
223  for(int i=cnt-1;i>=0; i--)
224  {
225  queue->enqueue(i); // the tasks: construct poly p[i]
226  }
227  for(int i=cpus;i>=0;i--)
228  {
229  queue->enqueue(-1); // stop sign, one for each child
230  }
231  // Create a queue of polys
232  VRef<Queue<VRef<VString> > > rqueue = vnew<Queue<VRef<VString> > >();
233  for (int i=0;i<cpus;i++)
234  {
235  int pid = fork_process();
236  if (pid==0) break; //child
237  }
238  if (parent_pid!=getpid()) // child ------------------------------------------
239  {
240  number *x=(number *)omAlloc(rl*sizeof(number));
241  poly *p=(poly *)omAlloc(rl*sizeof(poly));
242  CFArray inv_cache(rl);
244  n_SwitchChinRem=1;
245  loop
246  {
247  int ind=queue->dequeue();
248  if (ind== -1)
249  {
250  exit(0);
251  }
252 
253  for(int j=rl-1;j>=0;j--)
254  {
255  if(ind>=IDELEMS(xx[j])*xx[j]->nrows) // out of range of this ideal
256  p[j]=NULL;
257  else
258  p[j]=xx[j]->m[ind];
259  }
260  poly res=p_ChineseRemainder(p,x,q,rl,inv_cache,r);
261  long l=size_poly(res,r);
262  //printf("size: %ld kB\n",(l+1023)/1024);
263  VRef<VString> msg = vstring(l+1);
264  char *s=(char*)msg->str();
265  send_poly(s,ind,res,r);
266  rqueue->enqueue(msg);
267  if (TEST_OPT_PROT) printf(".");
268  }
269  }
270  else // parent ---------------------------------------------------
271  {
272  if (TEST_OPT_PROT) printf("%d childs created\n",cpus);
273  VRef<VString> msg;
274  while(cnt>0)
275  {
276  msg=rqueue->dequeue();
277  char *s=(char*)msg->str();
278  int ind;
279  poly p=NULL;
280  get_poly(s,ind,&p,r);
281  //printf("got res[%d]\n",ind);
282  result->m[ind]=p;
283  msg.free();
284  cnt--;
285  }
286  // removes queues
287  queue.free();
288  rqueue.free();
289  vmem_deinit();
290  }
291  return result;
292 }
293 
294 ideal id_Farey_0(ideal x, number N, const ring r)
295 {
296  int cnt=IDELEMS(x)*x->nrows;
297  int cpus=(int)(long)feOptValue(FE_OPT_CPUS);
298  if (2*cpus>=cnt) /* at least 2 polys for each process,
299  or switch to seriell version */
300  return id_Farey(x,N,r);
301  ideal result=idInit(cnt,x->rank);
302  result->nrows=x->nrows; // for lifting matrices
303  result->ncols=x->ncols; // for lifting matrices
304 
305  int parent_pid=getpid();
306  using namespace vspace;
307  vmem_init();
308  // Create a queue of int
309  VRef<Queue<int> > queue = vnew<Queue<int> >();
310  for(int i=cnt-1;i>=0; i--)
311  {
312  queue->enqueue(i); // the tasks: construct poly p[i]
313  }
314  for(int i=cpus;i>=0;i--)
315  {
316  queue->enqueue(-1); // stop sign, one for each child
317  }
318  // Create a queue of polys
319  VRef<Queue<VRef<VString> > > rqueue = vnew<Queue<VRef<VString> > >();
320  for (int i=0;i<cpus;i++)
321  {
322  int pid = fork_process();
323  if (pid==0) break; //child
324  }
325  if (parent_pid!=getpid()) // child ------------------------------------------
326  {
327  loop
328  {
329  int ind=queue->dequeue();
330  if (ind== -1)
331  {
332  exit(0);
333  }
334 
335  poly res=p_Farey(x->m[ind],N,r);
336  long l=size_poly(res,r);
337  VRef<VString> msg = vstring(l+1);
338  char *s=(char*)msg->str();
339  send_poly(s,ind,res,r);
340  rqueue->enqueue(msg);
341  if (TEST_OPT_PROT) printf(".");
342  }
343  }
344  else // parent ---------------------------------------------------
345  {
346  if (TEST_OPT_PROT) printf("%d childs created\n",cpus);
347  VRef<VString> msg;
348  while(cnt>0)
349  {
350  msg=rqueue->dequeue();
351  char *s=(char*)msg->str();
352  int ind;
353  poly p=NULL;
354  get_poly(s,ind,&p,r);
355  //printf("got res[%d]\n",ind);
356  result->m[ind]=p;
357  msg.free();
358  cnt--;
359  }
360  // removes queues
361  queue.free();
362  rqueue.free();
363  vmem_deinit();
364  }
365  return result;
366 }
367 
368 void test_n(poly n)
369 {
370  p_Write(n,currRing);
371  char *buf=(char*)omAlloc0(2048*1000);
372  int ll=size_poly(n,currRing);
373  printf("size: %d\n",ll);
374  char *s=send_poly(buf,12345,n,currRing);
375  printf("send len: %d\n",(int)(s-buf));
376  long *d=(long*)buf;
377  for(int i=0;i<=ll/SIZEOF_LONG;i++) printf("%ld ",d[i]);
378  printf("\n");
379  n=NULL;
380  s=get_poly(buf,ll,&n,currRing);
381  printf("read len: %d\n",(int)(s-buf));
382  Print(":index: %d\n",ll);
383  p_Write(n,currRing);
384  PrintLn();
385  omFree(buf);
386 }
387 #endif
All the auxiliary stuff.
int BOOLEAN
Definition: auxiliary.h:87
const CanonicalForm CFMap CFMap & N
Definition: cfEzgcd.cc:56
int l
Definition: cfEzgcd.cc:100
int m
Definition: cfEzgcd.cc:128
int i
Definition: cfEzgcd.cc:132
Variable x
Definition: cfModGcd.cc:4084
int p
Definition: cfModGcd.cc:4080
cl
Definition: cfModGcd.cc:4102
int int ncols
Definition: cf_linsys.cc:32
int nrows
Definition: cf_linsys.cc:32
#define Print
Definition: emacs.cc:80
return result
Definition: facAbsBiFact.cc:75
const CanonicalForm int s
Definition: facAbsFact.cc:51
CanonicalForm res
Definition: facAbsFact.cc:60
int j
Definition: facHensel.cc:110
void WerrorS(const char *s)
Definition: feFopen.cc:24
static void * feOptValue(feOptIndex opt)
Definition: feOpt.h:40
#define EXTERN_VAR
Definition: globaldefs.h:6
ideal id_Farey(ideal x, number N, const ring r)
Definition: ideals.cc:2836
static char * get_number(char *s, number *n)
Definition: kChinese.cc:55
#define mpz_isNeg(A)
Definition: kChinese.cc:15
static char * get_poly(char *s, int &ind, poly *p, const ring r)
Definition: kChinese.cc:165
static long size_poly(poly p, const ring r)
Definition: kChinese.cc:184
static char * send_poly(char *s, int ind, poly p, const ring r)
Definition: kChinese.cc:146
ideal id_Farey_0(ideal x, number N, const ring r)
Definition: kChinese.cc:294
void test_n(poly n)
Definition: kChinese.cc:368
static long size_number(number n)
Definition: kChinese.cc:94
number nlRInit(long i)
Definition: longrat.cc:2490
static char * send_mon(char *s, poly m, const ring r)
Definition: kChinese.cc:121
ideal id_ChineseRemainder_0(ideal *xx, number *q, int rl, const ring r)
Definition: kChinese.cc:195
static char * send_number(char *s, number n)
Definition: kChinese.cc:24
static char * get_mon(char *s, poly *m, const ring r)
Definition: kChinese.cc:130
static long size_mon(poly m, const ring r)
Definition: kChinese.cc:139
poly p_ChineseRemainder(poly *xx, mpz_ptr *x, mpz_ptr *q, int rl, mpz_ptr *C, const ring R)
VAR int n_SwitchChinRem
Definition: longrat.cc:3052
#define SR_INT
Definition: longrat.h:67
#define pIter(p)
Definition: monomials.h:37
#define pNext(p)
Definition: monomials.h:36
#define p_GetCoeff(p, r)
Definition: monomials.h:50
Definition: vspace.cc:9
pid_t fork_process()
Definition: vspace.cc:419
static void vmem_deinit()
Definition: vspace.h:407
static VRef< VString > vstring(const char *s)
Definition: vspace.h:766
static Status vmem_init()
Definition: vspace.h:403
#define omAlloc(size)
Definition: omAllocDecl.h:210
#define omFree(addr)
Definition: omAllocDecl.h:261
#define omAlloc0(size)
Definition: omAllocDecl.h:211
#define NULL
Definition: omList.c:12
#define TEST_OPT_PROT
Definition: options.h:103
poly p_Farey(poly p, number N, const ring r)
Definition: p_polys.cc:54
void p_Write(poly p, ring lmRing, ring tailRing)
Definition: polys0.cc:342
static poly pReverse(poly p)
Definition: p_polys.h:335
static unsigned pLength(poly a)
Definition: p_polys.h:191
static poly p_Init(const ring r, omBin bin)
Definition: p_polys.h:1280
VAR ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
Definition: polys.cc:13
Compatiblity layer for legacy polynomial operations (over currRing)
void PrintLn()
Definition: reporter.cc:310
#define mpz_size1(A)
Definition: si_gmp.h:12
int status int void * buf
Definition: si_signals.h:59
ideal idInit(int idsize, int rank)
initialise an ideal / module
Definition: simpleideals.cc:35
ideal id_ChineseRemainder(ideal *xx, number *q, int rl, const ring r)
#define IDELEMS(i)
Definition: simpleideals.h:23
#define loop
Definition: structs.h:80
void free()
Definition: vspace.h:470
#define SR_HDL(A)
Definition: tgb.cc:35