Changeset 4df33b in git
 Timestamp:
 Sep 28, 2010, 12:19:26 PM (13 years ago)
 Branches:
 (u'spielwiese', '0d6b7fcd9813a1ca1ed4220cfa2b104b97a0a003')
 Children:
 c94e60f986f0fb36f24beb94af60519b01cb313c
 Parents:
 9a24da7e54cf925089f9c1420b3c4b79dc39eea6
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

Singular/misc_ip.cc
r9a24da r4df33b 53 53 } 54 54 55 void divTimes_ui(mpz_t n, long d, int* times)55 void divTimes_ui(mpz_t n, unsigned long d, int* times) 56 56 { 57 57 *times = 0; … … 67 67 mpz_clear(r); 68 68 mpz_clear(q); 69 } 70 71 static inline void divTimes_ui_ui(unsigned long *n, unsigned long d, int* times) 72 { 73 *times = 0; 74 unsigned long q=(*n) / d; 75 unsigned long r=(*n) % d; 76 while (r==0) 77 { 78 (*times)++; 79 (*n)=q; 80 q=(*n)/d; r=(*n)%d; 81 } 69 82 } 70 83 … … 180 193 lists primes = (lists)omAllocBin(slists_bin); primes>Init(1000); 181 194 int* multiplicities = new int[1000]; 182 183 divTimes_ui(nn, 2, &tt); 184 if (tt > 0) 185 { 186 setListEntry_ui(primes, index, 2); 187 multiplicities[index++] = tt; 188 } 189 190 divTimes_ui(nn, 3, &tt); 191 if (tt > 0) 192 { 193 setListEntry_ui(primes, index, 3); 194 multiplicities[index++] = tt; 195 } 196 197 unsigned long p_ui=5; add = 2; 198 mpz_sqrt(sr, nn); 199 if ((mpz_cmp_ui(b, 0) == 0)  (mpz_cmp(pb, sr) > 0)) mpz_set(pb, sr); 200 int limit=mpz_get_ui(pb); 201 if ((limit==0)(mpz_cmp_ui(pb,limit)!=0)) limit=1<<31; 202 while (p_ui <=limit) 203 { 204 divTimes_ui(nn, p_ui, &tt); 195 int positive=1; 196 197 if (!nlIsZero(n)) 198 { 199 if (!nlGreaterZero(n)) 200 { 201 positive=1; 202 mpz_neg(nn,nn); 203 } 204 divTimes_ui(nn, 2, &tt); 205 205 if (tt > 0) 206 206 { 207 setListEntry_ui(primes, index, p_ui);207 setListEntry_ui(primes, index, 2); 208 208 multiplicities[index++] = tt; 209 //mpz_sqrt(sr, nn); 210 //if ((mpz_cmp_ui(b, 0) == 0)  (mpz_cmp(pb, sr) > 0)) mpz_set(pb, sr); 211 if (mpz_size1(nn)<=2) 212 { 209 } 210 211 divTimes_ui(nn, 3, &tt); 212 if (tt > 0) 213 { 214 setListEntry_ui(primes, index, 3); 215 multiplicities[index++] = tt; 216 } 217 218 unsigned long p_ui=5; add = 2; 219 mpz_sqrt(sr, nn); 220 if ((mpz_cmp_ui(b, 0) == 0)  (mpz_cmp(pb, sr) > 0)) mpz_set(pb, sr); 221 int limit=mpz_get_ui(pb); 222 if ((limit==0)(mpz_cmp_ui(pb,limit)!=0)) limit=1<<31; 223 while (p_ui <=limit) 224 { 225 divTimes_ui(nn, p_ui, &tt); 226 if (tt > 0) 227 { 228 setListEntry_ui(primes, index, p_ui); 229 multiplicities[index++] = tt; 230 //mpz_sqrt(sr, nn); 231 //if ((mpz_cmp_ui(b, 0) == 0)  (mpz_cmp(pb, sr) > 0)) mpz_set(pb, sr); 232 if (mpz_size1(nn)<=2) 233 { 234 mpz_sqrt(sr, nn); 235 if ((mpz_cmp_ui(b, 0) == 0)  (mpz_cmp(pb, sr) > 0)) mpz_set(pb, sr); 236 unsigned long l=mpz_get_ui(sr); 237 if (l<limit) limit=l; 238 if (mpz_size1(nn)<=1) 239 { 240 unsigned long nn_ui=mpz_get_ui(nn); 241 while ((p_ui <=limit)&&(nn_ui>1)) 242 { 243 divTimes_ui_ui(&nn_ui, p_ui, &tt); 244 if (tt > 0) 245 { 246 setListEntry_ui(primes, index, p_ui); 247 multiplicities[index++] = tt; 248 if (nn_ui>(limit/6)) limit=nn_ui/6; 249 } 250 p_ui +=add; 251 add += 2; if (add == 6) add = 2; 252 } 253 mpz_set_ui(nn,nn_ui); 254 break; 255 } 256 } 257 } 258 p_ui +=add; 259 add += 2; if (add == 6) add = 2; 260 } 261 mpz_set_ui(p, p_ui); 262 mpz_sqrt(sr, nn); 263 if ((mpz_cmp_ui(b, 0) == 0)  (mpz_cmp(pb, sr) > 0)) mpz_set(pb, sr); 264 while (mpz_cmp(pb, p) >= 0) 265 { 266 divTimes(nn, p, &tt); 267 if (tt > 0) 268 { 269 setListEntry(primes, index, p); 270 multiplicities[index++] = tt; 213 271 mpz_sqrt(sr, nn); 214 272 if ((mpz_cmp_ui(b, 0) == 0)  (mpz_cmp(pb, sr) > 0)) mpz_set(pb, sr); 215 unsigned long l=mpz_get_ui(sr); 216 if (l<limit) limit=l; 217 } 218 } 219 p_ui +=add; 220 add += 2; if (add == 6) add = 2; 221 } 222 mpz_set_ui(p, p_ui); 223 mpz_sqrt(sr, nn); 224 if ((mpz_cmp_ui(b, 0) == 0)  (mpz_cmp(pb, sr) > 0)) mpz_set(pb, sr); 225 while (mpz_cmp(pb, p) >= 0) 226 { 227 divTimes(nn, p, &tt); 228 if (tt > 0) 229 { 230 setListEntry(primes, index, p); 231 multiplicities[index++] = tt; 232 mpz_sqrt(sr, nn); 233 if ((mpz_cmp_ui(b, 0) == 0)  (mpz_cmp(pb, sr) > 0)) mpz_set(pb, sr); 234 } 235 mpz_add_ui(p, p, add); 236 add += 2; if (add == 6) add = 2; 237 } 238 if ((mpz_cmp_ui(nn, 1) > 0) && 239 ((mpz_cmp_ui(b, 0) == 0)  (mpz_cmp(nn, b) <= 0))) 240 { 241 setListEntry(primes, index, nn); 242 multiplicities[index++] = 1; 243 mpz_set_ui(nn, 1); 273 } 274 mpz_add_ui(p, p, add); 275 add += 2; if (add == 6) add = 2; 276 } 277 if ((mpz_cmp_ui(nn, 1) > 0) && 278 ((mpz_cmp_ui(b, 0) == 0)  (mpz_cmp(nn, b) <= 0))) 279 { 280 setListEntry(primes, index, nn); 281 multiplicities[index++] = 1; 282 mpz_set_ui(nn, 1); 283 } 244 284 } 245 285 … … 265 305 lists L=(lists)omAllocBin(slists_bin); 266 306 L>Init(4); 307 if (positive==1) mpz_neg(nn,nn); 267 308 setListEntry(L, 0, nn); 268 309 L>m[1].rtyp = LIST_CMD; L>m[1].data = (void*)primesL;
Note: See TracChangeset
for help on using the changeset viewer.