Changeset fb4075b in git
- Timestamp:
- Nov 8, 2010, 5:49:37 PM (13 years ago)
- Branches:
- (u'spielwiese', 'd1ba061a762c62d3a25159d8da8b6e17332291fa')
- Children:
- a7ee69c27ecef9708a4ea7677de1715f1ff6fa86
- Parents:
- f34215cd2606c6fc99bffc25e1a276b416786dca
- git-author:
- Hans Schoenemann <hannes@mathematik.uni-kl.de>2010-11-08 17:49:37+01:00
- git-committer:
- Mohamed Barakat <mohamed.barakat@rwth-aachen.de>2011-11-09 11:55:36+01:00
- Location:
- polys
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
polys/monomials/p_polys.cc
rf34215 rfb4075b 1272 1272 } 1273 1273 } 1274 /*2 1275 * assumes that the head term of b is a multiple of the head term of a 1276 * and return the multiplicant *m 1277 * Frank's observation: If LM(b) = LM(a)*m, then we may actually set 1278 * negative(!) exponents in the below loop. I suspect that the correct 1279 * comment should be "assumes that LM(a) = LM(b)*m, for some monomial m..." 1280 */ 1281 poly p_Divide(poly a, poly b, const ring r) 1282 { 1283 assume((p_GetComp(a,r)==p_GetComp(b,r)) || (p_GetComp(b,r)==0)); 1284 int i; 1285 poly result = pInit(); 1286 1287 for(i=(int)r->N; i; i--) 1288 p_SetExp(result,i, p_GetExp(a,i,r)- p_GetExp(b,i,r),r); 1289 p_SetComp(result, p_GetComp(a,r) - p_GetComp(b,r),r); 1290 p_Setm(result,r); 1291 return result; 1292 } 1293 1294 /*2 1295 * divides a by the monomial b, ignores monomials which are not divisible 1296 * assumes that b is not NULL 1297 */ 1298 poly p_DivideM(poly a, poly b, const ring r) 1299 { 1300 if (a==NULL) return NULL; 1301 poly result=a; 1302 poly prev=NULL; 1303 int i; 1304 #ifdef HAVE_RINGS 1305 number inv=pGetCoeff(b); 1306 #else 1307 number inv=n_Invers(pGetCoeff(b),r->cf); 1308 #endif 1309 1310 while (a!=NULL) 1311 { 1312 if (p_DivisibleBy(b,a,r)) 1313 { 1314 for(i=(int)r->N; i; i--) 1315 p_SubExp(a,i, p_GetExp(b,i,r),r); 1316 p_SubComp(a, p_GetComp(b,r),r); 1317 p_Setm(a,r); 1318 prev=a; 1319 pIter(a); 1320 } 1321 else 1322 { 1323 if (prev==NULL) 1324 { 1325 p_DeleteLm(&result,r); 1326 a=result; 1327 } 1328 else 1329 { 1330 p_DeleteLm(&pNext(prev),r); 1331 a=pNext(prev); 1332 } 1333 } 1334 } 1335 #ifdef HAVE_RINGS 1336 if (n_IsUnit(inv,r->cf)) 1337 { 1338 inv = n_Invers(inv,r->cf); 1339 p_Mult_nn(result,inv,r); 1340 n_Delete(&inv, r->cf); 1341 } 1342 else 1343 { 1344 p_Div_nn(result,inv,r); 1345 } 1346 #else 1347 p_Mult_nn(result,inv,r); 1348 n_Delete(&inv, r->cf); 1349 #endif 1350 p_Delete(&b, r); 1351 return result; 1352 } 1274 1353 1275 1354 /*************************************************************** -
polys/monomials/p_polys.h
rf34215 rfb4075b 1755 1755 poly p_mInit(const char *s, BOOLEAN &ok, const ring r); /* monom s -> poly, interpreter */ 1756 1756 const char * p_Read(const char *s, poly &p,const ring r); /* monom -> poly */ 1757 poly p_Divide(poly a, poly b, const ring r); 1758 poly p_DivideM(poly a, poly b, const ring r); 1757 1759 #endif // P_POLYS_H 1758 1760 -
polys/polys.cc
rf34215 rfb4075b 91 91 } 92 92 93 /*294 * assumes that the head term of b is a multiple of the head term of a95 * and return the multiplicant *m96 * Frank's observation: If LM(b) = LM(a)*m, then we may actually set97 * negative(!) exponents in the below loop. I suspect that the correct98 * comment should be "assumes that LM(a) = LM(b)*m, for some monomial m..."99 */100 #define pDivide(a,b) p_Divide(a,b,currRing)101 poly p_Divide(poly a, poly b, cont ring r)102 {103 assume((pGetComp(a)==pGetComp(b)) || (pGetComp(b)==0));104 int i;105 poly result = pInit();106 107 for(i=(int)r->N; i; i--)108 p_SetExp(result,i, p_GetExp(a,i,r)- p_GetExp(b,i,r),r);109 p_SetComp(result, p_GetComp(a,r) - p_GetComp(b,r),r);110 p_Setm(result,r);111 return result;112 }113 114 93 #ifdef HAVE_RINGS //TODO Oliver 115 94 #define pDiv_nn(p, n) p_Div_nn(p, n, currRing) … … 146 125 } 147 126 #endif 148 149 /*2150 * divides a by the monomial b, ignores monomials which are not divisible151 * assumes that b is not NULL, destroys b152 */153 poly p_DivideM(poly a, poly b, const ring r)154 {155 if (a==NULL) { pDelete(&b); return NULL; }156 poly result=a;157 poly prev=NULL;158 int i;159 #ifdef HAVE_RINGS160 number inv=pGetCoeff(b);161 #else162 number inv=nInvers(pGetCoeff(b));163 #endif164 165 while (a!=NULL)166 {167 if (p_DivisibleBy(b,a,r))168 {169 assume((p_GetComp(a,r)==p_GetComp(b,r)) || (p_GetComp(b,r)==0));170 for(i=(int)r->N; i; i--)171 p_SubExp(a,i, p_GetExp(b,i,r),r);172 p_SubComp(a, p_GetComp(b,r),r);173 p_Setm(a,r);174 prev=a;175 pIter(a);176 }177 else178 {179 if (prev==NULL)180 {181 p_DeleteLm(&result,r);182 a=result;183 }184 else185 {186 p_DeleteLm(&pNext(prev),r);187 a=pNext(prev);188 }189 }190 }191 #ifdef HAVE_RINGS192 if (n_IsUnit(inv,r->cf))193 {194 inv = n_Invers(inv,r->cf);195 p_Mult_nn(result,inv,r);196 n_Delete(&inv, r->cf);197 }198 else199 {200 p_Div_nn(result,inv,r);201 }202 #else203 p_Mult_nn(result,inv,r);204 n_Delete(&inv, r->cf);205 #endif206 p_Delete(&b, r);207 return result;208 }209 127 210 128 /*2 -
polys/polys.h
rf34215 rfb4075b 289 289 // ----------------- define to enable new p_procs -----*/ 290 290 291 poly p_Divide(poly a, poly b, const ring r); 292 poly p_DivideM(poly a, poly b, const ring r); 291 #define pDivide(a,b) p_Divide(a,b,currRing) 292 293 293 void pLcm(poly a, poly b, poly m); 294 294 poly pDiff(poly a, int k);
Note: See TracChangeset
for help on using the changeset viewer.