Changeset f2a4f3f in git
- Timestamp:
- Jul 25, 2008, 6:06:18 PM (16 years ago)
- Branches:
- (u'spielwiese', '5b153614cbc72bfa198d75b1e9e33dab2645d9fe')
- Children:
- 0a9e45327ccc99fcc83b828f8bd61c187d8e46ee
- Parents:
- 356960f1a09d63bc25ab2b04269bb1f6ce67fd28
- Location:
- kernel
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
kernel/gring.cc
r356960 rf2a4f3f 7 7 * Author: levandov (Viktor Levandovsky) 8 8 * Created: 8/00 - 11/00 9 * Version: $Id: gring.cc,v 1.6 6 2008-07-24 09:45:53motsak Exp $9 * Version: $Id: gring.cc,v 1.67 2008-07-25 16:06:18 motsak Exp $ 10 10 *******************************************************************/ 11 11 … … 53 53 static const bool bNoCache = false; // only formula whenever possible, only make sanse if bNoFormula is false! 54 54 55 56 // false, true, false == old "good" Plural 57 // false, false ==>> Plural + Cache + Direct Formula - not much 58 // false, false, true ==>> Plural Mult + Direct Formula (no ~cache) 59 // true, *, * == new OOP multiplication! 55 60 56 61 bool bUseExtensions = true; … … 495 500 cpower = cpower + F[i]; 496 501 } 497 cpower = cpower*G[j]; 502 cpower = cpower*G[j]; // bug! here may happen an arithmetic overflow!!! 498 503 tpower = tpower + cpower; 499 504 } … … 516 521 if (F[i]!=0) 517 522 { 518 cpower = F[i]*G[j]; 523 cpower = F[i]*G[j]; // bug! overflow danger!!! 519 524 cff = n_Copy(p_GetCoeff(MATELEM(r->GetNC()->COM,j,i),r),r); 520 525 nPower(cff,cpower,&tmp_num); … … 766 771 int *Nxt=(int*)omAlloc0((rN+1)*sizeof(int)); 767 772 int *lF=(int *)omAlloc0((rN+1)*sizeof(int)); 773 768 774 int cnt=0; int cnf=0; 769 775 /* splitting F wrt jG */ … … 773 779 if (F[i]!=0) cnf++; 774 780 } 775 if (cnf==0) freeT(Prv,rN); 781 782 if (cnf==0) 783 { 784 freeT(Prv,rN); Prv = NULL; 785 } 786 776 787 for (i=jG+1;i<=rN;i++) 777 788 { … … 790 801 poly Rout=pOne(); 791 802 out=gnc_uu_Mult_ww(q,Nxt[q],jG,bG,r); 792 freeT(Nxt,rN); 803 804 freeT(Nxt,rN); Nxt = NULL; 793 805 794 806 if (cnf!=0) … … 797 809 p_SetExpV(Rout,Prv,r); 798 810 p_Setm(Rout,r); 811 799 812 #ifdef PDEBUG 800 813 p_Test(Rout,r); 801 814 #endif 815 802 816 freeT(Prv,rN); 817 Prv = NULL; 818 803 819 out=gnc_mm_Mult_p(Rout,out,r); /* getting the final result */ 804 820 } 805 821 806 omFreeSize((ADDRESS)lF,(rN+1)*sizeof(int)); 822 freeT(lF,rN); 823 lF = NULL; 824 807 825 p_Delete(&Rout,r); 826 827 assume(Nxt == NULL); 828 assume(lF == NULL); 829 assume(Prv == NULL); 830 808 831 return (out); 809 832 } … … 816 839 i=cnt+2; /* later in freeN */ 817 840 int *Op=Nxt; 841 818 842 int *On=(int *)omAlloc0((rN+1)*sizeof(int)); 819 843 int *U=(int *)omAlloc0((rN+1)*sizeof(int)); … … 904 928 /* make leadterm */ 905 929 /* ??????????? we have done it already :-0 */ 930 906 931 Rout=pOne(); 907 932 p_SetExpV(Rout,U,r); 908 933 p_Setm(Rout,r); /* use again this name */ 909 934 p_SetCoeff(Rout,c[cnt+1],r); /* last computed coef */ 935 910 936 out=p_Add_q(out,Rout,r); 937 911 938 Rout=NULL; 912 freeT(U,rN); 913 freeN(c,i); 914 omFreeSize((ADDRESS)lF,(rN+1)*sizeof(int)); 939 940 freeT(U, rN); 941 freeN(c, i); 942 freeT(lF, rN); 915 943 916 944 if (cnf!=0) … … 919 947 p_SetExpV(Rout,Prv,r); 920 948 p_Setm(Rout,r); 921 freeT(Prv, rN);949 freeT(Prv, rN); 922 950 out=gnc_mm_Mult_p(Rout,out,r); /* getting the final result */ 923 951 p_Delete(&Rout,r); 924 952 } 953 925 954 return (out); 926 955 } -
kernel/ncSAMult.cc
r356960 rf2a4f3f 7 7 * Author: motsak 8 8 * Created: 9 * Version: $Id: ncSAMult.cc,v 1. 7 2008-07-23 07:09:46motsak Exp $9 * Version: $Id: ncSAMult.cc,v 1.8 2008-07-25 16:06:18 motsak Exp $ 10 10 *******************************************************************/ 11 11 … … 46 46 #endif 47 47 48 49 CGlobalMultiplier* const pMultiplier = r->GetNC()->GetGlobalMultiplier(); 50 assume( pMultiplier != NULL ); 51 52 poly pMonom = pMultiplier->LM(m, r); 53 poly pResult = pMultiplier->MultiplyPE(p, pMonom); 54 p_Delete(&pMonom, r); 48 poly pResult; 49 50 if (p_IsConstant(m, r)) 51 pResult = pp_Mult_nn(p, p_GetCoeff(m,r),r); 52 else 53 { 54 CGlobalMultiplier* const pMultiplier = r->GetNC()->GetGlobalMultiplier(); 55 assume( pMultiplier != NULL ); 56 57 poly pMonom = pMultiplier->LM(m, r); 58 pResult = pMultiplier->MultiplyPE(p, pMonom); 59 p_Delete(&pMonom, r); 60 p_Test(pResult, r); 61 pResult = p_Mult_nn(pResult, p_GetCoeff(m, r), r); 62 } 63 64 #if OUTPUT 55 65 p_Test(pResult, r); 56 pResult = p_Mult_nn(pResult, p_GetCoeff(m, r), r); 57 p_Test(pResult, r); 58 59 #if OUTPUT 66 60 67 Print("gnc_pp_Mult_mm(p, m) => "); p_Write(pResult, r); 61 68 PrintS("p: "); p_Write(p, r); … … 81 88 p_Write(m, r); 82 89 #endif 83 84 CGlobalMultiplier* const pMultiplier = r->GetNC()->GetGlobalMultiplier(); 85 assume( pMultiplier != NULL ); 86 87 poly pMonom = pMultiplier->LM(m, r); 88 poly pResult = pMultiplier->MultiplyPEDestroy(p, pMonom); 89 p_Delete(&pMonom, r); 90 91 poly pResult; 92 93 if (p_IsConstant(m, r)) 94 pResult = p_Mult_nn(p, p_GetCoeff(m,r),r); 95 else 96 { 97 CGlobalMultiplier* const pMultiplier = r->GetNC()->GetGlobalMultiplier(); 98 assume( pMultiplier != NULL ); 99 100 poly pMonom = pMultiplier->LM(m, r); 101 pResult = pMultiplier->MultiplyPEDestroy(p, pMonom); 102 p_Delete(&pMonom, r); 103 p_Test(pResult, r); 104 pResult = p_Mult_nn(pResult, p_GetCoeff(m, r), r); 105 } 106 107 #if OUTPUT 90 108 p_Test(pResult, r); 91 pResult = p_Mult_nn(pResult, p_GetCoeff(m, r), r); 92 p_Test(pResult, r); 93 94 #if OUTPUT 109 95 110 Print("gnc_p_Mult_mm(p, m) => "); p_Write(pResult, r); 96 111 // PrintS("p: "); p_Write(p, r); … … 117 132 PrintS("p: "); p_Write(p, r); 118 133 #endif 119 CGlobalMultiplier* const pMultiplier = r->GetNC()->GetGlobalMultiplier(); 120 assume( pMultiplier != NULL ); 121 122 poly pMonom = pMultiplier->LM(m, r); 123 poly pResult = pMultiplier->MultiplyEPDestroy(pMonom, p); 124 p_Delete(&pMonom, r); 134 135 poly pResult; 136 137 if (p_IsConstant(m, r)) 138 pResult = p_Mult_nn(p, p_GetCoeff(m,r),r); 139 else 140 { 141 CGlobalMultiplier* const pMultiplier = r->GetNC()->GetGlobalMultiplier(); 142 assume( pMultiplier != NULL ); 143 144 poly pMonom = pMultiplier->LM(m, r); 145 pResult = pMultiplier->MultiplyEPDestroy(pMonom, p); 146 p_Delete(&pMonom, r); 147 p_Test(pResult, r); 148 pResult = p_Mult_nn(pResult, p_GetCoeff(m, r), r); 149 } 150 151 #if OUTPUT 125 152 p_Test(pResult, r); 126 pResult = p_Mult_nn(pResult, p_GetCoeff(m, r), r); 127 p_Test(pResult, r); 128 129 130 #if OUTPUT 153 131 154 Print("gnc_mm_Mult_p(m, p) => "); p_Write(pResult, r); 132 155 // PrintS("p: "); p_Write(p, r); … … 153 176 #endif 154 177 155 156 CGlobalMultiplier* const pMultiplier = r->GetNC()->GetGlobalMultiplier(); 157 assume( pMultiplier != NULL ); 158 159 poly pMonom = pMultiplier->LM(m, r); 160 poly pResult = pMultiplier->MultiplyEP(pMonom, p); 161 p_Delete(&pMonom, r); 178 poly pResult; 179 180 if (p_IsConstant(m, r)) 181 pResult = pp_Mult_nn(p, p_GetCoeff(m,r),r); 182 else 183 { 184 CGlobalMultiplier* const pMultiplier = r->GetNC()->GetGlobalMultiplier(); 185 assume( pMultiplier != NULL ); 186 187 poly pMonom = pMultiplier->LM(m, r); 188 pResult = pMultiplier->MultiplyEP(pMonom, p); 189 p_Delete(&pMonom, r); 190 p_Test(pResult, r); 191 pResult = p_Mult_nn(pResult, p_GetCoeff(m, r), r); 192 } 193 194 #if OUTPUT 162 195 p_Test(pResult, r); 163 pResult = p_Mult_nn(pResult, p_GetCoeff(m, r), r); 164 p_Test(pResult, r); 165 166 #if OUTPUT 196 167 197 Print("gnc_mm_Mult_pp(m, p) => "); p_Write(pResult, r); 168 198 PrintS("p: "); p_Write(p, r); -
kernel/ncSAMult.h
r356960 rf2a4f3f 4 4 * Computer Algebra System SINGULAR * 5 5 *****************************************/ 6 /* $Id: ncSAMult.h,v 1. 8 2008-07-23 07:09:46motsak Exp $ */6 /* $Id: ncSAMult.h,v 1.9 2008-07-25 16:06:18 motsak Exp $ */ 7 7 #ifdef HAVE_PLURAL 8 8 … … 74 74 75 75 } 76 77 // Main templates:78 79 // Poly * Exponent80 inline poly MultiplyPE(const poly pPoly, const CExponent expRight)81 {82 bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET);83 CPolynomialSummator sum(GetBasering(), bUsePolynomial);84 85 for( poly q = pPoly; q !=NULL; q = pNext(q) )86 sum += MultiplyTE(q, expRight);87 88 return sum;89 }90 91 // Exponent * Poly92 inline poly MultiplyEP(const CExponent expLeft, const poly pPoly)93 {94 bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET);95 CPolynomialSummator sum(GetBasering(), bUsePolynomial);96 97 for( poly q = pPoly; q !=NULL; q = pNext(q) )98 sum += MultiplyET(expLeft, q);99 100 return sum;101 }102 103 // Poly * Exponent104 inline poly MultiplyPEDestroy(poly pPoly, const CExponent expRight)105 {106 bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET);107 CPolynomialSummator sum(GetBasering(), bUsePolynomial);108 109 for( ; pPoly!=NULL; pPoly = p_LmDeleteAndNext(pPoly, GetBasering()) )110 sum += MultiplyTE(pPoly, expRight);111 112 return sum;113 }114 115 // Exponent * Poly116 inline poly MultiplyEPDestroy(const CExponent expLeft, poly pPoly)117 {118 bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET);119 CPolynomialSummator sum(GetBasering(), bUsePolynomial);120 121 for( ; pPoly!=NULL; pPoly = p_LmDeleteAndNext(pPoly, GetBasering()) )122 sum += MultiplyET(expLeft, pPoly);123 124 return sum;125 }126 127 76 128 77 // protected: … … 251 200 virtual poly MultiplyEM(const CExponent expLeft, const poly pMonom); 252 201 253 202 // Main templates: 203 204 // Poly * Exponent 205 inline poly MultiplyPE(const poly pPoly, const CExponent expRight) 206 { 207 bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET); 208 CPolynomialSummator sum(GetBasering(), bUsePolynomial); 209 210 for( poly q = pPoly; q !=NULL; q = pNext(q) ) 211 sum += MultiplyTE(q, expRight); 212 213 return sum; 214 } 215 216 // Exponent * Poly 217 inline poly MultiplyEP(const CExponent expLeft, const poly pPoly) 218 { 219 bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET); 220 CPolynomialSummator sum(GetBasering(), bUsePolynomial); 221 222 for( poly q = pPoly; q !=NULL; q = pNext(q) ) 223 sum += MultiplyET(expLeft, q); 224 225 return sum; 226 } 227 228 // Poly * Exponent 229 inline poly MultiplyPEDestroy(poly pPoly, const CExponent expRight) 230 { 231 bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET); 232 CPolynomialSummator sum(GetBasering(), bUsePolynomial); 233 234 for( ; pPoly!=NULL; pPoly = p_LmDeleteAndNext(pPoly, GetBasering()) ) 235 sum += MultiplyTE(pPoly, expRight); 236 237 return sum; 238 } 239 240 // Exponent * Poly 241 inline poly MultiplyEPDestroy(const CExponent expLeft, poly pPoly) 242 { 243 bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET); 244 CPolynomialSummator sum(GetBasering(), bUsePolynomial); 245 246 for( ; pPoly!=NULL; pPoly = p_LmDeleteAndNext(pPoly, GetBasering()) ) 247 sum += MultiplyET(expLeft, pPoly); 248 249 return sum; 250 } 251 252 254 253 }; 255 254 … … 268 267 virtual ~CGlobalMultiplier(); 269 268 270 271 269 272 270 // protected: 273 271 typedef poly CExponent; … … 284 282 // Exponent * Monom 285 283 virtual poly MultiplyEM(const CExponent expLeft, const poly pMonom); 284 285 286 // Main templates: 287 288 // Poly * Exponent 289 inline poly MultiplyPE(const poly pPoly, const CExponent expRight) 290 { 291 assume( pPoly != NULL ); assume( expRight != NULL ); 292 const int iComponentMonom = p_GetComp(expRight, GetBasering()); 293 294 bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET); 295 CPolynomialSummator sum(GetBasering(), bUsePolynomial); 296 297 298 if( iComponentMonom!=0 ) 299 { 300 for( poly q = pPoly; q !=NULL; q = pNext(q) ) 301 { 302 #ifdef PDEBUG 303 { 304 const int iComponent = p_GetComp(q, GetBasering()); 305 assume(iComponent == 0); 306 if( iComponent!=0 ) 307 { 308 Werror("MultiplyPE: both sides have non-zero components: %d and %d!\n", iComponent, iComponentMonom); 309 // what should we do further?!? 310 return NULL; 311 } 312 313 } 314 #endif 315 sum += MultiplyTE(q, expRight); // NO Component!!! 316 } 317 poly t = sum; p_SetCompP(t, iComponentMonom, GetBasering()); 318 return t; 319 } // iComponentMonom != 0! 320 else 321 { // iComponentMonom == 0! 322 for( poly q = pPoly; q !=NULL; q = pNext(q) ) 323 { 324 const int iComponent = p_GetComp(q, GetBasering()); 325 326 #ifdef PDEBUG 327 if( iComponent!=0 ) 328 { 329 Warn("MultiplyPE: Multiplication in the left module from the right by component %d!\n", iComponent); 330 // what should we do further?!? 331 } 332 #endif 333 poly t = MultiplyTE(q, expRight); // NO Component!!! 334 p_SetCompP(t, iComponent, GetBasering()); 335 sum += t; 336 } 337 return sum; 338 } // iComponentMonom == 0! 339 } 340 341 // Exponent * Poly 342 inline poly MultiplyEP(const CExponent expLeft, const poly pPoly) 343 { 344 assume( pPoly != NULL ); assume( expLeft != NULL ); 345 const int iComponentMonom = p_GetComp(expLeft, GetBasering()); 346 347 bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET); 348 CPolynomialSummator sum(GetBasering(), bUsePolynomial); 349 350 if( iComponentMonom!=0 ) 351 { 352 for( poly q = pPoly; q !=NULL; q = pNext(q) ) 353 { 354 #ifdef PDEBUG 355 { 356 const int iComponent = p_GetComp(q, GetBasering()); 357 assume(iComponent == 0); 358 if( iComponent!=0 ) 359 { 360 Werror("MultiplyEP: both sides have non-zero components: %d and %d!\n", iComponent, iComponentMonom); 361 // what should we do further?!? 362 return NULL; 363 } 364 } 365 #endif 366 sum += MultiplyET(expLeft, q); 367 } 368 poly t = sum; p_SetCompP(t, iComponentMonom, GetBasering()); 369 return t; 370 } // iComponentMonom != 0! 371 else 372 { // iComponentMonom == 0! 373 for( poly q = pPoly; q !=NULL; q = pNext(q) ) 374 { 375 const int iComponent = p_GetComp(q, GetBasering()); 376 377 poly t = MultiplyET(expLeft, q); // NO Component!!! 378 p_SetCompP(t, iComponent, GetBasering()); 379 sum += t; 380 } 381 return sum; 382 } // iComponentMonom == 0! 383 } 384 385 // Poly * Exponent 386 inline poly MultiplyPEDestroy(poly pPoly, const CExponent expRight) 387 { 388 assume( pPoly != NULL ); assume( expRight != NULL ); 389 const int iComponentMonom = p_GetComp(expRight, GetBasering()); 390 391 bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET); 392 CPolynomialSummator sum(GetBasering(), bUsePolynomial); 393 394 395 if( iComponentMonom!=0 ) 396 { 397 for(poly q = pPoly ; q!=NULL; q = p_LmDeleteAndNext(q, GetBasering()) ) 398 { 399 #ifdef PDEBUG 400 { 401 const int iComponent = p_GetComp(q, GetBasering()); 402 assume(iComponent == 0); 403 if( iComponent!=0 ) 404 { 405 Werror("MultiplyPEDestroy: both sides have non-zero components: %d and %d!\n", iComponent, iComponentMonom); 406 // what should we do further?!? 407 return NULL; 408 } 409 410 } 411 #endif 412 sum += MultiplyTE(q, expRight); // NO Component!!! 413 } 414 poly t = sum; p_SetCompP(t, iComponentMonom, GetBasering()); 415 return t; 416 } // iComponentMonom != 0! 417 else 418 { // iComponentMonom == 0! 419 for(poly q = pPoly ; q!=NULL; q = p_LmDeleteAndNext(q, GetBasering()) ) 420 { 421 const int iComponent = p_GetComp(q, GetBasering()); 422 423 #ifdef PDEBUG 424 if( iComponent!=0 ) 425 { 426 Warn("MultiplyPEDestroy: Multiplication in the left module from the right by component %d!\n", iComponent); 427 // what should we do further?!? 428 } 429 #endif 430 poly t = MultiplyTE(q, expRight); // NO Component!!! 431 p_SetCompP(t, iComponent, GetBasering()); 432 sum += t; 433 } 434 return sum; 435 } // iComponentMonom == 0! 436 437 } 438 439 // Exponent * Poly 440 inline poly MultiplyEPDestroy(const CExponent expLeft, poly pPoly) 441 { 442 443 assume( pPoly != NULL ); assume( expLeft != NULL ); 444 const int iComponentMonom = p_GetComp(expLeft, GetBasering()); 445 446 bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET); 447 CPolynomialSummator sum(GetBasering(), bUsePolynomial); 448 449 if( iComponentMonom!=0 ) 450 { 451 for(poly q = pPoly ; q!=NULL; q = p_LmDeleteAndNext(q, GetBasering()) ) 452 { 453 #ifdef PDEBUG 454 { 455 const int iComponent = p_GetComp(q, GetBasering()); 456 assume(iComponent == 0); 457 if( iComponent!=0 ) 458 { 459 Werror("MultiplyEPDestroy: both sides have non-zero components: %d and %d!\n", iComponent, iComponentMonom); 460 // what should we do further?!? 461 return NULL; 462 } 463 } 464 #endif 465 sum += MultiplyET(expLeft, q); 466 } 467 poly t = sum; p_SetCompP(t, iComponentMonom, GetBasering()); 468 return t; 469 } // iComponentMonom != 0! 470 else 471 { // iComponentMonom == 0! 472 for(poly q = pPoly ; q!=NULL; q = p_LmDeleteAndNext(q, GetBasering()) ) 473 { 474 const int iComponent = p_GetComp(q, GetBasering()); 475 476 poly t = MultiplyET(expLeft, q); // NO Component!!! 477 p_SetCompP(t, iComponent, GetBasering()); 478 sum += t; 479 } 480 return sum; 481 } // iComponentMonom == 0! 482 483 } 484 485 486 286 487 287 488 }; -
kernel/summator.cc
r356960 rf2a4f3f 7 7 * Author: motsak 8 8 * Created: 9 * Version: $Id: summator.cc,v 1. 3 2008-07-08 13:30:54 SingularExp $9 * Version: $Id: summator.cc,v 1.4 2008-07-25 16:06:18 motsak Exp $ 10 10 *******************************************************************/ 11 11 … … 68 68 if(!m_bUsePolynomial) 69 69 { 70 poly out; 71 int pLength; 72 73 sBucketClearAdd(m_temp.m_bucket, &out, &pLength); 70 74 sBucketDestroy(&m_temp.m_bucket); 75 76 if(out != NULL) 77 p_Delete(&out, m_basering); 71 78 // m_temp.m_bucket = NULL; 72 79 }
Note: See TracChangeset
for help on using the changeset viewer.