Changeset 0001f9 in git
- Timestamp:
- Jan 5, 2006, 7:32:19 PM (18 years ago)
- Branches:
- (u'spielwiese', 'a719bcf0b8dbc648b128303a49777a094b57592c')
- Children:
- 44f7bcdefb6eab48841c3b44943126dc814d480e
- Parents:
- f764969054eb825a57286a74a605b9a43ee1a84c
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
Singular/walk.cc
rf76496 r0001f9 2 2 * Computer Algebra System SINGULAR * 3 3 *****************************************/ 4 /* $Id: walk.cc,v 1. 8 2005-05-31 07:42:20 brickenExp $ */4 /* $Id: walk.cc,v 1.9 2006-01-05 18:32:19 Singular Exp $ */ 5 5 /* 6 6 * ABSTRACT: Implementation of the Groebner walk … … 13 13 14 14 //#define UPPER_BOUND //for the original "Tran" algorithm 15 //#define REPRESENTATION_OF_SIGMA //if one perturbs sigma in Tran 15 //#define REPRESENTATION_OF_SIGMA //if one perturbs sigma in Tran 16 16 17 17 //#define TEST_OVERFLOW … … 108 108 clock_t xtif, xtstd, xtlift, xtred, xtnw; 109 109 clock_t xftostd, xtextra, xftinput, to; 110 110 111 111 /*2 112 112 *utilities for TSet, LSet … … 307 307 strat->sevT = initsevT(); 308 308 if (pOrdSgn == -1) strat->honey = TRUE; 309 310 309 310 311 311 //initSCC(F,Q,strat); 312 312 initS(F,Q,strat); … … 354 354 } 355 355 356 static void TimeString(clock_t tinput, clock_t tostd, clock_t tif,clock_t tstd, 357 358 { 359 double totm = ((double) (clock() - tinput))/1000000; 356 static void TimeString(clock_t tinput, clock_t tostd, clock_t tif,clock_t tstd, 357 clock_t tlf,clock_t tred, clock_t tnw, int step) 358 { 359 double totm = ((double) (clock() - tinput))/1000000; 360 360 double ostd,mostd, mif, mstd, mextra, mlf, mred, mnw, mxif,mxstd,mxlf,mxred,mxnw,tot; 361 361 … … 388 388 double res = (double) 100 - tot; 389 389 Print("\n// &%d&%.2f&%.2f&%.2f&%.2f&%.2f&%.2f&%.2f&%.2f&%.2f&%.2f&%.2f&%.2f&%.2f&%.2f&%.2f(%.2f)\\ \\", 390 391 392 } 393 394 static void TimeStringFractal(clock_t tinput, clock_t tostd, clock_t tif,clock_t tstd, 395 396 { 397 390 step, ostd, totm, mostd,mif,mstd,mlf,mred,mnw,mxif,mxstd,mxlf,mxred,mxnw,tot,res, 391 ((((double) xtextra)/1000000)/totm)*100); 392 } 393 394 static void TimeStringFractal(clock_t tinput, clock_t tostd, clock_t tif,clock_t tstd, 395 clock_t textra, clock_t tlf,clock_t tred, clock_t tnw) 396 { 397 398 398 double totm = ((double) (clock() - tinput))/1000000; 399 399 double ostd, mostd, mif, mstd, mextra, mlf, mred, mnw, tot, res; … … 415 415 tot = mostd+mif+mstd+mextra+mlf+mred+mnw; 416 416 res = (double) 100.00-tot; 417 Print("\n// &%.2f &%.2f&%.2f &%.2f &%.2f &%.2f &%.2f &%.2f &%.2f&%.2f&%.2f\\ \\ ", 418 417 Print("\n// &%.2f &%.2f&%.2f &%.2f &%.2f &%.2f &%.2f &%.2f &%.2f&%.2f&%.2f\\ \\ ", 418 ostd,totm,mostd,mif,mstd,mextra,mlf,mred,mnw,tot,res); 419 419 } 420 420 421 421 static void idString(ideal L, char* st) 422 422 { 423 int i, nL = IDELEMS(L); 423 int i, nL = IDELEMS(L); 424 424 425 425 Print("\n// ideal %s = ", st); 426 426 for(i=0; i<nL-1; i++) 427 427 Print(" %s, ", pString(L->m[i])); 428 428 429 429 Print(" %s;", pString(L->m[nL-1])); 430 430 } … … 432 432 static void headidString(ideal L, char* st) 433 433 { 434 int i, nL = IDELEMS(L); 434 int i, nL = IDELEMS(L); 435 435 436 436 Print("\n// ideal %s = ", st); 437 437 for(i=0; i<nL-1; i++) 438 438 Print(" %s, ", pString(pHead(L->m[i]))); 439 439 440 440 Print(" %s;", pString(pHead(L->m[nL-1]))); 441 441 } … … 443 443 static void idElements(ideal L, char* st) 444 444 { 445 int i, nL = IDELEMS(L); 446 int K[nL];445 int i, nL = IDELEMS(L); 446 int *K=(int *)omAlloc(nL*sizeof(int)); 447 447 448 448 Print("\n// #monoms of %s = ", st); … … 457 457 nsame = 1; 458 458 for(j=i+1; j<nL; j++){ 459 460 461 K[j]=0; 462 459 if(K[j]==K[i]){ 460 nsame ++; 461 K[j]=0; 462 } 463 463 } 464 464 if(nsame == 1) 465 465 Print("%d, ",K[i]); 466 466 else 467 Print("%d[%d], ", K[i], nsame); 468 } 469 } 467 Print("%d[%d], ", K[i], nsame); 468 } 469 } 470 omFree(K); 470 471 } 471 472 … … 495 496 Print("%d, ", (*ivb)[i]); 496 497 Print("%d) := (", (*ivb)[nV]); 497 498 498 499 for(i=0; i<nV; i++) 499 500 Print("%d, ", (*ivc)[i]); … … 558 559 559 560 int i, wgrad; 560 561 561 562 mpz_t zmul; 562 563 mpz_init(zmul); … … 565 566 mpz_t zsum; 566 567 mpz_init(zsum); 567 568 568 569 for (i=pVariables; i>0; i--) 569 570 { … … 585 586 } 586 587 } 587 588 588 589 return wgrad; 589 590 } … … 602 603 pIter(p); 603 604 604 if (maxtemp > max) 605 if (maxtemp > max) 605 606 max = maxtemp; 606 607 } … … 619 620 620 621 int i, wgrad; 621 622 622 623 mpz_t zmul; 623 624 mpz_init(zmul); … … 634 635 } 635 636 mpz_init_set(result, ztmp); 637 mpz_clear(ztmp); 638 mpz_clear(sing_int); 639 mpz_clear(zvec); 640 mpz_clear(zmul); 636 641 } 637 642 … … 647 652 mpz_t max; mpz_init(max); 648 653 mpz_t maxtmp; mpz_init(maxtmp); 649 650 poly hg, in_w_g = NULL; 654 655 poly hg, in_w_g = NULL; 651 656 652 657 while(g != NULL) … … 656 661 MLmWeightedDegree_gmp(maxtmp, hg, curr_weight); 657 662 658 if(mpz_cmp(maxtmp, max)>0) 663 if(mpz_cmp(maxtmp, max)>0) 659 664 { 660 665 mpz_init_set(max, maxtmp); 661 666 pDelete(&in_w_g); 662 667 in_w_g = pHead(hg); 663 } 668 } 664 669 else { 665 670 if(mpz_cmp(maxtmp, max)==0) … … 713 718 mi = MpolyInitialForm(G->m[i], iv); 714 719 gi = G->m[i]; 715 720 716 721 if(mi == NULL) 717 722 { 718 723 pDelete(&mi); 719 724 720 725 if(Overflow_Error == FALSE) 721 726 Overflow_Error = nError; 722 727 723 728 return 0; 724 } 729 } 725 730 if(!pLmEqual(mi, gi)) 726 731 { 727 732 pDelete(&mi); 728 733 729 734 if(Overflow_Error == FALSE) 730 735 Overflow_Error = nError; 731 736 732 737 return 0; 733 738 } 734 739 735 740 pDelete(&mi); 736 741 } … … 739 744 Overflow_Error = nError; 740 745 741 return 1; 746 return 1; 742 747 } 743 748 … … 746 751 static inline long Mlcm(long &i1, long &i2) 747 752 { 748 long temp = gcd(i1, i2); 753 long temp = gcd(i1, i2); 749 754 return ((i1 / temp)* i2); 750 755 } … … 759 764 int i, n = a->length(); 760 765 long result = 0; 761 766 762 767 for(i=n-1; i>=0; i--) 763 768 result += (*a)[i] * (*b)[i]; … … 770 775 { 771 776 assume( a->length() == b->length()); 772 int i, n = a->length(); 773 intvec* result = new intvec(n); 774 777 int i, n = a->length(); 778 intvec* result = new intvec(n); 779 775 780 for(i=n-1; i>=0; i--) 776 781 (*result)[i] = (*a)[i] - (*b)[i]; … … 786 791 int i, nR = currRing->N; 787 792 intvec* result = new intvec(nR); 788 793 789 794 for(i=nR-1; i>=0; i--) 790 795 (*result)[i] = pGetExp(f,i+1); … … 795 800 /* return 1, if two given intvecs are the same, otherwise 0*/ 796 801 int MivSame(intvec* u , intvec* v) 797 { 802 { 798 803 assume(u->length() == v->length()); 799 804 800 805 int i, niv = u->length(); 801 806 802 807 for (i=0; i<niv; i++) 803 808 if ((*u)[i] != (*v)[i]) 804 809 return 0; 805 810 806 return 1; 811 return 1; 807 812 } 808 813 809 814 int M3ivSame(intvec* temp, intvec* u , intvec* v) 810 { 815 { 811 816 assume(temp->length() == u->length() && u->length() == v->length()); 812 817 813 818 if((MivSame(temp, u)) == 1) 814 return 0; 819 return 0; 815 820 816 821 if((MivSame(temp, v)) == 1) 817 822 return 1; 818 823 819 return 2; 824 return 2; 820 825 } 821 826 … … 828 833 ideal G1 = kStd(G, NULL, testHomog, NULL); 829 834 test=save_test; 830 835 831 836 idSkipZeroes(G1); 832 837 return G1; … … 842 847 test=save_test; 843 848 844 idSkipZeroes(G1); 849 idSkipZeroes(G1); 845 850 return G1; 846 851 } … … 852 857 intvec* MivMatrixOrder(intvec* iv) 853 858 { 854 int i,j, nR = iv->length(); 859 int i,j, nR = iv->length(); 855 860 intvec* ivm = new intvec(nR*nR); 856 861 … … 867 872 intvec* Mivdp(int nR) 868 873 { 869 int i; 874 int i; 870 875 intvec* ivm = new intvec(nR); 871 876 … … 879 884 intvec* Mivlp(int nR) 880 885 { 881 int i; 886 int i; 882 887 intvec* ivm = new intvec(nR); 883 888 (*ivm)[0] = 1; 884 889 885 return ivm; 890 return ivm; 886 891 } 887 892 … … 891 896 int nV = currRing->N; 892 897 int nG = IDELEMS(G); 893 intvec* ivUnit = Mivdp(nV);//19.02 898 intvec* ivUnit = Mivdp(nV);//19.02 894 899 int i,j, tmpdeg, maxdeg=0; 895 900 number tmpcoeff , maxcoeff=nNULL; … … 898 903 { 899 904 tmpdeg = MwalkWeightDegree(G->m[i], ivUnit); 900 if (tmpdeg > maxdeg ) 905 if (tmpdeg > maxdeg ) 901 906 maxdeg = tmpdeg; 902 907 } … … 910 915 tmpcoeff = pGetCoeff(p); 911 916 if(nGreater(tmpcoeff,maxcoeff)) 912 917 maxcoeff = nCopy(tmpcoeff); 913 918 pIter(p); 914 919 } … … 934 939 * basering. * 935 940 * This programm computes a perturbated vector with a p_deg perturbation * 936 * degree which smaller than the numbers of variables * 941 * degree which smaller than the numbers of variables * 937 942 ******************************************************************************/ 938 943 /* ivtarget is a matrix order of a degree reverse lex. order */ … … 945 950 intvec* v_null = new intvec(nV); 946 951 947 952 948 953 //Checking that the perturbed degree is valid 949 954 if(pdeg > nV || pdeg <= 0) 950 { 955 { 951 956 WerrorS("//** The perturbed degree is wrong!!"); 952 957 return v_null; 953 958 } 954 959 delete v_null; 955 960 956 961 if(pdeg == 1) 957 962 return ivtarget; 958 963 959 mpz_t pert_vector[nV];960 964 mpz_t *pert_vector=(mpz_t*)omAlloc(nV*sizeof(mpz_t)); 965 961 966 for(i=0; i<nV; i++) 962 967 mpz_init_set_si(pert_vector[i], (*ivtarget)[i]); 963 964 968 969 965 970 // Calculate max1 = Max(A2)+Max(A3)+...+Max(Apdeg), 966 971 // where the Ai are the i-te rows of the matrix target_ord. 967 972 968 973 int ntemp, maxAi, maxA=0; 969 for(i=1; i<pdeg; i++) 974 for(i=1; i<pdeg; i++) 970 975 { 971 976 maxAi = (*ivtarget)[i*nV]; 972 977 if(maxAi<0) maxAi = -maxAi; 973 978 974 979 for(j=i*nV+1; j<(i+1)*nV; j++) 975 980 { 976 981 ntemp = (*ivtarget)[j]; 977 982 if(ntemp < 0) ntemp = -ntemp; 978 983 979 984 if(ntemp > maxAi) 980 985 maxAi = ntemp; 981 986 } 982 maxA += maxAi; 987 maxA += maxAi; 983 988 } 984 989 … … 986 991 987 992 intvec* ivUnit = Mivdp(nV); 988 993 989 994 mpz_t tot_deg; mpz_init(tot_deg); 990 995 mpz_t maxdeg; mpz_init(maxdeg); 991 996 mpz_t inveps; mpz_init(inveps); 992 997 993 998 994 999 for(i=nG-1; i>=0; i--) 995 1000 { 996 1001 mpz_set_ui(maxdeg, MwalkWeightDegree(G->m[i], ivUnit)); 997 if (mpz_cmp(maxdeg, tot_deg) > 0 ) 1002 if (mpz_cmp(maxdeg, tot_deg) > 0 ) 998 1003 mpz_set(tot_deg, maxdeg); 999 1004 } 1000 1001 delete ivUnit; 1005 1006 delete ivUnit; 1002 1007 mpz_mul_ui(inveps, tot_deg, maxA); 1003 1008 mpz_add_ui(inveps, inveps, 1); … … 1005 1010 1006 1011 //xx1.06.02 takes "small" inveps 1007 #ifdef INVEPS_SMALL_IN_MPERTVECTOR 1008 if(mpz_cmp_ui(inveps, pdeg)>0 && pdeg > 3) 1012 #ifdef INVEPS_SMALL_IN_MPERTVECTOR 1013 if(mpz_cmp_ui(inveps, pdeg)>0 && pdeg > 3) 1009 1014 { 1010 1015 /* 1011 Print("\n// choose the\"small\" inverse epsilon := %d / %d = ", 1012 mpz_get_si(inveps), pdeg); 1016 Print("\n// choose the\"small\" inverse epsilon := %d / %d = ", 1017 mpz_get_si(inveps), pdeg); 1013 1018 */ 1014 mpz_fdiv_q_ui(inveps, inveps, pdeg); 1019 mpz_fdiv_q_ui(inveps, inveps, pdeg); 1015 1020 //mpz_out_str(stdout, 10, inveps); 1016 1021 } 1017 1022 #else 1018 //PrintS("\n// the \"big\" inverse epsilon: "); 1023 //PrintS("\n// the \"big\" inverse epsilon: "); 1019 1024 mpz_out_str(stdout, 10, inveps); 1020 #endif 1025 #endif 1021 1026 1022 1027 // pert(A1) = inveps^(pdeg-1)*A1 + inveps^(pdeg-2)*A2+...+A_pdeg, … … 1031 1036 mpz_add_ui(pert_vector[j], pert_vector[j],(*ivtarget)[i*nV+j]); 1032 1037 } 1033 1038 1034 1039 mpz_t ztemp; 1035 1040 mpz_init(ztemp); … … 1054 1059 { 1055 1060 (*result)[i] = mpz_get_si(pert_vector[i]); 1056 1057 if(mpz_cmp(pert_vector[i], sing_int)>=0) 1058 { 1061 1062 if(mpz_cmp(pert_vector[i], sing_int)>=0) 1063 { 1059 1064 ntrue++; 1060 1065 if(Overflow_Error == FALSE) … … 1062 1067 Overflow_Error = TRUE; 1063 1068 PrintS("\n// ** OVERFLOW in \"MPertvectors\": "); 1064 1069 mpz_out_str( stdout, 10, pert_vector[i]); 1065 1070 PrintS(" is greater than 2147483647 (max. integer representation)"); 1066 1071 Print("\n// So vector[%d] := %d is wrong!!", i+1, (*result)[i]); 1067 } 1068 } 1069 } 1070 1071 if(Overflow_Error == TRUE) 1072 } 1073 } 1074 } 1075 1076 if(Overflow_Error == TRUE) 1072 1077 { 1073 1078 ivString(result, "pert_vector"); 1074 1079 Print("\n// %d element(s) of it is overflow!!", ntrue); 1075 1080 } 1076 1081 1082 mpz_clear(ztemp); 1083 mpz_clear(sing_int); 1084 omFree(pert_vector); 1077 1085 return result; 1078 1086 } … … 1087 1095 int i, j, nG = IDELEMS(G); 1088 1096 intvec* pert_vector = new intvec(nV); 1089 1097 1090 1098 //Checking that the perturbated degree is valid 1091 1099 if(pdeg > nV || pdeg <= 0) 1092 { 1100 { 1093 1101 WerrorS("//** The perturbed degree is wrong!!"); 1094 1102 return pert_vector; … … 1099 1107 if(pdeg == 1) 1100 1108 return pert_vector; 1101 1109 1102 1110 // Calculate max1 = Max(A2)+Max(A3)+...+Max(Apdeg), 1103 1111 // where the Ai are the i-te rows of the matrix target_ord. … … 1105 1113 for(i=1; i<pdeg; i++) 1106 1114 { 1107 maxAi = (*ivtarget)[i*nV]; 1115 maxAi = (*ivtarget)[i*nV]; 1108 1116 for(j=i*nV+1; j<(i+1)*nV; j++) 1109 1117 { … … 1112 1120 maxAi = ntemp; 1113 1121 } 1114 maxA += maxAi; 1115 } 1116 1122 maxA += maxAi; 1123 } 1124 1117 1125 // Calculate inveps := 1/eps, where 1/eps > deg(p)*max1 for all p in G. 1118 1126 int inveps, tot_deg = 0, maxdeg; 1119 1127 1120 intvec* ivUnit = Mivdp(nV);//19.02 1128 intvec* ivUnit = Mivdp(nV);//19.02 1121 1129 for(i=nG-1; i>=0; i--) 1122 1130 { 1123 1131 //maxdeg = pTotaldegree(G->m[i], currRing); //it's wrong for ex1,2,rose 1124 1132 maxdeg = MwalkWeightDegree(G->m[i], ivUnit); 1125 if (maxdeg > tot_deg ) 1133 if (maxdeg > tot_deg ) 1126 1134 tot_deg = maxdeg; 1127 1135 } … … 1129 1137 1130 1138 inveps = (tot_deg * maxA) + 1; 1131 1139 1132 1140 //9.10.01 1133 #ifdef INVEPS_SMALL_IN_FRACTAL 1141 #ifdef INVEPS_SMALL_IN_FRACTAL 1134 1142 /* 1135 Print("\n// choose the\"small\" inverse epsilon := %d / %d = ", 1136 inveps, pdeg); 1143 Print("\n// choose the\"small\" inverse epsilon := %d / %d = ", 1144 inveps, pdeg); 1137 1145 */ 1138 if(inveps > pdeg && pdeg > 3) 1146 if(inveps > pdeg && pdeg > 3) 1139 1147 inveps = inveps / pdeg; 1140 1148 1141 1149 //Print(" %d", inveps); 1142 1150 #else 1143 PrintS("\n// the \"big\" inverse epsilon %d", inveps); 1151 PrintS("\n// the \"big\" inverse epsilon %d", inveps); 1144 1152 #endif 1145 1153 1146 1154 // Pert(A1) = inveps^(pdeg-1)*A1 + inveps^(pdeg-2)*A2+...+A_pdeg, 1147 1155 for ( i=1; i < pdeg; i++ ) 1148 1156 for(j=0; j<nV; j++) 1149 1157 (*pert_vector)[j] = inveps*((*pert_vector)[j]) + (*ivtarget)[i*nV+j]; 1150 1158 1151 1159 int temp = (*pert_vector)[0]; 1152 1160 for(i=1; i<nV; i++) … … 1171 1179 int i; 1172 1180 intvec* ivM = new intvec(nV*nV); 1173 1181 1174 1182 for(i=0; i<nV; i++) 1175 1183 (*ivM)[i*nV + i] = 1; 1176 1184 1177 1185 return(ivM); 1178 1186 } … … 1183 1191 int i; 1184 1192 intvec* ivM = new intvec(nV*nV); 1185 1193 1186 1194 for(i=0; i<nV; i++) 1187 1195 (*ivM)[i] = 1; … … 1189 1197 for(i=1; i<nV; i++) 1190 1198 (*ivM)[(i+1)*nV - i] = -1; 1191 1199 1192 1200 return(ivM); 1193 1201 } … … 1199 1207 int nV = ivstart->length(); 1200 1208 intvec* ivM = new intvec(nV*nV); 1201 1209 1202 1210 for(i=0; i<nV; i++) 1203 1211 (*ivM)[i] = (*ivstart)[i]; … … 1205 1213 for(i=1; i<nV; i++) 1206 1214 (*ivM)[i*nV + i-1] = 1; 1207 1215 1208 1216 return(ivM); 1209 1217 } … … 1214 1222 int nV = ivstart->length(); 1215 1223 intvec* ivM = new intvec(nV*nV); 1216 1224 1217 1225 for(i=0; i<nV; i++) 1218 1226 (*ivM)[i] = (*ivstart)[i]; … … 1223 1231 for(i=2; i<nV; i++) 1224 1232 (*ivM)[(i+1)*nV - i] = -1; 1225 1233 1226 1234 return(ivM); 1227 1235 } … … 1231 1239 int i; 1232 1240 intvec* ivM = new intvec(nV*nV); 1233 1241 1234 1242 for(i=0; i<nV; i++) 1235 1243 (*ivM)[i] = 1; … … 1237 1245 for(i=1; i<nV; i++) 1238 1246 (*ivM)[(i+1)*nV - i] = -1; 1239 1247 1240 1248 return(ivM); 1241 1249 } … … 1245 1253 int i; 1246 1254 intvec* ivM = new intvec(nV); 1247 1255 1248 1256 for(i=nV-1; i>=0; i--) 1249 1257 (*ivM)[i] = 1; 1250 1258 1251 1259 return(ivM); 1252 1260 } … … 1267 1275 // where the Ai are the i-te rows of the matrix 'targer_ord'. 1268 1276 int ntemp, maxAi, maxA=0; 1269 for(i=1; i<nV; i++) 1277 for(i=1; i<nV; i++) 1270 1278 { 1271 1279 maxAi = (*ivtarget)[i*nV]; … … 1280 1288 maxAi = ntemp; 1281 1289 } 1282 maxA = maxA + maxAi; 1290 maxA = maxA + maxAi; 1283 1291 } 1284 1292 intvec* ivUnit = Mivdp(nV); 1285 1293 1286 1294 // Calculate inveps = 1/eps, where 1/eps > deg(p)*max1 for all p in G. 1287 1295 mpz_t tot_deg; mpz_init(tot_deg); 1288 1296 mpz_t maxdeg; mpz_init(maxdeg); 1289 1297 mpz_t inveps; mpz_init(inveps); 1290 1298 1291 1299 1292 1300 for(i=nG-1; i>=0; i--) 1293 1301 { 1294 1302 mpz_set_ui(maxdeg, MwalkWeightDegree(G->m[i], ivUnit)); 1295 if (mpz_cmp(maxdeg, tot_deg) > 0 ) 1303 if (mpz_cmp(maxdeg, tot_deg) > 0 ) 1296 1304 mpz_set(tot_deg, maxdeg); 1297 1305 } 1298 1299 delete ivUnit; 1300 //inveps = (tot_deg * maxA) + 1; 1306 1307 delete ivUnit; 1308 //inveps = (tot_deg * maxA) + 1; 1301 1309 mpz_mul_ui(inveps, tot_deg, maxA); 1302 1310 mpz_add_ui(inveps, inveps, 1); 1303 1311 1304 1312 //xx1.06.02 takes "small" inveps 1305 #ifdef INVEPS_SMALL_IN_FRACTAL 1313 #ifdef INVEPS_SMALL_IN_FRACTAL 1306 1314 if(mpz_cmp_ui(inveps, nV)>0 && nV > 3) 1307 1315 mpz_cdiv_q_ui(inveps, inveps, nV); 1308 1316 1309 1317 //PrintS("\n// choose the \"small\" inverse epsilon!"); 1310 #endif 1318 #endif 1311 1319 1312 1320 // PrintLn(); mpz_out_str(stdout, 10, inveps); 1313 1321 1314 1322 // Calculate the perturbed target orders: 1315 mpz_t ivtemp[nV];1316 mpz_t pert_vector[niv];1317 1323 mpz_t *ivtemp=(mpz_t *)omAlloc(nV*sizeof(mpz_t)); 1324 mpz_t *pert_vector=(mpz_t *)omAlloc(niv*sizeof(mpz_t)); 1325 1318 1326 for(i=0; i<nV; i++) 1319 1327 { … … 1321 1329 mpz_init_set_si(pert_vector[i], (*ivtarget)[i]); 1322 1330 } 1323 1331 1324 1332 mpz_t ztmp; mpz_init(ztmp); 1325 1333 BOOLEAN isneg = FALSE; … … 1330 1338 { 1331 1339 mpz_mul(ztmp, inveps, ivtemp[j]); 1332 1340 1333 1341 if((*ivtarget)[i*nV+j]<0) 1334 mpz_sub_ui(ivtemp[j], ztmp, -(*ivtarget)[i*nV+j]); 1335 else 1342 mpz_sub_ui(ivtemp[j], ztmp, -(*ivtarget)[i*nV+j]); 1343 else 1336 1344 mpz_add_ui(ivtemp[j], ztmp,(*ivtarget)[i*nV+j]); 1337 1345 } 1338 1346 1339 for(j=0; j<nV; j++) 1340 mpz_init_set(pert_vector[i*nV+j],ivtemp[j]); 1341 } 1342 1347 for(j=0; j<nV; j++) 1348 mpz_init_set(pert_vector[i*nV+j],ivtemp[j]); 1349 } 1350 1343 1351 /* 2147483647 is max. integer representation in SINGULAR */ 1344 1352 mpz_t sing_int; … … 1347 1355 intvec* result = new intvec(niv); 1348 1356 BOOLEAN nflow = FALSE; 1349 1350 // computes gcd 1351 mpz_set(ztmp, pert_vector[0]); 1357 1358 // computes gcd 1359 mpz_set(ztmp, pert_vector[0]); 1352 1360 for(i=0; i<niv; i++) 1353 1361 { … … 1358 1366 1359 1367 for(i=0; i<niv; i++) 1360 { 1368 { 1361 1369 mpz_divexact(pert_vector[i], pert_vector[i], ztmp); 1362 1370 (* result)[i] = mpz_get_si(pert_vector[i]); 1363 1371 1364 1372 if(mpz_cmp(pert_vector[i], sing_int)>0) 1365 1373 if(nflow == FALSE) … … 1379 1387 if(Overflow_Error == TRUE) 1380 1388 ivString(result, "new_vector"); 1381 1389 1390 omFree(pert_vector); 1391 omFree(ivtemp); 1392 mpz_clear(ztmp); 1393 1382 1394 return result; 1383 1395 } … … 1395 1407 return NULL; 1396 1408 1397 if(mB < mA) 1409 if(mB < mA) 1398 1410 mA = mB; 1399 1411 … … 1402 1414 int i, k=0; 1403 1415 for(i=0; i<mA; i++) 1404 { 1416 { 1405 1417 result->m[k] = pMult(A->m[i], pCopy(B->m[i])); 1406 1418 A->m[i]=NULL; 1407 1419 if (result->m[k]!=NULL) k++; 1408 1420 } 1409 1421 1410 1422 idDelete(&A); 1411 1423 idSkipZeroes(result); 1412 return result; 1424 return result; 1413 1425 } 1414 1426 … … 1422 1434 ********************************************************************/ 1423 1435 static ideal MLifttwoIdeal(ideal Gw, ideal M, ideal G) 1424 { 1436 { 1425 1437 ideal Mtmp = idLift(Gw, M, NULL, FALSE, TRUE, TRUE, NULL); 1426 1438 … … 1428 1440 //So, it is better, if one tests whether Gw is a GB 1429 1441 //in ideals.cc: 1430 //idLift (ideal mod, ideal submod,ideal * rest, BOOLEAN goodShape, 1442 //idLift (ideal mod, ideal submod,ideal * rest, BOOLEAN goodShape, 1431 1443 // BOOLEAN isSB,BOOLEAN divide,matrix * unit) 1432 1444 1433 1445 /* Let be Mtmp = {m1,...,ms}, where mi=sum hij.in_gj, for all i=1,...,s 1434 1446 We compute F = {f1,...,fs}, where fi=sum hij.gj */ … … 1436 1448 ideal idpol, idLG; 1437 1449 ideal F = idInit(nM, 1); 1438 1450 1439 1451 for(i=0; i<nM; i++) 1440 1452 { 1441 1453 idpol = idVec2Ideal(Mtmp->m[i]); 1442 1454 idLG = MidMult(idpol, G); 1443 idpol = NULL; 1455 idpol = NULL; 1444 1456 F->m[i] = NULL; 1445 1457 for(j=IDELEMS(idLG)-1; j>=0; j--) … … 1461 1473 int n = nG-1; 1462 1474 Print("\n//** Ideal %s besteht aus %d Polynomen mit ", Ch, nG); 1463 1475 1464 1476 for(i=0; i<nG; i++) 1465 1477 { … … 1477 1489 } 1478 1490 1479 static void HeadidString(ideal L, char* st) 1480 { 1481 int i, nL = IDELEMS(L)-1; 1491 static void HeadidString(ideal L, char* st) 1492 { 1493 int i, nL = IDELEMS(L)-1; 1482 1494 1483 1495 Print("// The head terms of the ideal %s = ", st); 1484 1496 for(i=0; i<nL; i++) 1485 1497 Print(" %s, ", pString(pHead(L->m[i]))); 1486 1498 1487 1499 Print(" %s;\n", pString(pHead(L->m[nL]))); 1488 1500 } … … 1505 1517 with respect to an ideal G. 1506 1518 */ 1507 static intvec* MwalkNextWeightCC(intvec* curr_weight, intvec* target_weight, 1519 static intvec* MwalkNextWeightCC(intvec* curr_weight, intvec* target_weight, 1508 1520 ideal G) 1509 1521 { … … 1511 1523 Overflow_Error = FALSE; 1512 1524 1513 assume(currRing != NULL && curr_weight != NULL && 1525 assume(currRing != NULL && curr_weight != NULL && 1514 1526 target_weight != NULL && G != NULL); 1515 1527 … … 1522 1534 mpz_init(t_nenner); 1523 1535 1524 mpz_t s_zaehler, s_nenner, temp, MwWd; 1536 mpz_t s_zaehler, s_nenner, temp, MwWd; 1525 1537 mpz_init(s_zaehler); 1526 1538 mpz_init(s_nenner); … … 1529 1541 1530 1542 1531 mpz_t deg_w0_p1, deg_d0_p1; 1543 mpz_t deg_w0_p1, deg_d0_p1; 1532 1544 mpz_init(deg_w0_p1); 1533 1545 mpz_init(deg_d0_p1); … … 1540 1552 1541 1553 mpz_t ggt; 1542 1554 1543 1555 int tn0, tn1, tz1, ncmp, gcd_tmp, ntmp; 1544 1556 intvec* diff_weight = MivSub(target_weight, curr_weight); 1545 1557 1546 1558 poly g, gw; 1547 1559 for (j=0; j<nG; j++) 1548 1560 { 1549 g = G->m[j]; 1550 if (g != NULL) 1561 g = G->m[j]; 1562 if (g != NULL) 1551 1563 { 1552 1564 ivtemp = MExpPol(g); … … 1554 1566 mpz_set_si(deg_d0_p1, MivDotProduct(ivtemp, diff_weight)); 1555 1567 delete ivtemp; 1556 1568 1557 1569 pIter(g); 1558 1570 while (g != NULL) … … 1566 1578 mpz_set_si(MwWd, MivDotProduct(ivtemp, diff_weight)); 1567 1579 mpz_sub(s_nenner, MwWd, deg_d0_p1); 1568 1580 1569 1581 // check for 0 < s <= 1 1570 if( (mpz_cmp(s_zaehler,t_null) > 0 && 1582 if( (mpz_cmp(s_zaehler,t_null) > 0 && 1571 1583 mpz_cmp(s_nenner, s_zaehler)>=0) || 1572 (mpz_cmp(s_zaehler, t_null) < 0 && 1584 (mpz_cmp(s_zaehler, t_null) < 0 && 1573 1585 mpz_cmp(s_nenner, s_zaehler)<=0)) 1574 1586 { … … 1579 1591 mpz_neg(s_nenner, s_nenner); 1580 1592 } 1581 1593 1582 1594 //compute a simply fraction of s 1583 1595 cancel(s_zaehler, s_nenner); 1584 1596 1585 1597 if(mpz_cmp(t_nenner, t_null) != 0) 1586 { 1598 { 1587 1599 mpz_mul(sztn, s_zaehler, t_nenner); 1588 1600 mpz_mul(sntz, s_nenner, t_zaehler); 1589 1601 1590 1602 if(mpz_cmp(sztn,sntz) < 0) 1591 { 1603 { 1592 1604 mpz_add(t_nenner, t_null, s_nenner); 1593 1605 mpz_add(t_zaehler,t_null, s_zaehler); 1594 } 1606 } 1595 1607 } 1596 1608 else … … 1607 1619 } 1608 1620 1609 mpz_t vec[nRing];1610 1611 /* there is no 0<t<1 and define the next weight vector that is equal to 1621 mpz_t *vec=(mpz_t*)omAlloc(nRing*sizeof(mpz_t)); 1622 1623 /* there is no 0<t<1 and define the next weight vector that is equal to 1612 1624 the current weight vector */ 1613 1625 if(mpz_cmp(t_nenner, t_null) == 0) 1614 { 1615 delete diff_weight; 1626 { 1627 delete diff_weight; 1616 1628 diff_weight = ivCopy(curr_weight);//take memory 1617 1629 goto FINISH; 1618 1630 } 1619 1631 1620 /* define the target vector as the next weight vector, if t = 1 */ 1632 /* define the target vector as the next weight vector, if t = 1 */ 1621 1633 if(mpz_cmp_si(t_nenner, 1)==0 && mpz_cmp_si(t_zaehler,1)==0) 1622 1634 { 1623 delete diff_weight; 1635 delete diff_weight; 1624 1636 diff_weight = ivCopy(target_weight); //this takes memory 1625 1637 goto FINISH; … … 1629 1641 //14.08.03 simplify the both vectors curr_weight and diff_weight (C-int) 1630 1642 gcd_tmp = (*curr_weight)[0]; 1631 1643 1632 1644 for (j=1; j<nRing; j++) 1633 1645 { … … 1636 1648 break; 1637 1649 } 1638 1650 1639 1651 if(gcd_tmp != 1) 1640 1652 for (j=0; j<nRing; j++) … … 1658 1670 PrintS("\n// t_zaehler: "); mpz_out_str( stdout, 10, t_zaehler); 1659 1671 PrintS(", t_nenner: "); mpz_out_str( stdout, 10, t_nenner); 1660 #endif 1672 #endif 1661 1673 1662 1674 mpz_t ddf; mpz_init(ddf); 1663 1675 mpz_t dcw; mpz_init(dcw); 1664 1676 BOOLEAN isdwpos; 1665 1677 1666 1678 // construct a new weight vector 1667 1679 for (j=0; j<nRing; j++) 1668 { 1680 { 1669 1681 mpz_set_si(dcw, (*curr_weight)[j]); 1670 1682 mpz_mul(s_nenner, t_nenner, dcw); 1671 1683 1672 1684 if( (*diff_weight)[j]>0) 1673 mpz_mul_ui(s_zaehler, t_zaehler, (*diff_weight)[j]); 1685 mpz_mul_ui(s_zaehler, t_zaehler, (*diff_weight)[j]); 1674 1686 else 1675 1687 { 1676 mpz_mul_ui(s_zaehler, t_zaehler, -(*diff_weight)[j]); 1677 mpz_neg(s_zaehler, s_zaehler); 1678 } 1679 1680 mpz_add(sntz, s_nenner, s_zaehler); 1688 mpz_mul_ui(s_zaehler, t_zaehler, -(*diff_weight)[j]); 1689 mpz_neg(s_zaehler, s_zaehler); 1690 } 1691 1692 mpz_add(sntz, s_nenner, s_zaehler); 1681 1693 1682 1694 mpz_init_set(vec[j], sntz); … … 1685 1697 Print("\n// j = %d ==> ", j); 1686 1698 PrintS("("); 1687 mpz_out_str( stdout, 10, t_nenner); 1688 Print(" * %d)", (*curr_weight)[j]); 1689 Print(" + ("); mpz_out_str( stdout, 10, t_zaehler); 1699 mpz_out_str( stdout, 10, t_nenner); 1700 Print(" * %d)", (*curr_weight)[j]); 1701 Print(" + ("); mpz_out_str( stdout, 10, t_zaehler); 1690 1702 Print(" * %d) = ", (*diff_weight)[j]); 1691 1703 mpz_out_str( stdout, 10, s_nenner); 1692 1704 PrintS(" + "); 1693 1705 mpz_out_str( stdout, 10, s_zaehler); 1694 PrintS(" = "); mpz_out_str( stdout, 10, sntz); 1706 PrintS(" = "); mpz_out_str( stdout, 10, sntz); 1695 1707 Print(" ==> vector[%d]: ", j); mpz_out_str(stdout, 10, vec[j]); 1696 #endif 1697 1698 if(j==0) 1708 #endif 1709 1710 if(j==0) 1699 1711 mpz_init_set(ggt, sntz); 1700 else 1712 else 1701 1713 if(mpz_cmp_si(ggt,1) != 0) 1702 mpz_gcd(ggt, ggt, sntz); 1703 1714 mpz_gcd(ggt, ggt, sntz); 1715 1704 1716 } 1705 1717 1706 1718 #ifdef NEXT_VECTORS_CC 1707 1719 PrintS("\n// gcd of elements of the vector: "); 1708 mpz_out_str( stdout, 10, ggt); 1709 #endif 1710 1720 mpz_out_str( stdout, 10, ggt); 1721 #endif 1722 1711 1723 mpz_t omega; 1712 1724 mpz_t sing_int; 1713 1725 mpz_init_set_ui(sing_int, 2147483647); 1714 1715 /* construct a new weight vector and check whether vec[j] is overflow!! 1726 1727 /* construct a new weight vector and check whether vec[j] is overflow!! 1716 1728 i.e. vec[j] > 2^31. 1717 1729 If vec[j] doesn't overflow, define a weight vector 1718 otherwise, report that overflow appears. 1719 In the second case test whether the defined new vector correct is 1730 otherwise, report that overflow appears. 1731 In the second case test whether the defined new vector correct is 1720 1732 plays an important rolle */ 1721 1733 … … 1728 1740 mpz_divexact(vec[j], vec[j], ggt); 1729 1741 (*diff_weight)[j] = mpz_get_si(vec[j]); 1730 } 1731 1732 if(mpz_cmp(vec[j], sing_int)>=0) 1742 } 1743 1744 if(mpz_cmp(vec[j], sing_int)>=0) 1733 1745 if(Overflow_Error == FALSE) 1734 1746 { 1735 1736 1737 1738 1739 1740 Print("\n// So vector[%d] := %d is wrong!!\n",j+1, (*diff_weight)[j]); 1741 } 1742 } 1743 1747 Overflow_Error = TRUE; 1748 1749 PrintS("\n// ** OVERFLOW in \"NextVector\": "); 1750 mpz_out_str( stdout, 10, vec[j]); 1751 PrintS(" is greater than 2147483647 (max. integer representation)"); 1752 Print("\n// So vector[%d] := %d is wrong!!\n",j+1, (*diff_weight)[j]); 1753 } 1754 } 1755 1744 1756 FINISH: 1745 1757 … … 1752 1764 mpz_clear(deg_w0_p1); 1753 1765 mpz_clear(deg_d0_p1); 1766 omFree(vec); 1754 1767 1755 1768 if(Overflow_Error == FALSE) … … 1759 1772 } 1760 1773 1761 /* 1762 compute an intermediate weight vector from iva to ivb w.r.t. 1774 /* 1775 compute an intermediate weight vector from iva to ivb w.r.t. 1763 1776 the reduced Groebner basis G. 1764 Return NULL, if it is equal to iva or iva = avb. 1777 Return NULL, if it is equal to iva or iva = avb. 1765 1778 */ 1766 1779 intvec* MkInterRedNextWeight(intvec* iva, intvec* ivb, ideal G) … … 1782 1795 return tmp; 1783 1796 } 1784 1797 1785 1798 delete tmp; 1786 return result; 1799 return result; 1787 1800 } 1788 1801 … … 1794 1807 //3.11.01 1795 1808 1796 if (ppNoether!=NULL) 1797 pDelete(&ppNoether); 1798 1799 if (((sLastPrinted.rtyp>BEGIN_RING) && (sLastPrinted.rtyp<END_RING)) || 1800 ((sLastPrinted.rtyp==LIST_CMD)&&(lRingDependend((lists)sLastPrinted.data)))) 1801 1802 { 1803 sLastPrinted.CleanUp(); 1804 memset(&sLastPrinted,0,sizeof(sleftv)); 1805 } 1809 if (ppNoether!=NULL) 1810 pDelete(&ppNoether); 1811 1812 if (((sLastPrinted.rtyp>BEGIN_RING) && (sLastPrinted.rtyp<END_RING)) || 1813 ((sLastPrinted.rtyp==LIST_CMD)&&(lRingDependend((lists)sLastPrinted.data)))) 1814 1815 { 1816 sLastPrinted.CleanUp(); 1817 memset(&sLastPrinted,0,sizeof(sleftv)); 1818 } 1806 1819 1807 1820 ring r = IDRING(tmp); … … 1834 1847 1835 1848 /* ringorder a for the first block: var 1..nv */ 1836 r->order[0] = ringorder_a; 1849 r->order[0] = ringorder_a; 1837 1850 r->block0[0] = 1; 1838 1851 r->block1[0] = nv; 1839 1852 1840 1853 /* ringorder lp for the second block: var 1..nv */ 1841 1854 r->order[1] = ringorder_lp; 1842 1855 r->block0[1] = 1; 1843 1856 r->block1[1] = nv; 1844 1845 /* ringorder C for the third block */ 1846 // it is very important within "idLift", 1857 1858 /* ringorder C for the third block */ 1859 // it is very important within "idLift", 1847 1860 // especially, by ring syz_ring=rCurrRingAssure_SyzComp(); 1848 1861 // therefore, nb must be (nBlocks(currRing) + 1) 1849 r->order[2] = ringorder_C; 1862 r->order[2] = ringorder_C; 1850 1863 1851 1864 /* the last block: everything is 0 */ … … 1857 1870 /* complete ring intializations */ 1858 1871 rComplete(r); 1859 1872 1860 1873 rChangeCurrRing(r); 1861 1874 currRingHdl = tmp; … … 1867 1880 { 1868 1881 idhdl tmp = enterid(IDID(currRingHdl),IDLEV(currRingHdl)+1,RING_CMD,&IDROOT,TRUE); 1869 1870 1871 if (ppNoether!=NULL) 1872 pDelete(&ppNoether); 1873 1874 if (((sLastPrinted.rtyp>BEGIN_RING) && (sLastPrinted.rtyp<END_RING)) || 1875 ((sLastPrinted.rtyp==LIST_CMD)&&(lRingDependend((lists)sLastPrinted.data)))) 1876 1877 { 1878 sLastPrinted.CleanUp(); 1879 memset(&sLastPrinted,0,sizeof(sleftv)); 1880 } 1882 1883 1884 if (ppNoether!=NULL) 1885 pDelete(&ppNoether); 1886 1887 if (((sLastPrinted.rtyp>BEGIN_RING) && (sLastPrinted.rtyp<END_RING)) || 1888 ((sLastPrinted.rtyp==LIST_CMD)&&(lRingDependend((lists)sLastPrinted.data)))) 1889 1890 { 1891 sLastPrinted.CleanUp(); 1892 memset(&sLastPrinted,0,sizeof(sleftv)); 1893 } 1881 1894 1882 1895 ring r = IDRING(tmp); … … 1897 1910 1898 1911 /*weights: entries for 3 blocks: NULL Made:???*/ 1899 1912 1900 1913 r->wvhdl = (int **)omAlloc0(nb * sizeof(int_ptr)); 1901 1914 … … 1909 1922 r->block0[0] = 1; 1910 1923 r->block1[0] = nv; 1911 1912 /* ringorder C for the second block */ 1913 r->order[1] = ringorder_C; 1924 1925 /* ringorder C for the second block */ 1926 r->order[1] = ringorder_C; 1914 1927 1915 1928 /* the last block: everything is 0 */ … … 1921 1934 /* complete ring intializations */ 1922 1935 rComplete(r); 1923 1936 1924 1937 //rSetHdl(tmp); 1925 1938 … … 1935 1948 int i, nv = currRing->N; 1936 1949 int nb = rBlocks(currRing) + 1; 1937 1950 1938 1951 ring res=(ring)omAllocBin(ip_sring_bin); 1939 1952 1940 1953 memcpy4(res,currRing,sizeof(ip_sring)); 1941 1954 1942 1955 res->VarOffset = NULL; 1943 1956 res->ref=0; … … 1950 1963 int l=rPar(currRing); 1951 1964 res->parameter=(char **)omAlloc(l*sizeof(char_ptr)); 1952 1965 1953 1966 for(i=l-1;i>=0;i--) 1954 1967 res->parameter[i]=omStrDup(currRing->parameter[i]); 1955 1968 } 1956 1969 1957 intvec* iva = va; 1970 intvec* iva = va; 1958 1971 1959 1972 /*weights: entries for 3 blocks: NULL Made:???*/ … … 1967 1980 res->block0 = (int *)omAlloc0(nb * sizeof(int *)); 1968 1981 res->block1 = (int *)omAlloc0(nb * sizeof(int *)); 1969 1982 1970 1983 /* ringorder a for the first block: var 1..nv */ 1971 res->order[0] = ringorder_a; 1984 res->order[0] = ringorder_a; 1972 1985 res->block0[0] = 1; 1973 1986 res->block1[0] = nv; 1974 1987 1975 1988 /* ringorder lp for the second block: var 1..nv */ 1976 1989 res->order[1] = ringorder_lp; 1977 1990 res->block0[1] = 1; 1978 1991 res->block1[1] = nv; 1979 1980 /* ringorder C for the third block */ 1981 // it is very important within "idLift", 1992 1993 /* ringorder C for the third block */ 1994 // it is very important within "idLift", 1982 1995 // especially, by ring syz_ring=rCurrRingAssure_SyzComp(); 1983 1996 // therefore, nb must be (nBlocks(currRing) + 1) 1984 res->order[2] = ringorder_C; 1997 res->order[2] = ringorder_C; 1985 1998 1986 1999 /* the last block: everything is 0 */ … … 1990 2003 res->OrdSgn = 1; 1991 2004 1992 2005 1993 2006 res->names = (char **)omAlloc0(nv * sizeof(char_ptr)); 1994 2007 for (i=nv-1; i>=0; i--) … … 1997 2010 /* complete ring intializations */ 1998 2011 rComplete(res); 1999 2000 2012 2013 2001 2014 // clean up history 2002 2015 if (sLastPrinted.RingDependend()) … … 2005 2018 memset(&sLastPrinted,0,sizeof(sleftv)); 2006 2019 } 2007 2020 2008 2021 2009 2022 /* execute the created ring */ … … 2017 2030 2018 2031 ring r=(ring)omAllocBin(ip_sring_bin); 2019 2032 2020 2033 memcpy4(r,currRing,sizeof(ip_sring)); 2021 2034 2022 2035 r->VarOffset = NULL; 2023 2036 r->ref=0; … … 2030 2043 int l=rPar(currRing); 2031 2044 r->parameter=(char **)omAlloc(l*sizeof(char_ptr)); 2032 2045 2033 2046 for(i=l-1;i>=0;i--) 2034 2047 r->parameter[i]=omStrDup(currRing->parameter[i]); … … 2050 2063 2051 2064 /*weights: entries for 3 blocks: NULL Made:???*/ 2052 2065 2053 2066 r->wvhdl = (int **)omAlloc0(nb * sizeof(int_ptr)); 2054 2067 … … 2062 2075 r->block0[0] = 1; 2063 2076 r->block1[0] = nv; 2064 2065 /* ringorder C for the second block */ 2066 r->order[1] = ringorder_C; 2077 2078 /* ringorder C for the second block */ 2079 r->order[1] = ringorder_C; 2067 2080 2068 2081 /* the last block: everything is 0 */ … … 2078 2091 int l=rPar(currRing); 2079 2092 r->parameter=(char **)omAlloc(l*sizeof(char_ptr)); 2080 2093 2081 2094 for(i=l-1;i>=0;i--) 2082 2095 r->parameter[i]=omStrDup(currRing->parameter[i]); … … 2085 2098 /* complete ring intializations */ 2086 2099 rComplete(r); 2087 2100 2088 2101 // clean up history 2089 2102 if (sLastPrinted.RingDependend()) … … 2104 2117 if((* hilb)[i]==0) 2105 2118 return 1; 2106 2119 2107 2120 return 0; 2108 2121 } … … 2119 2132 **************************************************************************/ 2120 2133 2121 static ideal LastGB(ideal G, intvec* curr_weight,int tp_deg) 2134 static ideal LastGB(ideal G, intvec* curr_weight,int tp_deg) 2122 2135 { 2123 2136 BOOLEAN nError = Overflow_Error; … … 2132 2145 intvec* target_weight; 2133 2146 intvec* iv_lp = Mivlp(nV); //define (1,0,...,0) 2134 intvec* pert_target_vector; 2147 intvec* pert_target_vector; 2135 2148 intvec* ivNull = new intvec(nV); 2136 2149 intvec* extra_curr_weight = new intvec(nV); … … 2172 2185 2173 2186 while(1) 2174 { 2187 { 2175 2188 nwalk++; 2176 2189 nstep++; … … 2179 2192 next_weight = MkInterRedNextWeight(curr_weight,target_weight, G); 2180 2193 xtnw=xtnw+clock()-to; 2181 #ifdef PRINT_VECTORS 2194 #ifdef PRINT_VECTORS 2182 2195 MivString(curr_weight, target_weight, next_weight); 2183 2196 #endif … … 2189 2202 nlast = 1; 2190 2203 delete next_weight; 2191 2204 2192 2205 //idElements(G, "G"); 2193 2206 //Print("\n// ring r%d_%d = %s;\n", tp_deg, nwalk, rString(currRing)); 2194 2207 2195 2208 break; 2196 2209 } … … 2219 2232 xtif=xtif+clock()-to; 2220 2233 2221 #ifdef ENDWALKS 2222 if(endwalks == 1){ 2223 Print("\n// ring r%d_%d = %s;\n", tp_deg, nwalk, rString(currRing)); 2224 idElements(Gomega, "Gw"); 2225 headidString(Gomega, "Gw"); 2226 } 2227 #endif 2228 2234 #ifdef ENDWALKS 2235 if(endwalks == 1){ 2236 Print("\n// ring r%d_%d = %s;\n", tp_deg, nwalk, rString(currRing)); 2237 idElements(Gomega, "Gw"); 2238 headidString(Gomega, "Gw"); 2239 } 2240 #endif 2241 2229 2242 #ifndef BUCHBERGER_ALG 2230 2243 if(isNolVector(curr_weight) == 0) … … 2233 2246 hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing); 2234 2247 #endif // BUCHBERGER_ALG 2235 2248 2236 2249 /* define a new ring that its ordering is "(a(curr_weight),lp) */ 2237 2250 //..25.03.03 VMrDefault(curr_weight); … … 2241 2254 VMrDefault(curr_weight); 2242 2255 2243 newRing = currRing; 2256 newRing = currRing; 2244 2257 Gomega1 = idrMoveR(Gomega, oldRing); 2245 2258 … … 2250 2263 #else 2251 2264 M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight); 2252 delete hilb_func; 2265 delete hilb_func; 2253 2266 #endif // BUCHBERGER_ALG 2254 2267 xtstd=xtstd+clock()-to; 2255 /* change the ring to oldRing */ 2268 /* change the ring to oldRing */ 2256 2269 rChangeCurrRing(oldRing); 2257 2270 M1 = idrMoveR(M, newRing); 2258 Gomega2 = idrMoveR(Gomega1, newRing); 2259 2271 Gomega2 = idrMoveR(Gomega1, newRing); 2272 2260 2273 to=clock(); 2261 2274 /* compute a reduced Groebner basis of <G> w.r.t. "newRing" */ … … 2263 2276 xtlift=xtlift+clock()-to; 2264 2277 2265 idDelete(&M1); 2278 idDelete(&M1); 2266 2279 idDelete(&G); 2267 2280 2268 /* change the ring to newRing */ 2281 /* change the ring to newRing */ 2269 2282 rChangeCurrRing(newRing); 2270 2283 F1 = idrMoveR(F, oldRing); 2271 2284 2272 2285 to=clock(); 2273 /* reduce the Groebner basis <G> w.r.t. new ring */ 2286 /* reduce the Groebner basis <G> w.r.t. new ring */ 2274 2287 G = kInterRedCC(F1, NULL); 2275 2288 xtred=xtred+clock()-to; 2276 2289 idDelete(&F1); 2277 2290 2278 2291 if(endwalks == 1){ 2279 2292 //Print("\n// ring r%d_%d = %s;\n", tp_deg, nwalk, rString(currRing)); … … 2295 2308 2296 2309 F1 = idrMoveR(G, newRing); 2297 2310 2298 2311 if(nnwinC == 0 || test_w_in_ConeCC(F1, pert_target_vector) != 1) 2299 2312 { … … 2301 2314 rChangeCurrRing(newRing); 2302 2315 G = idrMoveR(F1, oldRing); 2303 Print("\n// takes %d steps and calls the recursion of level %d:", 2304 2316 Print("\n// takes %d steps and calls the recursion of level %d:", 2317 nwalk, tp_deg-1); 2305 2318 2306 2319 F1 = LastGB(G,curr_weight, tp_deg-1); 2307 2320 } 2308 2321 2309 2322 TargetRing = currRing; 2310 2323 rChangeCurrRing(EXXRing); … … 2314 2327 { 2315 2328 if(nlast == 1) 2316 { 2329 { 2317 2330 //OMEGA_OVERFLOW_LASTGB: 2318 2331 /* 2319 if(MivSame(curr_weight, iv_lp) == 1) 2320 if (currRing->parameter != NULL) 2321 DefRingParlp(); 2322 else 2332 if(MivSame(curr_weight, iv_lp) == 1) 2333 if (currRing->parameter != NULL) 2334 DefRingParlp(); 2335 else 2323 2336 VMrDefaultlp(); 2324 else 2325 if (currRing->parameter != NULL) 2326 DefRingPar(curr_weight); 2327 else 2328 VMrDefault(curr_weight); 2337 else 2338 if (currRing->parameter != NULL) 2339 DefRingPar(curr_weight); 2340 else 2341 VMrDefault(curr_weight); 2329 2342 */ 2330 2331 2332 if (currRing->parameter != NULL) 2333 DefRingParlp(); 2334 else 2335 VMrDefaultlp(); 2336 2337 2343 2344 //..25.03.03 VMrDefaultlp();//define and execute the ring "lp" 2345 if (currRing->parameter != NULL) 2346 DefRingParlp(); 2347 else 2348 VMrDefaultlp(); 2349 2350 2338 2351 F1 = idrMoveR(G, newRing); 2339 2352 //Print("\n// Apply \"std\" in ring r%d_%d = %s;\n", tp_deg, nwalk, rString(currRing)); … … 2346 2359 rChangeCurrRing(EXXRing); 2347 2360 result = idrMoveR(G, newRing); 2348 } 2361 } 2349 2362 delete target_weight; 2350 2363 delete last_omega; … … 2354 2367 Overflow_Error = nError; 2355 2368 2356 return(result); 2369 return(result); 2357 2370 } 2358 2371 … … 2363 2376 int i; 2364 2377 for(i=IDELEMS(G)-1; i>=0; i--) 2365 #if 0 2378 #if 0 2366 2379 if(pLength(G->m[i])>2) 2367 2380 return 1; … … 2373 2386 //&& (G->m[i]->next->next->next->next!=NULL) /* len >=4 */ 2374 2387 ) return 1; 2375 #endif 2388 #endif 2376 2389 return 0; 2377 2390 } … … 2392 2405 2393 2406 2394 /* Implementation of the improved Groebner walk algorithm which is written 2407 /* Implementation of the improved Groebner walk algorithm which is written 2395 2408 by Quoc-Nam Tran (2000). 2396 One perturbs the original target weight vector, only if 2397 the next intermediate weight vector is equal to the current target weight 2409 One perturbs the original target weight vector, only if 2410 the next intermediate weight vector is equal to the current target weight 2398 2411 vector. This must be repeated until the wanted reduced Groebner basis 2399 2412 to reach. 2400 If the numbers of variables is big enough, the representation of the origin 2401 weight vector may be very big. Therefore, it is possible the intermediate 2402 weight vector doesn't stay in the correct Groebner cone. 2403 In this case we have just a reduced Groebner basis of the given ideal 2404 with respect to another monomial order. Then we have to compute 2413 If the numbers of variables is big enough, the representation of the origin 2414 weight vector may be very big. Therefore, it is possible the intermediate 2415 weight vector doesn't stay in the correct Groebner cone. 2416 In this case we have just a reduced Groebner basis of the given ideal 2417 with respect to another monomial order. Then we have to compute 2405 2418 a wanted reduced Groebner basis of it with respect to the given order. 2406 2419 At the following subroutine we use the improved Buchberger algorithm or … … 2414 2427 { 2415 2428 int i, nH =IDELEMS(h); 2416 2429 2417 2430 ideal m = idInit(nH,h->rank); 2418 2431 2419 2432 for (i=nH-1;i>=0; i--) 2420 2433 { 2421 if (h->m[i]!=NULL) 2434 if (h->m[i]!=NULL) 2422 2435 m->m[i]=pHead(h->m[i]); 2423 2436 } … … 2432 2445 if(nG != IDELEMS(H1)) 2433 2446 return 0; 2434 2447 2435 2448 for(i=nG-1; i>=0; i--) 2436 { 2437 #if 0 2449 { 2450 #if 0 2438 2451 poly t; 2439 2452 if((t=pSub(pCopy(H0->m[i]), pCopy(H1->m[i]))) != NULL) … … 2446 2459 if(!pEqualPolys(H0->m[i],H1->m[i])) 2447 2460 return 0; 2448 #endif 2461 #endif 2449 2462 } 2450 2463 … … 2460 2473 int degtmp, result = 0; 2461 2474 intvec* ivUnit = Mivdp(nV); 2462 2475 2463 2476 for(i=nG-1; i>=0; i--) 2464 2477 { … … 2472 2485 } 2473 2486 2474 /* perturb the weight vector iva w.r.t. the ideal G. 2487 /* perturb the weight vector iva w.r.t. the ideal G. 2475 2488 the monomial order of the current ring is the w_1 weight lex. order. 2476 2489 define w := d^(n-1)w_1+ d^(n-2)w_2, ...+ dw_(n-1)+ w_n … … 2490 2503 // define the sequence which expresses the current monomial ordering 2491 2504 // w_1 = iva; w_2 = (1,0,..,0); w_n = (0,...,0,1,0) 2492 intvec* ivMat = MivMatrixOrder(iva); 2505 intvec* ivMat = MivMatrixOrder(iva); 2493 2506 2494 2507 int mtmp, m=(*iva)[0]; … … 2497 2510 { 2498 2511 mtmp = (*ivMat)[i]; 2499 2512 2500 2513 if(mtmp <0) mtmp = -mtmp; 2501 2514 … … 2503 2516 m = mtmp; 2504 2517 } 2505 2518 2506 2519 /* define the maximal total degree of polynomials of G */ 2507 2520 mpz_t ndeg; … … 2520 2533 //ndeg = (2*maxdeg*maxdeg + (nV+1)*maxdeg)*m;//Kalkbrenner (1999) 2521 2534 mpz_pow_ui(ztmp, maxdeg, 2); 2522 mpz_mul_ui(ztmp, ztmp, 2); 2535 mpz_mul_ui(ztmp, ztmp, 2); 2523 2536 mpz_mul_ui(maxdeg, maxdeg, nV+1); 2524 2537 mpz_add(ndeg, ztmp, maxdeg); 2525 2538 mpz_mul_ui(ndeg, ndeg, m); 2526 2539 2527 2540 //PrintS("\n// with the new upper degree bound (2d^2+(n+1)d)*m "); 2528 2541 //Print("\n// where d = %d, n = %d and bound = %d", maxdeg, nV, ndeg); … … 2530 2543 2531 2544 /* 29.08.03*/ 2532 #ifdef INVEPS_SMALL_IN_TRAN 2545 #ifdef INVEPS_SMALL_IN_TRAN 2533 2546 if(mpz_cmp_ui(ndeg, nV)>0 && nV > 3) 2534 2547 mpz_cdiv_q_ui(ndeg, ndeg, nV); … … 2539 2552 mpz_t deg_tmp; 2540 2553 mpz_init_set(deg_tmp, ndeg); 2541 2542 mpz_t ivres[nV];2554 2555 mpz_t *ivres=( mpz_t *) omAlloc(nV*sizeof(mpz_t)); 2543 2556 mpz_init_set_si(ivres[nV-1],1); 2544 2557 2545 2558 for(i=nV-2; i>=0; i--) 2546 2559 { … … 2549 2562 } 2550 2563 2551 mpz_t ivtmp[nV];2564 mpz_t *ivtmp=(mpz_t *)omAlloc(nV*sizeof(mpz_t)); 2552 2565 for(i=0; i<nV; i++) 2553 2566 mpz_init(ivtmp[i]); 2554 2567 2555 2568 mpz_t sing_int; 2556 2569 mpz_init_set_ui(sing_int, 2147483647); 2557 2570 2558 2571 intvec* repr_vector = new intvec(nV); 2559 2572 2560 2573 /* define ivtmp := ndeg^(n-1).w_1 + ndeg^(n-2).w_2 + ... + w_n */ 2561 2574 for(i=0; i<nV; i++) … … 2584 2597 { 2585 2598 Overflow_Error = TRUE; 2586 2599 2587 2600 PrintS("\n// ** OVERFLOW in \"Repr.Vector\": "); 2588 2601 mpz_out_str( stdout, 10, ivtmp[i]); … … 2592 2605 } 2593 2606 } 2594 if(Overflow_Error == TRUE) 2607 if(Overflow_Error == TRUE) 2595 2608 { 2596 2609 ivString(repr_vector, "repvector"); … … 2598 2611 } 2599 2612 2600 if(Overflow_Error == FALSE) 2613 if(Overflow_Error == FALSE) 2601 2614 Overflow_Error=nError; 2602 2615 2616 omFree(ivres); 2617 omFree(ivtmp); 2603 2618 return repr_vector; 2604 2619 } … … 2630 2645 //ndeg = (2*maxdeg*maxdeg + (nV+1)*maxdeg);//Kalkbrenner (1999) 2631 2646 mpz_pow_ui(ztmp, maxdeg, 2); 2632 mpz_mul_ui(ztmp, ztmp, 2); 2647 mpz_mul_ui(ztmp, ztmp, 2); 2633 2648 mpz_mul_ui(maxdeg, maxdeg, nV+1); 2634 2649 mpz_add(ndeg, ztmp, maxdeg); 2635 2650 /* 2636 2651 PrintS("\n// with the new upper degree bound (2d^2+(n+1)d)*m "); 2637 Print("\n// where d = %d, n = %d and bound = %d", 2652 Print("\n// where d = %d, n = %d and bound = %d", 2638 2653 mpz_get_si(maxdeg), nV, mpz_get_si(ndeg)); 2639 2654 */ 2640 2655 #endif //UPPER_BOUND 2641 2656 2642 #ifdef INVEPS_SMALL_IN_TRAN 2657 #ifdef INVEPS_SMALL_IN_TRAN 2643 2658 if(mpz_cmp_ui(ndeg, nV)>0 && nV > 3) 2644 2659 mpz_cdiv_q_ui(ndeg, ndeg, nV); … … 2650 2665 mpz_t deg_tmp; 2651 2666 mpz_init_set(deg_tmp, ndeg); 2652 2653 mpz_t ivres[nV];2667 2668 mpz_t *ivres=(mpz_t *)omAlloc(nV*sizeof(mpz_t)); 2654 2669 mpz_init_set_si(ivres[nV-1], 1); 2655 2670 2656 2671 for(i=nV-2; i>=0; i--) 2657 2672 { … … 2680 2695 Print("\n// So vector[%d] := %d is wrong!!\n",i+1,(*repr_vector)[i]); 2681 2696 } 2682 } 2683 } 2684 if(Overflow_Error == TRUE) 2697 } 2698 } 2699 if(Overflow_Error == TRUE) 2685 2700 { 2686 2701 ivString(repr_vector, "repvector"); 2687 2702 Print("\n// %d element(s) of it are overflow!!", ntrue); 2688 2703 } 2689 if(Overflow_Error == FALSE) 2704 if(Overflow_Error == FALSE) 2690 2705 Overflow_Error = nError; 2691 2706 2707 omFree(ivres); 2692 2708 return repr_vector; 2693 2709 } … … 2705 2721 intvec* ivUnit = Mivdp(nV); 2706 2722 int degtmp, maxdeg = 0; 2707 2723 2708 2724 for(i=IDELEMS(G)-1; i>=0; i--) 2709 2725 { … … 2716 2732 mpz_t ztmp; 2717 2733 mpz_init_set_si(ztmp, maxdeg); 2718 mpz_t ivres[nV];2734 mpz_t *ivres=(mpz_t *)omAlloc(nV*sizeof(mpz_t)); 2719 2735 mpz_init_set_si(ivres[nV-1], 1); // (*ivres)[nV-1] = 1; 2720 2736 … … 2725 2741 } 2726 2742 2727 mpz_t ivtmp[nV];2743 mpz_t *ivtmp=(mpz_t*)omAlloc(nV*sizeof(mpz_t)); 2728 2744 for(i=0; i<nV; i++) 2729 2745 mpz_init(ivtmp[i]); … … 2738 2754 mpz_neg(ztmp, ztmp); 2739 2755 } 2740 else 2756 else 2741 2757 mpz_mul_ui(ztmp, ivres[i], (*M)[i*nV+j]); 2742 2758 2743 mpz_add(ivtmp[j], ivtmp[j], ztmp); 2759 mpz_add(ivtmp[j], ivtmp[j], ztmp); 2744 2760 } 2745 2761 … … 2762 2778 PrintS(" is greater than 2147483647 (max. integer representation)"); 2763 2779 Print("\n// So vector[%d] := %d is wrong!!\n",i+1,(*repvector)[i]); 2764 } 2765 } 2766 } 2767 if(Overflow_Error == TRUE) 2780 } 2781 } 2782 } 2783 if(Overflow_Error == TRUE) 2768 2784 { 2769 2785 ivString(repvector, "repvector"); … … 2771 2787 } 2772 2788 2773 if(Overflow_Error == FALSE) 2789 if(Overflow_Error == FALSE) 2774 2790 Overflow_Error = nError; 2775 2791 2792 mpz_clear(ztmp); 2793 omFree(ivtmp); 2794 omFree(ivres); 2776 2795 return repvector; 2777 2796 } … … 2781 2800 2782 2801 2783 /* The following subroutine is the implementation of our first improved 2802 /* The following subroutine is the implementation of our first improved 2784 2803 Groebner walk algorithm, i.e. the first altervative algorithm. 2785 2804 First we use the Grobner walk algorithm and then we call the changed 2786 2805 perturbation walk algorithm with decreased degree, if an intermediate 2787 weight vector is equal to the current target weight vector. 2788 This call will be only repeated until we get the wanted reduced Groebner 2806 weight vector is equal to the current target weight vector. 2807 This call will be only repeated until we get the wanted reduced Groebner 2789 2808 basis or n times, where n is the numbers of variables. 2790 2809 */ … … 2805 2824 /* 7 Februar 2002 */ 2806 2825 /* npwinc = 0, if curr_weight doesn't stay in the correct Groebner cone */ 2807 static ideal Rec_LastGB(ideal G, intvec* curr_weight, 2808 intvec* orig_target_weight, int tp_deg, int npwinc) 2826 static ideal Rec_LastGB(ideal G, intvec* curr_weight, 2827 intvec* orig_target_weight, int tp_deg, int npwinc) 2809 2828 { 2810 2829 BOOLEAN nError = Overflow_Error; … … 2829 2848 2830 2849 //08 Juli 03 2831 intvec* hilb_func; 2850 intvec* hilb_func; 2832 2851 /* to avoid (1,0,...,0) as the target vector */ 2833 2852 intvec* last_omega = new intvec(nV); … … 2856 2875 cone(k-1) is equal to cone(k) */ 2857 2876 if(test_G_GB_walk(H0_tmp,H1)==1) { 2858 idDelete(&H0_tmp); 2877 idDelete(&H0_tmp); 2859 2878 idDelete(&H1); 2860 2879 G = ssG; … … 2862 2881 newRing = currRing; 2863 2882 delete ivNull; 2864 2883 2865 2884 if(npwinc != 0) 2866 2885 goto LastGB_Finish; … … 2870 2889 } 2871 2890 } 2872 idDelete(&H0_tmp); 2891 idDelete(&H0_tmp); 2873 2892 idDelete(&H1); 2874 2893 2875 2894 iv_M_lp = MivMatrixOrderlp(nV); 2876 2895 target_weight = MPertVectors(ssG, iv_M_lp, tp_deg); … … 2881 2900 2882 2901 if(Overflow_Error == TRUE) { 2883 nOverflow_Error = Overflow_Error; 2902 nOverflow_Error = Overflow_Error; 2884 2903 NEG = 1; 2885 2904 newRing = currRing; … … 2889 2908 2890 2909 while(1) 2891 { 2910 { 2892 2911 nwalk ++; 2893 2912 nstep++; 2894 2913 2895 if(nwalk==1) 2914 if(nwalk==1) 2896 2915 goto FIRST_STEP; 2897 2916 … … 2909 2928 2910 2929 oldRing = currRing; 2911 2930 2912 2931 /* defiNe a new ring that its ordering is "(a(curr_weight),lp) */ 2913 2932 if (currRing->parameter != NULL) … … 2916 2935 VMrDefault(curr_weight); 2917 2936 2918 newRing = currRing; 2937 newRing = currRing; 2919 2938 Gomega1 = idrMoveR(Gomega, oldRing); 2920 2939 to=clock(); … … 2924 2943 #else 2925 2944 M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight); 2926 delete hilb_func; 2945 delete hilb_func; 2927 2946 #endif // BUCHBERGER_ALG 2928 2947 xtstd=xtstd+clock()-to; 2929 /* change the ring to oldRing */ 2948 /* change the ring to oldRing */ 2930 2949 rChangeCurrRing(oldRing); 2931 2950 M1 = idrMoveR(M, newRing); 2932 Gomega2 = idrMoveR(Gomega1, newRing); 2933 2951 Gomega2 = idrMoveR(Gomega1, newRing); 2952 2934 2953 to=clock(); 2935 /* compute a reduced Groebner basis of <G> w.r.t. "newRing" by the 2954 /* compute a reduced Groebner basis of <G> w.r.t. "newRing" by the 2936 2955 lifting process*/ 2937 2956 F = MLifttwoIdeal(Gomega2, M1, G); 2938 2957 xtlift=xtlift+clock()-to; 2939 idDelete(&M1); 2940 idDelete(&Gomega2); 2958 idDelete(&M1); 2959 idDelete(&Gomega2); 2941 2960 idDelete(&G); 2942 2943 /* change the ring to newRing */ 2961 2962 /* change the ring to newRing */ 2944 2963 rChangeCurrRing(newRing); 2945 2964 F1 = idrMoveR(F, oldRing); 2946 2965 2947 2966 to=clock(); 2948 /* reduce the Groebner basis <G> w.r.t. new ring */ 2967 /* reduce the Groebner basis <G> w.r.t. new ring */ 2949 2968 G = kInterRedCC(F1, NULL); 2950 2969 xtred=xtred+clock()-to; 2951 2970 idDelete(&F1); 2952 2971 2953 2972 if(endwalks == 1) 2954 2973 break; … … 2960 2979 next_weight = MkInterRedNextWeight(curr_weight,target_weight, G); 2961 2980 xtnw=xtnw+clock()-to; 2962 #ifdef PRINT_VECTORS 2981 #ifdef PRINT_VECTORS 2963 2982 MivString(curr_weight, target_weight, next_weight); 2964 2983 #endif 2965 2984 2966 2985 if(Overflow_Error == TRUE) { 2967 2986 //PrintS("\n// ** The next vector does NOT stay in Cone!!\n"); … … 2973 2992 if(tp_deg == nV) 2974 2993 nlast = 1; 2975 2994 2976 2995 delete next_weight; 2977 2996 break; … … 2983 3002 break; 2984 3003 } 2985 3004 2986 3005 if(MivComp(next_weight, target_weight) == 1) { 2987 3006 if(tp_deg == nV) … … 2989 3008 else { 2990 3009 REC_LAST_GB_ALT2: 2991 2992 2993 2994 Print("\n// takes %d steps and calls \"Rec_LastGB\" (%d):", 2995 2996 3010 nOverflow_Error = Overflow_Error; 3011 tproc=tproc+clock()-tinput; 3012 /* 3013 Print("\n// takes %d steps and calls \"Rec_LastGB\" (%d):", 3014 nwalk, tp_deg+1); 3015 */ 2997 3016 G = Rec_LastGB(G,curr_weight, orig_target_weight, tp_deg+1,nnwinC); 2998 3017 newRing = currRing; … … 3001 3020 } 3002 3021 } 3003 3022 3004 3023 for(i=nV-1; i>=0; i--) { 3005 3024 //(*extra_curr_weight)[i] = (*curr_weight)[i]; 3006 3025 (*curr_weight)[i] = (*next_weight)[i]; 3007 } 3008 delete next_weight; 3009 }//while 3026 } 3027 delete next_weight; 3028 }//while 3010 3029 3011 3030 delete ivNull; … … 3020 3039 3021 3040 F1 = idrMoveR(G, newRing); 3022 3041 3023 3042 if(nnwinC == 0 || test_w_in_ConeCC(F1, target_weight) != 1 ) { 3024 3043 nOverflow_Error = Overflow_Error; 3025 3044 //Print("\n// takes %d steps and calls \"Rec_LastGB (%d):", tp_deg+1); 3026 3045 tproc=tproc+clock()-tinput; 3027 F1 = Rec_LastGB(F1,curr_weight, orig_target_weight, tp_deg+1,nnwinC); 3028 } 3046 F1 = Rec_LastGB(F1,curr_weight, orig_target_weight, tp_deg+1,nnwinC); 3047 } 3029 3048 delete target_weight; 3030 3049 … … 3034 3053 } 3035 3054 else { 3036 if(nlast == 1) { 3055 if(nlast == 1) { 3037 3056 JUNI_STD: 3038 3057 3039 3058 newRing = currRing; 3040 3059 if (currRing->parameter != NULL) … … 3054 3073 xtextra=xtextra+clock()-to; 3055 3074 3056 3075 3057 3076 idDelete(&F1); 3058 3077 newRing = currRing; 3059 3078 } 3060 3079 3061 3080 LastGB_Finish: 3062 3081 rChangeCurrRing(EXXRing); 3063 result = idrMoveR(G, newRing); 3064 } 3065 3066 if(Overflow_Error == FALSE) 3082 result = idrMoveR(G, newRing); 3083 } 3084 3085 if(Overflow_Error == FALSE) 3067 3086 Overflow_Error=nError; 3068 3087 /* 3069 3088 Print("\n// \"Rec_LastGB\" (%d) took %d steps and %.2f sec.Overflow_Error (%d)", tp_deg, 3070 3089 nwalk, ((double) tproc)/1000000, nOverflow_Error); 3071 3090 */ 3072 return(result); 3073 } 3074 3075 /* The following subroutine is the implementation of our second improved 3091 return(result); 3092 } 3093 3094 /* The following subroutine is the implementation of our second improved 3076 3095 Groebner walk algorithm, i.e. the second altervative algorithm. 3077 3096 First we use the Grobner walk algorithm and then we call the changed 3078 3097 perturbation walk algorithm with increased degree, if an intermediate 3079 weight vector is equal to the current target weight vector. 3080 This call will be only repeated until we get the wanted reduced Groebner 3098 weight vector is equal to the current target weight vector. 3099 This call will be only repeated until we get the wanted reduced Groebner 3081 3100 basis or n times, where n is the numbers of variables. 3082 3101 */ … … 3096 3115 int i, nV = currRing->N; 3097 3116 int nwalk=0, endwalks=0, nhilb=1; 3098 3117 3099 3118 ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G, G1; 3100 3119 ring endRing, newRing, oldRing; … … 3106 3125 3107 3126 ring XXRing = currRing; 3108 3127 3109 3128 //Print("\n// ring r_input = %s;", rString(currRing)); 3110 to = clock(); 3111 /* compute the reduced Groebner basis of the given ideal w.r.t. 3112 a "fast" monomial order, e.g. degree reverse lex. order (dp) */ 3129 to = clock(); 3130 /* compute the reduced Groebner basis of the given ideal w.r.t. 3131 a "fast" monomial order, e.g. degree reverse lex. order (dp) */ 3113 3132 G = MstdCC(Go); 3114 3133 tostd=clock()-to; 3115 3134 3116 3135 /* 3117 Print("\n// Computation of the first std took = %.2f sec", 3136 Print("\n// Computation of the first std took = %.2f sec", 3118 3137 ((double) tostd)/1000000); 3119 3138 */ 3120 3139 if(currRing->order[0] == ringorder_a) 3121 3140 goto NEXT_VECTOR; 3122 3141 3123 3142 while(1) 3124 { 3143 { 3125 3144 nwalk ++; 3126 3145 nstep ++; 3127 to = clock(); 3146 to = clock(); 3128 3147 /* compute an initial form ideal of <G> w.r.t. "curr_vector" */ 3129 3148 Gomega = MwalkInitialForm(G, curr_weight); … … 3137 3156 goto LAST_GB_ALT2; 3138 3157 } 3139 #endif 3158 #endif 3140 3159 oldRing = currRing; 3141 3160 3142 3161 /* define a new ring that its ordering is "(a(curr_weight),lp) */ 3143 3162 if (currRing->parameter != NULL) … … 3146 3165 VMrDefault(curr_weight); 3147 3166 3148 newRing = currRing; 3167 newRing = currRing; 3149 3168 Gomega1 = idrMoveR(Gomega, oldRing); 3150 to = clock(); 3169 to = clock(); 3151 3170 /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */ 3152 3171 M = MstdhomCC(Gomega1); 3153 3172 xtstd=xtstd+clock()-to; 3154 /* change the ring to oldRing */ 3173 /* change the ring to oldRing */ 3155 3174 rChangeCurrRing(oldRing); 3156 3175 M1 = idrMoveR(M, newRing); 3157 3176 Gomega2 = idrMoveR(Gomega1, newRing); 3158 3177 3159 to = clock(); 3160 /* compute the reduced Groebner basis of <G> w.r.t. "newRing" 3178 to = clock(); 3179 /* compute the reduced Groebner basis of <G> w.r.t. "newRing" 3161 3180 by the liftig process */ 3162 3181 F = MLifttwoIdeal(Gomega2, M1, G); 3163 3182 xtlift=xtlift+clock()-to; 3164 idDelete(&M1); 3165 idDelete(&Gomega2); 3166 idDelete(&G); 3167 3168 /* change the ring to newRing */ 3183 idDelete(&M1); 3184 idDelete(&Gomega2); 3185 idDelete(&G); 3186 3187 /* change the ring to newRing */ 3169 3188 rChangeCurrRing(newRing); 3170 F1 = idrMoveR(F, oldRing); 3171 3172 to = clock(); 3173 /* reduce the Groebner basis <G> w.r.t. newRing */ 3189 F1 = idrMoveR(F, oldRing); 3190 3191 to = clock(); 3192 /* reduce the Groebner basis <G> w.r.t. newRing */ 3174 3193 G = kInterRedCC(F1, NULL); 3175 3194 xtred=xtred+clock()-to; 3176 3195 idDelete(&F1); 3177 3196 3178 3197 if(endwalks == 1) 3179 3198 break; 3180 3199 3181 3200 NEXT_VECTOR: 3182 to = clock(); 3201 to = clock(); 3183 3202 /* compute a next weight vector */ 3184 3203 next_weight = MkInterRedNextWeight(curr_weight,target_weight, G); 3185 3204 xtnw=xtnw+clock()-to; 3186 #ifdef PRINT_VECTORS 3205 #ifdef PRINT_VECTORS 3187 3206 MivString(curr_weight, target_weight, next_weight); 3188 3207 #endif 3189 3208 3190 3209 if(Overflow_Error == TRUE) 3191 3210 { 3192 3211 /* 3193 3194 3212 ivString(next_weight, "omega"); 3213 PrintS("\n// ** The weight vector does NOT stay in Cone!!\n"); 3195 3214 */ 3196 3215 #ifdef TEST_OVERFLOW … … 3198 3217 #endif 3199 3218 3200 3219 3201 3220 newRing = currRing; 3202 3221 if (currRing->parameter != NULL) … … 3205 3224 VMrDefault(target_weight); 3206 3225 3207 F1 = idrMoveR(G, newRing); 3226 F1 = idrMoveR(G, newRing); 3208 3227 G = MstdCC(F1); 3209 3228 idDelete(&F1); … … 3218 3237 break; 3219 3238 } 3220 3239 3221 3240 if(MivComp(next_weight, target_weight) == 1) 3222 3241 { 3223 3242 if(MivSame(target_weight, exivlp)==1) 3224 3243 { 3225 LAST_GB_ALT2: 3226 3227 3244 LAST_GB_ALT2: 3245 nOverflow_Error = Overflow_Error; 3246 tproc = clock()-xftinput; 3228 3247 //Print("\n// takes %d steps and calls the recursion of level 2:", nwalk); 3229 3248 /* call the changed perturbation walk algorithm with degree 2 */ … … 3235 3254 endwalks = 1; 3236 3255 } 3237 3256 3238 3257 for(i=nV-1; i>=0; i--) 3239 3258 { … … 3251 3270 #ifdef TIME_TEST 3252 3271 Print("\n// \"Main procedure\" took %d steps dnd %.2f sec. Overflow_Error (%d)", nwalk, 3253 3272 ((double) tproc)/1000000, nOverflow_Error); 3254 3273 3255 3274 TimeStringFractal(xftinput, tostd, xtif, xtstd, xtextra,xtlift, xtred,xtnw); … … 3259 3278 Print("\n// Awalk2 took %d steps!!", nstep); 3260 3279 #endif 3261 3280 3262 3281 return(G); 3263 3282 } … … 3267 3286 /* The implementation of the fractal walk algorithmus */ 3268 3287 3269 /* The main procedur Mfwalk calls the recursive Subroutine 3288 /* The main procedur Mfwalk calls the recursive Subroutine 3270 3289 rec_fractal_call to compute the wanted Gröbner basis. 3271 3290 At the main procedur we compute the reduced Gröbner basis w.r.t. a "fast" 3272 order, e.g. "dp" and a sequence of weight vectors which are row vectors 3291 order, e.g. "dp" and a sequence of weight vectors which are row vectors 3273 3292 of a matrix. This matrix defines the given monomial order, e.g. "lp" 3274 3293 */ … … 3297 3316 3298 3317 /*********************************************************************** 3299 The procedur REC_GB_Mwalk computes a GB for <G> w.r.t. the weight order 3318 The procedur REC_GB_Mwalk computes a GB for <G> w.r.t. the weight order 3300 3319 otw, where G is a reduced GB w.r.t. the weight order cw. 3301 3320 The new procedur Mwalk calls REC_GB. 3302 3321 ************************************************************************/ 3303 static ideal REC_GB_Mwalk(ideal G, intvec* curr_weight, intvec* orig_target_weight, 3304 int tp_deg, int npwinc) 3322 static ideal REC_GB_Mwalk(ideal G, intvec* curr_weight, intvec* orig_target_weight, 3323 int tp_deg, int npwinc) 3305 3324 { 3306 3325 BOOLEAN nError = Overflow_Error; 3307 3326 Overflow_Error = FALSE; 3308 3327 3309 3328 int i, nV = currRing->N, ntwC, npwinC; 3310 3329 int nwalk=0, endwalks=0, nnwinC=1, nlast = 0; … … 3313 3332 intvec* iv_M_lp; 3314 3333 intvec* target_weight; 3315 intvec* ivNull = new intvec(nV); 3334 intvec* ivNull = new intvec(nV); 3316 3335 intvec* hilb_func; 3317 3336 BOOLEAN isGB = FALSE; … … 3336 3355 TargetRing = currRing; 3337 3356 ssG = idrMoveR(G,EXXRing); 3338 3357 3339 3358 ideal H0_tmp = idrMoveR(H0,EXXRing); 3340 3359 ideal H1 = idHeadCC(ssG); … … 3356 3375 else 3357 3376 goto LastGB_Finish; 3358 } 3377 } 3359 3378 idDelete(&H0_tmp); idDelete(&H1); 3360 3379 3361 3380 intvec* ivlp = Mivlp(nV); 3362 3381 if( MivSame(orig_target_weight, ivlp)==1 ) … … 3367 3386 //target_weight = MPertVectorslp(ssG, iv_M_lp, tp_deg); 3368 3387 target_weight = MPertVectors(ssG, iv_M_lp, tp_deg); 3369 3388 3370 3389 delete ivlp; 3371 3390 delete iv_M_lp; … … 3374 3393 G = idrMoveR(ssG, TargetRing); 3375 3394 } 3376 3395 3377 3396 while(1) 3378 { 3397 { 3379 3398 nwalk ++; 3380 3399 nstep++; … … 3394 3413 hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing); 3395 3414 #endif // BUCHBERGER_ALG 3396 3415 3397 3416 oldRing = currRing; 3398 3417 … … 3403 3422 VMrDefault(curr_weight); 3404 3423 3405 newRing = currRing; 3424 newRing = currRing; 3406 3425 Gomega1 = idrMoveR(Gomega, oldRing); 3407 3426 … … 3412 3431 #else 3413 3432 M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight); 3414 delete hilb_func; 3433 delete hilb_func; 3415 3434 #endif // BUCHBERGER_ALG 3416 3435 xtstd = xtstd + clock() - to; 3417 3436 3418 /* change the ring to oldRing */ 3437 /* change the ring to oldRing */ 3419 3438 rChangeCurrRing(oldRing); 3420 3439 3421 3440 M1 = idrMoveR(M, newRing); 3422 Gomega2 = idrMoveR(Gomega1, newRing); 3423 3441 Gomega2 = idrMoveR(Gomega1, newRing); 3442 3424 3443 to = clock(); 3425 3444 F = MLifttwoIdeal(Gomega2, M1, G); 3426 3445 xtlift = xtlift + clock() -to; 3427 3446 3428 idDelete(&M1); 3429 idDelete(&Gomega2); 3447 idDelete(&M1); 3448 idDelete(&Gomega2); 3430 3449 idDelete(&G); 3431 3450 3432 3451 3433 /* change the ring to newRing */ 3452 /* change the ring to newRing */ 3434 3453 rChangeCurrRing(newRing); 3435 3454 F1 = idrMoveR(F, oldRing); 3436 3455 3437 to = clock(); 3438 /* reduce the Groebner basis <G> w.r.t. new ring */ 3456 to = clock(); 3457 /* reduce the Groebner basis <G> w.r.t. new ring */ 3439 3458 G = kInterRedCC(F1, NULL); 3440 3459 xtred = xtred + clock() -to; … … 3451 3470 xtnw = xtnw + clock() - to; 3452 3471 3453 #ifdef PRINT_VECTORS 3472 #ifdef PRINT_VECTORS 3454 3473 MivString(curr_weight, target_weight, next_weight); 3455 3474 #endif 3456 3475 3457 /*check whether the computed vector does in the correct cone */ 3476 /*check whether the computed vector does in the correct cone */ 3458 3477 //ntwC = test_w_in_ConeCC(G, next_weight); 3459 3478 //if(ntwC != 1) … … 3473 3492 break; 3474 3493 } 3475 3494 3476 3495 if(MivComp(next_weight, target_weight) == 1) 3477 3496 { … … 3490 3509 (*curr_weight)[i] = (*next_weight)[i]; 3491 3510 3492 delete next_weight; 3493 }//while 3511 delete next_weight; 3512 }//while 3494 3513 3495 3514 delete ivNull; … … 3499 3518 //28.07.03 3500 3519 newRing = currRing; 3501 3520 3502 3521 if (currRing->parameter != NULL) 3503 3522 // DefRingParlp(); // … … 3508 3527 3509 3528 F1 = idrMoveR(G, newRing); 3510 3529 3511 3530 if(nnwinC == 0) 3512 3531 F1 = REC_GB_Mwalk(F1,curr_weight, orig_target_weight, tp_deg+1,nnwinC); … … 3514 3533 if(test_w_in_ConeCC(F1, target_weight) != 1) 3515 3534 F1 = REC_GB_Mwalk(F1,curr_weight, orig_target_weight,tp_deg+1,nnwinC); 3516 3535 3517 3536 delete target_weight; 3518 3537 … … 3524 3543 { 3525 3544 if(nlast == 1) 3526 { 3545 { 3527 3546 if (currRing->parameter != NULL) 3528 3547 DefRingPar(orig_target_weight); … … 3543 3562 newRing = currRing; 3544 3563 } 3545 3564 3546 3565 LastGB_Finish: 3547 3566 rChangeCurrRing(EXXRing); 3548 3567 result = idrMoveR(G, newRing); 3549 } 3550 3568 } 3569 3551 3570 if(Overflow_Error == FALSE) 3552 3571 Overflow_Error = nError; 3553 3554 return(result); 3572 3573 return(result); 3555 3574 } 3556 3575 … … 3558 3577 /* 08.09.02 */ 3559 3578 /******** THE NEW GRÖBNER WALK ALGORITHM **********/ 3560 /* Gröbnerwalk with a recursive "second" alternative GW, REC_GB_Mwalk 3579 /* Gröbnerwalk with a recursive "second" alternative GW, REC_GB_Mwalk 3561 3580 that only computes the last reduced GB */ 3562 3581 ideal Mwalk(ideal Go, intvec* curr_weight, intvec* target_weight) … … 3589 3608 (*last_omega)[i] = 1; 3590 3609 (*last_omega)[0] = 10000; 3591 3610 3592 3611 ring XXRing = currRing; 3593 3612 … … 3599 3618 if(currRing->order[0] == ringorder_a) 3600 3619 goto NEXT_VECTOR; 3601 3620 3602 3621 while(1) 3603 { 3622 { 3604 3623 nwalk ++; 3605 3624 nstep ++; … … 3616 3635 tim = clock(); 3617 3636 /* 3618 3619 3620 3637 Print("\n// **** Gröbnerwalk took %d steps and ", nwalk); 3638 PrintS("\n// **** call the rec. Pert. Walk to compute a red GB of:"); 3639 idElements(Gomega, "G_omega"); 3621 3640 */ 3622 3641 3623 3642 if(MivSame(exivlp, target_weight)==1) 3624 M = REC_GB_Mwalk(idCopy(Gomega), tmp_weight, curr_weight, 2,1); 3643 M = REC_GB_Mwalk(idCopy(Gomega), tmp_weight, curr_weight, 2,1); 3625 3644 else 3626 3645 goto NORMAL_GW; 3627 3646 /* 3628 Print("\n// time for the last std(Gw) = %.2f sec", 3629 3630 3647 Print("\n// time for the last std(Gw) = %.2f sec", 3648 ((double) (clock()-tim)/1000000)); 3649 PrintS("\n// ***************************************************\n"); 3631 3650 */ 3632 #ifdef CHECK_IDEAL_MWALK 3651 #ifdef CHECK_IDEAL_MWALK 3633 3652 idElements(Gomega, "G_omega"); 3634 3653 headidString(Gomega, "Gw"); 3635 3654 idElements(M, "M"); 3636 3655 //headidString(M, "M"); 3637 #endif 3656 #endif 3638 3657 to = clock(); 3639 3658 F = MLifttwoIdeal(Gomega, M, G); … … 3646 3665 oldRing = currRing; 3647 3666 3648 /* create a new ring newRing */ 3667 /* create a new ring newRing */ 3649 3668 if (currRing->parameter != NULL) 3650 3669 DefRingPar(curr_weight); … … 3652 3671 VMrDefault(curr_weight); 3653 3672 3654 newRing = currRing; 3655 F1 = idrMoveR(F, oldRing); 3673 newRing = currRing; 3674 F1 = idrMoveR(F, oldRing); 3656 3675 } 3657 3676 else … … 3660 3679 #ifndef BUCHBERGER_ALG 3661 3680 if(isNolVector(curr_weight) == 0) 3662 3681 hilb_func = hFirstSeries(Gomega,NULL,NULL,curr_weight,currRing); 3663 3682 else 3664 3683 hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing); 3665 3684 #endif // BUCHBERGER_ALG 3666 3685 3667 3686 /* define a new ring that its ordering is "(a(curr_weight),lp) */ 3668 3687 if (currRing->parameter != NULL) 3669 3688 DefRingPar(curr_weight); 3670 3689 else 3671 3672 3673 newRing = currRing; 3690 VMrDefault(curr_weight); 3691 3692 newRing = currRing; 3674 3693 Gomega1 = idrMoveR(Gomega, oldRing); 3675 3694 3676 3695 to = clock(); 3677 3696 /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */ … … 3680 3699 #else 3681 3700 M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight); 3682 delete hilb_func; 3701 delete hilb_func; 3683 3702 #endif // BUCHBERGER_ALG 3684 3703 tstd = tstd + clock() - to; 3685 3704 3686 /* change the ring to oldRing */ 3705 /* change the ring to oldRing */ 3687 3706 rChangeCurrRing(oldRing); 3688 3707 M1 = idrMoveR(M, newRing); … … 3691 3710 to = clock(); 3692 3711 /* compute a representation of the generators of submod (M) 3693 with respect to those of mod (Gomega). 3712 with respect to those of mod (Gomega). 3694 3713 Gomega is a reduced Groebner basis w.r.t. the current ring */ 3695 3714 F = MLifttwoIdeal(Gomega2, M1, G); 3696 3715 tlift = tlift + clock() - to; 3697 3716 3698 idDelete(&M1); 3699 idDelete(&Gomega2); 3700 idDelete(&G); 3701 3702 /* change the ring to newRing */ 3717 idDelete(&M1); 3718 idDelete(&Gomega2); 3719 idDelete(&G); 3720 3721 /* change the ring to newRing */ 3703 3722 rChangeCurrRing(newRing); 3704 F1 = idrMoveR(F, oldRing); 3723 F1 = idrMoveR(F, oldRing); 3705 3724 } 3706 3725 3707 3726 to = clock(); 3708 /* reduce the Groebner basis <G> w.r.t. new ring */ 3727 /* reduce the Groebner basis <G> w.r.t. new ring */ 3709 3728 G = kInterRedCC(F1, NULL); 3710 3729 if(endwalks != 1) … … 3716 3735 if(endwalks == 1) 3717 3736 break; 3718 3737 3719 3738 NEXT_VECTOR: 3720 3739 to = clock(); … … 3722 3741 intvec* next_weight = MkInterRedNextWeight(curr_weight,target_weight,G); 3723 3742 tnw = tnw + clock() - to; 3724 #ifdef PRINT_VECTORS 3743 #ifdef PRINT_VECTORS 3725 3744 MivString(curr_weight, target_weight, next_weight); 3726 3745 #endif … … 3736 3755 else 3737 3756 VMrDefault(target_weight); 3738 3739 F1 = idrMoveR(G, newRing); 3757 3758 F1 = idrMoveR(G, newRing); 3740 3759 G = MstdCC(F1); 3741 3760 idDelete(&F1); 3742 3761 3743 3762 newRing = currRing; 3744 3763 break; … … 3765 3784 delete tmp_weight; 3766 3785 delete ivNull; 3767 delete exivlp; 3768 3786 delete exivlp; 3787 3769 3788 #ifdef TIME_TEST 3770 3789 TimeString(tinput, tostd, tif, tstd, tlift, tred, tnw, nstep); … … 3783 3802 int i, nV = currRing->N; 3784 3803 int nwalk=0, endwalks=0; 3785 3804 3786 3805 ideal Gomega, M, F, Gomega1, Gomega2, M1, F1, G, G1; 3787 3806 ring endRing, newRing, oldRing; … … 3805 3824 3806 3825 while(1) 3807 { 3826 { 3808 3827 nwalk ++; 3809 3828 //Print("\n// Entering the %d-th step:", nwalk); 3810 3829 //Print("\n// ring r[%d] = %s;", nwalk, rString(currRing)); 3811 idString(G,"G"); 3830 idString(G,"G"); 3812 3831 /* compute an initial form ideal of <G> w.r.t. "curr_vector" */ 3813 3832 Gomega = MwalkInitialForm(G, curr_weight); 3814 3833 //ivString(curr_weight, "omega"); 3815 idString(Gomega,"Gw"); 3834 idString(Gomega,"Gw"); 3816 3835 3817 3836 #ifndef BUCHBERGER_ALG … … 3824 3843 3825 3844 oldRing = currRing; 3826 3845 3827 3846 /* define a new ring that its ordering is "(a(curr_weight),lp) */ 3828 3847 VMrDefault(curr_weight); 3829 newRing = currRing; 3830 3848 newRing = currRing; 3849 3831 3850 Gomega1 = idrMoveR(Gomega, oldRing); 3832 3851 3833 3852 /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */ 3834 3853 #ifdef BUCHBERGER_ALG … … 3836 3855 #else 3837 3856 M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight); 3838 delete hilb_func; 3857 delete hilb_func; 3839 3858 #endif // BUCHBERGER_ALG 3840 3859 3841 3860 idString(M,"M"); 3842 3861 3843 /* change the ring to oldRing */ 3862 /* change the ring to oldRing */ 3844 3863 rChangeCurrRing(oldRing); 3845 3864 M1 = idrMoveR(M, newRing); … … 3847 3866 3848 3867 /* compute a representation of the generators of submod (M) 3849 with respect to those of mod (Gomega). 3850 3868 with respect to those of mod (Gomega). 3869 Gomega is a reduced Groebner basis w.r.t. the current ring */ 3851 3870 F = MLifttwoIdeal(Gomega2, M1, G); 3852 idDelete(&M1); 3853 idDelete(&Gomega2); 3854 idDelete(&G); 3871 idDelete(&M1); 3872 idDelete(&Gomega2); 3873 idDelete(&G); 3855 3874 idString(F,"F"); 3856 3857 /* change the ring to newRing */ 3875 3876 /* change the ring to newRing */ 3858 3877 rChangeCurrRing(newRing); 3859 F1 = idrMoveR(F, oldRing); 3860 3861 /* reduce the Groebner basis <G> w.r.t. new ring */ 3878 F1 = idrMoveR(F, oldRing); 3879 3880 /* reduce the Groebner basis <G> w.r.t. new ring */ 3862 3881 G = kInterRedCC(F1, NULL); 3863 3882 //idSkipZeroes(G);//done by kInterRed … … 3866 3885 if(endwalks == 1) 3867 3886 break; 3868 3887 3869 3888 /* compute a next weight vector */ 3870 3889 intvec* next_weight = MkInterRedNextWeight(curr_weight,target_weight,G); 3871 #ifdef PRINT_VECTORS 3890 #ifdef PRINT_VECTORS 3872 3891 MivString(curr_weight, target_weight, next_weight); 3873 3892 #endif … … 3904 3923 /**************************************************************/ 3905 3924 /* If the perturbed target weight vector or an intermediate weight vector 3906 doesn't stay in the correct Groebner cone, we have only 3907 a reduced Groebner basis for the given ideal with respect to 3925 doesn't stay in the correct Groebner cone, we have only 3926 a reduced Groebner basis for the given ideal with respect to 3908 3927 a monomial order which differs to the given order. 3909 3928 Then we have to compute the wanted reduced Groebner basis for it. … … 3914 3933 /* use kStd, if nP = 0, else call LastGB */ 3915 3934 ideal Mpwalk(ideal Go, int op_deg, int tp_deg,intvec* curr_weight, 3916 3917 { 3935 intvec* target_weight, int nP) 3936 { 3918 3937 Set_Error(FALSE ); 3919 3938 Overflow_Error = FALSE; … … 3930 3949 int i, ntwC=1, ntestw=1, nV = currRing->N, op_tmp = op_deg; 3931 3950 int endwalks=0, nhilb=0, ntestomega=0; 3932 3951 3933 3952 ideal Gomega, M, F, G, Gomega1, Gomega2, M1,F1,Eresult,ssG; 3934 3953 ring newRing, oldRing, TargetRing; … … 3938 3957 intvec* orig_target = target_weight; 3939 3958 intvec* pert_target_vector = target_weight; 3940 intvec* ivNull = new intvec(nV); 3959 intvec* ivNull = new intvec(nV); 3941 3960 intvec* iv_dp = MivUnit(nV);// define (1,1,...,1) 3942 3961 intvec* cw_tmp = curr_weight; … … 3944 3963 intvec* next_weight; 3945 3964 intvec* extra_curr_weight = new intvec(nV); 3946 3965 3947 3966 /* to avoid (1,0,...,0) as the target vector */ 3948 3967 intvec* last_omega = new intvec(nV); … … 3950 3969 (*last_omega)[i] = 1; 3951 3970 (*last_omega)[0] = 10000; 3952 3971 3953 3972 ring XXRing = currRing; 3954 3973 … … 3965 3984 curr_weight = MPertVectors(G, iv_M_dp, op_deg); 3966 3985 } 3967 } 3986 } 3968 3987 else 3969 3988 { 3970 //ring order := (a(curr_weight),lp); 3989 //ring order := (a(curr_weight),lp); 3971 3990 if (currRing->parameter != NULL) 3972 3991 DefRingPar(curr_weight); 3973 3992 else 3974 3993 VMrDefault(curr_weight); 3975 3994 3976 3995 G = idrMoveR(Go, XXRing); 3977 3996 G = MstdCC(G); … … 3984 4003 delete iv_dp; 3985 4004 if(op_deg != 1) delete iv_M_dp; 3986 4005 3987 4006 ring HelpRing = currRing; 3988 4007 … … 4008 4027 iv_M_lp = MivMatrixOrder(target_weight); 4009 4028 //target_weight = MPertVectorslp(ssG, iv_M_lp, tp_deg); 4010 target_weight = MPertVectors(ssG, iv_M_lp, tp_deg); 4011 } 4029 target_weight = MPertVectors(ssG, iv_M_lp, tp_deg); 4030 } 4012 4031 delete iv_M_lp; 4013 4032 pert_target_vector = target_weight; //vor 19. mai 2003//test 19 Junu 03 … … 4021 4040 */ 4022 4041 while(1) 4023 { 4042 { 4024 4043 nstep ++; 4025 4044 to = clock(); 4026 /* compute an initial form ideal of <G> w.r.t. the weight vector 4045 /* compute an initial form ideal of <G> w.r.t. the weight vector 4027 4046 "curr_weight" */ 4028 4047 Gomega = MwalkInitialForm(G, curr_weight); 4029 4048 4030 4049 4031 4050 #ifdef ENDWALKS 4032 if(endwalks == 1){ 4033 Print("\n// ring r%d = %s;\n", nstep, rString(currRing)); 4051 if(endwalks == 1){ 4052 Print("\n// ring r%d = %s;\n", nstep, rString(currRing)); 4034 4053 idElements(G, "G"); 4035 4054 // idElements(Gomega, "Gw"); … … 4056 4075 VMrDefault(curr_weight); 4057 4076 4058 newRing = currRing; 4077 newRing = currRing; 4059 4078 Gomega1 = idrMoveR(Gomega, oldRing); 4060 4079 4061 #ifdef ENDWALKS 4080 #ifdef ENDWALKS 4062 4081 if(endwalks==1) 4063 4082 { 4064 4083 Print("\n// ring r%d = %s;\n", nstep, rString(currRing)); 4065 4084 idElements(Gomega1, "Gw"); 4066 headidString(Gomega1, "headGw"); 4085 headidString(Gomega1, "headGw"); 4067 4086 PrintS("\n// compute a rGB of Gw:\n"); 4068 4087 4069 #ifndef BUCHBERGER_ALG 4088 #ifndef BUCHBERGER_ALG 4070 4089 ivString(hilb_func, "w"); 4071 #endif4072 }4073 4090 #endif 4074 4091 } 4092 #endif 4093 4075 4094 tim = clock(); 4076 4095 to = clock(); … … 4080 4099 #else 4081 4100 M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight); 4082 delete hilb_func; 4101 delete hilb_func; 4083 4102 #endif // BUCHBERGER_ALG 4084 4103 4085 4104 if(endwalks == 1){ 4086 4105 xtstd = xtstd+clock()-to; 4087 #ifdef ENDWALKS 4088 Print("\n// time for the last std(Gw) = %.2f sec\n", 4106 #ifdef ENDWALKS 4107 Print("\n// time for the last std(Gw) = %.2f sec\n", 4089 4108 ((double) clock())/1000000 -((double)tim) /1000000); 4090 4109 #endif … … 4093 4112 tstd=tstd+clock()-to; 4094 4113 4095 /* change the ring to oldRing */ 4114 /* change the ring to oldRing */ 4096 4115 rChangeCurrRing(oldRing); 4097 4116 M1 = idrMoveR(M, newRing); 4098 Gomega2 = idrMoveR(Gomega1, newRing); 4099 4117 Gomega2 = idrMoveR(Gomega1, newRing); 4118 4100 4119 //if(endwalks==1) PrintS("\n// Lifting is working:.."); 4101 4120 4102 4121 to=clock(); 4103 4122 /* compute a representation of the generators of submod (M) 4104 with respect to those of mod (Gomega). 4123 with respect to those of mod (Gomega). 4105 4124 Gomega is a reduced Groebner basis w.r.t. the current ring */ 4106 4125 F = MLifttwoIdeal(Gomega2, M1, G); … … 4110 4129 xtlift=clock()-to; 4111 4130 4112 idDelete(&M1); 4113 idDelete(&Gomega2); 4131 idDelete(&M1); 4132 idDelete(&Gomega2); 4114 4133 idDelete(&G); 4115 4134 4116 /* change the ring to newRing */ 4135 /* change the ring to newRing */ 4117 4136 rChangeCurrRing(newRing); 4118 4137 F1 = idrMoveR(F, oldRing); … … 4121 4140 4122 4141 to=clock(); 4123 /* reduce the Groebner basis <G> w.r.t. new ring */ 4142 /* reduce the Groebner basis <G> w.r.t. new ring */ 4124 4143 G = kInterRedCC(F1, NULL); 4125 4144 if(endwalks != 1) … … 4129 4148 4130 4149 idDelete(&F1); 4131 4150 4132 4151 if(endwalks == 1) 4133 4152 break; 4134 4153 4135 4154 to=clock(); 4136 4155 /* compute a next weight vector */ 4137 4156 next_weight = MkInterRedNextWeight(curr_weight,target_weight, G); 4138 4157 tnw=tnw+clock()-to; 4139 #ifdef PRINT_VECTORS 4158 #ifdef PRINT_VECTORS 4140 4159 MivString(curr_weight, target_weight, next_weight); 4141 4160 #endif … … 4145 4164 ntwC = 0; 4146 4165 ntestomega = 1; 4147 //Print("\n// ring r%d = %s;\n", nstep, rString(currRing)); 4166 //Print("\n// ring r%d = %s;\n", nstep, rString(currRing)); 4148 4167 //idElements(G, "G"); 4149 4168 delete next_weight; 4150 goto FINISH_160302; 4169 goto FINISH_160302; 4151 4170 } 4152 4171 if(MivComp(next_weight, ivNull) == 1){ … … 4159 4178 endwalks = 1; 4160 4179 4161 for(i=nV-1; i>=0; i--) 4180 for(i=nV-1; i>=0; i--) 4162 4181 (*curr_weight)[i] = (*next_weight)[i]; 4163 4182 4164 4183 delete next_weight; 4165 4184 }//while … … 4168 4187 { 4169 4188 FINISH_160302://16.03.02 4170 if(MivSame(orig_target, exivlp) == 1) 4189 if(MivSame(orig_target, exivlp) == 1) 4171 4190 if (currRing->parameter != NULL) 4172 4191 DefRingParlp(); … … 4178 4197 else 4179 4198 VMrDefault(orig_target); 4180 4199 4181 4200 TargetRing=currRing; 4182 4201 F1 = idrMoveR(G, newRing); … … 4184 4203 headidString(G, "G"); 4185 4204 #endif 4186 4205 4187 4206 4188 4207 // check whether the pertubed target vector stays in the correct cone … … 4194 4213 { 4195 4214 /* 4196 if(ntestw != 1){ 4197 ivString(pert_target_vector, "tau"); 4198 PrintS("\n// ** perturbed target vector doesn't stay in cone!!"); 4199 Print("\n// ring r%d = %s;\n", nstep, rString(currRing)); 4200 idElements(F1, "G"); 4215 if(ntestw != 1){ 4216 ivString(pert_target_vector, "tau"); 4217 PrintS("\n// ** perturbed target vector doesn't stay in cone!!"); 4218 Print("\n// ring r%d = %s;\n", nstep, rString(currRing)); 4219 idElements(F1, "G"); 4201 4220 } 4202 4221 */ … … 4205 4224 ideal eF1; 4206 4225 if(nP == 0 || tp_deg == 1 || MivSame(orig_target, exivlp) != 1){ 4207 // PrintS("\n// ** calls \"std\" to compute a GB"); 4226 // PrintS("\n// ** calls \"std\" to compute a GB"); 4208 4227 eF1 = MstdCC(F1); 4209 4228 idDelete(&F1); 4210 4229 } 4211 4230 else { 4212 // PrintS("\n// ** calls \"LastGB\" to compute a GB"); 4213 4214 4215 4216 4231 // PrintS("\n// ** calls \"LastGB\" to compute a GB"); 4232 rChangeCurrRing(newRing); 4233 ideal F2 = idrMoveR(F1, TargetRing); 4234 eF1 = LastGB(F2, curr_weight, tp_deg-1); 4235 F2=NULL; 4217 4236 } 4218 4237 xtextra=clock()-to; … … 4221 4240 rChangeCurrRing(XXRing); 4222 4241 Eresult = idrMoveR(eF1, exTargetRing); 4223 } 4242 } 4224 4243 else{ 4225 4244 rChangeCurrRing(XXRing); … … 4230 4249 rChangeCurrRing(XXRing); 4231 4250 Eresult = idrMoveR(G, newRing); 4232 } 4251 } 4233 4252 delete ivNull; 4234 4253 if(tp_deg != 1) … … 4242 4261 4243 4262 #ifdef TIME_TEST 4244 TimeStringFractal(tinput, tostd, tif+xtif, tstd+xtstd,0, tlift+xtlift, tred+xtred, 4245 4263 TimeStringFractal(tinput, tostd, tif+xtif, tstd+xtstd,0, tlift+xtlift, tred+xtred, 4264 tnw+xtnw); 4246 4265 4247 4266 Print("\n// pSetm_Error = (%d)", ErrorCheck()); 4248 4267 Print("\n// It took %d steps and Overflow_Error? (%d)\n", nstep, Overflow_Error); 4249 4268 #endif 4250 return(Eresult); 4269 return(Eresult); 4251 4270 } 4252 4271 … … 4258 4277 int i,j; 4259 4278 intvec* ivM = new intvec(nV*nV); 4260 4279 4261 4280 for(i=0; i<nV; i++) 4262 4281 for(j=0; j<nV; j++) 4263 4282 (*ivM)[i*nV + j] = 1; 4264 4283 4265 4284 return(ivM); 4266 4285 } … … 4303 4322 (*omega)[i] = (*Xsigma)[(nV*(nlev-1))+i]; 4304 4323 4305 (*omega2)[i] = (*Xtau)[(nlev-1)*nV+i]; 4306 } 4307 4324 (*omega2)[i] = (*Xtau)[(nlev-1)*nV+i]; 4325 } 4326 4308 4327 if(nlev == 1) Xcall = 1; 4309 4328 else Xcall = 0; 4310 4329 4311 4330 ring oRing = currRing; 4312 4331 4313 4332 while(1) 4314 4333 { 4315 4334 #ifdef FIRST_STEP_FRACTAL 4316 // perturb the current weight vector only on the top level or 4335 // perturb the current weight vector only on the top level or 4317 4336 // after perturbation of the both vectors, nlev = 2 as the top level 4318 4337 if((nlev == 1 && Xcall == 0) || (nlev == 2 && Xngleich == 1)) 4319 4338 if(islengthpoly2(G) == 1) 4320 4339 { 4321 4322 4323 4324 4340 Mwlp = MivWeightOrderlp(omega); 4341 Xsigma = Mfpertvector(G, Mwlp); 4342 delete Mwlp; 4343 Overflow_Error = FALSE; 4325 4344 } 4326 4345 #endif … … 4331 4350 next_vect = MkInterRedNextWeight(omega,omega2,G); 4332 4351 xtnw=xtnw+clock()-to; 4333 #ifdef PRINT_VECTORS 4352 #ifdef PRINT_VECTORS 4334 4353 MivString(omega, omega2, next_vect); 4335 4354 #endif 4336 4355 oRing = currRing; 4337 4356 4338 4357 /* We only perturb the current target vector at the recursion level 1 */ 4339 if(Xngleich == 0 && nlev == 1) //(ngleich == 0) important, e.g. ex2, ex3 4340 if (MivComp(next_vect, omega2) == 1) 4358 if(Xngleich == 0 && nlev == 1) //(ngleich == 0) important, e.g. ex2, ex3 4359 if (MivComp(next_vect, omega2) == 1) 4341 4360 { 4342 4361 /* to dispense with taking initial (and lifting/interreducing 4343 4362 after the call of recursion */ 4344 4363 //Print("\n\n// ** Perturb the both vectors with degree %d with",nlev); 4345 4364 //idElements(G, "G"); 4346 4365 4347 4366 Xngleich = 1; … … 4350 4369 if (currRing->parameter != NULL) 4351 4370 DefRingPar(omtmp); 4352 else 4371 else 4353 4372 VMrDefault(omtmp); 4354 4373 4355 4374 testring = currRing; 4356 4375 Gt = idrMoveR(G, oRing); 4357 4376 4358 4377 /* perturb the original target vector w.r.t. the current GB */ 4359 delete Xtau; 4378 delete Xtau; 4360 4379 Xtau = NewVectorlp(Gt); 4361 4362 rChangeCurrRing(oRing); 4380 4381 rChangeCurrRing(oRing); 4363 4382 G = idrMoveR(Gt, testring); 4364 4383 … … 4367 4386 Xsigma = Mfpertvector(G, Mwlp); 4368 4387 delete Mwlp; 4369 4388 4370 4389 for(i=nV-1; i>=0; i--) { 4371 4390 (*omega2)[i] = (*Xtau)[nV+i]; 4372 4391 (*omega)[i] = (*Xsigma)[nV+i]; 4373 4392 } 4374 4393 4375 4394 delete next_vect; 4376 4377 4378 4379 4380 4381 next_vect = MkInterRedNextWeight(omega,omega2,G); 4382 4383 4384 #ifdef PRINT_VECTORS 4395 to=clock(); 4396 4397 /* to avoid the value of Overflow_Error that occur in Mfpertvector*/ 4398 Overflow_Error = FALSE; 4399 4400 next_vect = MkInterRedNextWeight(omega,omega2,G); 4401 xtnw=xtnw+clock()-to; 4402 4403 #ifdef PRINT_VECTORS 4385 4404 MivString(omega, omega2, next_vect); 4386 4405 #endif 4387 4406 } 4388 4389 4407 4408 4390 4409 /* check whether the the computed vector is in the correct cone */ 4391 /* If no, the reduced GB of an omega-homogeneous ideal will be 4410 /* If no, the reduced GB of an omega-homogeneous ideal will be 4392 4411 computed by Buchberger algorithm and stop this recursion step*/ 4393 4412 //if(test_w_in_ConeCC(G, next_vect) != 1) //e.g. Example s7, cyc6 … … 4414 4433 Gt = NULL; 4415 4434 4416 delete omega2; 4417 delete altomega; 4435 delete omega2; 4436 delete altomega; 4418 4437 4419 4438 //Print("\n// Leaving the %d-th recursion with %d steps", nlev, nwalks); 4420 4439 //Print(" ** Overflow_Error? (%d)", Overflow_Error); 4421 4440 nnflow ++; 4422 4441 4423 4442 Overflow_Error = FALSE; 4424 4443 return (G1); 4425 4444 } 4426 4445 4427 4428 /* If the perturbed target vector stays in the correct cone, 4446 4447 /* If the perturbed target vector stays in the correct cone, 4429 4448 return the current GB, 4430 otherwise, return the computed GB by the Buchberger-algorithm. 4431 Then we update the perturbed target vectors w.r.t. this GB. */ 4432 4449 otherwise, return the computed GB by the Buchberger-algorithm. 4450 Then we update the perturbed target vectors w.r.t. this GB. */ 4451 4433 4452 /* the computed vector is equal to the origin vector, since 4434 4453 t is not defined */ 4435 if (MivComp(next_vect, XivNull) == 1) 4454 if (MivComp(next_vect, XivNull) == 1) 4436 4455 { 4437 4456 if (currRing->parameter != NULL) … … 4444 4463 4445 4464 if(test_w_in_ConeCC(Gt, omega2) == 1) { 4446 delete omega2; 4465 delete omega2; 4447 4466 delete next_vect; 4448 4467 delete altomega; 4449 4468 //Print("\n// Leaving the %d-th recursion with %d steps ",nlev, nwalks); 4450 4469 //Print(" ** Overflow_Error? (%d)", Overflow_Error); 4451 4470 4452 4471 return (Gt); 4453 4472 } … … 4456 4475 //ivString(omega2, "tau'"); 4457 4476 //Print("\n// tau' doesn't stay in the correct cone!!"); 4458 4477 4459 4478 #ifndef MSTDCC_FRACTAL 4460 4479 //07.08.03 … … 4475 4494 Xtau = Xtautmp; 4476 4495 Xtautmp = NULL; 4477 //ivString(Xtau, "new Xtau"); 4478 4479 for(i=nV-1; i>=0; i--) 4496 //ivString(Xtau, "new Xtau"); 4497 4498 for(i=nV-1; i>=0; i--) 4480 4499 (*omega2)[i] = (*Xtau)[(nlev-1)*nV+i]; 4481 4500 4482 4501 //Print("\n// ring tau = %s;", rString(currRing)); 4483 4502 rChangeCurrRing(oRing); … … 4491 4510 to=clock(); 4492 4511 G = MstdCC(Gt); 4493 4494 4495 oRing = currRing; 4496 4497 // update the original target vector w.r.t. the current GB 4512 xtextra=xtextra+clock()-to; 4513 4514 oRing = currRing; 4515 4516 // update the original target vector w.r.t. the current GB 4498 4517 if(MivSame(Xivinput, Xivlp) == 1) 4499 4518 if (currRing->parameter != NULL) … … 4515 4534 rChangeCurrRing(oRing); 4516 4535 G = idrMoveR(Gt, testring); 4517 4518 delete omega2; 4536 4537 delete omega2; 4519 4538 delete next_vect; 4520 delete altomega; 4521 /* 4522 4523 4524 4539 delete altomega; 4540 /* 4541 Print("\n// Leaving the %d-th recursion with %d steps,", nlev,nwalks); 4542 Print(" ** Overflow_Error? (%d)", Overflow_Error); 4543 */ 4525 4544 if(Overflow_Error == TRUE) 4526 4545 nnflow ++; … … 4529 4548 return(G); 4530 4549 } 4531 } 4532 4533 for(i=nV-1; i>=0; i--) { 4550 } 4551 4552 for(i=nV-1; i>=0; i--) { 4534 4553 (*altomega)[i] = (*omega)[i]; 4535 4554 (*omega)[i] = (*next_vect)[i]; … … 4541 4560 Gomega = MwalkInitialForm(G, omega); 4542 4561 xtif=xtif+clock()-to; 4543 4562 4544 4563 #ifndef BUCHBERGER_ALG 4545 4564 if(isNolVector(omega) == 0) … … 4548 4567 hilb_func = hFirstSeries(Gomega,NULL,NULL,last_omega,currRing); 4549 4568 #endif // BUCHBERGER_ALG 4550 4569 4551 4570 if (currRing->parameter != NULL) 4552 4571 DefRingPar(omega); … … 4555 4574 4556 4575 Gomega1 = idrMoveR(Gomega, oRing); 4557 4558 /* Maximal recursion depth, to compute a red. GB */ 4576 4577 /* Maximal recursion depth, to compute a red. GB */ 4559 4578 /* Fractal walk with the alternative recursion */ 4560 4579 /* alternative recursion */ … … 4564 4583 { 4565 4584 /* 4566 if(Xnlev != nV) 4567 { 4568 Print("\n// ** Xnlev = %d", Xnlev); 4569 ivString(Xtau, "Xtau"); 4570 } 4585 if(Xnlev != nV) 4586 { 4587 Print("\n// ** Xnlev = %d", Xnlev); 4588 ivString(Xtau, "Xtau"); 4589 } 4571 4590 */ 4572 4591 to=clock(); … … 4575 4594 #else 4576 4595 Gresult =kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,omega); 4577 delete hilb_func; 4596 delete hilb_func; 4578 4597 #endif // BUCHBERGER_ALG 4579 4598 xtstd=xtstd+clock()-to; 4580 } 4599 } 4581 4600 else { 4582 4601 rChangeCurrRing(oRing); … … 4585 4604 } 4586 4605 4587 //convert a Groebner basis from a ring to another ring, 4606 //convert a Groebner basis from a ring to another ring, 4588 4607 new_ring = currRing; 4589 4590 rChangeCurrRing(oRing); 4608 4609 rChangeCurrRing(oRing); 4591 4610 Gresult1 = idrMoveR(Gresult, new_ring); 4592 4611 Gomega2 = idrMoveR(Gomega1, new_ring); 4593 4612 4594 4613 to=clock(); 4595 4614 /* Lifting process */ 4596 4615 F = MLifttwoIdeal(Gomega2, Gresult1, G); 4597 4616 xtlift=xtlift+clock()-to; 4598 idDelete(&Gresult1); 4599 idDelete(&Gomega2); 4617 idDelete(&Gresult1); 4618 idDelete(&Gomega2); 4600 4619 idDelete(&G); 4601 4620 4602 4621 rChangeCurrRing(new_ring); 4603 4622 F1 = idrMoveR(F, oRing); … … 4608 4627 xtred=xtred+clock()-to; 4609 4628 idDelete(&F1); 4610 } 4629 } 4611 4630 } 4612 4631 … … 4617 4636 //Print("// pSetm_Error = (%d)", ErrorCheck()); 4618 4637 //Print("\n// ring ro = %s;", rString(currRing)); 4619 4638 4620 4639 nnflow = 0; 4621 4640 Xngleich = 0; … … 4623 4642 xtif=0; xtstd=0; xtlift=0; xtred=0; xtnw=0; xtextra=0; 4624 4643 xftinput = clock(); 4625 4644 4626 4645 ring oldRing = currRing; 4627 4646 int i, nV = currRing->N; … … 4644 4663 && (Gw->m[i]->next!=NULL) /* len >=1 */ 4645 4664 && (Gw->m[i]->next->next!=NULL)) /* len >=2 */ 4646 { 4665 { 4647 4666 intvec* iv_dp = MivUnit(nV);// define (1,1,...,1) 4648 4667 intvec* Mdp; 4649 4668 4650 4669 if(MivSame(ivstart, iv_dp) != 1) 4651 4652 else 4653 4670 Mdp = MivWeightOrderdp(ivstart); 4671 else 4672 Mdp = MivMatrixOrderdp(nV); 4654 4673 4655 4674 Xsigma = Mfpertvector(I, Mdp); … … 4668 4687 Xivlp = Mivlp(nV); 4669 4688 4670 if(MivComp(ivtarget, Xivlp) != 1) 4689 if(MivComp(ivtarget, Xivlp) != 1) 4671 4690 { 4672 4691 if (currRing->parameter != NULL) … … 4698 4717 id_Delete(&I, oldRing); 4699 4718 ring tRing = currRing; 4700 4719 4701 4720 if (currRing->parameter != NULL) 4702 4721 DefRingPar(ivstart); … … 4712 4731 ideal resF; 4713 4732 ring helpRing = currRing; 4714 4715 J = rec_fractal_call(J, 1, ivtarget); 4716 4717 rChangeCurrRing(oldRing); 4733 4734 J = rec_fractal_call(J, 1, ivtarget); 4735 4736 rChangeCurrRing(oldRing); 4718 4737 resF = idrMoveR(J, helpRing); 4719 4738 idSkipZeroes(resF); … … 4725 4744 4726 4745 #ifdef TIME_TEST 4727 TimeStringFractal(xftinput, xftostd, xtif, xtstd, xtextra, 4728 xtlift, xtred, xtnw); 4746 TimeStringFractal(xftinput, xftostd, xtif, xtstd, xtextra, 4747 xtlift, xtred, xtnw); 4729 4748 4730 4749 … … 4732 4751 Print("\n// Overflow_Error? (%d)\n", Overflow_Error); 4733 4752 Print("\n// the numbers of Overflow_Error (%d)", nnflow); 4734 #endif 4735 4753 #endif 4754 4736 4755 return(resF); 4737 4756 } … … 4751 4770 4752 4771 int nsteppert=0, i, nV = currRing->N, nwalk=0, npert_tmp=0; 4753 int npert[2*nV];4772 int *npert=(int*)omAlloc(2*nV*sizeof(int)); 4754 4773 ideal Gomega, M,F, G1, Gomega1, Gomega2, M1, F1; 4755 4774 ring endRing, newRing, oldRing, lpRing; … … 4774 4793 for(i=nV-1; i>=0; i--) 4775 4794 (*target_weight)[i] = (*target_tmp)[i]; 4776 4795 4777 4796 ring XXRing = currRing; 4778 4797 newRing = currRing; … … 4796 4815 #ifdef REPRESENTATION_OF_SIGMA 4797 4816 ideal Gw = MwalkInitialForm(G, curr_weight); 4798 4817 4799 4818 if(islengthpoly2(Gw)==1) 4800 4819 { 4801 4820 intvec* MDp; 4802 4821 if(MivComp(curr_weight, iv_dp) == 1) 4803 MDp = MatrixOrderdp(nV); //MivWeightOrderlp(iv_dp); 4822 MDp = MatrixOrderdp(nV); //MivWeightOrderlp(iv_dp); 4804 4823 else 4805 4824 MDp = MivWeightOrderlp(curr_weight); 4806 4825 4807 4826 curr_weight = RepresentationMatrix_Dp(G, MDp); 4808 4827 … … 4810 4829 4811 4830 ring exring = currRing; 4812 4831 4813 4832 if (currRing->parameter != NULL) 4814 4833 DefRingPar(curr_weight); 4815 else 4834 else 4816 4835 VMrDefault(curr_weight); 4817 4836 to=clock(); 4818 4837 Gw = idrMoveR(G, exring); 4819 G = MstdCC(Gw); 4838 G = MstdCC(Gw); 4820 4839 Gw = NULL; 4821 4840 tostd=tostd+clock()-to; … … 4823 4842 goto COMPUTE_NEW_VECTOR; 4824 4843 } 4825 4844 4826 4845 idDelete(&Gw); 4827 4846 delete iv_dp; … … 4830 4849 4831 4850 while(1) 4832 { 4851 { 4833 4852 to=clock(); 4834 4853 /* compute an initial form ideal of <G> w.r.t. "curr_vector" */ … … 4851 4870 VMrDefault(curr_weight); 4852 4871 4853 newRing = currRing; 4854 Gomega1 = idrMoveR(Gomega, oldRing); 4855 4872 newRing = currRing; 4873 Gomega1 = idrMoveR(Gomega, oldRing); 4874 4856 4875 to=clock(); 4857 4876 /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */ … … 4860 4879 #else 4861 4880 M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight); 4862 delete hilb_func; 4881 delete hilb_func; 4863 4882 #endif // BUCHBERGER_ALG 4864 4883 tstd=tstd+clock()-to; 4865 4884 4866 /* change the ring to oldRing */ 4885 /* change the ring to oldRing */ 4867 4886 rChangeCurrRing(oldRing); 4868 4887 M1 = idrMoveR(M, newRing); 4869 Gomega2 = idrMoveR(Gomega1, newRing); 4888 Gomega2 = idrMoveR(Gomega1, newRing); 4870 4889 4871 4890 to=clock(); 4872 4891 /* compute a representation of the generators of submod (M) 4873 with respect to those of mod (Gomega). 4892 with respect to those of mod (Gomega). 4874 4893 Gomega is a reduced Groebner basis w.r.t. the current ring */ 4875 4894 F = MLifttwoIdeal(Gomega2, M1, G); 4876 4895 tlift=tlift+clock()-to; 4877 4896 4878 idDelete(&M1); 4879 idDelete(&Gomega2); 4897 idDelete(&M1); 4898 idDelete(&Gomega2); 4880 4899 idDelete(&G); 4881 4882 /* change the ring to newRing */ 4900 4901 /* change the ring to newRing */ 4883 4902 rChangeCurrRing(newRing); 4884 4903 F1 = idrMoveR(F, oldRing); 4885 4904 4886 4905 to=clock(); 4887 /* reduce the Groebner basis <G> w.r.t. new ring */ 4906 /* reduce the Groebner basis <G> w.r.t. new ring */ 4888 4907 G = kInterRedCC(F1, NULL); 4889 4908 tred=tred+clock()-to; 4890 4909 idDelete(&F1); 4891 4910 4892 4911 4893 4912 COMPUTE_NEW_VECTOR: 4894 4913 newRing = currRing; … … 4899 4918 next_weight = MwalkNextWeightCC(curr_weight,target_weight, G); 4900 4919 tnw=tnw+clock()-to; 4901 #ifdef PRINT_VECTORS 4920 #ifdef PRINT_VECTORS 4902 4921 MivString(curr_weight, target_weight, next_weight); 4903 4922 #endif 4904 4923 4905 4906 /* check whether the computed intermediate weight vector in 4924 4925 /* check whether the computed intermediate weight vector in 4907 4926 the correct cone is, since sometimes it is very big e.g. s7, cyc7. 4908 If it is NOT in the cone, then one has directly to compute 4927 If it is NOT in the cone, then one has directly to compute 4909 4928 a reduced Groebner basis with respect to the lexicographic ordering 4910 4929 for the known Groebner basis that it is computed in the last step. … … 4942 4961 /*apply kStd or LastGB to compute a lex. red. Groebner basis of <G>*/ 4943 4962 if(nP == 0 || MivSame(target_tmp, iv_lp) == 0){ 4944 //Print("\n\n// calls \"std in ring r_%d = %s;", nwalk, rString(currRing)); 4963 //Print("\n\n// calls \"std in ring r_%d = %s;", nwalk, rString(currRing)); 4945 4964 G = MstdCC(G1);//no result for qnt1 4946 4965 } 4947 4966 else { 4948 4949 4950 4951 4952 4953 4954 4955 4967 rChangeCurrRing(newRing); 4968 G1 = idrMoveR(G1, lpRing); 4969 4970 //Print("\n\n// calls \"LastGB\" (%d) to compute a GB", nV-1); 4971 G = LastGB(G1, curr_weight, nV-1); //no result for kats7 4972 4973 rChangeCurrRing(lpRing); 4974 G = idrMoveR(G, newRing); 4956 4975 } 4957 4976 textra=clock()-to; … … 4959 4978 npert_tmp = nwalk; 4960 4979 endwalks ++; 4961 break; 4962 } 4963 4964 /* check whether the computed Groebner basis a really Groebner basis is. 4965 if no, we perturb the target vector with the maximal "perturbation" 4980 break; 4981 } 4982 4983 /* check whether the computed Groebner basis a really Groebner basis is. 4984 if no, we perturb the target vector with the maximal "perturbation" 4966 4985 degree.*/ 4967 if(MivComp(next_weight, target_weight) == 1 || 4986 if(MivComp(next_weight, target_weight) == 1 || 4968 4987 MivComp(next_weight, curr_weight) == 1 ) 4969 { 4988 { 4970 4989 //Print("\n//ring r_%d = %s;", nwalk, rString(currRing)); 4971 4990 … … 4974 4993 npert[endwalks]=nwalk-npert_tmp; 4975 4994 npert_tmp = nwalk; 4976 4995 4977 4996 endwalks ++; 4978 4997 4979 /*it is very important if the walk only uses one step, e.g. Fate, liu*/ 4998 /*it is very important if the walk only uses one step, e.g. Fate, liu*/ 4980 4999 if(endwalks == 1 && MivComp(next_weight, curr_weight) == 1){ 4981 5000 rChangeCurrRing(XXRing); … … 4983 5002 goto FINISH; 4984 5003 } 4985 H0 = idHead(G); 5004 H0 = idHead(G); 4986 5005 4987 5006 if(MivSame(target_tmp, iv_lp) == 1) … … 4995 5014 else 4996 5015 VMrDefault(target_tmp); 4997 4998 lpRing = currRing; 4999 Glp = idrMoveR(G, newRing); 5016 5017 lpRing = currRing; 5018 Glp = idrMoveR(G, newRing); 5000 5019 H2 = idrMoveR(H0, newRing); 5001 5020 5002 5021 /* Apply Lemma 2.2 in Collart et. al (1997) to check whether 5003 5022 cone(k-1) is equal to cone(k) */ … … 5007 5026 poly t; 5008 5027 if((t=pSub(pHead(Glp->m[i]), pCopy(H2->m[i]))) != NULL) 5009 { 5028 { 5010 5029 pDelete(&t); 5011 5030 idDelete(&H2);//5.5.02 … … 5019 5038 5020 5039 if(nGB == 1) 5021 { 5040 { 5022 5041 G = Glp; 5023 5042 Glp = NULL; … … 5037 5056 { 5038 5057 Overflow_Error = FALSE; 5039 5040 5041 5042 5043 5044 target_weight = MPertVectors(Glp, Mlp, nV); 5045 5046 5047 5048 5049 5050 5051 5052 5053 5054 5058 5059 for(i=0; i<nV; i++) 5060 (*vector_tmp)[i] = (*target_weight)[i]; 5061 5062 delete target_weight; 5063 target_weight = MPertVectors(Glp, Mlp, nV); 5064 5065 if(MivComp(vector_tmp, target_weight)==1) 5066 { 5067 //PrintS("\n// The old and new representaion vector are the same!!"); 5068 G = Glp; 5069 newRing = currRing; 5070 goto OMEGA_OVERFLOW_TRAN_NEW; 5071 } 5072 5073 if(Overflow_Error == TRUE) 5055 5074 { 5056 5075 rChangeCurrRing(newRing); 5057 5076 G = idrMoveR(Glp, lpRing); 5058 5077 goto OMEGA_OVERFLOW_TRAN_NEW; 5059 5078 } 5060 5079 5061 5080 plength3 = TRUE; 5062 5081 pDelete(&p); … … 5072 5091 goto TRAN_LIFTING; 5073 5092 } 5074 5093 5075 5094 5076 5095 npertstep = nwalk; 5077 5096 nwalkpert = 1; 5078 5097 nsteppert ++; 5079 5098 5080 5099 /* 5081 Print("\n// Subroutine needs (%d) steps.", nwalk); 5082 idElements(Glp, "last G in walk:"); 5100 Print("\n// Subroutine needs (%d) steps.", nwalk); 5101 idElements(Glp, "last G in walk:"); 5083 5102 PrintS("\n//****************************************"); 5084 5103 Print("\n// Perturb the original target vector (%d): ", nsteppert); 5085 ivString(target_weight, "new target"); 5086 PrintS("\n//****************************************\n"); 5104 ivString(target_weight, "new target"); 5105 PrintS("\n//****************************************\n"); 5087 5106 */ 5088 5107 rChangeCurrRing(newRing); … … 5090 5109 5091 5110 delete next_weight; 5092 5111 5093 5112 //Print("\n// ring rNEW = %s;", rString(currRing)); 5094 5113 goto COMPUTE_NEW_VECTOR; … … 5098 5117 for(i=nV-1; i>=0; i--) 5099 5118 (*curr_weight)[i] = (*next_weight)[i]; 5100 5119 5101 5120 delete next_weight; 5102 5121 }//while … … 5106 5125 G = idrMoveR(G, lpRing); 5107 5126 5108 FINISH: 5127 FINISH: 5109 5128 delete ivNull; 5110 5129 delete next_weight; 5111 5130 delete iv_lp; 5131 omFree(npert); 5112 5132 5113 5133 #ifdef TIME_TEST 5114 Print("\n// Computation took %d steps and %.2f sec", 5134 Print("\n// Computation took %d steps and %.2f sec", 5115 5135 nwalk, ((double) (clock()-mtim)/1000000)); 5116 5136 … … 5133 5153 clock_t tinput=clock(); 5134 5154 int i, nV = currRing->N; 5135 int nwalk=0, endwalks=0, ntestwinC=1; 5155 int nwalk=0, endwalks=0, ntestwinC=1; 5136 5156 int tp_deg_tmp = tp_deg; 5137 5157 ideal Gomega, M, F, G, M1, F1, Gomega1, Gomega2, G1; 5138 5158 ring endRing, newRing, oldRing, TargetRing; 5139 5159 intvec* next_weight; 5140 intvec* ivNull = new intvec(nV); 5160 intvec* ivNull = new intvec(nV); 5141 5161 intvec* extra_curr_weight = new intvec(nV); 5142 5162 5143 5163 ring YXXRing = currRing; 5144 5164 5145 5165 intvec* iv_M_dpp = MivMatrixOrderlp(nV); 5146 5166 intvec* target_weight;// = Mivlp(nV); … … 5170 5190 if(Overflow_Error == FALSE) 5171 5191 break; 5172 5192 5173 5193 Overflow_Error = TRUE; 5174 5194 tp_deg --; … … 5182 5202 // Print("\n// tp_deg = %d", tp_deg); 5183 5203 // ivString(target_weight, "pert target"); 5184 5185 delete iv_M_dpp; 5204 5205 delete iv_M_dpp; 5186 5206 intvec* hilb_func; 5187 5207 … … 5196 5216 5197 5217 while(1) 5198 { 5218 { 5199 5219 nwalk ++; 5200 5220 nstep ++; 5201 5221 5202 if(nwalk==1) 5222 if(nwalk==1) 5203 5223 goto FIRST_STEP; 5204 5224 … … 5231 5251 VMrDefault(curr_weight); 5232 5252 5233 newRing = currRing; 5253 newRing = currRing; 5234 5254 Gomega1 = idrMoveR(Gomega, oldRing); 5235 5236 #ifdef ENDWALKS 5237 if(endwalks == 1) 5238 { 5239 Print("\n// it is %d-th step!!", nwalk); 5240 idElements(Gomega1, "Gw"); 5241 PrintS("\n// compute a rGB of Gw:"); 5242 } 5255 5256 #ifdef ENDWALKS 5257 if(endwalks == 1) 5258 { 5259 Print("\n// it is %d-th step!!", nwalk); 5260 idElements(Gomega1, "Gw"); 5261 PrintS("\n// compute a rGB of Gw:"); 5262 } 5243 5263 #endif 5244 5264 … … 5249 5269 #else 5250 5270 M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight); 5251 delete hilb_func; 5271 delete hilb_func; 5252 5272 #endif // BUCHBERGER_ALG 5253 5273 xtstd=xtstd+clock()-to; … … 5256 5276 rChangeCurrRing(oldRing); 5257 5277 M1 = idrMoveR(M, newRing); 5258 Gomega2 = idrMoveR(Gomega1, newRing); 5259 to=clock(); 5278 Gomega2 = idrMoveR(Gomega1, newRing); 5279 to=clock(); 5260 5280 5261 5281 // if(endwalks == 1) PrintS("\n// Lifting is still working:"); … … 5266 5286 xtlift=xtlift+clock()-to; 5267 5287 5268 idDelete(&M1); 5269 idDelete(&Gomega2); 5288 idDelete(&M1); 5289 idDelete(&Gomega2); 5270 5290 idDelete(&G); 5271 5291 5272 /* change the ring to newRing */ 5292 /* change the ring to newRing */ 5273 5293 rChangeCurrRing(newRing); 5274 5294 F1 = idrMoveR(F, oldRing); … … 5278 5298 G = kInterRedCC(F1, NULL); 5279 5299 xtred=xtred+clock()-to; 5280 idDelete(&F1); 5300 idDelete(&F1); 5281 5301 5282 5302 if(endwalks == 1) 5283 5303 break; 5284 5304 5285 5305 FIRST_STEP: 5286 5306 Overflow_Error=FALSE; … … 5289 5309 next_weight = MkInterRedNextWeight(curr_weight,target_weight, G); 5290 5310 xtnw=xtnw+clock()-to; 5291 #ifdef PRINT_VECTORS 5311 #ifdef PRINT_VECTORS 5292 5312 MivString(curr_weight, target_weight, next_weight); 5293 5313 #endif … … 5297 5317 delete next_weight; 5298 5318 5299 LASTGB_ALT1: 5319 LASTGB_ALT1: 5300 5320 if(tp_deg > 1){ 5301 5302 5321 nOverflow_Error = Overflow_Error; 5322 tproc = tproc+clock()-tinput; 5303 5323 /* 5304 Print("\n// A subroutine takes %d steps and calls \"Mpwalk\" (1,%d):", 5305 5306 5324 Print("\n// A subroutine takes %d steps and calls \"Mpwalk\" (1,%d):", 5325 nwalk, tp_deg-1); 5326 */ 5307 5327 G1 = Mpwalk_MAltwalk1(G, curr_weight, tp_deg-1); 5308 5328 goto MPW_Finish; … … 5311 5331 newRing = currRing; 5312 5332 ntestwinC = 0; 5313 break; 5333 break; 5314 5334 } 5315 5335 } 5316 5336 5317 5337 if(MivComp(next_weight, ivNull) == 1) 5318 { 5338 { 5319 5339 newRing = currRing; 5320 5340 delete next_weight; … … 5345 5365 { 5346 5366 PrintS("\n// The perturbed target vector doesn't STAY in the correct cone!!"); 5347 if(tp_deg == 1){ 5367 if(tp_deg == 1){ 5348 5368 //Print("\n// subroutine takes %d steps and applys \"std\"", nwalk); 5349 5369 to=clock(); … … 5358 5378 tproc = tproc+clock()-tinput; 5359 5379 /* 5360 Print("\n// B subroutine takes %d steps and calls \"Mpwalk\" (1,%d) :", 5361 5380 Print("\n// B subroutine takes %d steps and calls \"Mpwalk\" (1,%d) :", 5381 nwalk, tp_deg-1); 5362 5382 */ 5363 5383 G1 = Mpwalk_MAltwalk1(G1, curr_weight, tp_deg-1); 5364 5384 } 5365 } 5385 } 5366 5386 5367 5387 MPW_Finish: … … 5369 5389 rChangeCurrRing(YXXRing); 5370 5390 ideal result = idrMoveR(G1, newRing); 5371 5391 5372 5392 delete ivNull; 5373 5393 delete target_weight; … … 5375 5395 /* 5376 5396 Print("\n// \"Mpwalk\" (1,%d) took %d steps and %.2f sec. Overflow_Error (%d)", tp_deg, 5377 5397 nwalk, ((double) clock()-tinput)/1000000, nOverflow_Error); 5378 5398 */ 5379 5399 … … 5397 5417 int i, nV = currRing->N; 5398 5418 int nwalk=0, endwalks=0; 5399 int op_tmp = op_deg; 5419 int op_tmp = op_deg; 5400 5420 ideal Gomega, M, F, G, G1, Gomega1, Gomega2, M1, F1; 5401 5421 ring endRing, newRing, oldRing, TargetRing; 5402 5422 intvec* next_weight; 5403 5423 intvec* iv_M_dp; 5404 intvec* ivNull = new intvec(nV); 5424 intvec* ivNull = new intvec(nV); 5405 5425 intvec* iv_dp = MivUnit(nV);// define (1,1,...,1) 5406 5426 intvec* exivlp = Mivlp(nV); … … 5415 5435 (*last_omega)[i] = 1; 5416 5436 (*last_omega)[0] = 10000; 5417 5437 5418 5438 ring XXRing = currRing; 5419 5439 5420 5440 to=clock(); 5421 5441 /* compute a pertubed weight vector of the original weight vector. 5422 5442 The perturbation degree is recursive decrease until that vector 5423 5443 stays inn the correct cone. */ 5424 while(1) 5444 while(1) 5425 5445 { 5426 5446 if(Overflow_Error == FALSE) … … 5441 5461 else 5442 5462 VMrDefault(cw_tmp); 5443 5463 5444 5464 G = idrMoveR(Go, XXRing); 5445 5465 G = MstdCC(G); … … 5457 5477 if(Overflow_Error == FALSE) 5458 5478 break; 5459 5479 5460 5480 Overflow_Error = TRUE; 5461 5481 op_deg --; … … 5471 5491 5472 5492 while(1) 5473 { 5493 { 5474 5494 nwalk ++; 5475 5495 nstep ++; … … 5485 5505 (*curr_weight)[i] = (*extra_curr_weight)[i]; 5486 5506 delete extra_curr_weight; 5487 5507 5488 5508 newRing = currRing; 5489 5509 goto MSTD_ALT1; … … 5505 5525 VMrDefault(curr_weight); 5506 5526 5507 newRing = currRing; 5508 Gomega1 = idrMoveR(Gomega, oldRing); 5509 5527 newRing = currRing; 5528 Gomega1 = idrMoveR(Gomega, oldRing); 5529 5510 5530 to=clock(); 5511 5531 /* compute a reduced Groebner basis of <Gomega> w.r.t. "newRing" */ … … 5514 5534 #else 5515 5535 M=kStd(Gomega1,NULL,isHomog,NULL,hilb_func,0,NULL,curr_weight); 5516 delete hilb_func; 5536 delete hilb_func; 5517 5537 #endif // BUCHBERGER_ALG 5518 5538 xtstd=xtstd+clock()-to; 5519 5539 5520 /* change the ring to oldRing */ 5540 /* change the ring to oldRing */ 5521 5541 rChangeCurrRing(oldRing); 5522 5542 M1 = idrMoveR(M, newRing); … … 5529 5549 xtlift=xtlift+clock()-to; 5530 5550 5531 idDelete(&M1); 5532 idDelete(&Gomega2); 5533 idDelete(&G); 5534 5535 /* change the ring to newRing */ 5551 idDelete(&M1); 5552 idDelete(&Gomega2); 5553 idDelete(&G); 5554 5555 /* change the ring to newRing */ 5536 5556 rChangeCurrRing(newRing); 5537 5557 F1 = idrMoveR(F, oldRing); 5538 5558 5539 5559 to=clock(); 5540 /* reduce the Groebner basis <G> w.r.t. new ring */ 5560 /* reduce the Groebner basis <G> w.r.t. new ring */ 5541 5561 G = kInterRedCC(F1, NULL); 5542 5562 xtred=xtred+clock()-to; … … 5551 5571 next_weight = MkInterRedNextWeight(curr_weight,target_weight, G); 5552 5572 xtnw=xtnw+clock()-to; 5553 #ifdef PRINT_VECTORS 5573 #ifdef PRINT_VECTORS 5554 5574 MivString(curr_weight, target_weight, next_weight); 5555 5575 #endif 5556 5576 5557 5577 if(Overflow_Error == TRUE) 5558 5578 { 5559 5579 newRing = currRing; 5560 5580 5561 5581 if (currRing->parameter != NULL) 5562 5582 DefRingPar(target_weight); … … 5564 5584 VMrDefault(target_weight); 5565 5585 5566 F1 = idrMoveR(G, newRing); 5586 F1 = idrMoveR(G, newRing); 5567 5587 G = MstdCC(F1); 5568 5588 idDelete(&F1); … … 5587 5607 { 5588 5608 MSTD_ALT1: 5589 5590 5591 5592 Print("\n// main routine takes %d steps and calls \"Mpwalk\" (1,%d):", 5593 nwalk, tp_deg); 5594 5595 5596 // the "recursive-modified" perturbation walk alg (1,tp_deg) 5609 nOverflow_Error = Overflow_Error; 5610 tproc = clock()-xftinput; 5611 /* 5612 Print("\n// main routine takes %d steps and calls \"Mpwalk\" (1,%d):", 5613 nwalk, tp_deg); 5614 */ 5615 // compute the red. GB of <G> w.r.t. the lex order by 5616 // the "recursive-modified" perturbation walk alg (1,tp_deg) 5597 5617 G = Mpwalk_MAltwalk1(G, curr_weight, tp_deg); 5598 5618 delete next_weight; … … 5621 5641 #ifdef TIME_TEST 5622 5642 5623 Print("\n// \"Main procedure\" took %d steps, %.2f sec. and Overflow_Error(%d)", 5624 5643 Print("\n// \"Main procedure\" took %d steps, %.2f sec. and Overflow_Error(%d)", 5644 nwalk, ((double) tproc)/1000000, nOverflow_Error); 5625 5645 5626 5646 TimeStringFractal(xftinput, tostd, xtif, xtstd,xtextra, xtlift, xtred, xtnw); … … 5629 5649 Print("\n// Overflow_Error? (%d)", Overflow_Error); 5630 5650 Print("\n// Awalk1 took %d steps.\n", nstep); 5631 #endif 5651 #endif 5632 5652 return(result); 5633 5653 }
Note: See TracChangeset
for help on using the changeset viewer.