source: git/Singular/maps_ip.cc @ d828d63

spielwiese
Last change on this file since d828d63 was 6ce030f, checked in by Oleksandr Motsak <motsak@…>, 12 years ago
removal of the $Id$ svn tag from everywhere NOTE: the git SHA1 may be used instead (only on special places) NOTE: the libraries Singular/LIB/*.lib still contain the marker due to our current use of svn
  • Property mode set to 100644
File size: 9.8 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/*
5* ABSTRACT - the mapping of polynomials to other rings
6*/
7#define TRANSEXT_PRIVATES
8
9#include "config.h"
10#include <kernel/mod2.h>
11#include <omalloc/omalloc.h>
12
13#include <coeffs/numbers.h>
14#include <coeffs/coeffs.h>
15
16#include <polys/monomials/ring.h>
17#include <polys/monomials/maps.h>
18#include <polys/matpol.h>
19#include <polys/prCopy.h>
20#include <polys/ext_fields/transext.h>
21
22//#include <libpolys/polys/ext_fields/longtrans.h>
23// #include <kernel/longalg.h>
24
25#include <kernel/febase.h>
26#include <kernel/kstd1.h>
27
28#include "maps_ip.h"
29#include "ipid.h"
30
31
32#include "lists.h"
33#include "tok.h"
34
35/* debug output: Tok2Cmdname in maApplyFetch*/
36#include "ipshell.h"
37
38// define this if you want to use the fast_map routine for mapping ideals
39//#define FAST_MAP
40
41#ifdef FAST_MAP
42#include <polys/monomials/maps.h>
43#endif
44
45
46/*2
47* maps the expression w to res,
48* switch what: MAP_CMD: use theMap for mapping, N for preimage ring
49*              //FETCH_CMD: use pOrdPoly for mapping
50*              IMAP_CMD: use perm for mapping, N for preimage ring
51*              default: map only poly-structures,
52*                       use perm and par_perm, N and P,
53*/
54BOOLEAN maApplyFetch(int what,map theMap,leftv res, leftv w, ring preimage_r,
55                     int *perm, int *par_perm, int P, nMapFunc nMap)
56{
57  int i;
58  int N = preimage_r->N;
59#if 0
60  Print("N=%d what=%s ",N,Tok2Cmdname(what));
61  if (perm!=NULL) for(i=1;i<=N;i++) Print("%d -> %d ",i,perm[i]);
62  PrintS("\n");
63  Print("P=%d ",P);
64  if (par_perm!=NULL) for(i=0;i<P;i++) Print("%d -> %d ",i,par_perm[i]);
65  PrintS("\n");
66#endif
67
68  void *data=w->Data();
69  res->rtyp = w->rtyp;
70  switch (w->rtyp)
71  {
72    case NUMBER_CMD:
73      if (P!=0)
74      {
75//        WerrorS("Sorry 'napPermNumber' was lost in the refactoring process (due to Frank): needs to be fixed");
76//        return TRUE;
77#if 1
78// poly n_PermNumber(const number z, const int *par_perm, const int OldPar, const ring src, const ring dst);
79        res->data= (void *) n_PermNumber((number)data, par_perm, P, preimage_r, currRing);
80#endif
81        res->rtyp=POLY_CMD;
82        if (nCoeff_is_Extension(currRing->cf))
83          res->data=(void *)p_MinPolyNormalize((poly)res->data, currRing);
84        pTest((poly) res->data);       
85      }
86      else
87      {
88        assume( nMap != NULL );
89         
90        number a = nMap((number)data, preimage_r->cf, currRing->cf);
91       
92 
93        if (nCoeff_is_Extension(currRing->cf))
94        {
95          n_Normalize(a, currRing->cf); // ???
96/*
97          number a = (number)res->data;
98          number one = nInit(1);
99          number product = nMult(a, one );
100          nDelete(&one);
101          nDelete(&a);
102          res->data=(void *)product;
103 */
104        }
105        #ifdef LDEBUG
106        n_Test(a, currRing->cf);
107        #endif
108        res->data=(void *)a;
109
110      }
111      break;
112    case POLY_CMD:
113    case VECTOR_CMD:
114      if ((what==FETCH_CMD)&& (preimage_r->cf==currRing->cf))
115        res->data=(void *)prCopyR( (poly)data, preimage_r, currRing);
116      else
117      if ((what==IMAP_CMD) || ((what==FETCH_CMD) /* && (nMap!=nCopy)*/))
118        res->data=(void *)p_PermPoly((poly)data,perm,preimage_r,currRing, nMap,par_perm,P);
119      else /*if (what==MAP_CMD)*/
120      {
121        matrix s=mpNew(N,maMaxDeg_P((poly)data, preimage_r));
122        res->data=(void *)maEval(theMap, (poly)data, preimage_r, nMap, (ideal)s, currRing);
123        idDelete((ideal *)&s);
124      }
125      if (nCoeff_is_Extension(currRing->cf))
126        res->data=(void *)p_MinPolyNormalize((poly)res->data, currRing);
127      pTest((poly)res->data);
128      break;
129    case MODUL_CMD:
130    case MATRIX_CMD:
131    case IDEAL_CMD:
132    case MAP_CMD:
133    {
134      int C=((matrix)data)->cols();
135      int R;
136      if (w->rtyp==MAP_CMD) R=1;
137      else R=((matrix)data)->rows();
138      matrix m=mpNew(R,C);
139      char *tmpR=NULL;
140      if(w->rtyp==MAP_CMD)
141      {
142        tmpR=((map)data)->preimage;
143        ((matrix)data)->rank=((matrix)data)->rows();
144      }
145      if ((what==FETCH_CMD)&& (preimage_r->cf == currRing->cf))
146      {
147        for (i=R*C-1;i>=0;i--)
148        {
149          m->m[i]=prCopyR(((ideal)data)->m[i], preimage_r, currRing);
150          pTest(m->m[i]);
151        }
152      }
153      else
154      if ((what==IMAP_CMD) || ((what==FETCH_CMD) /* && (nMap!=nCopy)*/))
155      {
156        for (i=R*C-1;i>=0;i--)
157        {
158          m->m[i]=p_PermPoly(((ideal)data)->m[i],perm,preimage_r,currRing,
159                          nMap,par_perm,P);
160          pTest(m->m[i]);
161        }
162      }
163      else /* if(what==MAP_CMD) */
164      {
165        matrix s=mpNew(N,maMaxDeg_Ma((ideal)data,preimage_r));
166        for (i=R*C-1;i>=0;i--)
167        {
168          m->m[i]=maEval(theMap, ((ideal)data)->m[i], preimage_r, nMap, (ideal)s, currRing);
169          pTest(m->m[i]);
170        }
171        idDelete((ideal *)&s);
172      }
173      if (nCoeff_is_Extension(currRing->cf))
174      {
175        for (i=R*C-1;i>=0;i--)
176        {
177          m->m[i]=p_MinPolyNormalize(m->m[i], currRing);
178          pTest(m->m[i]);
179        }
180      }
181      if(w->rtyp==MAP_CMD)
182      {
183        ((map)data)->preimage=tmpR;
184        ((map)m)->preimage=omStrDup(tmpR);
185      }
186      else
187      {
188        m->rank=((matrix)data)->rank;
189      }
190      res->data=(char *)m;
191      idTest((ideal) m);
192      break;
193    }
194
195    case LIST_CMD:
196    {
197      lists l=(lists)data;
198      lists ml=(lists)omAllocBin(slists_bin);
199      ml->Init(l->nr+1);
200      for(i=0;i<=l->nr;i++)
201      {
202        if (((l->m[i].rtyp>BEGIN_RING)&&(l->m[i].rtyp<END_RING))
203        ||(l->m[i].rtyp==LIST_CMD))
204        {
205          if (maApplyFetch(what,theMap,&ml->m[i],&l->m[i],
206                           preimage_r,perm,par_perm,P,nMap))
207          {
208            ml->Clean();
209            omFreeBin((ADDRESS)ml, slists_bin);
210            res->rtyp=0;
211            return TRUE;
212          }
213        }
214        else
215        {
216          ml->m[i].Copy(&l->m[i]);
217        }
218      }
219      res->data=(char *)ml;
220      break;
221    }
222    default:
223    {
224      return TRUE;
225    }
226  }
227  return FALSE;
228}
229
230/*2
231* substitutes the parameter par (from 1..N) by image,
232* does not destroy p and  image
233*/
234poly pSubstPar(poly p, int par, poly image)
235{
236  ideal theMapI = idInit(rPar(currRing),1);
237  nMapFunc nMap = n_SetMap(currRing->cf->extRing->cf, currRing->cf);
238
239  int i;
240  for(i = rPar(currRing);i>0;i--)
241  {
242    if (i != par)
243      theMapI->m[i-1]= p_NSet(n_Param(i, currRing), currRing);
244    else
245      theMapI->m[i-1] = p_Copy(image, currRing);
246  }
247 
248
249  map theMap=(map)theMapI;
250  theMap->preimage=NULL;
251
252  leftv v=(leftv)omAllocBin(sleftv_bin);
253  sleftv tmpW;
254  poly res=NULL;
255
256  while (p!=NULL)
257  {
258    memset(v,0,sizeof(sleftv));
259
260    number d = n_GetDenom(p_GetCoeff(p, currRing), currRing);
261    if (!n_IsOne (d, currRing->cf)) WarnS("ignoring denominators of coefficients...");
262    n_Delete(&d, currRing);
263
264    number num = n_GetNumerator(p_GetCoeff(p, currRing), currRing); 
265
266    memset(&tmpW,0,sizeof(sleftv));
267    tmpW.rtyp = POLY_CMD;
268    tmpW.data = NUM (num); // a copy of this poly will be used
269     
270    if (maApplyFetch(MAP_CMD,theMap,v,&tmpW,currRing->cf->extRing,NULL,NULL,0,nMap))
271    {
272      WerrorS("map failed");
273      v->data=NULL;
274    }
275    n_Delete(&num, currRing);
276
277    //TODO check for memory leaks
278    poly pp = pHead(p);
279    //PrintS("map:");pWrite(pp);
280    pSetCoeff(pp, nInit(1));
281    //PrintS("->");pWrite((poly)(v->data));
282    poly ppp = pMult((poly)(v->data),pp);
283    //PrintS("->");pWrite(ppp);
284    res=pAdd(res,ppp);
285    pIter(p);
286  }
287  idDelete((ideal *)(&theMap));
288  omFreeBin((ADDRESS)v, sleftv_bin);
289  return res;
290}
291
292/*2
293* substitute the n-th parameter by the poly e in id
294* does not destroy id and e
295*/
296ideal  idSubstPar(ideal id, int n, poly e)
297{
298  int k=MATROWS((matrix)id)*MATCOLS((matrix)id);
299  ideal res=(ideal)mpNew(MATROWS((matrix)id),MATCOLS((matrix)id));
300
301  res->rank = id->rank;
302  for(k--;k>=0;k--)
303  {
304    res->m[k]=pSubstPar(id->m[k],n,e);
305  }
306  return res;
307}
308
309/*2
310* substitutes the variable var (from 1..N) by image,
311* does not destroy p and  image
312*/
313poly pSubstPoly(poly p, int var, poly image)
314{
315  if (p==NULL) return NULL;
316#ifdef HAVE_PLURAL
317  if (rIsPluralRing(currRing))
318  {
319    return pSubst(pCopy(p),var,image);
320  }
321#endif
322  map theMap=(map)idMaxIdeal(1);
323  theMap->preimage=NULL;
324  pDelete(&(theMap->m[var-1]));
325  theMap->m[var-1]=pCopy(image);
326
327  poly res=NULL;
328#ifdef FAST_MAP
329  if (pGetComp(p)==0)
330  {
331    ideal src_id=idInit(1,1);
332    src_id->m[0]=p;
333    ideal res_id=fast_map(src_id,currRing,(ideal)theMap,currRing);
334    res=res_id->m[0];
335    res_id->m[0]=NULL; idDelete(&res_id);
336    src_id->m[0]=NULL; idDelete(&src_id);
337  }
338  else
339#endif
340  {
341    sleftv tmpW;
342    memset(&tmpW,0,sizeof(sleftv));
343    tmpW.rtyp=POLY_CMD;
344    tmpW.data=p;
345    leftv v=(leftv)omAlloc0Bin(sleftv_bin);
346    if (maApplyFetch(MAP_CMD,theMap,v,&tmpW,currRing,NULL,NULL,0,
347                            n_SetMap(currRing->cf, currRing->cf)))
348    {
349      WerrorS("map failed");
350      v->data=NULL;
351    }
352    res=(poly)(v->data);
353    omFreeBin((ADDRESS)v, sleftv_bin);
354  }
355  idDelete((ideal *)(&theMap));
356  return res;
357}
358
359/*2
360* substitute the n-th variable by the poly e in id
361* does not destroy id and e
362*/
363ideal  idSubstPoly(ideal id, int n, poly e)
364{
365
366#ifdef HAVE_PLURAL
367  if (rIsPluralRing(currRing))
368  {
369    int k=MATROWS((matrix)id)*MATCOLS((matrix)id);
370    ideal res=(ideal)mpNew(MATROWS((matrix)id),MATCOLS((matrix)id));
371    res->rank = id->rank;
372    for(k--;k>=0;k--)
373    {
374      res->m[k]=pSubst(pCopy(id->m[k]),n,e);
375    }
376    return res;
377  }
378#endif
379  map theMap=(map)idMaxIdeal(1);
380  theMap->preimage=NULL;
381  pDelete(&(theMap->m[n-1]));
382  theMap->m[n-1]=pCopy(e);
383
384  leftv v=(leftv)omAlloc0Bin(sleftv_bin);
385  sleftv tmpW;
386  memset(&tmpW,0,sizeof(sleftv));
387  tmpW.rtyp=IDEAL_CMD;
388  tmpW.data=id;
389  if (maApplyFetch(MAP_CMD,theMap,v,&tmpW,currRing,NULL,NULL,0,
390                          n_SetMap(currRing->cf, currRing->cf)))
391  {
392    WerrorS("map failed");
393    v->data=NULL;
394  }
395  ideal res=(ideal)(v->data);
396  idDelete((ideal *)(&theMap));
397  omFreeBin((ADDRESS)v, sleftv_bin);
398  return res;
399}
Note: See TracBrowser for help on using the repository browser.