[7bf145] | 1 | /*****************************************************************************\ |
---|
[806c18] | 2 | * Computer Algebra System SINGULAR |
---|
[7bf145] | 3 | \*****************************************************************************/ |
---|
| 4 | /** @file facFqBivarUtil.cc |
---|
[806c18] | 5 | * |
---|
| 6 | * This file provides utility functions for bivariate factorization |
---|
[7bf145] | 7 | * |
---|
| 8 | * @author Martin Lee |
---|
| 9 | * |
---|
| 10 | **/ |
---|
| 11 | /*****************************************************************************/ |
---|
| 12 | |
---|
[e4fe2b] | 13 | #include "config.h" |
---|
[7bf145] | 14 | |
---|
| 15 | #include "cf_map.h" |
---|
[81d96c] | 16 | #include "algext.h" |
---|
[7bf145] | 17 | #include "cf_map_ext.h" |
---|
[6db552] | 18 | #include "templates/ftmpl_functions.h" |
---|
[7bf145] | 19 | #include "ExtensionInfo.h" |
---|
[fecc08] | 20 | #include "cf_algorithm.h" |
---|
| 21 | #include "cf_factory.h" |
---|
[07718c3] | 22 | #include "cf_util.h" |
---|
[fecc08] | 23 | #include "imm.h" |
---|
| 24 | #include "cf_iter.h" |
---|
[7bf145] | 25 | #include "facFqBivarUtil.h" |
---|
[f876a66] | 26 | #include "cfNewtonPolygon.h" |
---|
| 27 | #include "facHensel.h" |
---|
[0e2e23] | 28 | #include "facMul.h" |
---|
[7bf145] | 29 | |
---|
| 30 | |
---|
[806c18] | 31 | void append (CFList& factors1, const CFList& factors2) |
---|
[7bf145] | 32 | { |
---|
[806c18] | 33 | for (CFListIterator i= factors2; i.hasItem(); i++) |
---|
[7bf145] | 34 | { |
---|
| 35 | if (!i.getItem().inCoeffDomain()) |
---|
| 36 | factors1.append (i.getItem()); |
---|
| 37 | } |
---|
| 38 | return; |
---|
| 39 | } |
---|
| 40 | |
---|
[806c18] | 41 | void decompress (CFList& factors, const CFMap& N) |
---|
[7bf145] | 42 | { |
---|
[806c18] | 43 | for (CFListIterator i= factors; i.hasItem(); i++) |
---|
[7bf145] | 44 | i.getItem()= N (i.getItem()); |
---|
[f876a66] | 45 | } |
---|
| 46 | |
---|
| 47 | void decompress (CFFList& factors, const CFMap& N) |
---|
| 48 | { |
---|
| 49 | for (CFFListIterator i= factors; i.hasItem(); i++) |
---|
| 50 | i.getItem()= CFFactor (N (i.getItem().factor()), i.getItem().exp()); |
---|
[7bf145] | 51 | } |
---|
| 52 | |
---|
[806c18] | 53 | void appendSwapDecompress (CFList& factors1, const CFList& factors2, |
---|
| 54 | const CFList& factors3, const bool swap1, |
---|
| 55 | const bool swap2, const CFMap& N) |
---|
[7bf145] | 56 | { |
---|
| 57 | Variable x= Variable (1); |
---|
| 58 | Variable y= Variable (2); |
---|
[806c18] | 59 | for (CFListIterator i= factors1; i.hasItem(); i++) |
---|
[7bf145] | 60 | { |
---|
[806c18] | 61 | if (swap1) |
---|
[7bf145] | 62 | { |
---|
[806c18] | 63 | if (!swap2) |
---|
[7bf145] | 64 | i.getItem()= swapvar (i.getItem(), x, y); |
---|
| 65 | } |
---|
[806c18] | 66 | else |
---|
[7bf145] | 67 | { |
---|
[806c18] | 68 | if (swap2) |
---|
[7bf145] | 69 | i.getItem()= swapvar (i.getItem(), y, x); |
---|
| 70 | } |
---|
| 71 | i.getItem()= N (i.getItem()); |
---|
| 72 | } |
---|
[806c18] | 73 | for (CFListIterator i= factors2; i.hasItem(); i++) |
---|
[7bf145] | 74 | factors1.append (N (i.getItem())); |
---|
[806c18] | 75 | for (CFListIterator i= factors3; i.hasItem(); i++) |
---|
[7bf145] | 76 | factors1.append (N (i.getItem())); |
---|
| 77 | return; |
---|
| 78 | } |
---|
| 79 | |
---|
[806c18] | 80 | void swapDecompress (CFList& factors, const bool swap, const CFMap& N) |
---|
[7bf145] | 81 | { |
---|
| 82 | Variable x= Variable (1); |
---|
| 83 | Variable y= Variable (2); |
---|
[806c18] | 84 | for (CFListIterator i= factors; i.hasItem(); i++) |
---|
[7bf145] | 85 | { |
---|
[806c18] | 86 | if (swap) |
---|
[7bf145] | 87 | i.getItem()= swapvar (i.getItem(), x, y); |
---|
| 88 | i.getItem()= N (i.getItem()); |
---|
| 89 | } |
---|
| 90 | return; |
---|
| 91 | } |
---|
| 92 | |
---|
[806c18] | 93 | static inline |
---|
| 94 | bool GFInExtensionHelper (const CanonicalForm& F, const int number) |
---|
[7bf145] | 95 | { |
---|
| 96 | if (F.isOne()) return false; |
---|
| 97 | InternalCF* buf; |
---|
| 98 | int exp; |
---|
| 99 | bool result= false; |
---|
[806c18] | 100 | if (F.inBaseDomain()) |
---|
[7bf145] | 101 | { |
---|
| 102 | buf= F.getval(); |
---|
| 103 | exp= imm2int (buf); |
---|
[806c18] | 104 | if (exp%number != 0) |
---|
[7bf145] | 105 | return true; |
---|
[806c18] | 106 | else |
---|
[7bf145] | 107 | return result; |
---|
| 108 | } |
---|
[806c18] | 109 | else |
---|
[7bf145] | 110 | { |
---|
[806c18] | 111 | for (CFIterator i= F; i.hasTerms(); i++) |
---|
[7bf145] | 112 | { |
---|
| 113 | result= GFInExtensionHelper (i.coeff(), number); |
---|
| 114 | if (result == true) |
---|
| 115 | return result; |
---|
| 116 | } |
---|
| 117 | } |
---|
| 118 | return result; |
---|
| 119 | } |
---|
| 120 | |
---|
[806c18] | 121 | static inline |
---|
[f876a66] | 122 | bool FqInExtensionHelper (const CanonicalForm& F, const CanonicalForm& gamma, |
---|
| 123 | const CanonicalForm& delta, CFList& source, |
---|
| 124 | CFList& dest) |
---|
[7bf145] | 125 | { |
---|
| 126 | bool result= false; |
---|
[806c18] | 127 | if (F.inBaseDomain()) |
---|
[7bf145] | 128 | return result; |
---|
[806c18] | 129 | else if (F.inCoeffDomain()) |
---|
[7bf145] | 130 | { |
---|
[806c18] | 131 | if (!fdivides (gamma, F)) |
---|
[7bf145] | 132 | return true; |
---|
[806c18] | 133 | else |
---|
[f876a66] | 134 | { |
---|
| 135 | int pos= findItem (source, F); |
---|
| 136 | if (pos > 0) |
---|
| 137 | return false; |
---|
| 138 | Variable a; |
---|
| 139 | hasFirstAlgVar (F, a); |
---|
| 140 | int bound= ipower (getCharacteristic(), degree (getMipo (a))); |
---|
| 141 | CanonicalForm buf= 1; |
---|
| 142 | for (int i= 1; i < bound; i++) |
---|
| 143 | { |
---|
| 144 | buf *= gamma; |
---|
| 145 | if (buf == F) |
---|
| 146 | { |
---|
| 147 | source.append (buf); |
---|
| 148 | dest.append (power (delta, i)); |
---|
| 149 | return false; |
---|
| 150 | } |
---|
| 151 | } |
---|
| 152 | return true; |
---|
| 153 | } |
---|
[7bf145] | 154 | } |
---|
[806c18] | 155 | else |
---|
[7bf145] | 156 | { |
---|
[806c18] | 157 | for (CFIterator i= F; i.hasTerms(); i++) |
---|
[7bf145] | 158 | { |
---|
[f876a66] | 159 | result= FqInExtensionHelper (i.coeff(), gamma, delta, source, dest); |
---|
[7bf145] | 160 | if (result == true) |
---|
[806c18] | 161 | return result; |
---|
[7bf145] | 162 | } |
---|
| 163 | } |
---|
| 164 | return result; |
---|
| 165 | } |
---|
| 166 | |
---|
| 167 | bool isInExtension (const CanonicalForm& F, const CanonicalForm& gamma, |
---|
[f876a66] | 168 | const int k, const CanonicalForm& delta, |
---|
| 169 | CFList& source, CFList& dest) |
---|
[7bf145] | 170 | { |
---|
| 171 | bool result; |
---|
[806c18] | 172 | if (CFFactory::gettype() == GaloisFieldDomain) |
---|
[7bf145] | 173 | { |
---|
| 174 | int p= getCharacteristic(); |
---|
| 175 | int orderFieldExtension= ipower (p, getGFDegree()) - 1; |
---|
| 176 | int order= ipower (p, k) - 1; |
---|
| 177 | int number= orderFieldExtension/order; |
---|
| 178 | result= GFInExtensionHelper (F, number); |
---|
| 179 | return result; |
---|
| 180 | } |
---|
[806c18] | 181 | else |
---|
[7bf145] | 182 | { |
---|
[f876a66] | 183 | result= FqInExtensionHelper (F, gamma, delta, source, dest); |
---|
[7bf145] | 184 | return result; |
---|
| 185 | } |
---|
| 186 | } |
---|
| 187 | |
---|
[806c18] | 188 | CanonicalForm |
---|
| 189 | mapDown (const CanonicalForm& F, const ExtensionInfo& info, CFList& source, |
---|
| 190 | CFList& dest) |
---|
[7bf145] | 191 | { |
---|
| 192 | int k= info.getGFDegree(); |
---|
| 193 | Variable beta= info.getAlpha(); |
---|
| 194 | CanonicalForm primElem= info.getGamma(); |
---|
| 195 | CanonicalForm imPrimElem= info.getDelta(); |
---|
| 196 | if (k > 1) |
---|
| 197 | return GFMapDown (F, k); |
---|
| 198 | else if (k == 1) |
---|
| 199 | return F; |
---|
[c1b9927] | 200 | if (/*k==0 &&*/ beta == Variable (1)) |
---|
[7bf145] | 201 | return F; |
---|
[c1b9927] | 202 | else /*if (k==0 && beta != Variable (1))*/ |
---|
[618da5] | 203 | return mapDown (F, imPrimElem, primElem, beta, source, dest); |
---|
[7bf145] | 204 | } |
---|
| 205 | |
---|
[806c18] | 206 | void appendTestMapDown (CFList& factors, const CanonicalForm& f, |
---|
| 207 | const ExtensionInfo& info, CFList& source, CFList& dest) |
---|
[7bf145] | 208 | { |
---|
| 209 | int k= info.getGFDegree(); |
---|
| 210 | Variable beta= info.getBeta(); |
---|
| 211 | Variable alpha= info.getAlpha(); |
---|
| 212 | CanonicalForm delta= info.getDelta(); |
---|
| 213 | CanonicalForm gamma= info.getGamma(); |
---|
| 214 | CanonicalForm g= f; |
---|
| 215 | int degMipoBeta; |
---|
[1372ae] | 216 | if (!k && beta.level() == 1) |
---|
[7bf145] | 217 | degMipoBeta= 1; |
---|
[1372ae] | 218 | else if (!k && beta.level() != 1) |
---|
[7bf145] | 219 | degMipoBeta= degree (getMipo (beta)); |
---|
[806c18] | 220 | if (k > 1) |
---|
[7bf145] | 221 | { |
---|
[f876a66] | 222 | if (!isInExtension (g, gamma, k, delta, source, dest)) |
---|
[7bf145] | 223 | { |
---|
| 224 | g= GFMapDown (g, k); |
---|
| 225 | factors.append (g); |
---|
| 226 | } |
---|
| 227 | } |
---|
[806c18] | 228 | else if (k == 1) |
---|
| 229 | { |
---|
[1709e1] | 230 | if (!isInExtension (g, gamma, k, delta, source, dest)) |
---|
[7bf145] | 231 | factors.append (g); |
---|
| 232 | } |
---|
[806c18] | 233 | else if (!k && beta == Variable (1)) |
---|
[7bf145] | 234 | { |
---|
| 235 | if (degree (g, alpha) < degMipoBeta) |
---|
| 236 | factors.append (g); |
---|
| 237 | } |
---|
[806c18] | 238 | else if (!k && beta != Variable (1)) |
---|
[7bf145] | 239 | { |
---|
[f876a66] | 240 | if (!isInExtension (g, gamma, k, delta, source, dest)) |
---|
[7bf145] | 241 | { |
---|
[618da5] | 242 | g= mapDown (g, delta, gamma, alpha, source, dest); |
---|
[7bf145] | 243 | factors.append (g); |
---|
[806c18] | 244 | } |
---|
[7bf145] | 245 | } |
---|
| 246 | return; |
---|
| 247 | } |
---|
| 248 | |
---|
[806c18] | 249 | void |
---|
| 250 | appendMapDown (CFList& factors, const CanonicalForm& g, |
---|
| 251 | const ExtensionInfo& info, CFList& source, CFList& dest) |
---|
[7bf145] | 252 | { |
---|
| 253 | int k= info.getGFDegree(); |
---|
| 254 | Variable beta= info.getBeta(); |
---|
| 255 | Variable alpha= info.getAlpha(); |
---|
| 256 | CanonicalForm gamma= info.getGamma(); |
---|
| 257 | CanonicalForm delta= info.getDelta(); |
---|
[806c18] | 258 | if (k > 1) |
---|
[7bf145] | 259 | factors.append (GFMapDown (g, k)); |
---|
[806c18] | 260 | else if (k == 1) |
---|
[7bf145] | 261 | factors.append (g); |
---|
[806c18] | 262 | else if (!k && beta == Variable (1)) |
---|
[7bf145] | 263 | factors.append (g); |
---|
[806c18] | 264 | else if (!k && beta != Variable (1)) |
---|
[618da5] | 265 | factors.append (mapDown (g, delta, gamma, alpha, source, dest)); |
---|
[7bf145] | 266 | return; |
---|
| 267 | } |
---|
| 268 | |
---|
[806c18] | 269 | void normalize (CFList& factors) |
---|
[7bf145] | 270 | { |
---|
[1682691] | 271 | CanonicalForm lcinv; |
---|
[806c18] | 272 | for (CFListIterator i= factors; i.hasItem(); i++) |
---|
[1682691] | 273 | { |
---|
| 274 | lcinv= 1/Lc (i.getItem()); |
---|
| 275 | i.getItem() *= lcinv; |
---|
| 276 | } |
---|
[7bf145] | 277 | return; |
---|
| 278 | } |
---|
| 279 | |
---|
[f876a66] | 280 | void normalize (CFFList& factors) |
---|
| 281 | { |
---|
[1682691] | 282 | CanonicalForm lcinv; |
---|
[c1b9927] | 283 | for (CFFListIterator i= factors; i.hasItem(); i++) |
---|
[1682691] | 284 | { |
---|
| 285 | lcinv= 1/ Lc (i.getItem().factor()); |
---|
| 286 | i.getItem()= CFFactor (i.getItem().factor()*lcinv, |
---|
[f876a66] | 287 | i.getItem().exp()); |
---|
[1682691] | 288 | } |
---|
[f876a66] | 289 | return; |
---|
| 290 | } |
---|
| 291 | |
---|
[806c18] | 292 | CFList subset (int index [], const int& s, const CFArray& elements, |
---|
| 293 | bool& noSubset) |
---|
[7bf145] | 294 | { |
---|
| 295 | int r= elements.size(); |
---|
| 296 | int i= 0; |
---|
| 297 | CFList result; |
---|
| 298 | noSubset= false; |
---|
[806c18] | 299 | if (index[s - 1] == 0) |
---|
[7bf145] | 300 | { |
---|
[806c18] | 301 | while (i < s) |
---|
[7bf145] | 302 | { |
---|
| 303 | index[i]= i + 1; |
---|
| 304 | result.append (elements[i]); |
---|
| 305 | i++; |
---|
[806c18] | 306 | } |
---|
[7bf145] | 307 | return result; |
---|
| 308 | } |
---|
| 309 | int buf; |
---|
| 310 | int k; |
---|
| 311 | bool found= false; |
---|
[806c18] | 312 | if (index[s - 1] == r) |
---|
[7bf145] | 313 | { |
---|
[806c18] | 314 | if (index[0] == r - s + 1) |
---|
[7bf145] | 315 | { |
---|
| 316 | noSubset= true; |
---|
| 317 | return result; |
---|
| 318 | } |
---|
| 319 | else { |
---|
[806c18] | 320 | while (found == false) |
---|
[7bf145] | 321 | { |
---|
| 322 | if (index[s - 2 - i] < r - i - 1) |
---|
| 323 | found= true; |
---|
[806c18] | 324 | i++; |
---|
[7bf145] | 325 | } |
---|
| 326 | buf= index[s - i - 1]; |
---|
| 327 | k= 0; |
---|
[806c18] | 328 | while (s - i - 1 + k < s) |
---|
[7bf145] | 329 | { |
---|
| 330 | index[s - i - 1 + k]= buf + k + 1; |
---|
| 331 | k++; |
---|
[806c18] | 332 | } |
---|
[7bf145] | 333 | } |
---|
[806c18] | 334 | for (int j= 0; j < s; j++) |
---|
[7bf145] | 335 | result.append (elements[index[j] - 1]); |
---|
| 336 | return result; |
---|
| 337 | } |
---|
[806c18] | 338 | else |
---|
[7bf145] | 339 | { |
---|
| 340 | index[s - 1] += 1; |
---|
| 341 | for (int j= 0; j < s; j++) |
---|
| 342 | result.append (elements[index[j] - 1]); |
---|
| 343 | return result; |
---|
| 344 | } |
---|
[806c18] | 345 | } |
---|
[7bf145] | 346 | |
---|
[806c18] | 347 | CFArray copy (const CFList& list) |
---|
[7bf145] | 348 | { |
---|
| 349 | CFArray array= CFArray (list.length()); |
---|
| 350 | int j= 0; |
---|
[806c18] | 351 | for (CFListIterator i= list; i.hasItem(); i++, j++) |
---|
[7bf145] | 352 | array[j]= i.getItem(); |
---|
| 353 | return array; |
---|
| 354 | } |
---|
| 355 | |
---|
[806c18] | 356 | void indexUpdate (int index [], const int& subsetSize, const int& setSize, |
---|
| 357 | bool& noSubset) |
---|
[7bf145] | 358 | { |
---|
| 359 | noSubset= false; |
---|
[806c18] | 360 | if (subsetSize > setSize) |
---|
[7bf145] | 361 | { |
---|
| 362 | noSubset= true; |
---|
| 363 | return; |
---|
| 364 | } |
---|
[1673386] | 365 | int * v= new int [setSize]; |
---|
[806c18] | 366 | for (int i= 0; i < setSize; i++) |
---|
[7bf145] | 367 | v[i]= index[i]; |
---|
[806c18] | 368 | if (subsetSize == 1) |
---|
[7bf145] | 369 | { |
---|
| 370 | v[0]= v[0] - 1; |
---|
[806c18] | 371 | if (v[0] >= setSize) |
---|
[7bf145] | 372 | { |
---|
| 373 | noSubset= true; |
---|
[1673386] | 374 | delete [] v; |
---|
[7bf145] | 375 | return; |
---|
| 376 | } |
---|
| 377 | } |
---|
[806c18] | 378 | else |
---|
[7bf145] | 379 | { |
---|
[806c18] | 380 | if (v[subsetSize - 1] - v[0] + 1 == subsetSize && v[0] > 1) |
---|
[7bf145] | 381 | { |
---|
[806c18] | 382 | if (v[0] + subsetSize - 1 > setSize) |
---|
| 383 | { |
---|
[7bf145] | 384 | noSubset= true; |
---|
[1673386] | 385 | delete [] v; |
---|
[7bf145] | 386 | return; |
---|
| 387 | } |
---|
| 388 | v[0]= v[0] - 1; |
---|
[806c18] | 389 | for (int i= 1; i < subsetSize - 1; i++) |
---|
[7bf145] | 390 | v[i]= v[i - 1] + 1; |
---|
| 391 | v[subsetSize - 1]= v[subsetSize - 2]; |
---|
| 392 | } |
---|
[806c18] | 393 | else |
---|
[7bf145] | 394 | { |
---|
[806c18] | 395 | if (v[0] + subsetSize - 1 > setSize) |
---|
| 396 | { |
---|
[7bf145] | 397 | noSubset= true; |
---|
[1673386] | 398 | delete [] v; |
---|
[7bf145] | 399 | return; |
---|
| 400 | } |
---|
| 401 | for (int i= 1; i < subsetSize - 1; i++) |
---|
| 402 | v[i]= v[i - 1] + 1; |
---|
| 403 | v[subsetSize - 1]= v[subsetSize - 2]; |
---|
| 404 | } |
---|
| 405 | } |
---|
| 406 | for (int i= 0; i < setSize; i++) |
---|
| 407 | index[i]= v[i]; |
---|
[1673386] | 408 | delete [] v; |
---|
[7bf145] | 409 | } |
---|
| 410 | |
---|
[806c18] | 411 | int subsetDegree (const CFList& S) |
---|
[7bf145] | 412 | { |
---|
| 413 | int result= 0; |
---|
| 414 | for (CFListIterator i= S; i.hasItem(); i++) |
---|
| 415 | result += degree (i.getItem(), Variable (1)); |
---|
| 416 | return result; |
---|
| 417 | } |
---|
| 418 | |
---|
| 419 | CFFList multiplicity (CanonicalForm& F, const CFList& factors) |
---|
| 420 | { |
---|
| 421 | if (F.inCoeffDomain()) |
---|
| 422 | return CFFList (CFFactor (F, 1)); |
---|
| 423 | CFFList result; |
---|
| 424 | int multi= 0; |
---|
[21b8f4] | 425 | CanonicalForm quot; |
---|
[7bf145] | 426 | for (CFListIterator i= factors; i.hasItem(); i++) |
---|
| 427 | { |
---|
[21b8f4] | 428 | while (fdivides (i.getItem(), F, quot)) |
---|
[7bf145] | 429 | { |
---|
| 430 | multi++; |
---|
[21b8f4] | 431 | F= quot; |
---|
[7bf145] | 432 | } |
---|
[575e1d] | 433 | if (multi > 0) |
---|
| 434 | result.append (CFFactor (i.getItem(), multi)); |
---|
[7bf145] | 435 | multi= 0; |
---|
| 436 | } |
---|
| 437 | return result; |
---|
| 438 | } |
---|
| 439 | |
---|
[615ca8] | 440 | #ifdef HAVE_NTL |
---|
[f876a66] | 441 | CFArray |
---|
| 442 | logarithmicDerivative (const CanonicalForm& F, const CanonicalForm& G, int l, |
---|
| 443 | CanonicalForm& Q |
---|
| 444 | ) |
---|
| 445 | { |
---|
| 446 | Variable x= Variable (2); |
---|
| 447 | Variable y= Variable (1); |
---|
| 448 | CanonicalForm xToL= power (x, l); |
---|
| 449 | CanonicalForm q,r; |
---|
| 450 | CanonicalForm logDeriv; |
---|
| 451 | |
---|
| 452 | q= newtonDiv (F, G, xToL); |
---|
| 453 | |
---|
| 454 | logDeriv= mulMod2 (q, deriv (G, y), xToL); |
---|
| 455 | |
---|
[9614bb] | 456 | int j= degree (logDeriv, y) + 1; |
---|
[f876a66] | 457 | CFArray result= CFArray (j); |
---|
[9614bb] | 458 | CFIterator ii; |
---|
[f876a66] | 459 | for (CFIterator i= logDeriv; i.hasTerms() && !logDeriv.isZero(); i++) |
---|
| 460 | { |
---|
[9614bb] | 461 | for (ii= i.coeff(); ii.hasTerms(); ii++) |
---|
| 462 | result[ii.exp()] += ii.coeff()*power (x,i.exp()); |
---|
[f876a66] | 463 | } |
---|
| 464 | Q= q; |
---|
| 465 | return result; |
---|
| 466 | } |
---|
| 467 | |
---|
| 468 | CFArray |
---|
| 469 | logarithmicDerivative (const CanonicalForm& F, const CanonicalForm& G, int l, |
---|
| 470 | int oldL, const CanonicalForm& oldQ, CanonicalForm& Q |
---|
| 471 | ) |
---|
| 472 | { |
---|
| 473 | Variable x= Variable (2); |
---|
| 474 | Variable y= Variable (1); |
---|
| 475 | CanonicalForm xToL= power (x, l); |
---|
| 476 | CanonicalForm xToOldL= power (x, oldL); |
---|
| 477 | CanonicalForm xToLOldL= power (x, l-oldL); |
---|
| 478 | CanonicalForm q,r; |
---|
| 479 | CanonicalForm logDeriv; |
---|
| 480 | |
---|
[ac3fcca] | 481 | CanonicalForm bufF; |
---|
| 482 | if ((oldL > 100 && l - oldL < 50) || (oldL < 100 && l - oldL < 30)) |
---|
| 483 | { |
---|
| 484 | bufF= F; |
---|
| 485 | CanonicalForm oldF= mulMod2 (G, oldQ, xToL); |
---|
| 486 | bufF -= oldF; |
---|
| 487 | bufF= div (bufF, xToOldL); |
---|
| 488 | } |
---|
| 489 | else |
---|
| 490 | { |
---|
| 491 | //middle product style computation of [G*oldQ]^{l}_{oldL} |
---|
| 492 | CanonicalForm G3= div (G, xToOldL); |
---|
| 493 | CanonicalForm Up= mulMod2 (G3, oldQ, xToLOldL); |
---|
[fd80670] | 494 | CanonicalForm xToOldL2= power (x, (oldL+1)/2); |
---|
[ac3fcca] | 495 | CanonicalForm G2= mod (G, xToOldL); |
---|
| 496 | CanonicalForm G1= div (G2, xToOldL2); |
---|
| 497 | CanonicalForm G0= mod (G2, xToOldL2); |
---|
| 498 | CanonicalForm oldQ1= div (oldQ, xToOldL2); |
---|
| 499 | CanonicalForm oldQ0= mod (oldQ, xToOldL2); |
---|
[fd80670] | 500 | CanonicalForm Mid; |
---|
| 501 | if (oldL % 2 == 1) |
---|
| 502 | Mid= mulMod2 (G1, oldQ1*x, xToLOldL); |
---|
| 503 | else |
---|
| 504 | Mid= mulMod2 (G1, oldQ1, xToLOldL); |
---|
[ac3fcca] | 505 | //computation of Low might be faster using a real middle product? |
---|
| 506 | CanonicalForm Low= mulMod2 (G0, oldQ1, xToOldL)+mulMod2 (G1, oldQ0, xToOldL); |
---|
[fd80670] | 507 | Low= div (Low, power (x, oldL/2)); |
---|
| 508 | Low= mod (Low, xToLOldL); |
---|
[ac3fcca] | 509 | Up += Mid + Low; |
---|
| 510 | bufF= div (F, xToOldL); |
---|
| 511 | bufF -= Up; |
---|
| 512 | } |
---|
[f876a66] | 513 | |
---|
[72f1e4b] | 514 | if (l-oldL > 0) |
---|
| 515 | q= newtonDiv (bufF, G, xToLOldL); |
---|
| 516 | else |
---|
| 517 | q= 0; |
---|
[f876a66] | 518 | q *= xToOldL; |
---|
| 519 | q += oldQ; |
---|
| 520 | |
---|
| 521 | logDeriv= mulMod2 (q, deriv (G, y), xToL); |
---|
| 522 | |
---|
[9614bb] | 523 | int j= degree (logDeriv,y) + 1; |
---|
[f876a66] | 524 | CFArray result= CFArray (j); |
---|
[9614bb] | 525 | CFIterator ii; |
---|
[f876a66] | 526 | for (CFIterator i= logDeriv; i.hasTerms() && !logDeriv.isZero(); i++) |
---|
| 527 | { |
---|
[9614bb] | 528 | for (ii= i.coeff(); ii.hasTerms(); ii++) |
---|
| 529 | result[ii.exp()] += ii.coeff()*power (x,i.exp()); |
---|
[f876a66] | 530 | } |
---|
| 531 | Q= q; |
---|
| 532 | return result; |
---|
| 533 | } |
---|
[615ca8] | 534 | #endif |
---|
[f876a66] | 535 | |
---|
| 536 | void |
---|
| 537 | writeInMatrix (CFMatrix& M, const CFArray& A, const int column, |
---|
| 538 | const int startIndex |
---|
| 539 | ) |
---|
| 540 | { |
---|
[050d1b] | 541 | ASSERT (A.size () - startIndex >= 0, "wrong starting index"); |
---|
| 542 | ASSERT (A.size () - startIndex <= M.rows(), "wrong starting index"); |
---|
[f876a66] | 543 | ASSERT (column > 0 && column <= M.columns(), "wrong column"); |
---|
| 544 | if (A.size() - startIndex <= 0) return; |
---|
| 545 | int j= 1; |
---|
| 546 | for (int i= startIndex; i < A.size(); i++, j++) |
---|
| 547 | M (j, column)= A [i]; |
---|
| 548 | } |
---|
| 549 | |
---|
| 550 | CFArray getCoeffs (const CanonicalForm& F, const int k) |
---|
| 551 | { |
---|
[050d1b] | 552 | ASSERT (F.isUnivariate() || F.inCoeffDomain(), "univariate input expected"); |
---|
[f876a66] | 553 | if (degree (F, 2) < k) |
---|
| 554 | return CFArray(); |
---|
| 555 | |
---|
| 556 | CFArray result= CFArray (degree (F) - k + 1); |
---|
| 557 | CFIterator j= F; |
---|
| 558 | for (int i= degree (F); i >= k; i--) |
---|
| 559 | { |
---|
| 560 | if (j.exp() == i) |
---|
| 561 | { |
---|
| 562 | result [i - k]= j.coeff(); |
---|
[c1b9927] | 563 | j++; |
---|
[f876a66] | 564 | if (!j.hasTerms()) |
---|
| 565 | return result; |
---|
| 566 | } |
---|
| 567 | else |
---|
| 568 | result[i - k]= 0; |
---|
| 569 | } |
---|
| 570 | return result; |
---|
| 571 | } |
---|
| 572 | |
---|
| 573 | CFArray getCoeffs (const CanonicalForm& F, const int k, const Variable& alpha) |
---|
| 574 | { |
---|
[050d1b] | 575 | ASSERT (F.isUnivariate() || F.inCoeffDomain(), "univariate input expected"); |
---|
[f876a66] | 576 | if (degree (F, 2) < k) |
---|
| 577 | return CFArray (); |
---|
| 578 | |
---|
| 579 | int d= degree (getMipo (alpha)); |
---|
| 580 | CFArray result= CFArray ((degree (F) - k + 1)*d); |
---|
| 581 | CFIterator j= F; |
---|
| 582 | CanonicalForm buf; |
---|
| 583 | CFIterator iter; |
---|
| 584 | for (int i= degree (F); i >= k; i--) |
---|
| 585 | { |
---|
| 586 | if (j.exp() == i) |
---|
| 587 | { |
---|
| 588 | iter= j.coeff(); |
---|
| 589 | for (int l= degree (j.coeff(), alpha); l >= 0; l--) |
---|
| 590 | { |
---|
| 591 | if (iter.exp() == l) |
---|
| 592 | { |
---|
| 593 | result [(i - k)*d + l]= iter.coeff(); |
---|
| 594 | iter++; |
---|
| 595 | if (!iter.hasTerms()) |
---|
| 596 | break; |
---|
| 597 | } |
---|
| 598 | } |
---|
| 599 | j++; |
---|
| 600 | if (!j.hasTerms()) |
---|
| 601 | return result; |
---|
| 602 | } |
---|
| 603 | else |
---|
| 604 | { |
---|
| 605 | for (int l= 0; l < d; l++) |
---|
| 606 | result[(i - k)*d + l]= 0; |
---|
| 607 | } |
---|
| 608 | } |
---|
| 609 | return result; |
---|
| 610 | } |
---|
| 611 | |
---|
| 612 | #ifdef HAVE_NTL |
---|
| 613 | CFArray |
---|
| 614 | getCoeffs (const CanonicalForm& G, const int k, const int l, const int degMipo, |
---|
| 615 | const Variable& alpha, const CanonicalForm& evaluation, |
---|
| 616 | const mat_zz_p& M) |
---|
| 617 | { |
---|
[050d1b] | 618 | ASSERT (G.isUnivariate() || G.inCoeffDomain(), "univariate input expected"); |
---|
[f876a66] | 619 | CanonicalForm F= G (G.mvar() - evaluation, G.mvar()); |
---|
| 620 | if (F.isZero()) |
---|
| 621 | return CFArray (); |
---|
| 622 | |
---|
| 623 | Variable y= Variable (2); |
---|
| 624 | F= F (power (y, degMipo), y); |
---|
| 625 | F= F (y, alpha); |
---|
| 626 | zz_pX NTLF= convertFacCF2NTLzzpX (F); |
---|
| 627 | NTLF.rep.SetLength (l*degMipo); |
---|
| 628 | NTLF.rep= M*NTLF.rep; |
---|
| 629 | NTLF.normalize(); |
---|
| 630 | F= convertNTLzzpX2CF (NTLF, y); |
---|
| 631 | |
---|
| 632 | if (degree (F, 2) < k) |
---|
| 633 | return CFArray(); |
---|
| 634 | |
---|
| 635 | CFArray result= CFArray (degree (F) - k + 1); |
---|
| 636 | |
---|
| 637 | CFIterator j= F; |
---|
| 638 | for (int i= degree (F); i >= k; i--) |
---|
| 639 | { |
---|
| 640 | if (j.exp() == i) |
---|
| 641 | { |
---|
| 642 | result [i - k]= j.coeff(); |
---|
[c1b9927] | 643 | j++; |
---|
[f876a66] | 644 | if (!j.hasTerms()) |
---|
| 645 | return result; |
---|
| 646 | } |
---|
| 647 | else |
---|
| 648 | result[i - k]= 0; |
---|
| 649 | } |
---|
| 650 | return result; |
---|
| 651 | } |
---|
| 652 | #endif |
---|
| 653 | |
---|
[542864] | 654 | int * computeBounds (const CanonicalForm& F, int& n, bool& isIrreducible) |
---|
[f876a66] | 655 | { |
---|
| 656 | n= degree (F, 1); |
---|
| 657 | int* result= new int [n]; |
---|
| 658 | int sizeOfNewtonPolygon; |
---|
| 659 | int** newtonPolyg= newtonPolygon (F, sizeOfNewtonPolygon); |
---|
| 660 | |
---|
[542864] | 661 | isIrreducible= false; |
---|
| 662 | if (sizeOfNewtonPolygon == 3) |
---|
| 663 | { |
---|
| 664 | bool check1= |
---|
| 665 | (newtonPolyg[0][0]==0 || newtonPolyg[1][0]==0 || newtonPolyg[2][0]==0); |
---|
| 666 | if (check1) |
---|
| 667 | { |
---|
| 668 | bool check2= |
---|
| 669 | (newtonPolyg[0][1]==0 || newtonPolyg[1][1]==0 || newtonPolyg[2][0]==0); |
---|
| 670 | if (check2) |
---|
| 671 | { |
---|
| 672 | int p=getCharacteristic(); |
---|
| 673 | int d=1; |
---|
| 674 | char bufGFName='Z'; |
---|
| 675 | bool GF= (CFFactory::gettype()==GaloisFieldDomain); |
---|
| 676 | if (GF) |
---|
| 677 | { |
---|
| 678 | d= getGFDegree(); |
---|
| 679 | bufGFName=gf_name; |
---|
| 680 | } |
---|
| 681 | setCharacteristic(0); |
---|
| 682 | CanonicalForm tmp= gcd (newtonPolyg[0][0],newtonPolyg[0][1]); |
---|
| 683 | tmp= gcd (tmp, newtonPolyg[1][0]); |
---|
| 684 | tmp= gcd (tmp, newtonPolyg[1][1]); |
---|
| 685 | tmp= gcd (tmp, newtonPolyg[2][0]); |
---|
| 686 | tmp= gcd (tmp, newtonPolyg[2][1]); |
---|
| 687 | isIrreducible= (tmp==1); |
---|
| 688 | if (GF) |
---|
| 689 | setCharacteristic (p, d, bufGFName); |
---|
| 690 | else |
---|
| 691 | setCharacteristic(p); |
---|
| 692 | } |
---|
| 693 | } |
---|
| 694 | } |
---|
| 695 | |
---|
[f876a66] | 696 | int minX, minY, maxX, maxY; |
---|
| 697 | minX= newtonPolyg [0] [0]; |
---|
| 698 | minY= newtonPolyg [0] [1]; |
---|
| 699 | maxX= minX; |
---|
| 700 | maxY= minY; |
---|
| 701 | for (int i= 1; i < sizeOfNewtonPolygon; i++) |
---|
| 702 | { |
---|
| 703 | if (minX > newtonPolyg [i] [0]) |
---|
| 704 | minX= newtonPolyg [i] [0]; |
---|
| 705 | if (maxX < newtonPolyg [i] [0]) |
---|
| 706 | maxX= newtonPolyg [i] [0]; |
---|
| 707 | if (minY > newtonPolyg [i] [1]) |
---|
| 708 | minY= newtonPolyg [i] [1]; |
---|
| 709 | if (maxY < newtonPolyg [i] [1]) |
---|
| 710 | maxY= newtonPolyg [i] [1]; |
---|
| 711 | } |
---|
| 712 | |
---|
| 713 | int k= maxX; |
---|
| 714 | for (int i= 0; i < n; i++) |
---|
| 715 | { |
---|
| 716 | if (i + 1 > maxY || i + 1 < minY) |
---|
| 717 | { |
---|
| 718 | result [i]= 0; |
---|
| 719 | continue; |
---|
| 720 | } |
---|
| 721 | int* point= new int [2]; |
---|
| 722 | point [0]= k; |
---|
| 723 | point [1]= i + 1; |
---|
| 724 | while (!isInPolygon (newtonPolyg, sizeOfNewtonPolygon, point) && k > 0) |
---|
| 725 | { |
---|
| 726 | k--; |
---|
| 727 | point [0]= k; |
---|
| 728 | } |
---|
| 729 | result [i]= k; |
---|
| 730 | k= maxX; |
---|
| 731 | delete [] point; |
---|
| 732 | } |
---|
| 733 | |
---|
| 734 | return result; |
---|
| 735 | } |
---|
| 736 | |
---|
[dd8047] | 737 | int |
---|
| 738 | substituteCheck (const CanonicalForm& F, const Variable& x) |
---|
| 739 | { |
---|
| 740 | if (F.inCoeffDomain()) |
---|
| 741 | return 0; |
---|
| 742 | if (degree (F, x) < 0) |
---|
| 743 | return 0; |
---|
| 744 | CanonicalForm f= swapvar (F, F.mvar(), x); |
---|
| 745 | int sizef= 0; |
---|
| 746 | for (CFIterator i= f; i.hasTerms(); i++, sizef++) |
---|
| 747 | { |
---|
| 748 | if (i.exp() == 1) |
---|
| 749 | return 0; |
---|
| 750 | } |
---|
| 751 | int * expf= new int [sizef]; |
---|
| 752 | int j= 0; |
---|
| 753 | for (CFIterator i= f; i.hasTerms(); i++, j++) |
---|
| 754 | expf [j]= i.exp(); |
---|
| 755 | |
---|
| 756 | int indf= sizef - 1; |
---|
| 757 | if (expf[indf] == 0) |
---|
| 758 | indf--; |
---|
| 759 | |
---|
| 760 | int result= expf[indf]; |
---|
| 761 | for (int i= indf - 1; i >= 0; i--) |
---|
| 762 | { |
---|
| 763 | if (expf [i]%result != 0) |
---|
| 764 | { |
---|
| 765 | delete [] expf; |
---|
| 766 | return 0; |
---|
| 767 | } |
---|
| 768 | } |
---|
| 769 | |
---|
| 770 | delete [] expf; |
---|
| 771 | return result; |
---|
| 772 | } |
---|
| 773 | |
---|
[686ce3] | 774 | static int |
---|
| 775 | substituteCheck (const CanonicalForm& F, const CanonicalForm& G) |
---|
| 776 | { |
---|
| 777 | if (F.inCoeffDomain() || G.inCoeffDomain()) |
---|
| 778 | return 0; |
---|
| 779 | Variable x= Variable (1); |
---|
| 780 | if (degree (F, x) <= 1 || degree (G, x) <= 1) |
---|
| 781 | return 0; |
---|
| 782 | CanonicalForm f= swapvar (F, F.mvar(), x); |
---|
| 783 | CanonicalForm g= swapvar (G, G.mvar(), x); |
---|
| 784 | int sizef= 0; |
---|
| 785 | int sizeg= 0; |
---|
| 786 | for (CFIterator i= f; i.hasTerms(); i++, sizef++) |
---|
| 787 | { |
---|
| 788 | if (i.exp() == 1) |
---|
| 789 | return 0; |
---|
| 790 | } |
---|
| 791 | for (CFIterator i= g; i.hasTerms(); i++, sizeg++) |
---|
| 792 | { |
---|
| 793 | if (i.exp() == 1) |
---|
| 794 | return 0; |
---|
| 795 | } |
---|
| 796 | int * expf= new int [sizef]; |
---|
| 797 | int * expg= new int [sizeg]; |
---|
| 798 | int j= 0; |
---|
| 799 | for (CFIterator i= f; i.hasTerms(); i++, j++) |
---|
| 800 | { |
---|
| 801 | expf [j]= i.exp(); |
---|
| 802 | } |
---|
| 803 | j= 0; |
---|
| 804 | for (CFIterator i= g; i.hasTerms(); i++, j++) |
---|
| 805 | { |
---|
| 806 | expg [j]= i.exp(); |
---|
| 807 | } |
---|
| 808 | |
---|
| 809 | int indf= sizef - 1; |
---|
| 810 | int indg= sizeg - 1; |
---|
| 811 | if (expf[indf] == 0) |
---|
| 812 | indf--; |
---|
| 813 | if (expg[indg] == 0) |
---|
| 814 | indg--; |
---|
| 815 | |
---|
| 816 | if ((expg[indg]%expf [indf] != 0 && expf[indf]%expg[indg] != 0) || |
---|
| 817 | (expg[indg] == 1 && expf[indf] == 1)) |
---|
| 818 | { |
---|
| 819 | delete [] expg; |
---|
| 820 | delete [] expf; |
---|
| 821 | return 0; |
---|
| 822 | } |
---|
| 823 | |
---|
| 824 | int result; |
---|
| 825 | if (expg [indg]%expf [indf] == 0) |
---|
| 826 | result= expf[indf]; |
---|
| 827 | else |
---|
| 828 | result= expg[indg]; |
---|
| 829 | for (int i= indf - 1; i >= 0; i--) |
---|
| 830 | { |
---|
| 831 | if (expf [i]%result != 0) |
---|
| 832 | { |
---|
| 833 | delete [] expf; |
---|
| 834 | delete [] expg; |
---|
| 835 | return 0; |
---|
| 836 | } |
---|
| 837 | } |
---|
| 838 | |
---|
| 839 | for (int i= indg - 1; i >= 0; i--) |
---|
| 840 | { |
---|
| 841 | if (expg [i]%result != 0) |
---|
| 842 | { |
---|
| 843 | delete [] expf; |
---|
| 844 | delete [] expg; |
---|
| 845 | return 0; |
---|
| 846 | } |
---|
| 847 | } |
---|
| 848 | |
---|
| 849 | delete [] expg; |
---|
| 850 | delete [] expf; |
---|
| 851 | return result; |
---|
| 852 | } |
---|
| 853 | |
---|
| 854 | int recSubstituteCheck (const CanonicalForm& F, const int d) |
---|
| 855 | { |
---|
| 856 | if (F.inCoeffDomain()) |
---|
| 857 | return 0; |
---|
| 858 | Variable x= Variable (1); |
---|
| 859 | if (degree (F, x) <= 1) |
---|
| 860 | return 0; |
---|
| 861 | CanonicalForm f= swapvar (F, F.mvar(), x); |
---|
| 862 | int sizef= 0; |
---|
| 863 | for (CFIterator i= f; i.hasTerms(); i++, sizef++) |
---|
| 864 | { |
---|
| 865 | if (i.exp() == 1) |
---|
| 866 | return 0; |
---|
| 867 | } |
---|
| 868 | int * expf= new int [sizef]; |
---|
| 869 | int j= 0; |
---|
| 870 | for (CFIterator i= f; i.hasTerms(); i++, j++) |
---|
| 871 | { |
---|
| 872 | expf [j]= i.exp(); |
---|
| 873 | } |
---|
| 874 | |
---|
| 875 | int indf= sizef - 1; |
---|
| 876 | if (expf[indf] == 0) |
---|
| 877 | indf--; |
---|
| 878 | |
---|
| 879 | if ((d%expf [indf] != 0 && expf[indf]%d != 0) || (expf[indf] == 1)) |
---|
| 880 | { |
---|
| 881 | delete [] expf; |
---|
| 882 | return 0; |
---|
| 883 | } |
---|
| 884 | |
---|
| 885 | int result; |
---|
| 886 | if (d%expf [indf] == 0) |
---|
| 887 | result= expf[indf]; |
---|
| 888 | else |
---|
| 889 | result= d; |
---|
| 890 | for (int i= indf - 1; i >= 0; i--) |
---|
| 891 | { |
---|
| 892 | if (expf [i]%result != 0) |
---|
| 893 | { |
---|
| 894 | delete [] expf; |
---|
| 895 | return 0; |
---|
| 896 | } |
---|
| 897 | } |
---|
| 898 | |
---|
| 899 | delete [] expf; |
---|
| 900 | return result; |
---|
| 901 | } |
---|
| 902 | |
---|
| 903 | int substituteCheck (const CFList& L) |
---|
| 904 | { |
---|
| 905 | ASSERT (L.length() > 1, "expected a list of at least two elements"); |
---|
| 906 | if (L.length() < 2) |
---|
| 907 | return 0; |
---|
| 908 | CFListIterator i= L; |
---|
| 909 | i++; |
---|
| 910 | int result= substituteCheck (L.getFirst(), i.getItem()); |
---|
| 911 | if (result <= 1) |
---|
| 912 | return result; |
---|
| 913 | i++; |
---|
| 914 | for (;i.hasItem(); i++) |
---|
| 915 | { |
---|
| 916 | result= recSubstituteCheck (i.getItem(), result); |
---|
| 917 | if (result <= 1) |
---|
| 918 | return result; |
---|
| 919 | } |
---|
| 920 | return result; |
---|
| 921 | } |
---|
| 922 | |
---|
| 923 | void |
---|
| 924 | subst (const CanonicalForm& F, CanonicalForm& A, const int d, const Variable& x) |
---|
| 925 | { |
---|
| 926 | if (d <= 1) |
---|
| 927 | { |
---|
| 928 | A= F; |
---|
| 929 | return; |
---|
| 930 | } |
---|
| 931 | if (degree (F, x) <= 0) |
---|
| 932 | { |
---|
| 933 | A= F; |
---|
| 934 | return; |
---|
| 935 | } |
---|
| 936 | CanonicalForm C= 0; |
---|
| 937 | CanonicalForm f= swapvar (F, x, F.mvar()); |
---|
| 938 | for (CFIterator i= f; i.hasTerms(); i++) |
---|
| 939 | C += i.coeff()*power (f.mvar(), i.exp()/ d); |
---|
| 940 | A= swapvar (C, x, F.mvar()); |
---|
| 941 | } |
---|
| 942 | |
---|
| 943 | CanonicalForm |
---|
| 944 | reverseSubst (const CanonicalForm& F, const int d, const Variable& x) |
---|
| 945 | { |
---|
| 946 | if (d <= 1) |
---|
| 947 | return F; |
---|
| 948 | if (degree (F, x) <= 0) |
---|
| 949 | return F; |
---|
| 950 | CanonicalForm f= swapvar (F, x, F.mvar()); |
---|
| 951 | CanonicalForm result= 0; |
---|
| 952 | for (CFIterator i= f; i.hasTerms(); i++) |
---|
| 953 | result += i.coeff()*power (f.mvar(), d*i.exp()); |
---|
| 954 | return swapvar (result, x, F.mvar()); |
---|
| 955 | } |
---|
| 956 | |
---|
| 957 | void |
---|
| 958 | reverseSubst (CFList& L, const int d, const Variable& x) |
---|
| 959 | { |
---|
| 960 | for (CFListIterator i= L; i.hasItem(); i++) |
---|
| 961 | i.getItem()= reverseSubst (i.getItem(), d, x); |
---|
| 962 | } |
---|
| 963 | |
---|