Changeset 880433 in git
- Timestamp:
- Oct 31, 2012, 5:50:34 PM (11 years ago)
- Branches:
- (u'spielwiese', 'd1b01e9d51ade4b46b745d3bada5c5f3696be3a8')
- Children:
- 8e3dbbf7777bf1a5f6acd0133065b3542d71b4fc
- Parents:
- f7286a23b786af032ab46ebb6351640245468492
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
Singular/interpolation.cc
rf7286a r880433 37 37 #define exponent int 38 38 39 modp_number myp; // all modp computation done mod myp40 int myp_index; // index of small prime in Singular table of primes39 static modp_number myp; // all modp computation done mod myp 40 static int myp_index; // index of small prime in Singular table of primes 41 41 42 42 static inline modp_number modp_mul (modp_number x,modp_number y) … … 87 87 typedef bool *coord_exist_table; 88 88 89 int final_base_dim; // dimension of the quotient space, known from the beginning90 int last_solve_column; // last non-zero column in "solve" part of matrix, used for speed up91 int n_points; // modp_number of ideals (points)92 int *multiplicity; // multiplicities of points93 int variables; // modp_number of variables94 int max_coord; // maximal possible coordinate product used during Evaluation95 bool only_modp; // perform only one modp computations96 97 modp_coordinates *modp_points; // coordinates of points for modp problem - used by Evaluate (this is also initial data for only modp)98 q_coordinates *q_points; // coordinates of points for rational data (not used for modp)99 int_coordinates *int_points; // coordinates of points for integer data - used to check generators (not used for modp)100 coord_exist_table *coord_exist; // checks whether all coordinates has been initialized101 mon_list_entry *check_list; // monomials to be checked in next stages102 coordinates *points; // power products of coordinates of points used in modp cycles103 condition_type *condition_list; // conditions stored in an array104 mon_list_entry *lt_list; // leading terms found so far105 mon_list_entry *base_list; // standard monomials found so far106 row_list_entry *row_list; // rows of the matrix (both parts)107 modp_number *my_row; // one special row to perform operations108 modp_number *my_solve_row; // one special row to find the linear dependence ("solve" part)109 mono_type *column_name; // monomials assigned to columns in solve_row110 111 int n_results; // modp_number of performed modp computations (not discarded)112 modp_number modp_denom; // denominator of mod p computations113 modp_result_entry *modp_result; // list of results for various mod p calculations (used for modp - first result is the desired one)114 modp_result_entry *cur_result; // pointer to current result (as before)115 modp_number *congr; // primes used in computations (chinese remainder theorem) (not used for modp)116 modp_number *in_gamma; // inverts used in chinese remainder theorem (not used for modp)117 mpz_t bigcongr; // result, in fact, is given in mod bigcongr (not used for modp)118 119 mpz_t *polycoef; // polynomial integercoefficients (not used for modp)120 mono_type *polyexp; // polynomial exponents89 static int final_base_dim; // dimension of the quotient space, known from the beginning 90 static int last_solve_column; // last non-zero column in "solve" part of matrix, used for speed up 91 static int n_points; // modp_number of ideals (points) 92 static int *multiplicity; // multiplicities of points 93 static int variables; // modp_number of variables 94 static int max_coord; // maximal possible coordinate product used during Evaluation 95 static bool only_modp; // perform only one modp computations 96 97 static modp_coordinates *modp_points; // coordinates of points for modp problem - used by Evaluate (this is also initial data for only modp) 98 static q_coordinates *q_points; // coordinates of points for rational data (not used for modp) 99 static int_coordinates *int_points; // coordinates of points for integer data - used to check generators (not used for modp) 100 static coord_exist_table *coord_exist; // checks whether all coordinates has been initialized 101 static mon_list_entry *check_list; // monomials to be checked in next stages 102 static coordinates *points; // power products of coordinates of points used in modp cycles 103 static condition_type *condition_list; // conditions stored in an array 104 static mon_list_entry *lt_list; // leading terms found so far 105 static mon_list_entry *base_list; // standard monomials found so far 106 static row_list_entry *row_list; // rows of the matrix (both parts) 107 static modp_number *my_row; // one special row to perform operations 108 static modp_number *my_solve_row; // one special row to find the linear dependence ("solve" part) 109 static mono_type *column_name; // monomials assigned to columns in solve_row 110 111 static int n_results; // modp_number of performed modp computations (not discarded) 112 static modp_number modp_denom; // denominator of mod p computations 113 static modp_result_entry *modp_result; // list of results for various mod p calculations (used for modp - first result is the desired one) 114 static modp_result_entry *cur_result; // pointer to current result (as before) 115 static modp_number *congr; // primes used in computations (chinese remainder theorem) (not used for modp) 116 static modp_number *in_gamma; // inverts used in chinese remainder theorem (not used for modp) 117 static mpz_t bigcongr; // result, in fact, is given in mod bigcongr (not used for modp) 118 119 static mpz_t *polycoef; // polynomial integercoefficients (not used for modp) 120 static mono_type *polyexp; // polynomial exponents 121 121 122 122 struct gen_list_struct {mpz_t *polycoef; … … 125 125 typedef struct gen_list_struct gen_list_entry; 126 126 127 gen_list_entry *gen_list=NULL; // list of resulting generators - output data (integer version)128 129 int generic_n_generators; // modp_number of generators - should be the same for all modp comp (not used for modp)130 mono_type *generic_column_name; // monomials assigned to columns in solve_row - should be the same for all modp comp (!!! used for modp)131 mon_list_entry *generic_lt=NULL; // leading terms for ordered generators - should be the same for all modp comp (not used for modp)132 int good_primes; // modp_number of good primes so far;133 int bad_primes; // modp_number of bad primes so far;134 mpz_t common_denom; // common denominator used to force points coordinates to Z (not used for modp)135 bool denom_divisible; // common denominator is divisible by p (not used for modp)136 137 poly comparizon_p1; //polynomials used to do comparizons by Singular138 poly comparizon_p2;139 140 modp_number *modp_Reverse; // reverses in mod p141 142 bool protocol; // true to show the protocol127 static gen_list_entry *gen_list=NULL; // list of resulting generators - output data (integer version) 128 129 static int generic_n_generators; // modp_number of generators - should be the same for all modp comp (not used for modp) 130 static mono_type *generic_column_name; // monomials assigned to columns in solve_row - should be the same for all modp comp (!!! used for modp) 131 static mon_list_entry *generic_lt=NULL; // leading terms for ordered generators - should be the same for all modp comp (not used for modp) 132 static int good_primes; // modp_number of good primes so far; 133 static int bad_primes; // modp_number of bad primes so far; 134 static mpz_t common_denom; // common denominator used to force points coordinates to Z (not used for modp) 135 static bool denom_divisible; // common denominator is divisible by p (not used for modp) 136 137 static poly comparizon_p1; //polynomials used to do comparizons by Singular 138 static poly comparizon_p2; 139 140 static modp_number *modp_Reverse; // reverses in mod p 141 142 static bool protocol; // true to show the protocol 143 143 144 144 #ifdef checksize 145 int maximal_size=0;145 static int maximal_size=0; 146 146 #endif 147 147 … … 225 225 #endif 226 226 227 modp_number OneInverse(modp_number a,modp_number p) // computes inverse of d mod p without using tables227 static modp_number OneInverse(modp_number a,modp_number p) // computes inverse of d mod p without using tables 228 228 { 229 229 long u, v, u0, u1, u2, q, r; … … 254 254 } 255 255 256 int CalcBaseDim () // returns the maximal (expected) dimension of quotient space256 static int CalcBaseDim () // returns the maximal (expected) dimension of quotient space 257 257 { 258 258 int c; … … 273 273 } 274 274 275 bool EqualMon (mono_type m1, mono_type m2) // compares two monomials, true if equal, false otherwise275 static bool EqualMon (mono_type m1, mono_type m2) // compares two monomials, true if equal, false otherwise 276 276 { 277 277 int i; … … 281 281 } 282 282 283 exponent MonDegree (mono_type mon) // computes the degree of a monomial283 static exponent MonDegree (mono_type mon) // computes the degree of a monomial 284 284 { 285 285 exponent v=0; … … 289 289 } 290 290 291 bool Greater (mono_type m1, mono_type m2) // true if m1 > m2, false otherwise. uses Singular comparing function291 static bool Greater (mono_type m1, mono_type m2) // true if m1 > m2, false otherwise. uses Singular comparing function 292 292 { 293 293 for (int j = variables; j; j--) … … 302 302 } 303 303 304 mon_list_entry* MonListAdd (mon_list_entry *list, mono_type mon) // adds one entry to the list of monomials structure304 static mon_list_entry* MonListAdd (mon_list_entry *list, mono_type mon) // adds one entry to the list of monomials structure 305 305 { 306 306 mon_list_entry *curptr=list; … … 327 327 } 328 328 329 mono_type MonListElement (mon_list_entry *list, int n) // returns ith element from list of monomials - no error checking!329 static mono_type MonListElement (mon_list_entry *list, int n) // returns ith element from list of monomials - no error checking! 330 330 { 331 331 mon_list_entry *cur=list; … … 335 335 } 336 336 337 mono_type ZeroMonomial () // allocates memory for new monomial with all enries equal 0337 static mono_type ZeroMonomial () // allocates memory for new monomial with all enries equal 0 338 338 { 339 339 mono_type m; … … 342 342 } 343 343 344 void GeneralInit () // general initialization345 { 346 int i,j ,k;344 static void GeneralInit () // general initialization 345 { 346 int i,j; 347 347 points=(coordinates*)omAlloc(sizeof(coordinates)*n_points); 348 348 for (i=0;i<n_points;i++) … … 403 403 } 404 404 405 void InitProcData () // initialization for procedures which do computations mod p406 { 407 int i ,j;405 static void InitProcData () // initialization for procedures which do computations mod p 406 { 407 int i; 408 408 check_list=MonListAdd (check_list,ZeroMonomial ()); 409 409 my_row=(modp_number*)omAlloc0(sizeof(modp_number)*final_base_dim); … … 439 439 } 440 440 441 mon_list_entry* FreeMonList (mon_list_entry *list) // destroys a list of monomials, returns NULL pointer441 static mon_list_entry* FreeMonList (mon_list_entry *list) // destroys a list of monomials, returns NULL pointer 442 442 { 443 443 mon_list_entry *ptr; … … 453 453 } 454 454 455 void GeneralDone () // to be called before exit to free memory456 { 457 int i,j ,k;455 static void GeneralDone () // to be called before exit to free memory 456 { 457 int i,j; 458 458 for (i=0;i<n_points;i++) 459 459 { … … 504 504 } 505 505 506 void FreeProcData () // to be called after one modp computation to free memory506 static void FreeProcData () // to be called after one modp computation to free memory 507 507 { 508 508 int i; … … 531 531 } 532 532 533 void modp_Evaluate(modp_number *ev, mono_type mon, condition_type con) // evaluates monomial on condition (modp)533 static void modp_Evaluate(modp_number *ev, mono_type mon, condition_type con) // evaluates monomial on condition (modp) 534 534 { 535 535 int i; … … 554 554 } 555 555 556 void int_Evaluate(mpz_t ev, mono_type mon, condition_type con) // (***) evaluates monomial on condition for integer numbers556 static void int_Evaluate(mpz_t ev, mono_type mon, condition_type con) // (***) evaluates monomial on condition for integer numbers 557 557 { 558 558 int i; … … 581 581 } 582 582 583 void ProduceRow(mono_type mon) // produces a row for monomial - first part is an evaluated part, second 0 to obtain the coefs of dependence583 static void ProduceRow(mono_type mon) // produces a row for monomial - first part is an evaluated part, second 0 to obtain the coefs of dependence 584 584 { 585 585 modp_number *row; … … 602 602 } 603 603 604 void IntegerPoints () // produces integer points from rationals by scaling the coordinate system604 static void IntegerPoints () // produces integer points from rationals by scaling the coordinate system 605 605 { 606 606 int i,j; … … 630 630 } 631 631 632 void int_PrepareProducts () // prepares coordinates of points and products for modp case (from integer coefs) 633 { 634 bool ok=true; 632 static void int_PrepareProducts () // prepares coordinates of points and products for modp case (from integer coefs) 633 { 635 634 int i,j,k; 636 635 mpz_t big_myp; … … 655 654 } 656 655 657 void modp_PrepareProducts () // prepares products for modp computation from modp data656 static void modp_PrepareProducts () // prepares products for modp computation from modp data 658 657 { 659 658 int i,j,k; … … 669 668 } 670 669 671 void MakeConditions () // prepares a list of conditions from list of multiplicities670 static void MakeConditions () // prepares a list of conditions from list of multiplicities 672 671 { 673 672 condition_type *con; … … 700 699 } 701 700 702 void ReduceRow () // reduces my_row by previous rows, does the same for solve row701 static void ReduceRow () // reduces my_row by previous rows, does the same for solve row 703 702 { 704 703 if (row_list==NULL) return ; … … 776 775 } 777 776 778 bool RowIsZero () // check whether a row is zero777 static bool RowIsZero () // check whether a row is zero 779 778 { 780 779 bool zero=1; … … 790 789 } 791 790 792 bool DivisibleMon (mono_type m1, mono_type m2) // checks whether m1 is divisible by m2791 static bool DivisibleMon (mono_type m1, mono_type m2) // checks whether m1 is divisible by m2 793 792 { 794 793 int i; … … 798 797 } 799 798 800 void ReduceCheckListByMon (mono_type m) // from check_list monomials divisible by m are thrown out799 static void ReduceCheckListByMon (mono_type m) // from check_list monomials divisible by m are thrown out 801 800 { 802 801 mon_list_entry *c_ptr; … … 826 825 } 827 826 828 void TakeNextMonomial (mono_type mon) // reads first monomial from check_list, then it is deleted827 static void TakeNextMonomial (mono_type mon) // reads first monomial from check_list, then it is deleted 829 828 { 830 829 mon_list_entry *n_check_list; … … 839 838 } 840 839 841 void UpdateCheckList (mono_type m) // adds all monomials which are "next" to m (divisible by m and degree +1)840 static void UpdateCheckList (mono_type m) // adds all monomials which are "next" to m (divisible by m and degree +1) 842 841 { 843 842 int i; … … 850 849 } 851 850 852 void ReduceCheckListByLTs () // deletes all monomials from check_list which are divisible by one of the leading terms851 static void ReduceCheckListByLTs () // deletes all monomials from check_list which are divisible by one of the leading terms 853 852 { 854 853 mon_list_entry *ptr; … … 861 860 } 862 861 863 void RowListAdd (int first_col, mono_type mon) // puts a row into matrix862 static void RowListAdd (int first_col, mono_type mon) // puts a row into matrix 864 863 { 865 864 row_list_entry *ptr; … … 888 887 889 888 890 void PrepareRow (mono_type mon) // prepares a row and puts it into matrix889 static void PrepareRow (mono_type mon) // prepares a row and puts it into matrix 891 890 { 892 891 modp_number *row; … … 941 940 } 942 941 943 void NewResultEntry () // creates an entry for new modp result942 static void NewResultEntry () // creates an entry for new modp result 944 943 { 945 944 modp_result_entry *temp; 946 int i;947 945 temp=(modp_result_entry*)omAlloc0(sizeof(modp_result_entry)); 948 946 if (cur_result==NULL) … … 964 962 } 965 963 966 void FreeResultEntry (modp_result_entry *e) // destroys the result entry, without worrying about where it is964 static void FreeResultEntry (modp_result_entry *e) // destroys the result entry, without worrying about where it is 967 965 { 968 966 generator_entry *cur_gen; … … 981 979 982 980 983 void NewGenerator (mono_type mon) // new generator in modp comp found, shoul be stored on the list 984 { 985 int i; 981 static void NewGenerator (mono_type mon) // new generator in modp comp found, shoul be stored on the list 982 { 986 983 generator_entry *cur_ptr; 987 984 generator_entry *prev_ptr; … … 1006 1003 } 1007 1004 1008 void MultGenerators () // before reconstructing, all denominators must be cancelled1005 static void MultGenerators () // before reconstructing, all denominators must be cancelled 1009 1006 { 1010 1007 #ifndef integerstrategy … … 1044 1041 #endif 1045 1042 1046 modp_number TakePrime (modp_number p) // takes "previous" (smaller) prime1043 static modp_number TakePrime (modp_number p) // takes "previous" (smaller) prime 1047 1044 { 1048 1045 #ifdef HAVE_FACTORY … … 1054 1051 } 1055 1052 1056 void PrepareChinese (int n) // initialization for CRA1053 static void PrepareChinese (int n) // initialization for CRA 1057 1054 { 1058 1055 int i,k; … … 1081 1078 } 1082 1079 1083 void CloseChinese (int n) // after CRA1080 static void CloseChinese () // after CRA 1084 1081 { 1085 1082 omFree(in_gamma); … … 1088 1085 } 1089 1086 1090 void ClearGCD () // divides polynomials in basis by gcd of coefficients1087 static void ClearGCD () // divides polynomials in basis by gcd of coefficients 1091 1088 { 1092 1089 bool first_gcd=true; … … 1111 1108 } 1112 1109 1113 void ReconstructGenerator (int ngen,int n,bool show) // recostruction of generator from various modp comp 1114 { 1115 int size; 1110 static void ReconstructGenerator (int ngen,int n) // recostruction of generator from various modp comp 1111 { 1116 1112 int i,j,k; 1117 1113 int coef; … … 1177 1173 } 1178 1174 1179 void Discard () // some unlucky prime occures1175 static void Discard () // some unlucky prime occures 1180 1176 { 1181 1177 modp_result_entry *temp; … … 1227 1223 } 1228 1224 1229 void modp_SetColumnNames () // used by modp - sets column names, the old table will be destroyed1225 static void modp_SetColumnNames () // used by modp - sets column names, the old table will be destroyed 1230 1226 { 1231 1227 int i; … … 1233 1229 } 1234 1230 1235 void CheckColumnSequence () // checks if scheme of computations is as generic one1231 static void CheckColumnSequence () // checks if scheme of computations is as generic one 1236 1232 { 1237 1233 int i; … … 1309 1305 #endif 1310 1306 1311 bool CheckGenerator () // evaluates generator to check whether it is good1312 { 1313 mpz_t val,sum ,temp;1307 static bool CheckGenerator () // evaluates generator to check whether it is good 1308 { 1309 mpz_t val,sum; 1314 1310 int con,i; 1315 1311 mpz_init(val); … … 1336 1332 } 1337 1333 1338 void ClearGenList ()1334 static void ClearGenList () 1339 1335 { 1340 1336 gen_list_entry *temp; … … 1355 1351 } 1356 1352 1357 void UpdateGenList ()1353 static void UpdateGenList () 1358 1354 { 1359 1355 gen_list_entry *temp,*prev; … … 1416 1412 1417 1413 1418 void modp_Main ()1414 static void modp_Main () 1419 1415 { 1420 1416 mono_type cur_mon; … … 1473 1469 } 1474 1470 1475 void ResolveCoeff (mpq_t c, number m)1476 { 1477 1478 1479 1480 1481 1482 1483 1484 1485 1486 1487 1488 1489 1490 1491 1492 1493 1494 1471 static void ResolveCoeff (mpq_t c, number m) 1472 { 1473 if ((long)m & SR_INT) 1474 { 1475 long m_val=SR_TO_INT(m); 1476 mpq_set_si(c,m_val,1); 1477 } 1478 else 1479 { 1480 if (m->s<2) 1481 { 1482 mpz_set(mpq_numref(c),m->z); 1483 mpz_set(mpq_denref(c),m->n); 1484 mpq_canonicalize(c); 1485 } 1486 else 1487 { 1488 mpq_set_z(c,m->z); 1489 } 1490 } 1495 1491 } 1496 1492 … … 1673 1669 for (i=0;i<generic_n_generators;i++) 1674 1670 { 1675 ReconstructGenerator (i,modp_cycles ,false);1671 ReconstructGenerator (i,modp_cycles); 1676 1672 correct_gen=CheckGenerator (); 1677 1673 if (!correct_gen) … … 1693 1689 Print("maximal size of output: %d, precision bound: %d.\n",maximalsize,mpz_sizeinbase(bigcongr,10)); 1694 1690 #endif 1695 CloseChinese ( modp_cycles);1691 CloseChinese (); 1696 1692 modp_cycles=modp_cycles+10; 1697 1693 }
Note: See TracChangeset
for help on using the changeset viewer.