Coverage for /usr/lib/python3/dist-packages/sympy/polys/domains/gaussiandomains.py: 44%

1"""Domains of Gaussian type."""

3from sympy.core.numbers import I

4from sympy.polys.polyerrors import CoercionFailed

5from sympy.polys.domains.integerring import ZZ

6from sympy.polys.domains.rationalfield import QQ

7from sympy.polys.domains.algebraicfield import AlgebraicField

8from sympy.polys.domains.domain import Domain

9from sympy.polys.domains.domainelement import DomainElement

10from sympy.polys.domains.field import Field

11from sympy.polys.domains.ring import Ring

14class GaussianElement(DomainElement):

15 """Base class for elements of Gaussian type domains."""

16 base: Domain

17 _parent: Domain

19 __slots__ = ('x', 'y')

21 def __new__(cls, x, y=0):

22 conv = cls.base.convert

23 return cls.new(conv(x), conv(y))

25 @classmethod

26 def new(cls, x, y):

27 """Create a new GaussianElement of the same domain."""

28 obj = super().__new__(cls)

29 obj.x = x

30 obj.y = y

31 return obj

33 def parent(self):

34 """The domain that this is an element of (ZZ_I or QQ_I)"""

35 return self._parent

37 def __hash__(self):

38 return hash((self.x, self.y))

40 def __eq__(self, other):

41 if isinstance(other, self.__class__):

42 return self.x == other.x and self.y == other.y

43 else:

44 return NotImplemented

46 def __lt__(self, other):

47 if not isinstance(other, GaussianElement):

48 return NotImplemented

49 return [self.y, self.x] < [other.y, other.x]

51 def __pos__(self):

52 return self

54 def __neg__(self):

55 return self.new(-self.x, -self.y)

57 def __repr__(self):

58 return "%s(%s, %s)" % (self._parent.rep, self.x, self.y)

60 def __str__(self):

61 return str(self._parent.to_sympy(self))

63 @classmethod

64 def _get_xy(cls, other):

65 if not isinstance(other, cls):

66 try:

67 other = cls._parent.convert(other)

68 except CoercionFailed:

69 return None, None

70 return other.x, other.y

72 def __add__(self, other):

73 x, y = self._get_xy(other)

74 if x is not None:

75 return self.new(self.x + x, self.y + y)

76 else:

77 return NotImplemented

79 __radd__ = __add__

81 def __sub__(self, other):

82 x, y = self._get_xy(other)

83 if x is not None:

84 return self.new(self.x - x, self.y - y)

85 else:

86 return NotImplemented

88 def __rsub__(self, other):

89 x, y = self._get_xy(other)

90 if x is not None:

91 return self.new(x - self.x, y - self.y)

92 else:

93 return NotImplemented

95 def __mul__(self, other):

96 x, y = self._get_xy(other)

97 if x is not None:

98 return self.new(self.x*x - self.y*y, self.x*y + self.y*x)

99 else:

100 return NotImplemented

101

102 __rmul__ = __mul__

103

104 def __pow__(self, exp):

105 if exp == 0:

106 return self.new(1, 0)

107 if exp < 0:

108 self, exp = 1/self, -exp

109 if exp == 1:

110 return self

111 pow2 = self

112 prod = self if exp % 2 else self._parent.one

113 exp //= 2

114 while exp:

115 pow2 *= pow2

116 if exp % 2:

117 prod *= pow2

118 exp //= 2

119 return prod

120

121 def __bool__(self):

122 return bool(self.x) or bool(self.y)

123

124 def quadrant(self):

125 """Return quadrant index 0-3.

126

127 0 is included in quadrant 0.

128 """

129 if self.y > 0:

130 return 0 if self.x > 0 else 1

131 elif self.y < 0:

132 return 2 if self.x < 0 else 3

133 else:

134 return 0 if self.x >= 0 else 2

135

136 def __rdivmod__(self, other):

137 try:

138 other = self._parent.convert(other)

139 except CoercionFailed:

140 return NotImplemented

141 else:

142 return other.__divmod__(self)

143

144 def __rtruediv__(self, other):

145 try:

146 other = QQ_I.convert(other)

147 except CoercionFailed:

148 return NotImplemented

149 else:

150 return other.__truediv__(self)

151

152 def __floordiv__(self, other):

153 qr = self.__divmod__(other)

154 return qr if qr is NotImplemented else qr[0]

155

156 def __rfloordiv__(self, other):

157 qr = self.__rdivmod__(other)

158 return qr if qr is NotImplemented else qr[0]

159

160 def __mod__(self, other):

161 qr = self.__divmod__(other)

162 return qr if qr is NotImplemented else qr[1]

163

164 def __rmod__(self, other):

165 qr = self.__rdivmod__(other)

166 return qr if qr is NotImplemented else qr[1]

167

168

169class GaussianInteger(GaussianElement):

170 """Gaussian integer: domain element for :ref:`ZZ_I`

171

172 >>> from sympy import ZZ_I

173 >>> z = ZZ_I(2, 3)

174 >>> z

175 (2 + 3*I)

176 >>> type(z)

177 <class 'sympy.polys.domains.gaussiandomains.GaussianInteger'>

178 """

179 base = ZZ

180

181 def __truediv__(self, other):

182 """Return a Gaussian rational."""

183 return QQ_I.convert(self)/other

184

185 def __divmod__(self, other):

186 if not other:

187 raise ZeroDivisionError('divmod({}, 0)'.format(self))

188 x, y = self._get_xy(other)

189 if x is None:

190 return NotImplemented

191

192 # multiply self and other by x - I*y

193 # self/other == (a + I*b)/c

194 a, b = self.x*x + self.y*y, -self.x*y + self.y*x

195 c = x*x + y*y

196

197 # find integers qx and qy such that

198 # |a - qx*c| <= c/2 and |b - qy*c| <= c/2

199 qx = (2*a + c) // (2*c) # -c <= 2*a - qx*2*c < c

200 qy = (2*b + c) // (2*c)

201

202 q = GaussianInteger(qx, qy)

203 # |self/other - q| < 1 since

204 # |a/c - qx|**2 + |b/c - qy|**2 <= 1/4 + 1/4 < 1

205

206 return q, self - q*other # |r| < |other|

207

208

209class GaussianRational(GaussianElement):

210 """Gaussian rational: domain element for :ref:`QQ_I`

211

212 >>> from sympy import QQ_I, QQ

213 >>> z = QQ_I(QQ(2, 3), QQ(4, 5))

214 >>> z

215 (2/3 + 4/5*I)

216 >>> type(z)

217 <class 'sympy.polys.domains.gaussiandomains.GaussianRational'>

218 """

219 base = QQ

220

221 def __truediv__(self, other):

222 """Return a Gaussian rational."""

223 if not other:

224 raise ZeroDivisionError('{} / 0'.format(self))

225 x, y = self._get_xy(other)

226 if x is None:

227 return NotImplemented

228 c = x*x + y*y

229

230 return GaussianRational((self.x*x + self.y*y)/c,

231 (-self.x*y + self.y*x)/c)

232

233 def __divmod__(self, other):

234 try:

235 other = self._parent.convert(other)

236 except CoercionFailed:

237 return NotImplemented

238 if not other:

239 raise ZeroDivisionError('{} % 0'.format(self))

240 else:

241 return self/other, QQ_I.zero

242

243

244class GaussianDomain():

245 """Base class for Gaussian domains."""

246 dom = None # type: Domain

247

248 is_Numerical = True

249 is_Exact = True

250

251 has_assoc_Ring = True

252 has_assoc_Field = True

253

254 def to_sympy(self, a):

255 """Convert ``a`` to a SymPy object. """

256 conv = self.dom.to_sympy

257 return conv(a.x) + I*conv(a.y)

258

259 def from_sympy(self, a):

260 """Convert a SymPy object to ``self.dtype``."""

261 r, b = a.as_coeff_Add()

262 x = self.dom.from_sympy(r) # may raise CoercionFailed

263 if not b:

264 return self.new(x, 0)

265 r, b = b.as_coeff_Mul()

266 y = self.dom.from_sympy(r)

267 if b is I:

268 return self.new(x, y)

269 else:

270 raise CoercionFailed("{} is not Gaussian".format(a))

271

272 def inject(self, *gens):

273 """Inject generators into this domain. """

274 return self.poly_ring(*gens)

275

276 def canonical_unit(self, d):

277 unit = self.units[-d.quadrant()] # - for inverse power

278 return unit

279

280 def is_negative(self, element):

281 """Returns ``False`` for any ``GaussianElement``. """

282 return False

283

284 def is_positive(self, element):

285 """Returns ``False`` for any ``GaussianElement``. """

286 return False

287

288 def is_nonnegative(self, element):

289 """Returns ``False`` for any ``GaussianElement``. """

290 return False

291

292 def is_nonpositive(self, element):

293 """Returns ``False`` for any ``GaussianElement``. """

294 return False

295

296 def from_ZZ_gmpy(K1, a, K0):

297 """Convert a GMPY mpz to ``self.dtype``."""

298 return K1(a)

299

300 def from_ZZ(K1, a, K0):

301 """Convert a ZZ_python element to ``self.dtype``."""

302 return K1(a)

303

304 def from_ZZ_python(K1, a, K0):

305 """Convert a ZZ_python element to ``self.dtype``."""

306 return K1(a)

307

308 def from_QQ(K1, a, K0):

309 """Convert a GMPY mpq to ``self.dtype``."""

310 return K1(a)

311

312 def from_QQ_gmpy(K1, a, K0):

313 """Convert a GMPY mpq to ``self.dtype``."""

314 return K1(a)

315

316 def from_QQ_python(K1, a, K0):

317 """Convert a QQ_python element to ``self.dtype``."""

318 return K1(a)

319

320 def from_AlgebraicField(K1, a, K0):

321 """Convert an element from ZZ<I> or QQ<I> to ``self.dtype``."""

322 if K0.ext.args[0] == I:

323 return K1.from_sympy(K0.to_sympy(a))

324

325

326class GaussianIntegerRing(GaussianDomain, Ring):

327 r"""Ring of Gaussian integers ``ZZ_I``

328

329 The :ref:`ZZ_I` domain represents the `Gaussian integers`_ `\mathbb{Z}[i]`

330 as a :py:class:`~.Domain` in the domain system (see

331 :ref:`polys-domainsintro`).

332

333 By default a :py:class:`~.Poly` created from an expression with

334 coefficients that are combinations of integers and ``I`` (`\sqrt{-1}`)

335 will have the domain :ref:`ZZ_I`.

336

337 >>> from sympy import Poly, Symbol, I

338 >>> x = Symbol('x')

339 >>> p = Poly(x**2 + I)

340 >>> p

341 Poly(x**2 + I, x, domain='ZZ_I')

342 >>> p.domain

343 ZZ_I

344

345 The :ref:`ZZ_I` domain can be used to factorise polynomials that are

346 reducible over the Gaussian integers.

347

348 >>> from sympy import factor

349 >>> factor(x**2 + 1)

350 x**2 + 1

351 >>> factor(x**2 + 1, domain='ZZ_I')

352 (x - I)*(x + I)

353

354 The corresponding `field of fractions`_ is the domain of the Gaussian

355 rationals :ref:`QQ_I`. Conversely :ref:`ZZ_I` is the `ring of integers`_

356 of :ref:`QQ_I`.

357

358 >>> from sympy import ZZ_I, QQ_I

359 >>> ZZ_I.get_field()

360 QQ_I

361 >>> QQ_I.get_ring()

362 ZZ_I

363

364 When using the domain directly :ref:`ZZ_I` can be used as a constructor.

365

366 >>> ZZ_I(3, 4)

367 (3 + 4*I)

368 >>> ZZ_I(5)

369 (5 + 0*I)

370

371 The domain elements of :ref:`ZZ_I` are instances of

372 :py:class:`~.GaussianInteger` which support the rings operations

373 ``+,-,*,**``.

374

375 >>> z1 = ZZ_I(5, 1)

376 >>> z2 = ZZ_I(2, 3)

377 >>> z1

378 (5 + 1*I)

379 >>> z2

380 (2 + 3*I)

381 >>> z1 + z2

382 (7 + 4*I)

383 >>> z1 * z2

384 (7 + 17*I)

385 >>> z1 ** 2

386 (24 + 10*I)

387

388 Both floor (``//``) and modulo (``%``) division work with

389 :py:class:`~.GaussianInteger` (see the :py:meth:`~.Domain.div` method).

390

391 >>> z3, z4 = ZZ_I(5), ZZ_I(1, 3)

392 >>> z3 // z4 # floor division

393 (1 + -1*I)

394 >>> z3 % z4 # modulo division (remainder)

395 (1 + -2*I)

396 >>> (z3//z4)*z4 + z3%z4 == z3

397 True

398

399 True division (``/``) in :ref:`ZZ_I` gives an element of :ref:`QQ_I`. The

400 :py:meth:`~.Domain.exquo` method can be used to divide in :ref:`ZZ_I` when

401 exact division is possible.

402

403 >>> z1 / z2

404 (1 + -1*I)

405 >>> ZZ_I.exquo(z1, z2)

406 (1 + -1*I)

407 >>> z3 / z4

408 (1/2 + -3/2*I)

409 >>> ZZ_I.exquo(z3, z4)

410 Traceback (most recent call last):

411 ...

412 ExactQuotientFailed: (1 + 3*I) does not divide (5 + 0*I) in ZZ_I

413

414 The :py:meth:`~.Domain.gcd` method can be used to compute the `gcd`_ of any

415 two elements.

416

417 >>> ZZ_I.gcd(ZZ_I(10), ZZ_I(2))

418 (2 + 0*I)

419 >>> ZZ_I.gcd(ZZ_I(5), ZZ_I(2, 1))

420 (2 + 1*I)

421

422 .. _Gaussian integers: https://en.wikipedia.org/wiki/Gaussian_integer

423 .. _gcd: https://en.wikipedia.org/wiki/Greatest_common_divisor

424

425 """

426 dom = ZZ

427 dtype = GaussianInteger

428 zero = dtype(ZZ(0), ZZ(0))

429 one = dtype(ZZ(1), ZZ(0))

430 imag_unit = dtype(ZZ(0), ZZ(1))

431 units = (one, imag_unit, -one, -imag_unit) # powers of i

432

433 rep = 'ZZ_I'

434

435 is_GaussianRing = True

436 is_ZZ_I = True

437

438 def __init__(self): # override Domain.__init__

439 """For constructing ZZ_I."""

440

441 def get_ring(self):

442 """Returns a ring associated with ``self``. """

443 return self

444

445 def get_field(self):

446 """Returns a field associated with ``self``. """

447 return QQ_I

448

449 def normalize(self, d, *args):

450 """Return first quadrant element associated with ``d``.

451

452 Also multiply the other arguments by the same power of i.

453 """

454 unit = self.canonical_unit(d)

455 d *= unit

456 args = tuple(a*unit for a in args)

457 return (d,) + args if args else d

458

459 def gcd(self, a, b):

460 """Greatest common divisor of a and b over ZZ_I."""

461 while b:

462 a, b = b, a % b

463 return self.normalize(a)

464

465 def lcm(self, a, b):

466 """Least common multiple of a and b over ZZ_I."""

467 return (a * b) // self.gcd(a, b)

468

469 def from_GaussianIntegerRing(K1, a, K0):

470 """Convert a ZZ_I element to ZZ_I."""

471 return a

472

473 def from_GaussianRationalField(K1, a, K0):

474 """Convert a QQ_I element to ZZ_I."""

475 return K1.new(ZZ.convert(a.x), ZZ.convert(a.y))

476

477ZZ_I = GaussianInteger._parent = GaussianIntegerRing()

478

479

480class GaussianRationalField(GaussianDomain, Field):

481 r"""Field of Gaussian rationals ``QQ_I``

482

483 The :ref:`QQ_I` domain represents the `Gaussian rationals`_ `\mathbb{Q}(i)`

484 as a :py:class:`~.Domain` in the domain system (see

485 :ref:`polys-domainsintro`).

486

487 By default a :py:class:`~.Poly` created from an expression with

488 coefficients that are combinations of rationals and ``I`` (`\sqrt{-1}`)

489 will have the domain :ref:`QQ_I`.

490

491 >>> from sympy import Poly, Symbol, I

492 >>> x = Symbol('x')

493 >>> p = Poly(x**2 + I/2)

494 >>> p

495 Poly(x**2 + I/2, x, domain='QQ_I')

496 >>> p.domain

497 QQ_I

498

499 The polys option ``gaussian=True`` can be used to specify that the domain

500 should be :ref:`QQ_I` even if the coefficients do not contain ``I`` or are

501 all integers.

502

503 >>> Poly(x**2)

504 Poly(x**2, x, domain='ZZ')

505 >>> Poly(x**2 + I)

506 Poly(x**2 + I, x, domain='ZZ_I')

507 >>> Poly(x**2/2)

508 Poly(1/2*x**2, x, domain='QQ')

509 >>> Poly(x**2, gaussian=True)

510 Poly(x**2, x, domain='QQ_I')

511 >>> Poly(x**2 + I, gaussian=True)

512 Poly(x**2 + I, x, domain='QQ_I')

513 >>> Poly(x**2/2, gaussian=True)

514 Poly(1/2*x**2, x, domain='QQ_I')

515

516 The :ref:`QQ_I` domain can be used to factorise polynomials that are

517 reducible over the Gaussian rationals.

518

519 >>> from sympy import factor, QQ_I

520 >>> factor(x**2/4 + 1)

521 (x**2 + 4)/4

522 >>> factor(x**2/4 + 1, domain='QQ_I')

523 (x - 2*I)*(x + 2*I)/4

524 >>> factor(x**2/4 + 1, domain=QQ_I)

525 (x - 2*I)*(x + 2*I)/4

526

527 It is also possible to specify the :ref:`QQ_I` domain explicitly with

528 polys functions like :py:func:`~.apart`.

529

530 >>> from sympy import apart

531 >>> apart(1/(1 + x**2))

532 1/(x**2 + 1)

533 >>> apart(1/(1 + x**2), domain=QQ_I)

534 I/(2*(x + I)) - I/(2*(x - I))

535

536 The corresponding `ring of integers`_ is the domain of the Gaussian

537 integers :ref:`ZZ_I`. Conversely :ref:`QQ_I` is the `field of fractions`_

538 of :ref:`ZZ_I`.

539

540 >>> from sympy import ZZ_I, QQ_I, QQ

541 >>> ZZ_I.get_field()

542 QQ_I

543 >>> QQ_I.get_ring()

544 ZZ_I

545

546 When using the domain directly :ref:`QQ_I` can be used as a constructor.

547

548 >>> QQ_I(3, 4)

549 (3 + 4*I)

550 >>> QQ_I(5)

551 (5 + 0*I)

552 >>> QQ_I(QQ(2, 3), QQ(4, 5))

553 (2/3 + 4/5*I)

554

555 The domain elements of :ref:`QQ_I` are instances of

556 :py:class:`~.GaussianRational` which support the field operations

557 ``+,-,*,**,/``.

558

559 >>> z1 = QQ_I(5, 1)

560 >>> z2 = QQ_I(2, QQ(1, 2))

561 >>> z1

562 (5 + 1*I)

563 >>> z2

564 (2 + 1/2*I)

565 >>> z1 + z2

566 (7 + 3/2*I)

567 >>> z1 * z2

568 (19/2 + 9/2*I)

569 >>> z2 ** 2

570 (15/4 + 2*I)

571

572 True division (``/``) in :ref:`QQ_I` gives an element of :ref:`QQ_I` and

573 is always exact.

574

575 >>> z1 / z2

576 (42/17 + -2/17*I)

577 >>> QQ_I.exquo(z1, z2)

578 (42/17 + -2/17*I)

579 >>> z1 == (z1/z2)*z2

580 True

581

582 Both floor (``//``) and modulo (``%``) division can be used with

583 :py:class:`~.GaussianRational` (see :py:meth:`~.Domain.div`)

584 but division is always exact so there is no remainder.

585

586 >>> z1 // z2

587 (42/17 + -2/17*I)

588 >>> z1 % z2

589 (0 + 0*I)

590 >>> QQ_I.div(z1, z2)

591 ((42/17 + -2/17*I), (0 + 0*I))

592 >>> (z1//z2)*z2 + z1%z2 == z1

593 True

594

595 .. _Gaussian rationals: https://en.wikipedia.org/wiki/Gaussian_rational

596 """

597 dom = QQ

598 dtype = GaussianRational

599 zero = dtype(QQ(0), QQ(0))

600 one = dtype(QQ(1), QQ(0))

601 imag_unit = dtype(QQ(0), QQ(1))

602 units = (one, imag_unit, -one, -imag_unit) # powers of i

603

604 rep = 'QQ_I'

605

606 is_GaussianField = True

607 is_QQ_I = True

608

609 def __init__(self): # override Domain.__init__

610 """For constructing QQ_I."""

611

612 def get_ring(self):

613 """Returns a ring associated with ``self``. """

614 return ZZ_I

615

616 def get_field(self):

617 """Returns a field associated with ``self``. """

618 return self

619

620 def as_AlgebraicField(self):

621 """Get equivalent domain as an ``AlgebraicField``. """

622 return AlgebraicField(self.dom, I)

623

624 def numer(self, a):

625 """Get the numerator of ``a``."""

626 ZZ_I = self.get_ring()

627 return ZZ_I.convert(a * self.denom(a))

628

629 def denom(self, a):

630 """Get the denominator of ``a``."""

631 ZZ = self.dom.get_ring()

632 QQ = self.dom

633 ZZ_I = self.get_ring()

634 denom_ZZ = ZZ.lcm(QQ.denom(a.x), QQ.denom(a.y))

635 return ZZ_I(denom_ZZ, ZZ.zero)

636

637 def from_GaussianIntegerRing(K1, a, K0):

638 """Convert a ZZ_I element to QQ_I."""

639 return K1.new(a.x, a.y)

640

641 def from_GaussianRationalField(K1, a, K0):

642 """Convert a QQ_I element to QQ_I."""

643 return a

644

645QQ_I = GaussianRational._parent = GaussianRationalField()