My Project
Loading...
Searching...
No Matches
maps.cc
Go to the documentation of this file.
1/****************************************
2* Computer Algebra System SINGULAR *
3****************************************/
4/*
5* ABSTRACT - the mapping of polynomials to other rings
6*/
7
8#include "misc/auxiliary.h"
9#include "misc/options.h"
10
11#include "coeffs/coeffs.h"
12#include "coeffs/numbers.h"
13
16#include "polys/simpleideals.h"
17#include "polys/prCopy.h"
19
20#ifdef HAVE_PLURAL
21#include "polys/nc/nc.h"
22#endif
23
24// This is a very dirty way to "normalize" numbers w.r.t. a
25// MinPoly
26
27#define MAX_MAP_DEG 128
28
29/*2
30* copy a map
31*/
32map maCopy(map theMap, const ring r)
33{
34 int i;
35 map m=(map)idInit(IDELEMS(theMap),0);
36 for (i=IDELEMS(theMap)-1; i>=0; i--)
37 m->m[i] = p_Copy(theMap->m[i],r);
38 m->preimage=omStrDup(theMap->preimage);
39 return m;
40}
41
42
43/*2
44* return the image of var(v)^pExp, where var(v) maps to p
45*/
46poly maEvalVariable(poly p, int v,int pExp, ideal s, const ring dst_r)
47{
48 if (pExp==1)
49 return p_Copy(p,dst_r);
50
51 poly res;
52
53 if((s!=NULL)&&(pExp<MAX_MAP_DEG))
54 {
55 int j=2;
56 poly p0=p;
57 // find starting point
58 if(MATELEM(s,v,1)==NULL)
59 {
60 MATELEM(s,v,1)=p_Copy(p/*theMap->m[v-1]*/,dst_r);
61 }
62 else
63 {
64 while((j<=pExp)&&(MATELEM(s,v,j)!=NULL))
65 {
66 j++;
67 }
68 p0=MATELEM(s,v,j-1);
69 }
70 // multiply
71 for(;j<=pExp;j++)
72 {
73 p0=MATELEM(s,v,j)=pp_Mult_qq(p0, p,dst_r);
74 p_Normalize(p0, dst_r);
75 }
76 res=p_Copy(p0/*MATELEM(s,v,pExp)*/,dst_r);
77 }
78 else //if ((p->next!=NULL)&&(p->next->next==NULL))
79 {
80 res=p_Power(p_Copy(p,dst_r),pExp,dst_r);
81 }
82 return res;
83}
84
85static poly maEvalMonom(map theMap, poly p,ring preimage_r, ideal s,
86 nMapFunc nMap, const ring dst_r)
87{
88 p_Test(p,preimage_r);
89 poly q=p_NSet(nMap(pGetCoeff(p),preimage_r->cf,dst_r->cf),dst_r);
90
91 int i;
92 for(i=1;i<=preimage_r->N; i++)
93 {
94 int pExp=p_GetExp( p,i,preimage_r);
95 if (pExp != 0)
96 {
97 if (theMap->m[i-1]!=NULL)
98 {
99 poly p1=theMap->m[i-1];
100 poly pp=maEvalVariable(p1,i,pExp,s,dst_r);
101 q = p_Mult_q(q,pp,dst_r);
102 }
103 else
104 {
105 p_Delete(&q,dst_r);
106 break;
107 }
108 }
109 }
110 int modulComp = p_GetComp( p,preimage_r);
111 if (q!=NULL) p_SetCompP(q,modulComp,dst_r);
112 return q;
113}
114
115poly maEval(map theMap, poly p,ring preimage_r,nMapFunc nMap, ideal s, const ring dst_r)
116{
117 poly result = NULL;
118 int i;
119
120// for(i=1; i<=preimage_r->N; i++)
121// {
122// pTest(theMap->m[i-1]);
123// }
124// while (p!=NULL)
125// {
126// poly q=maEvalMonom(theMap,p,preimage_r,s);
127// result = pAdd(result,q);
128// pIter(p);
129// }
130 if (p!=NULL)
131 {
132 int l = pLength(p)-1;
133 poly* monoms;
134 if (l>0)
135 {
136 monoms = (poly*) omAlloc(l*sizeof(poly));
137
138 for (i=0; i<l; i++)
139 {
140 monoms[i]=maEvalMonom(theMap,p,preimage_r,s, nMap, dst_r);
141 pIter(p);
142 }
143 }
144 result=maEvalMonom(theMap,p,preimage_r,s, nMap, dst_r);
145 if (l>0)
146 {
147 for(i = l-1; i>=0; i--)
148 {
149 result=p_Add_q(result, monoms[i], dst_r);
150 }
151 omFreeSize((ADDRESS)monoms,l*sizeof(poly));
152 }
153
154 assume(dst_r != NULL);
155 assume(dst_r->cf != NULL);
156
157 if (nCoeff_is_algExt(dst_r->cf))
159 }
160 return result;
161}
162
163void maFindPerm(char const * const * const preim_names, int preim_n, char const * const * const preim_par, int preim_p,
164 char const * const * const names, int n, char const * const * const par, int nop,
165 int * perm, int *par_perm, n_coeffType ch)
166{
167 int i,j;
168 /* find correspondig vars */
169 for (i=0; i<preim_n; i++)
170 {
171 for(j=0; j<n; j++)
172 {
173 if (strcmp(preim_names[i],names[j])==0)
174 {
175 if (BVERBOSE(V_IMAP))
176 Print("// var %s: nr %d -> nr %d\n",preim_names[i],i+1,j+1);
177 /* var i+1 from preimage ring is var j+1 (index j+1) from image ring */
178 perm[i+1]=j+1;
179 break;
180 }
181 }
182 if ((perm[i+1]==0)&&(par!=NULL)
183 // do not consider par of Fq
184 && (ch!=n_GF))
185 {
186 for(j=0; j<nop; j++)
187 {
188 if (strcmp(preim_names[i],par[j])==0)
189 {
190 if (BVERBOSE(V_IMAP))
191 Print("// var %s: nr %d -> par %d\n",preim_names[i],i+1,j+1);
192 /* var i+1 from preimage ring is par j+1 (index j) from image ring */
193 perm[i+1]=-(j+1);
194 }
195 }
196 }
197 }
198 if (par_perm!=NULL)
199 {
200 for (i=0; i<preim_p; i++)
201 {
202 for(j=0; j<n; j++)
203 {
204 if (strcmp(preim_par[i],names[j])==0)
205 {
206 if (BVERBOSE(V_IMAP))
207 Print("// par %s: par %d -> nr %d\n",preim_par[i],i+1,j+1);
208 /*par i+1 from preimage ring is var j+1 (index j+1) from image ring*/
209 par_perm[i]=j+1;
210 break;
211 }
212 }
213 if ((par!=NULL) && (par_perm[i]==0))
214 {
215 for(j=0; j<nop; j++)
216 {
217 if (strcmp(preim_par[i],par[j])==0)
218 {
219 if (BVERBOSE(V_IMAP))
220 Print("// par %s: nr %d -> par %d\n",preim_par[i],i+1,j+1);
221 /*par i+1 from preimage ring is par j+1 (index j) from image ring */
222 par_perm[i]=-(j+1);
223 }
224 }
225 }
226 }
227 }
228}
229
230#ifdef HAVE_SHIFTBBA
231void maFindPermLP(char const * const * const preim_names, int preim_n, char const * const * const preim_par, int preim_p,
232 char const * const * const names, int n, char const * const * const par, int nop,
233 int * perm, int *par_perm, n_coeffType ch, int lV)
234{
235 int i,j,b;
236 /* find correspondig vars */
237 for (b=0;b<preim_n/lV;b++)
238 {
239 for (i=b*lV; i<(b+1)*lV; i++)
240 {
241 int cnt=0;
242 for(j=0; j<n; j++)
243 {
244 if (strcmp(preim_names[i],names[j])==0)
245 {
246 if (cnt==b)
247 {
248 if (BVERBOSE(V_IMAP))
249 Print("// var %s: nr %d -> nr %d\n",preim_names[i],i+1,j+1);
250 /* var i+1 from preimage ring is var j+1 (index j+1) from image ring */
251 perm[i+1]=j+1;
252 break;
253 }
254 else cnt++;
255 }
256 }
257 if ((perm[i+1]==0)&&(par!=NULL)
258 // do not consider par of Fq
259 && (ch!=n_GF))
260 {
261 for(j=0; j<nop; j++)
262 {
263 if (strcmp(preim_names[i],par[j])==0)
264 {
265 if (BVERBOSE(V_IMAP))
266 Print("// var %s: nr %d -> par %d\n",preim_names[i],i+1,j+1);
267 /* var i+1 from preimage ring is par j+1 (index j) from image ring */
268 perm[i+1]=-(j+1);
269 }
270 }
271 }
272 }
273 }
274 if (par_perm!=NULL)
275 {
276 for (i=0; i<preim_p; i++)
277 {
278 for(j=0; j<n; j++)
279 {
280 if (strcmp(preim_par[i],names[j])==0)
281 {
282 if (BVERBOSE(V_IMAP))
283 Print("// par %s: par %d -> nr %d\n",preim_par[i],i+1,j+1);
284 /*par i+1 from preimage ring is var j+1 (index j+1) from image ring*/
285 par_perm[i]=j+1;
286 break;
287 }
288 }
289 if ((par!=NULL) && (par_perm[i]==0))
290 {
291 for(j=0; j<nop; j++)
292 {
293 if (strcmp(preim_par[i],par[j])==0)
294 {
295 if (BVERBOSE(V_IMAP))
296 Print("// par %s: nr %d -> par %d\n",preim_par[i],i+1,j+1);
297 /*par i+1 from preimage ring is par j+1 (index j) from image ring */
298 par_perm[i]=-(j+1);
299 }
300 }
301 }
302 }
303 }
304}
305
306void maFetchPermLP(const ring preimage_r, const ring dst_r, int * perm)
307/* perm=(int *)omAlloc0(preimage_r->N+1 * sizeof(int)); */
308{
309 for (int i = 0; i < preimage_r->N + 1; i++) { perm[i] = 0; }
310
311 int preimage_lV = preimage_r->isLPring;
312 int r_lV = dst_r->isLPring;
313
314 int preimage_ncgens = preimage_r->LPncGenCount;
315 int r_ncges = dst_r->LPncGenCount;
316
317 int preimage_vars = preimage_lV - preimage_ncgens;
318 int r_vars = r_lV - r_ncges;
319
320 // for each block
321 for (int i = 0; i < si_min(preimage_r->N / preimage_lV, dst_r->N / r_lV); i++)
322 {
323 // align variables
324 for (int j = 1; j <= si_min(preimage_vars, r_vars); j++)
325 {
326 perm[(i * preimage_lV) + j] = (i * r_lV) + j;
327 }
328
329 // align ncgens
330 for (int j = 1; j <= si_min(preimage_ncgens, r_ncges); j++)
331 {
332 perm[(i * preimage_lV) + preimage_vars + j] = (i * r_lV) + r_vars + j;
333 }
334 }
335}
336#endif
337
338/*2
339* embeds poly p from the subring r into the current ring
340*/
341poly maIMap(ring r, poly p, const ring dst_r)
342{
343 /* the simplest case:*/
344 if(r==dst_r) return p_Copy(p,dst_r);
345 nMapFunc nMap=n_SetMap(r->cf,dst_r->cf);
346 int *perm=(int *)omAlloc0((r->N+1)*sizeof(int));
347 //int *par_perm=(int *)omAlloc0(rPar(r)*sizeof(int));
348 maFindPerm(r->names, rVar(r), rParameter(r), rPar(r),
349 dst_r->names, rVar(dst_r),rParameter(dst_r), rPar(dst_r),
350 perm,NULL, dst_r->cf->type);
351 poly res=p_PermPoly(p,perm,r,dst_r, nMap /*,par_perm,rPar(r)*/);
352 omFreeSize((ADDRESS)perm,(r->N+1)*sizeof(int));
353 //omFreeSize((ADDRESS)par_perm,rPar(r)*sizeof(int));
354 return res;
355}
356
357/*3
358* find the max. degree in one variable, but not larger than MAX_MAP_DEG
359*/
360int maMaxDeg_Ma(ideal a,ring preimage_r)
361{
362 int i,j;
363 int N = preimage_r->N;
364 poly p;
365 int *m=(int *)omAlloc0(N*sizeof(int));
366
367 for (i=MATROWS(a)*MATCOLS(a)-1;i>=0;i--)
368 {
369 p=a->m[i];
370 //pTest(p); // cannot test p because it is from another ring
371 while(p!=NULL)
372 {
373 for(j=N-1;j>=0;j--)
374 {
375 m[j]=si_max(m[j],(int)p_GetExp( p,j+1,preimage_r));
376 if (m[j]>=MAX_MAP_DEG)
377 {
379 goto max_deg_fertig_id;
380 }
381 }
382 pIter(p);
383 }
384 }
385 i=m[0];
386 for(j=N-1;j>0;j--)
387 {
388 i=si_max(i,m[j]);
389 }
390max_deg_fertig_id:
391 omFreeSize((ADDRESS)m,N*sizeof(int));
392 return i;
393}
394
395/*3
396* find the max. degree in one variable, but not larger than MAX_MAP_DEG
397*/
398int maMaxDeg_P(poly p,ring preimage_r)
399{
400 int i,j;
401 int N = preimage_r->N;
402 int *m=(int *)omAlloc0(N*sizeof(int));
403
404// pTest(p);
405 while(p!=NULL)
406 {
407 for(j=N-1;j>=0;j--)
408 {
409 m[j]=si_max(m[j],(int)p_GetExp(p,j+1,preimage_r));
410 if (m[j]>=MAX_MAP_DEG)
411 {
413 goto max_deg_fertig_p;
414 }
415 }
416 pIter(p);
417 }
418 i=m[0];
419 for(j=N-1;j>0;j--)
420 {
421 i=si_max(i,m[j]);
422 }
423max_deg_fertig_p:
424 omFreeSize((ADDRESS)m,N*sizeof(int));
425 return i;
426}
427
428// This is a very dirty way to cancel monoms whose number equals the
429// MinPoly
430poly p_MinPolyNormalize(poly p, const ring r)
431{
432 const coeffs C = r->cf;
433 number one = n_Init(1, C);
434 spolyrec rp;
435
436 poly q = &rp;
437
438 while (p != NULL)
439 {
440 // this returns 0, if p == MinPoly
441 number product = n_Mult(p_GetCoeff(p, r), one, C);
442 if ((product == NULL)||(n_IsZero(product, C)))
443 {
444 p_LmDelete(&p, r);
445 }
446 else
447 {
448 p_SetCoeff(p, product, r);
449 pNext(q) = p;
450 q = p;
451 p = pNext(p);
452 }
453 }
454 pNext(q) = NULL;
455 n_Delete(&one, C);
456 return rp.next;
457}
All the auxiliary stuff.
static int si_max(const int a, const int b)
Definition: auxiliary.h:124
void * ADDRESS
Definition: auxiliary.h:119
static int si_min(const int a, const int b)
Definition: auxiliary.h:125
CanonicalForm FACTORY_PUBLIC pp(const CanonicalForm &)
CanonicalForm pp ( const CanonicalForm & f )
Definition: cf_gcd.cc:676
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
int p
Definition: cfModGcd.cc:4078
CanonicalForm b
Definition: cfModGcd.cc:4103
CanonicalForm map(const CanonicalForm &primElem, const Variable &alpha, const CanonicalForm &F, const Variable &beta)
map from to such that is mapped onto
Definition: cf_map_ext.cc:504
Coefficient rings, fields and other domains suitable for Singular polynomials.
static FORCE_INLINE number n_Mult(number a, number b, const coeffs r)
return the product of 'a' and 'b', i.e., a*b
Definition: coeffs.h:633
n_coeffType
Definition: coeffs.h:27
@ n_GF
\GF{p^n < 2^16}
Definition: coeffs.h:32
static FORCE_INLINE nMapFunc n_SetMap(const coeffs src, const coeffs dst)
set the mapping function pointers for translating numbers from src to dst
Definition: coeffs.h:697
static FORCE_INLINE BOOLEAN n_IsZero(number n, const coeffs r)
TRUE iff 'n' represents the zero element.
Definition: coeffs.h:461
static FORCE_INLINE void n_Delete(number *p, const coeffs r)
delete 'p'
Definition: coeffs.h:452
static FORCE_INLINE number n_Init(long i, const coeffs r)
a number representing i in the given coeff field/ring r
Definition: coeffs.h:535
static FORCE_INLINE BOOLEAN nCoeff_is_algExt(const coeffs r)
TRUE iff r represents an algebraic extension field.
Definition: coeffs.h:907
number(* nMapFunc)(number a, const coeffs src, const coeffs dst)
maps "a", which lives in src, into dst
Definition: coeffs.h:73
#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
const Variable & v
< [in] a sqrfree bivariate poly
Definition: facBivar.h:39
int j
Definition: facHensel.cc:110
map maCopy(map theMap, const ring r)
Definition: maps.cc:32
poly maEval(map theMap, poly p, ring preimage_r, nMapFunc nMap, ideal s, const ring dst_r)
Definition: maps.cc:115
void maFetchPermLP(const ring preimage_r, const ring dst_r, int *perm)
Definition: maps.cc:306
poly maIMap(ring r, poly p, const ring dst_r)
Definition: maps.cc:341
int maMaxDeg_P(poly p, ring preimage_r)
Definition: maps.cc:398
void maFindPerm(char const *const *const preim_names, int preim_n, char const *const *const preim_par, int preim_p, char const *const *const names, int n, char const *const *const par, int nop, int *perm, int *par_perm, n_coeffType ch)
Definition: maps.cc:163
static poly maEvalMonom(map theMap, poly p, ring preimage_r, ideal s, nMapFunc nMap, const ring dst_r)
Definition: maps.cc:85
void maFindPermLP(char const *const *const preim_names, int preim_n, char const *const *const preim_par, int preim_p, char const *const *const names, int n, char const *const *const par, int nop, int *perm, int *par_perm, n_coeffType ch, int lV)
Definition: maps.cc:231
poly maEvalVariable(poly p, int v, int pExp, ideal s, const ring dst_r)
Definition: maps.cc:46
int maMaxDeg_Ma(ideal a, ring preimage_r)
Definition: maps.cc:360
poly p_MinPolyNormalize(poly p, const ring r)
Definition: maps.cc:430
#define MAX_MAP_DEG
Definition: maps.cc:27
#define MATELEM(mat, i, j)
1-based access to matrix
Definition: matpol.h:29
#define MATROWS(i)
Definition: matpol.h:26
#define MATCOLS(i)
Definition: matpol.h:27
#define assume(x)
Definition: mod2.h:389
#define p_GetComp(p, r)
Definition: monomials.h:64
poly next
Definition: monomials.h:24
#define pIter(p)
Definition: monomials.h:37
#define pNext(p)
Definition: monomials.h:36
#define p_GetCoeff(p, r)
Definition: monomials.h:50
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy
Definition: monomials.h:44
The main handler for Singular numbers which are suitable for Singular polynomials.
#define omStrDup(s)
Definition: omAllocDecl.h:263
#define omFreeSize(addr, size)
Definition: omAllocDecl.h:260
#define omAlloc(size)
Definition: omAllocDecl.h:210
#define omAlloc0(size)
Definition: omAllocDecl.h:211
#define NULL
Definition: omList.c:12
#define BVERBOSE(a)
Definition: options.h:35
#define V_IMAP
Definition: options.h:53
poly p_PermPoly(poly p, const int *perm, const ring oldRing, const ring dst, nMapFunc nMap, const int *par_perm, int OldPar, BOOLEAN use_mult)
Definition: p_polys.cc:4130
poly p_Power(poly p, int i, const ring r)
Definition: p_polys.cc:2197
void p_Normalize(poly p, const ring r)
Definition: p_polys.cc:3813
poly p_NSet(number n, const ring r)
returns the poly representing the number n, destroys n
Definition: p_polys.cc:1473
static int pLength(poly a)
Definition: p_polys.h:188
static poly p_Add_q(poly p, poly q, const ring r)
Definition: p_polys.h:934
static void p_LmDelete(poly p, const ring r)
Definition: p_polys.h:721
static poly p_Mult_q(poly p, poly q, const ring r)
Definition: p_polys.h:1112
static void p_SetCompP(poly p, int i, ring r)
Definition: p_polys.h:252
static number p_SetCoeff(poly p, number n, ring r)
Definition: p_polys.h:410
static long p_GetExp(const poly p, const unsigned long iBitmask, const int VarOffset)
get a single variable exponent @Note: the integer VarOffset encodes:
Definition: p_polys.h:467
static void p_Delete(poly *p, const ring r)
Definition: p_polys.h:899
static poly pp_Mult_qq(poly p, poly q, const ring r)
Definition: p_polys.h:1149
static poly p_Copy(poly p, const ring r)
returns a copy of p
Definition: p_polys.h:844
#define p_Test(p, r)
Definition: p_polys.h:159
static int rPar(const ring r)
(r->cf->P)
Definition: ring.h:599
static char const ** rParameter(const ring r)
(r->cf->parameter)
Definition: ring.h:625
static short rVar(const ring r)
#define rVar(r) (r->N)
Definition: ring.h:592
ideal idInit(int idsize, int rank)
initialise an ideal / module
Definition: simpleideals.cc:35
#define IDELEMS(i)
Definition: simpleideals.h:23