001/* 002 * $Id$ 003 */ 004 005package edu.jas.ufd; 006 007 008import org.apache.logging.log4j.LogManager; 009import org.apache.logging.log4j.Logger; 010 011import edu.jas.poly.GenPolynomial; 012import edu.jas.poly.PolyUtil; 013import edu.jas.structure.GcdRingElem; 014import edu.jas.structure.RingFactory; 015 016 017/** 018 * Greatest common divisor algorithms with monic polynomial remainder sequence. 019 * If C is a field, then the monic PRS (on coefficients) is computed otherwise 020 * no simplifications in the reduction are made. 021 * @author Heinz Kredel 022 */ 023 024public class GreatestCommonDivisorSimple<C extends GcdRingElem<C>> extends GreatestCommonDivisorAbstract<C> { 025 026 027 private static final Logger logger = LogManager.getLogger(GreatestCommonDivisorSimple.class); 028 029 030 private static final boolean debug = logger.isDebugEnabled(); 031 032 033 /** 034 * Univariate GenPolynomial greatest common divisor. Uses pseudoRemainder for 035 * remainder. 036 * @param P univariate GenPolynomial. 037 * @param S univariate GenPolynomial. 038 * @return gcd(P,S). 039 */ 040 @Override 041 public GenPolynomial<C> baseGcd(GenPolynomial<C> P, GenPolynomial<C> S) { 042 if (S == null || S.isZERO()) { 043 return P; 044 } 045 if (P == null || P.isZERO()) { 046 return S; 047 } 048 if (P.ring.nvar > 1) { 049 throw new IllegalArgumentException(this.getClass().getName() + " no univariate polynomial"); 050 } 051 boolean field = P.ring.coFac.isField(); 052 long e = P.degree(0); 053 long f = S.degree(0); 054 GenPolynomial<C> q; 055 GenPolynomial<C> r; 056 if (f > e) { 057 r = P; 058 q = S; 059 long g = f; 060 f = e; 061 e = g; 062 } else { 063 q = P; 064 r = S; 065 } 066 if (debug) { 067 logger.debug("degrees: e = {}, f = {}", e, f); 068 } 069 C c; 070 if (field) { 071 r = r.monic(); 072 q = q.monic(); 073 c = P.ring.getONECoefficient(); 074 } else { 075 r = r.abs(); 076 q = q.abs(); 077 C a = baseContent(r); 078 C b = baseContent(q); 079 c = gcd(a, b); // indirection 080 r = divide(r, a); // indirection 081 q = divide(q, b); // indirection 082 } 083 if (r.isONE()) { 084 return r.multiply(c); 085 } 086 if (q.isONE()) { 087 return q.multiply(c); 088 } 089 GenPolynomial<C> x; 090 //System.out.println("q = " + q); 091 //System.out.println("r = " + r); 092 while (!r.isZERO()) { 093 x = PolyUtil.<C> baseSparsePseudoRemainder(q, r); 094 q = r; 095 if (field) { 096 r = x.monic(); 097 } else { 098 r = x; 099 } 100 //System.out.println("q = " + q); 101 //System.out.println("r = " + r); 102 } 103 q = basePrimitivePart(q); 104 return (q.multiply(c)).abs(); 105 } 106 107 108 /** 109 * Univariate GenPolynomial recursive greatest common divisor. Uses 110 * pseudoRemainder for remainder. 111 * @param P univariate recursive GenPolynomial. 112 * @param S univariate recursive GenPolynomial. 113 * @return gcd(P,S). 114 */ 115 @Override 116 public GenPolynomial<GenPolynomial<C>> recursiveUnivariateGcd(GenPolynomial<GenPolynomial<C>> P, 117 GenPolynomial<GenPolynomial<C>> S) { 118 if (S == null || S.isZERO()) { 119 return P; 120 } 121 if (P == null || P.isZERO()) { 122 return S; 123 } 124 if (P.ring.nvar > 1) { 125 throw new IllegalArgumentException(this.getClass().getName() + " no univariate polynomial"); 126 } 127 boolean field = P.leadingBaseCoefficient().ring.coFac.isField(); 128 long e = P.degree(0); 129 long f = S.degree(0); 130 GenPolynomial<GenPolynomial<C>> q; 131 GenPolynomial<GenPolynomial<C>> r; 132 if (f > e) { 133 r = P; 134 q = S; 135 long g = f; 136 f = e; 137 e = g; 138 } else { 139 q = P; 140 r = S; 141 } 142 if (debug) { 143 logger.debug("degrees: e = {}, f = {}", e, f); 144 } 145 if (field) { 146 r = PolyUtil.<C> monic(r); 147 q = PolyUtil.<C> monic(q); 148 } else { 149 r = r.abs(); 150 q = q.abs(); 151 } 152 GenPolynomial<C> a = recursiveContent(r); 153 GenPolynomial<C> b = recursiveContent(q); 154 155 GenPolynomial<C> c = gcd(a, b); // go to recursion 156 //System.out.println("rgcd c = " + c); 157 r = PolyUtil.<C> recursiveDivide(r, a); 158 q = PolyUtil.<C> recursiveDivide(q, b); 159 if (r.isONE()) { 160 return r.multiply(c); 161 } 162 if (q.isONE()) { 163 return q.multiply(c); 164 } 165 GenPolynomial<GenPolynomial<C>> x; 166 while (!r.isZERO()) { 167 x = PolyUtil.<C> recursiveSparsePseudoRemainder(q, r); 168 if (logger.isDebugEnabled()) { 169 logger.info("recursiveSparsePseudoRemainder.bits = {}", x.bitLength()); 170 } 171 q = r; 172 if (field) { 173 r = PolyUtil.<C> monic(x); 174 } else { 175 r = x; 176 } 177 } 178 q = recursivePrimitivePart(q); 179 q = q.abs().multiply(c); 180 return q; 181 } 182 183 184 /** 185 * Univariate GenPolynomial resultant. 186 * @param P univariate GenPolynomial. 187 * @param S univariate GenPolynomial. 188 * @return res(P,S). 189 */ 190 @Override 191 public GenPolynomial<C> baseResultant(GenPolynomial<C> P, GenPolynomial<C> S) { 192 if (S == null || S.isZERO()) { 193 return S; 194 } 195 if (P == null || P.isZERO()) { 196 return P; 197 } 198 if (P.ring.nvar > 1 || P.ring.nvar == 0) { 199 throw new IllegalArgumentException("no univariate polynomial"); 200 } 201 long e = P.degree(0); 202 long f = S.degree(0); 203 if (f == 0 && e == 0) { 204 return P.ring.getONE(); 205 } 206 if (e == 0) { 207 return P.power(f); 208 } 209 if (f == 0) { 210 return S.power(e); 211 } 212 GenPolynomial<C> q; 213 GenPolynomial<C> r; 214 int s = 0; // sign is +, 1 for sign is - 215 if (e < f) { 216 r = P; 217 q = S; 218 long t = e; 219 e = f; 220 f = t; 221 if ((e % 2 != 0) && (f % 2 != 0)) { // odd(e) && odd(f) 222 s = 1; 223 } 224 } else { 225 q = P; 226 r = S; 227 } 228 RingFactory<C> cofac = P.ring.coFac; 229 boolean field = cofac.isField(); 230 C c = cofac.getONE(); 231 GenPolynomial<C> x; 232 long g; 233 do { 234 if (field) { 235 x = q.remainder(r); 236 } else { 237 x = PolyUtil.<C> baseSparsePseudoRemainder(q, r); 238 //System.out.println("x_s = " + x + ", lbcf(r) = " + r.leadingBaseCoefficient()); 239 } 240 if (x.isZERO()) { 241 return x; 242 } 243 //System.out.println("x = " + x); 244 e = q.degree(0); 245 f = r.degree(0); 246 if ((e % 2 != 0) && (f % 2 != 0)) { // odd(e) && odd(f) 247 s = 1 - s; 248 } 249 g = x.degree(0); 250 C c2 = r.leadingBaseCoefficient(); 251 for (int i = 0; i < (e - g); i++) { 252 c = c.multiply(c2); 253 } 254 q = r; 255 r = x; 256 } while (g != 0); 257 C c2 = r.leadingBaseCoefficient(); 258 for (int i = 0; i < f; i++) { 259 c = c.multiply(c2); 260 } 261 if (s == 1) { 262 c = c.negate(); 263 } 264 x = P.ring.getONE().multiply(c); 265 return x; 266 } 267 268 269 /** 270 * Univariate GenPolynomial recursive resultant. 271 * @param P univariate recursive GenPolynomial. 272 * @param S univariate recursive GenPolynomial. 273 * @return res(P,S). 274 */ 275 @Override 276 public GenPolynomial<GenPolynomial<C>> recursiveUnivariateResultant(GenPolynomial<GenPolynomial<C>> P, 277 GenPolynomial<GenPolynomial<C>> S) { 278 if (S == null || S.isZERO()) { 279 return S; 280 } 281 if (P == null || P.isZERO()) { 282 return P; 283 } 284 if (P.ring.nvar > 1 || P.ring.nvar == 0) { 285 throw new IllegalArgumentException("no recursive univariate polynomial"); 286 } 287 long e = P.degree(0); 288 long f = S.degree(0); 289 if (f == 0 && e == 0) { 290 // if coeffs are multivariate (and non constant) 291 // otherwise it would be 1 292 GenPolynomial<C> t = resultant(P.leadingBaseCoefficient(), S.leadingBaseCoefficient()); 293 return P.ring.getONE().multiply(t); 294 } 295 if (e == 0) { 296 return P.power(f); 297 } 298 if (f == 0) { 299 return S.power(e); 300 } 301 GenPolynomial<GenPolynomial<C>> q; 302 GenPolynomial<GenPolynomial<C>> r; 303 int s = 0; // sign is +, 1 for sign is - 304 if (f > e) { 305 r = P; 306 q = S; 307 long g = f; 308 f = e; 309 e = g; 310 if ((e % 2 != 0) && (f % 2 != 0)) { // odd(e) && odd(f) 311 s = 1; 312 } 313 } else { 314 q = P; 315 r = S; 316 } 317 GenPolynomial<GenPolynomial<C>> x; 318 RingFactory<GenPolynomial<C>> cofac = P.ring.coFac; 319 GenPolynomial<C> c = cofac.getONE(); 320 long g; 321 do { 322 x = PolyUtil.<C> recursiveSparsePseudoRemainder(q, r); 323 //x = PolyUtil.<C>recursiveDensePseudoRemainder(q,r); 324 if (x.isZERO()) { 325 return x; 326 } 327 //no: x = recursivePrimitivePart(x); 328 //System.out.println("x = " + x); 329 e = q.degree(0); 330 f = r.degree(0); 331 if ((e % 2 != 0) && (f % 2 != 0)) { // odd(e) && odd(f) 332 s = 1 - s; 333 } 334 g = x.degree(0); 335 GenPolynomial<C> c2 = r.leadingBaseCoefficient(); 336 for (int i = 0; i < (e - g); i++) { 337 c = c.multiply(c2); 338 } 339 q = r; 340 r = x; 341 } while (g != 0); 342 GenPolynomial<C> c2 = r.leadingBaseCoefficient(); 343 for (int i = 0; i < f; i++) { 344 c = c.multiply(c2); 345 } 346 if (s == 1) { 347 c = c.negate(); 348 } 349 x = P.ring.getONE().multiply(c); 350 return x; 351 } 352 353}