Coverage for /usr/lib/python3/dist-packages/sympy/geometry/ellipse.py: 21%

440 statements  

« prev     ^ index     » next       coverage.py v7.9.1, created at 2025-06-14 15:55 +0200

1"""Elliptical geometrical entities. 

2 

3Contains 

4* Ellipse 

5* Circle 

6 

7""" 

8 

9from sympy.core.expr import Expr 

10from sympy.core.relational import Eq 

11from sympy.core import S, pi, sympify 

12from sympy.core.evalf import N 

13from sympy.core.parameters import global_parameters 

14from sympy.core.logic import fuzzy_bool 

15from sympy.core.numbers import Rational, oo 

16from sympy.core.sorting import ordered 

17from sympy.core.symbol import Dummy, uniquely_named_symbol, _symbol 

18from sympy.simplify import simplify, trigsimp 

19from sympy.functions.elementary.miscellaneous import sqrt, Max 

20from sympy.functions.elementary.trigonometric import cos, sin 

21from sympy.functions.special.elliptic_integrals import elliptic_e 

22from .entity import GeometryEntity, GeometrySet 

23from .exceptions import GeometryError 

24from .line import Line, Segment, Ray2D, Segment2D, Line2D, LinearEntity3D 

25from .point import Point, Point2D, Point3D 

26from .util import idiff, find 

27from sympy.polys import DomainError, Poly, PolynomialError 

28from sympy.polys.polyutils import _not_a_coeff, _nsort 

29from sympy.solvers import solve 

30from sympy.solvers.solveset import linear_coeffs 

31from sympy.utilities.misc import filldedent, func_name 

32 

33from mpmath.libmp.libmpf import prec_to_dps 

34 

35import random 

36 

37x, y = [Dummy('ellipse_dummy', real=True) for i in range(2)] 

38 

39 

40class Ellipse(GeometrySet): 

41 """An elliptical GeometryEntity. 

42 

43 Parameters 

44 ========== 

45 

46 center : Point, optional 

47 Default value is Point(0, 0) 

48 hradius : number or SymPy expression, optional 

49 vradius : number or SymPy expression, optional 

50 eccentricity : number or SymPy expression, optional 

51 Two of `hradius`, `vradius` and `eccentricity` must be supplied to 

52 create an Ellipse. The third is derived from the two supplied. 

53 

54 Attributes 

55 ========== 

56 

57 center 

58 hradius 

59 vradius 

60 area 

61 circumference 

62 eccentricity 

63 periapsis 

64 apoapsis 

65 focus_distance 

66 foci 

67 

68 Raises 

69 ====== 

70 

71 GeometryError 

72 When `hradius`, `vradius` and `eccentricity` are incorrectly supplied 

73 as parameters. 

74 TypeError 

75 When `center` is not a Point. 

76 

77 See Also 

78 ======== 

79 

80 Circle 

81 

82 Notes 

83 ----- 

84 Constructed from a center and two radii, the first being the horizontal 

85 radius (along the x-axis) and the second being the vertical radius (along 

86 the y-axis). 

87 

88 When symbolic value for hradius and vradius are used, any calculation that 

89 refers to the foci or the major or minor axis will assume that the ellipse 

90 has its major radius on the x-axis. If this is not true then a manual 

91 rotation is necessary. 

92 

93 Examples 

94 ======== 

95 

96 >>> from sympy import Ellipse, Point, Rational 

97 >>> e1 = Ellipse(Point(0, 0), 5, 1) 

98 >>> e1.hradius, e1.vradius 

99 (5, 1) 

100 >>> e2 = Ellipse(Point(3, 1), hradius=3, eccentricity=Rational(4, 5)) 

101 >>> e2 

102 Ellipse(Point2D(3, 1), 3, 9/5) 

103 

104 """ 

105 

106 def __contains__(self, o): 

107 if isinstance(o, Point): 

108 res = self.equation(x, y).subs({x: o.x, y: o.y}) 

109 return trigsimp(simplify(res)) is S.Zero 

110 elif isinstance(o, Ellipse): 

111 return self == o 

112 return False 

113 

114 def __eq__(self, o): 

115 """Is the other GeometryEntity the same as this ellipse?""" 

116 return isinstance(o, Ellipse) and (self.center == o.center and 

117 self.hradius == o.hradius and 

118 self.vradius == o.vradius) 

119 

120 def __hash__(self): 

121 return super().__hash__() 

122 

123 def __new__( 

124 cls, center=None, hradius=None, vradius=None, eccentricity=None, **kwargs): 

125 

126 hradius = sympify(hradius) 

127 vradius = sympify(vradius) 

128 

129 if center is None: 

130 center = Point(0, 0) 

131 else: 

132 if len(center) != 2: 

133 raise ValueError('The center of "{}" must be a two dimensional point'.format(cls)) 

134 center = Point(center, dim=2) 

135 

136 if len(list(filter(lambda x: x is not None, (hradius, vradius, eccentricity)))) != 2: 

137 raise ValueError(filldedent(''' 

138 Exactly two arguments of "hradius", "vradius", and 

139 "eccentricity" must not be None.''')) 

140 

141 if eccentricity is not None: 

142 eccentricity = sympify(eccentricity) 

143 if eccentricity.is_negative: 

144 raise GeometryError("Eccentricity of ellipse/circle should lie between [0, 1)") 

145 elif hradius is None: 

146 hradius = vradius / sqrt(1 - eccentricity**2) 

147 elif vradius is None: 

148 vradius = hradius * sqrt(1 - eccentricity**2) 

149 

150 if hradius == vradius: 

151 return Circle(center, hradius, **kwargs) 

152 

153 if S.Zero in (hradius, vradius): 

154 return Segment(Point(center[0] - hradius, center[1] - vradius), Point(center[0] + hradius, center[1] + vradius)) 

155 

156 if hradius.is_real is False or vradius.is_real is False: 

157 raise GeometryError("Invalid value encountered when computing hradius / vradius.") 

158 

159 return GeometryEntity.__new__(cls, center, hradius, vradius, **kwargs) 

160 

161 def _svg(self, scale_factor=1., fill_color="#66cc99"): 

162 """Returns SVG ellipse element for the Ellipse. 

163 

164 Parameters 

165 ========== 

166 

167 scale_factor : float 

168 Multiplication factor for the SVG stroke-width. Default is 1. 

169 fill_color : str, optional 

170 Hex string for fill color. Default is "#66cc99". 

171 """ 

172 

173 c = N(self.center) 

174 h, v = N(self.hradius), N(self.vradius) 

175 return ( 

176 '<ellipse fill="{1}" stroke="#555555" ' 

177 'stroke-width="{0}" opacity="0.6" cx="{2}" cy="{3}" rx="{4}" ry="{5}"/>' 

178 ).format(2. * scale_factor, fill_color, c.x, c.y, h, v) 

179 

180 @property 

181 def ambient_dimension(self): 

182 return 2 

183 

184 @property 

185 def apoapsis(self): 

186 """The apoapsis of the ellipse. 

187 

188 The greatest distance between the focus and the contour. 

189 

190 Returns 

191 ======= 

192 

193 apoapsis : number 

194 

195 See Also 

196 ======== 

197 

198 periapsis : Returns shortest distance between foci and contour 

199 

200 Examples 

201 ======== 

202 

203 >>> from sympy import Point, Ellipse 

204 >>> p1 = Point(0, 0) 

205 >>> e1 = Ellipse(p1, 3, 1) 

206 >>> e1.apoapsis 

207 2*sqrt(2) + 3 

208 

209 """ 

210 return self.major * (1 + self.eccentricity) 

211 

212 def arbitrary_point(self, parameter='t'): 

213 """A parameterized point on the ellipse. 

214 

215 Parameters 

216 ========== 

217 

218 parameter : str, optional 

219 Default value is 't'. 

220 

221 Returns 

222 ======= 

223 

224 arbitrary_point : Point 

225 

226 Raises 

227 ====== 

228 

229 ValueError 

230 When `parameter` already appears in the functions. 

231 

232 See Also 

233 ======== 

234 

235 sympy.geometry.point.Point 

236 

237 Examples 

238 ======== 

239 

240 >>> from sympy import Point, Ellipse 

241 >>> e1 = Ellipse(Point(0, 0), 3, 2) 

242 >>> e1.arbitrary_point() 

243 Point2D(3*cos(t), 2*sin(t)) 

244 

245 """ 

246 t = _symbol(parameter, real=True) 

247 if t.name in (f.name for f in self.free_symbols): 

248 raise ValueError(filldedent('Symbol %s already appears in object ' 

249 'and cannot be used as a parameter.' % t.name)) 

250 return Point(self.center.x + self.hradius*cos(t), 

251 self.center.y + self.vradius*sin(t)) 

252 

253 @property 

254 def area(self): 

255 """The area of the ellipse. 

256 

257 Returns 

258 ======= 

259 

260 area : number 

261 

262 Examples 

263 ======== 

264 

265 >>> from sympy import Point, Ellipse 

266 >>> p1 = Point(0, 0) 

267 >>> e1 = Ellipse(p1, 3, 1) 

268 >>> e1.area 

269 3*pi 

270 

271 """ 

272 return simplify(S.Pi * self.hradius * self.vradius) 

273 

274 @property 

275 def bounds(self): 

276 """Return a tuple (xmin, ymin, xmax, ymax) representing the bounding 

277 rectangle for the geometric figure. 

278 

279 """ 

280 

281 h, v = self.hradius, self.vradius 

282 return (self.center.x - h, self.center.y - v, self.center.x + h, self.center.y + v) 

283 

284 @property 

285 def center(self): 

286 """The center of the ellipse. 

287 

288 Returns 

289 ======= 

290 

291 center : number 

292 

293 See Also 

294 ======== 

295 

296 sympy.geometry.point.Point 

297 

298 Examples 

299 ======== 

300 

301 >>> from sympy import Point, Ellipse 

302 >>> p1 = Point(0, 0) 

303 >>> e1 = Ellipse(p1, 3, 1) 

304 >>> e1.center 

305 Point2D(0, 0) 

306 

307 """ 

308 return self.args[0] 

309 

310 @property 

311 def circumference(self): 

312 """The circumference of the ellipse. 

313 

314 Examples 

315 ======== 

316 

317 >>> from sympy import Point, Ellipse 

318 >>> p1 = Point(0, 0) 

319 >>> e1 = Ellipse(p1, 3, 1) 

320 >>> e1.circumference 

321 12*elliptic_e(8/9) 

322 

323 """ 

324 if self.eccentricity == 1: 

325 # degenerate 

326 return 4*self.major 

327 elif self.eccentricity == 0: 

328 # circle 

329 return 2*pi*self.hradius 

330 else: 

331 return 4*self.major*elliptic_e(self.eccentricity**2) 

332 

333 @property 

334 def eccentricity(self): 

335 """The eccentricity of the ellipse. 

336 

337 Returns 

338 ======= 

339 

340 eccentricity : number 

341 

342 Examples 

343 ======== 

344 

345 >>> from sympy import Point, Ellipse, sqrt 

346 >>> p1 = Point(0, 0) 

347 >>> e1 = Ellipse(p1, 3, sqrt(2)) 

348 >>> e1.eccentricity 

349 sqrt(7)/3 

350 

351 """ 

352 return self.focus_distance / self.major 

353 

354 def encloses_point(self, p): 

355 """ 

356 Return True if p is enclosed by (is inside of) self. 

357 

358 Notes 

359 ----- 

360 Being on the border of self is considered False. 

361 

362 Parameters 

363 ========== 

364 

365 p : Point 

366 

367 Returns 

368 ======= 

369 

370 encloses_point : True, False or None 

371 

372 See Also 

373 ======== 

374 

375 sympy.geometry.point.Point 

376 

377 Examples 

378 ======== 

379 

380 >>> from sympy import Ellipse, S 

381 >>> from sympy.abc import t 

382 >>> e = Ellipse((0, 0), 3, 2) 

383 >>> e.encloses_point((0, 0)) 

384 True 

385 >>> e.encloses_point(e.arbitrary_point(t).subs(t, S.Half)) 

386 False 

387 >>> e.encloses_point((4, 0)) 

388 False 

389 

390 """ 

391 p = Point(p, dim=2) 

392 if p in self: 

393 return False 

394 

395 if len(self.foci) == 2: 

396 # if the combined distance from the foci to p (h1 + h2) is less 

397 # than the combined distance from the foci to the minor axis 

398 # (which is the same as the major axis length) then p is inside 

399 # the ellipse 

400 h1, h2 = [f.distance(p) for f in self.foci] 

401 test = 2*self.major - (h1 + h2) 

402 else: 

403 test = self.radius - self.center.distance(p) 

404 

405 return fuzzy_bool(test.is_positive) 

406 

407 def equation(self, x='x', y='y', _slope=None): 

408 """ 

409 Returns the equation of an ellipse aligned with the x and y axes; 

410 when slope is given, the equation returned corresponds to an ellipse 

411 with a major axis having that slope. 

412 

413 Parameters 

414 ========== 

415 

416 x : str, optional 

417 Label for the x-axis. Default value is 'x'. 

418 y : str, optional 

419 Label for the y-axis. Default value is 'y'. 

420 _slope : Expr, optional 

421 The slope of the major axis. Ignored when 'None'. 

422 

423 Returns 

424 ======= 

425 

426 equation : SymPy expression 

427 

428 See Also 

429 ======== 

430 

431 arbitrary_point : Returns parameterized point on ellipse 

432 

433 Examples 

434 ======== 

435 

436 >>> from sympy import Point, Ellipse, pi 

437 >>> from sympy.abc import x, y 

438 >>> e1 = Ellipse(Point(1, 0), 3, 2) 

439 >>> eq1 = e1.equation(x, y); eq1 

440 y**2/4 + (x/3 - 1/3)**2 - 1 

441 >>> eq2 = e1.equation(x, y, _slope=1); eq2 

442 (-x + y + 1)**2/8 + (x + y - 1)**2/18 - 1 

443 

444 A point on e1 satisfies eq1. Let's use one on the x-axis: 

445 

446 >>> p1 = e1.center + Point(e1.major, 0) 

447 >>> assert eq1.subs(x, p1.x).subs(y, p1.y) == 0 

448 

449 When rotated the same as the rotated ellipse, about the center 

450 point of the ellipse, it will satisfy the rotated ellipse's 

451 equation, too: 

452 

453 >>> r1 = p1.rotate(pi/4, e1.center) 

454 >>> assert eq2.subs(x, r1.x).subs(y, r1.y) == 0 

455 

456 References 

457 ========== 

458 

459 .. [1] https://math.stackexchange.com/questions/108270/what-is-the-equation-of-an-ellipse-that-is-not-aligned-with-the-axis 

460 .. [2] https://en.wikipedia.org/wiki/Ellipse#Shifted_ellipse 

461 

462 """ 

463 

464 x = _symbol(x, real=True) 

465 y = _symbol(y, real=True) 

466 

467 dx = x - self.center.x 

468 dy = y - self.center.y 

469 

470 if _slope is not None: 

471 L = (dy - _slope*dx)**2 

472 l = (_slope*dy + dx)**2 

473 h = 1 + _slope**2 

474 b = h*self.major**2 

475 a = h*self.minor**2 

476 return l/b + L/a - 1 

477 

478 else: 

479 t1 = (dx/self.hradius)**2 

480 t2 = (dy/self.vradius)**2 

481 return t1 + t2 - 1 

482 

483 def evolute(self, x='x', y='y'): 

484 """The equation of evolute of the ellipse. 

485 

486 Parameters 

487 ========== 

488 

489 x : str, optional 

490 Label for the x-axis. Default value is 'x'. 

491 y : str, optional 

492 Label for the y-axis. Default value is 'y'. 

493 

494 Returns 

495 ======= 

496 

497 equation : SymPy expression 

498 

499 Examples 

500 ======== 

501 

502 >>> from sympy import Point, Ellipse 

503 >>> e1 = Ellipse(Point(1, 0), 3, 2) 

504 >>> e1.evolute() 

505 2**(2/3)*y**(2/3) + (3*x - 3)**(2/3) - 5**(2/3) 

506 """ 

507 if len(self.args) != 3: 

508 raise NotImplementedError('Evolute of arbitrary Ellipse is not supported.') 

509 x = _symbol(x, real=True) 

510 y = _symbol(y, real=True) 

511 t1 = (self.hradius*(x - self.center.x))**Rational(2, 3) 

512 t2 = (self.vradius*(y - self.center.y))**Rational(2, 3) 

513 return t1 + t2 - (self.hradius**2 - self.vradius**2)**Rational(2, 3) 

514 

515 @property 

516 def foci(self): 

517 """The foci of the ellipse. 

518 

519 Notes 

520 ----- 

521 The foci can only be calculated if the major/minor axes are known. 

522 

523 Raises 

524 ====== 

525 

526 ValueError 

527 When the major and minor axis cannot be determined. 

528 

529 See Also 

530 ======== 

531 

532 sympy.geometry.point.Point 

533 focus_distance : Returns the distance between focus and center 

534 

535 Examples 

536 ======== 

537 

538 >>> from sympy import Point, Ellipse 

539 >>> p1 = Point(0, 0) 

540 >>> e1 = Ellipse(p1, 3, 1) 

541 >>> e1.foci 

542 (Point2D(-2*sqrt(2), 0), Point2D(2*sqrt(2), 0)) 

543 

544 """ 

545 c = self.center 

546 hr, vr = self.hradius, self.vradius 

547 if hr == vr: 

548 return (c, c) 

549 

550 # calculate focus distance manually, since focus_distance calls this 

551 # routine 

552 fd = sqrt(self.major**2 - self.minor**2) 

553 if hr == self.minor: 

554 # foci on the y-axis 

555 return (c + Point(0, -fd), c + Point(0, fd)) 

556 elif hr == self.major: 

557 # foci on the x-axis 

558 return (c + Point(-fd, 0), c + Point(fd, 0)) 

559 

560 @property 

561 def focus_distance(self): 

562 """The focal distance of the ellipse. 

563 

564 The distance between the center and one focus. 

565 

566 Returns 

567 ======= 

568 

569 focus_distance : number 

570 

571 See Also 

572 ======== 

573 

574 foci 

575 

576 Examples 

577 ======== 

578 

579 >>> from sympy import Point, Ellipse 

580 >>> p1 = Point(0, 0) 

581 >>> e1 = Ellipse(p1, 3, 1) 

582 >>> e1.focus_distance 

583 2*sqrt(2) 

584 

585 """ 

586 return Point.distance(self.center, self.foci[0]) 

587 

588 @property 

589 def hradius(self): 

590 """The horizontal radius of the ellipse. 

591 

592 Returns 

593 ======= 

594 

595 hradius : number 

596 

597 See Also 

598 ======== 

599 

600 vradius, major, minor 

601 

602 Examples 

603 ======== 

604 

605 >>> from sympy import Point, Ellipse 

606 >>> p1 = Point(0, 0) 

607 >>> e1 = Ellipse(p1, 3, 1) 

608 >>> e1.hradius 

609 3 

610 

611 """ 

612 return self.args[1] 

613 

614 def intersection(self, o): 

615 """The intersection of this ellipse and another geometrical entity 

616 `o`. 

617 

618 Parameters 

619 ========== 

620 

621 o : GeometryEntity 

622 

623 Returns 

624 ======= 

625 

626 intersection : list of GeometryEntity objects 

627 

628 Notes 

629 ----- 

630 Currently supports intersections with Point, Line, Segment, Ray, 

631 Circle and Ellipse types. 

632 

633 See Also 

634 ======== 

635 

636 sympy.geometry.entity.GeometryEntity 

637 

638 Examples 

639 ======== 

640 

641 >>> from sympy import Ellipse, Point, Line 

642 >>> e = Ellipse(Point(0, 0), 5, 7) 

643 >>> e.intersection(Point(0, 0)) 

644 [] 

645 >>> e.intersection(Point(5, 0)) 

646 [Point2D(5, 0)] 

647 >>> e.intersection(Line(Point(0,0), Point(0, 1))) 

648 [Point2D(0, -7), Point2D(0, 7)] 

649 >>> e.intersection(Line(Point(5,0), Point(5, 1))) 

650 [Point2D(5, 0)] 

651 >>> e.intersection(Line(Point(6,0), Point(6, 1))) 

652 [] 

653 >>> e = Ellipse(Point(-1, 0), 4, 3) 

654 >>> e.intersection(Ellipse(Point(1, 0), 4, 3)) 

655 [Point2D(0, -3*sqrt(15)/4), Point2D(0, 3*sqrt(15)/4)] 

656 >>> e.intersection(Ellipse(Point(5, 0), 4, 3)) 

657 [Point2D(2, -3*sqrt(7)/4), Point2D(2, 3*sqrt(7)/4)] 

658 >>> e.intersection(Ellipse(Point(100500, 0), 4, 3)) 

659 [] 

660 >>> e.intersection(Ellipse(Point(0, 0), 3, 4)) 

661 [Point2D(3, 0), Point2D(-363/175, -48*sqrt(111)/175), Point2D(-363/175, 48*sqrt(111)/175)] 

662 >>> e.intersection(Ellipse(Point(-1, 0), 3, 4)) 

663 [Point2D(-17/5, -12/5), Point2D(-17/5, 12/5), Point2D(7/5, -12/5), Point2D(7/5, 12/5)] 

664 """ 

665 # TODO: Replace solve with nonlinsolve, when nonlinsolve will be able to solve in real domain 

666 

667 if isinstance(o, Point): 

668 if o in self: 

669 return [o] 

670 else: 

671 return [] 

672 

673 elif isinstance(o, (Segment2D, Ray2D)): 

674 ellipse_equation = self.equation(x, y) 

675 result = solve([ellipse_equation, Line( 

676 o.points[0], o.points[1]).equation(x, y)], [x, y], 

677 set=True)[1] 

678 return list(ordered([Point(i) for i in result if i in o])) 

679 

680 elif isinstance(o, Polygon): 

681 return o.intersection(self) 

682 

683 elif isinstance(o, (Ellipse, Line2D)): 

684 if o == self: 

685 return self 

686 else: 

687 ellipse_equation = self.equation(x, y) 

688 return list(ordered([Point(i) for i in solve( 

689 [ellipse_equation, o.equation(x, y)], [x, y], 

690 set=True)[1]])) 

691 elif isinstance(o, LinearEntity3D): 

692 raise TypeError('Entity must be two dimensional, not three dimensional') 

693 else: 

694 raise TypeError('Intersection not handled for %s' % func_name(o)) 

695 

696 def is_tangent(self, o): 

697 """Is `o` tangent to the ellipse? 

698 

699 Parameters 

700 ========== 

701 

702 o : GeometryEntity 

703 An Ellipse, LinearEntity or Polygon 

704 

705 Raises 

706 ====== 

707 

708 NotImplementedError 

709 When the wrong type of argument is supplied. 

710 

711 Returns 

712 ======= 

713 

714 is_tangent: boolean 

715 True if o is tangent to the ellipse, False otherwise. 

716 

717 See Also 

718 ======== 

719 

720 tangent_lines 

721 

722 Examples 

723 ======== 

724 

725 >>> from sympy import Point, Ellipse, Line 

726 >>> p0, p1, p2 = Point(0, 0), Point(3, 0), Point(3, 3) 

727 >>> e1 = Ellipse(p0, 3, 2) 

728 >>> l1 = Line(p1, p2) 

729 >>> e1.is_tangent(l1) 

730 True 

731 

732 """ 

733 if isinstance(o, Point2D): 

734 return False 

735 elif isinstance(o, Ellipse): 

736 intersect = self.intersection(o) 

737 if isinstance(intersect, Ellipse): 

738 return True 

739 elif intersect: 

740 return all((self.tangent_lines(i)[0]).equals(o.tangent_lines(i)[0]) for i in intersect) 

741 else: 

742 return False 

743 elif isinstance(o, Line2D): 

744 hit = self.intersection(o) 

745 if not hit: 

746 return False 

747 if len(hit) == 1: 

748 return True 

749 # might return None if it can't decide 

750 return hit[0].equals(hit[1]) 

751 elif isinstance(o, Ray2D): 

752 intersect = self.intersection(o) 

753 if len(intersect) == 1: 

754 return intersect[0] != o.source and not self.encloses_point(o.source) 

755 else: 

756 return False 

757 elif isinstance(o, (Segment2D, Polygon)): 

758 all_tangents = False 

759 segments = o.sides if isinstance(o, Polygon) else [o] 

760 for segment in segments: 

761 intersect = self.intersection(segment) 

762 if len(intersect) == 1: 

763 if not any(intersect[0] in i for i in segment.points) \ 

764 and not any(self.encloses_point(i) for i in segment.points): 

765 all_tangents = True 

766 continue 

767 else: 

768 return False 

769 else: 

770 return all_tangents 

771 return all_tangents 

772 elif isinstance(o, (LinearEntity3D, Point3D)): 

773 raise TypeError('Entity must be two dimensional, not three dimensional') 

774 else: 

775 raise TypeError('Is_tangent not handled for %s' % func_name(o)) 

776 

777 @property 

778 def major(self): 

779 """Longer axis of the ellipse (if it can be determined) else hradius. 

780 

781 Returns 

782 ======= 

783 

784 major : number or expression 

785 

786 See Also 

787 ======== 

788 

789 hradius, vradius, minor 

790 

791 Examples 

792 ======== 

793 

794 >>> from sympy import Point, Ellipse, Symbol 

795 >>> p1 = Point(0, 0) 

796 >>> e1 = Ellipse(p1, 3, 1) 

797 >>> e1.major 

798 3 

799 

800 >>> a = Symbol('a') 

801 >>> b = Symbol('b') 

802 >>> Ellipse(p1, a, b).major 

803 a 

804 >>> Ellipse(p1, b, a).major 

805 b 

806 

807 >>> m = Symbol('m') 

808 >>> M = m + 1 

809 >>> Ellipse(p1, m, M).major 

810 m + 1 

811 

812 """ 

813 ab = self.args[1:3] 

814 if len(ab) == 1: 

815 return ab[0] 

816 a, b = ab 

817 o = b - a < 0 

818 if o == True: 

819 return a 

820 elif o == False: 

821 return b 

822 return self.hradius 

823 

824 @property 

825 def minor(self): 

826 """Shorter axis of the ellipse (if it can be determined) else vradius. 

827 

828 Returns 

829 ======= 

830 

831 minor : number or expression 

832 

833 See Also 

834 ======== 

835 

836 hradius, vradius, major 

837 

838 Examples 

839 ======== 

840 

841 >>> from sympy import Point, Ellipse, Symbol 

842 >>> p1 = Point(0, 0) 

843 >>> e1 = Ellipse(p1, 3, 1) 

844 >>> e1.minor 

845 1 

846 

847 >>> a = Symbol('a') 

848 >>> b = Symbol('b') 

849 >>> Ellipse(p1, a, b).minor 

850 b 

851 >>> Ellipse(p1, b, a).minor 

852 a 

853 

854 >>> m = Symbol('m') 

855 >>> M = m + 1 

856 >>> Ellipse(p1, m, M).minor 

857 m 

858 

859 """ 

860 ab = self.args[1:3] 

861 if len(ab) == 1: 

862 return ab[0] 

863 a, b = ab 

864 o = a - b < 0 

865 if o == True: 

866 return a 

867 elif o == False: 

868 return b 

869 return self.vradius 

870 

871 def normal_lines(self, p, prec=None): 

872 """Normal lines between `p` and the ellipse. 

873 

874 Parameters 

875 ========== 

876 

877 p : Point 

878 

879 Returns 

880 ======= 

881 

882 normal_lines : list with 1, 2 or 4 Lines 

883 

884 Examples 

885 ======== 

886 

887 >>> from sympy import Point, Ellipse 

888 >>> e = Ellipse((0, 0), 2, 3) 

889 >>> c = e.center 

890 >>> e.normal_lines(c + Point(1, 0)) 

891 [Line2D(Point2D(0, 0), Point2D(1, 0))] 

892 >>> e.normal_lines(c) 

893 [Line2D(Point2D(0, 0), Point2D(0, 1)), Line2D(Point2D(0, 0), Point2D(1, 0))] 

894 

895 Off-axis points require the solution of a quartic equation. This 

896 often leads to very large expressions that may be of little practical 

897 use. An approximate solution of `prec` digits can be obtained by 

898 passing in the desired value: 

899 

900 >>> e.normal_lines((3, 3), prec=2) 

901 [Line2D(Point2D(-0.81, -2.7), Point2D(0.19, -1.2)), 

902 Line2D(Point2D(1.5, -2.0), Point2D(2.5, -2.7))] 

903 

904 Whereas the above solution has an operation count of 12, the exact 

905 solution has an operation count of 2020. 

906 """ 

907 p = Point(p, dim=2) 

908 

909 # XXX change True to something like self.angle == 0 if the arbitrarily 

910 # rotated ellipse is introduced. 

911 # https://github.com/sympy/sympy/issues/2815) 

912 if True: 

913 rv = [] 

914 if p.x == self.center.x: 

915 rv.append(Line(self.center, slope=oo)) 

916 if p.y == self.center.y: 

917 rv.append(Line(self.center, slope=0)) 

918 if rv: 

919 # at these special orientations of p either 1 or 2 normals 

920 # exist and we are done 

921 return rv 

922 

923 # find the 4 normal points and construct lines through them with 

924 # the corresponding slope 

925 eq = self.equation(x, y) 

926 dydx = idiff(eq, y, x) 

927 norm = -1/dydx 

928 slope = Line(p, (x, y)).slope 

929 seq = slope - norm 

930 

931 # TODO: Replace solve with solveset, when this line is tested 

932 yis = solve(seq, y)[0] 

933 xeq = eq.subs(y, yis).as_numer_denom()[0].expand() 

934 if len(xeq.free_symbols) == 1: 

935 try: 

936 # this is so much faster, it's worth a try 

937 xsol = Poly(xeq, x).real_roots() 

938 except (DomainError, PolynomialError, NotImplementedError): 

939 # TODO: Replace solve with solveset, when these lines are tested 

940 xsol = _nsort(solve(xeq, x), separated=True)[0] 

941 points = [Point(i, solve(eq.subs(x, i), y)[0]) for i in xsol] 

942 else: 

943 raise NotImplementedError( 

944 'intersections for the general ellipse are not supported') 

945 slopes = [norm.subs(zip((x, y), pt.args)) for pt in points] 

946 if prec is not None: 

947 points = [pt.n(prec) for pt in points] 

948 slopes = [i if _not_a_coeff(i) else i.n(prec) for i in slopes] 

949 return [Line(pt, slope=s) for pt, s in zip(points, slopes)] 

950 

951 @property 

952 def periapsis(self): 

953 """The periapsis of the ellipse. 

954 

955 The shortest distance between the focus and the contour. 

956 

957 Returns 

958 ======= 

959 

960 periapsis : number 

961 

962 See Also 

963 ======== 

964 

965 apoapsis : Returns greatest distance between focus and contour 

966 

967 Examples 

968 ======== 

969 

970 >>> from sympy import Point, Ellipse 

971 >>> p1 = Point(0, 0) 

972 >>> e1 = Ellipse(p1, 3, 1) 

973 >>> e1.periapsis 

974 3 - 2*sqrt(2) 

975 

976 """ 

977 return self.major * (1 - self.eccentricity) 

978 

979 @property 

980 def semilatus_rectum(self): 

981 """ 

982 Calculates the semi-latus rectum of the Ellipse. 

983 

984 Semi-latus rectum is defined as one half of the chord through a 

985 focus parallel to the conic section directrix of a conic section. 

986 

987 Returns 

988 ======= 

989 

990 semilatus_rectum : number 

991 

992 See Also 

993 ======== 

994 

995 apoapsis : Returns greatest distance between focus and contour 

996 

997 periapsis : The shortest distance between the focus and the contour 

998 

999 Examples 

1000 ======== 

1001 

1002 >>> from sympy import Point, Ellipse 

1003 >>> p1 = Point(0, 0) 

1004 >>> e1 = Ellipse(p1, 3, 1) 

1005 >>> e1.semilatus_rectum 

1006 1/3 

1007 

1008 References 

1009 ========== 

1010 

1011 .. [1] https://mathworld.wolfram.com/SemilatusRectum.html 

1012 .. [2] https://en.wikipedia.org/wiki/Ellipse#Semi-latus_rectum 

1013 

1014 """ 

1015 return self.major * (1 - self.eccentricity ** 2) 

1016 

1017 def auxiliary_circle(self): 

1018 """Returns a Circle whose diameter is the major axis of the ellipse. 

1019 

1020 Examples 

1021 ======== 

1022 

1023 >>> from sympy import Ellipse, Point, symbols 

1024 >>> c = Point(1, 2) 

1025 >>> Ellipse(c, 8, 7).auxiliary_circle() 

1026 Circle(Point2D(1, 2), 8) 

1027 >>> a, b = symbols('a b') 

1028 >>> Ellipse(c, a, b).auxiliary_circle() 

1029 Circle(Point2D(1, 2), Max(a, b)) 

1030 """ 

1031 return Circle(self.center, Max(self.hradius, self.vradius)) 

1032 

1033 def director_circle(self): 

1034 """ 

1035 Returns a Circle consisting of all points where two perpendicular 

1036 tangent lines to the ellipse cross each other. 

1037 

1038 Returns 

1039 ======= 

1040 

1041 Circle 

1042 A director circle returned as a geometric object. 

1043 

1044 Examples 

1045 ======== 

1046 

1047 >>> from sympy import Ellipse, Point, symbols 

1048 >>> c = Point(3,8) 

1049 >>> Ellipse(c, 7, 9).director_circle() 

1050 Circle(Point2D(3, 8), sqrt(130)) 

1051 >>> a, b = symbols('a b') 

1052 >>> Ellipse(c, a, b).director_circle() 

1053 Circle(Point2D(3, 8), sqrt(a**2 + b**2)) 

1054 

1055 References 

1056 ========== 

1057 

1058 .. [1] https://en.wikipedia.org/wiki/Director_circle 

1059 

1060 """ 

1061 return Circle(self.center, sqrt(self.hradius**2 + self.vradius**2)) 

1062 

1063 def plot_interval(self, parameter='t'): 

1064 """The plot interval for the default geometric plot of the Ellipse. 

1065 

1066 Parameters 

1067 ========== 

1068 

1069 parameter : str, optional 

1070 Default value is 't'. 

1071 

1072 Returns 

1073 ======= 

1074 

1075 plot_interval : list 

1076 [parameter, lower_bound, upper_bound] 

1077 

1078 Examples 

1079 ======== 

1080 

1081 >>> from sympy import Point, Ellipse 

1082 >>> e1 = Ellipse(Point(0, 0), 3, 2) 

1083 >>> e1.plot_interval() 

1084 [t, -pi, pi] 

1085 

1086 """ 

1087 t = _symbol(parameter, real=True) 

1088 return [t, -S.Pi, S.Pi] 

1089 

1090 def random_point(self, seed=None): 

1091 """A random point on the ellipse. 

1092 

1093 Returns 

1094 ======= 

1095 

1096 point : Point 

1097 

1098 Examples 

1099 ======== 

1100 

1101 >>> from sympy import Point, Ellipse 

1102 >>> e1 = Ellipse(Point(0, 0), 3, 2) 

1103 >>> e1.random_point() # gives some random point 

1104 Point2D(...) 

1105 >>> p1 = e1.random_point(seed=0); p1.n(2) 

1106 Point2D(2.1, 1.4) 

1107 

1108 Notes 

1109 ===== 

1110 

1111 When creating a random point, one may simply replace the 

1112 parameter with a random number. When doing so, however, the 

1113 random number should be made a Rational or else the point 

1114 may not test as being in the ellipse: 

1115 

1116 >>> from sympy.abc import t 

1117 >>> from sympy import Rational 

1118 >>> arb = e1.arbitrary_point(t); arb 

1119 Point2D(3*cos(t), 2*sin(t)) 

1120 >>> arb.subs(t, .1) in e1 

1121 False 

1122 >>> arb.subs(t, Rational(.1)) in e1 

1123 True 

1124 >>> arb.subs(t, Rational('.1')) in e1 

1125 True 

1126 

1127 See Also 

1128 ======== 

1129 sympy.geometry.point.Point 

1130 arbitrary_point : Returns parameterized point on ellipse 

1131 """ 

1132 t = _symbol('t', real=True) 

1133 x, y = self.arbitrary_point(t).args 

1134 # get a random value in [-1, 1) corresponding to cos(t) 

1135 # and confirm that it will test as being in the ellipse 

1136 if seed is not None: 

1137 rng = random.Random(seed) 

1138 else: 

1139 rng = random 

1140 # simplify this now or else the Float will turn s into a Float 

1141 r = Rational(rng.random()) 

1142 c = 2*r - 1 

1143 s = sqrt(1 - c**2) 

1144 return Point(x.subs(cos(t), c), y.subs(sin(t), s)) 

1145 

1146 def reflect(self, line): 

1147 """Override GeometryEntity.reflect since the radius 

1148 is not a GeometryEntity. 

1149 

1150 Examples 

1151 ======== 

1152 

1153 >>> from sympy import Circle, Line 

1154 >>> Circle((0, 1), 1).reflect(Line((0, 0), (1, 1))) 

1155 Circle(Point2D(1, 0), -1) 

1156 >>> from sympy import Ellipse, Line, Point 

1157 >>> Ellipse(Point(3, 4), 1, 3).reflect(Line(Point(0, -4), Point(5, 0))) 

1158 Traceback (most recent call last): 

1159 ... 

1160 NotImplementedError: 

1161 General Ellipse is not supported but the equation of the reflected 

1162 Ellipse is given by the zeros of: f(x, y) = (9*x/41 + 40*y/41 + 

1163 37/41)**2 + (40*x/123 - 3*y/41 - 364/123)**2 - 1 

1164 

1165 Notes 

1166 ===== 

1167 

1168 Until the general ellipse (with no axis parallel to the x-axis) is 

1169 supported a NotImplemented error is raised and the equation whose 

1170 zeros define the rotated ellipse is given. 

1171 

1172 """ 

1173 

1174 if line.slope in (0, oo): 

1175 c = self.center 

1176 c = c.reflect(line) 

1177 return self.func(c, -self.hradius, self.vradius) 

1178 else: 

1179 x, y = [uniquely_named_symbol( 

1180 name, (self, line), modify=lambda s: '_' + s, real=True) 

1181 for name in 'xy'] 

1182 expr = self.equation(x, y) 

1183 p = Point(x, y).reflect(line) 

1184 result = expr.subs(zip((x, y), p.args 

1185 ), simultaneous=True) 

1186 raise NotImplementedError(filldedent( 

1187 'General Ellipse is not supported but the equation ' 

1188 'of the reflected Ellipse is given by the zeros of: ' + 

1189 "f(%s, %s) = %s" % (str(x), str(y), str(result)))) 

1190 

1191 def rotate(self, angle=0, pt=None): 

1192 """Rotate ``angle`` radians counterclockwise about Point ``pt``. 

1193 

1194 Note: since the general ellipse is not supported, only rotations that 

1195 are integer multiples of pi/2 are allowed. 

1196 

1197 Examples 

1198 ======== 

1199 

1200 >>> from sympy import Ellipse, pi 

1201 >>> Ellipse((1, 0), 2, 1).rotate(pi/2) 

1202 Ellipse(Point2D(0, 1), 1, 2) 

1203 >>> Ellipse((1, 0), 2, 1).rotate(pi) 

1204 Ellipse(Point2D(-1, 0), 2, 1) 

1205 """ 

1206 if self.hradius == self.vradius: 

1207 return self.func(self.center.rotate(angle, pt), self.hradius) 

1208 if (angle/S.Pi).is_integer: 

1209 return super().rotate(angle, pt) 

1210 if (2*angle/S.Pi).is_integer: 

1211 return self.func(self.center.rotate(angle, pt), self.vradius, self.hradius) 

1212 # XXX see https://github.com/sympy/sympy/issues/2815 for general ellipes 

1213 raise NotImplementedError('Only rotations of pi/2 are currently supported for Ellipse.') 

1214 

1215 def scale(self, x=1, y=1, pt=None): 

1216 """Override GeometryEntity.scale since it is the major and minor 

1217 axes which must be scaled and they are not GeometryEntities. 

1218 

1219 Examples 

1220 ======== 

1221 

1222 >>> from sympy import Ellipse 

1223 >>> Ellipse((0, 0), 2, 1).scale(2, 4) 

1224 Circle(Point2D(0, 0), 4) 

1225 >>> Ellipse((0, 0), 2, 1).scale(2) 

1226 Ellipse(Point2D(0, 0), 4, 1) 

1227 """ 

1228 c = self.center 

1229 if pt: 

1230 pt = Point(pt, dim=2) 

1231 return self.translate(*(-pt).args).scale(x, y).translate(*pt.args) 

1232 h = self.hradius 

1233 v = self.vradius 

1234 return self.func(c.scale(x, y), hradius=h*x, vradius=v*y) 

1235 

1236 def tangent_lines(self, p): 

1237 """Tangent lines between `p` and the ellipse. 

1238 

1239 If `p` is on the ellipse, returns the tangent line through point `p`. 

1240 Otherwise, returns the tangent line(s) from `p` to the ellipse, or 

1241 None if no tangent line is possible (e.g., `p` inside ellipse). 

1242 

1243 Parameters 

1244 ========== 

1245 

1246 p : Point 

1247 

1248 Returns 

1249 ======= 

1250 

1251 tangent_lines : list with 1 or 2 Lines 

1252 

1253 Raises 

1254 ====== 

1255 

1256 NotImplementedError 

1257 Can only find tangent lines for a point, `p`, on the ellipse. 

1258 

1259 See Also 

1260 ======== 

1261 

1262 sympy.geometry.point.Point, sympy.geometry.line.Line 

1263 

1264 Examples 

1265 ======== 

1266 

1267 >>> from sympy import Point, Ellipse 

1268 >>> e1 = Ellipse(Point(0, 0), 3, 2) 

1269 >>> e1.tangent_lines(Point(3, 0)) 

1270 [Line2D(Point2D(3, 0), Point2D(3, -12))] 

1271 

1272 """ 

1273 p = Point(p, dim=2) 

1274 if self.encloses_point(p): 

1275 return [] 

1276 

1277 if p in self: 

1278 delta = self.center - p 

1279 rise = (self.vradius**2)*delta.x 

1280 run = -(self.hradius**2)*delta.y 

1281 p2 = Point(simplify(p.x + run), 

1282 simplify(p.y + rise)) 

1283 return [Line(p, p2)] 

1284 else: 

1285 if len(self.foci) == 2: 

1286 f1, f2 = self.foci 

1287 maj = self.hradius 

1288 test = (2*maj - 

1289 Point.distance(f1, p) - 

1290 Point.distance(f2, p)) 

1291 else: 

1292 test = self.radius - Point.distance(self.center, p) 

1293 if test.is_number and test.is_positive: 

1294 return [] 

1295 # else p is outside the ellipse or we can't tell. In case of the 

1296 # latter, the solutions returned will only be valid if 

1297 # the point is not inside the ellipse; if it is, nan will result. 

1298 eq = self.equation(x, y) 

1299 dydx = idiff(eq, y, x) 

1300 slope = Line(p, Point(x, y)).slope 

1301 

1302 # TODO: Replace solve with solveset, when this line is tested 

1303 tangent_points = solve([slope - dydx, eq], [x, y]) 

1304 

1305 # handle horizontal and vertical tangent lines 

1306 if len(tangent_points) == 1: 

1307 if tangent_points[0][ 

1308 0] == p.x or tangent_points[0][1] == p.y: 

1309 return [Line(p, p + Point(1, 0)), Line(p, p + Point(0, 1))] 

1310 else: 

1311 return [Line(p, p + Point(0, 1)), Line(p, tangent_points[0])] 

1312 

1313 # others 

1314 return [Line(p, tangent_points[0]), Line(p, tangent_points[1])] 

1315 

1316 @property 

1317 def vradius(self): 

1318 """The vertical radius of the ellipse. 

1319 

1320 Returns 

1321 ======= 

1322 

1323 vradius : number 

1324 

1325 See Also 

1326 ======== 

1327 

1328 hradius, major, minor 

1329 

1330 Examples 

1331 ======== 

1332 

1333 >>> from sympy import Point, Ellipse 

1334 >>> p1 = Point(0, 0) 

1335 >>> e1 = Ellipse(p1, 3, 1) 

1336 >>> e1.vradius 

1337 1 

1338 

1339 """ 

1340 return self.args[2] 

1341 

1342 

1343 def second_moment_of_area(self, point=None): 

1344 """Returns the second moment and product moment area of an ellipse. 

1345 

1346 Parameters 

1347 ========== 

1348 

1349 point : Point, two-tuple of sympifiable objects, or None(default=None) 

1350 point is the point about which second moment of area is to be found. 

1351 If "point=None" it will be calculated about the axis passing through the 

1352 centroid of the ellipse. 

1353 

1354 Returns 

1355 ======= 

1356 

1357 I_xx, I_yy, I_xy : number or SymPy expression 

1358 I_xx, I_yy are second moment of area of an ellise. 

1359 I_xy is product moment of area of an ellipse. 

1360 

1361 Examples 

1362 ======== 

1363 

1364 >>> from sympy import Point, Ellipse 

1365 >>> p1 = Point(0, 0) 

1366 >>> e1 = Ellipse(p1, 3, 1) 

1367 >>> e1.second_moment_of_area() 

1368 (3*pi/4, 27*pi/4, 0) 

1369 

1370 References 

1371 ========== 

1372 

1373 .. [1] https://en.wikipedia.org/wiki/List_of_second_moments_of_area 

1374 

1375 """ 

1376 

1377 I_xx = (S.Pi*(self.hradius)*(self.vradius**3))/4 

1378 I_yy = (S.Pi*(self.hradius**3)*(self.vradius))/4 

1379 I_xy = 0 

1380 

1381 if point is None: 

1382 return I_xx, I_yy, I_xy 

1383 

1384 # parallel axis theorem 

1385 I_xx = I_xx + self.area*((point[1] - self.center.y)**2) 

1386 I_yy = I_yy + self.area*((point[0] - self.center.x)**2) 

1387 I_xy = I_xy + self.area*(point[0] - self.center.x)*(point[1] - self.center.y) 

1388 

1389 return I_xx, I_yy, I_xy 

1390 

1391 

1392 def polar_second_moment_of_area(self): 

1393 """Returns the polar second moment of area of an Ellipse 

1394 

1395 It is a constituent of the second moment of area, linked through 

1396 the perpendicular axis theorem. While the planar second moment of 

1397 area describes an object's resistance to deflection (bending) when 

1398 subjected to a force applied to a plane parallel to the central 

1399 axis, the polar second moment of area describes an object's 

1400 resistance to deflection when subjected to a moment applied in a 

1401 plane perpendicular to the object's central axis (i.e. parallel to 

1402 the cross-section) 

1403 

1404 Examples 

1405 ======== 

1406 

1407 >>> from sympy import symbols, Circle, Ellipse 

1408 >>> c = Circle((5, 5), 4) 

1409 >>> c.polar_second_moment_of_area() 

1410 128*pi 

1411 >>> a, b = symbols('a, b') 

1412 >>> e = Ellipse((0, 0), a, b) 

1413 >>> e.polar_second_moment_of_area() 

1414 pi*a**3*b/4 + pi*a*b**3/4 

1415 

1416 References 

1417 ========== 

1418 

1419 .. [1] https://en.wikipedia.org/wiki/Polar_moment_of_inertia 

1420 

1421 """ 

1422 second_moment = self.second_moment_of_area() 

1423 return second_moment[0] + second_moment[1] 

1424 

1425 

1426 def section_modulus(self, point=None): 

1427 """Returns a tuple with the section modulus of an ellipse 

1428 

1429 Section modulus is a geometric property of an ellipse defined as the 

1430 ratio of second moment of area to the distance of the extreme end of 

1431 the ellipse from the centroidal axis. 

1432 

1433 Parameters 

1434 ========== 

1435 

1436 point : Point, two-tuple of sympifyable objects, or None(default=None) 

1437 point is the point at which section modulus is to be found. 

1438 If "point=None" section modulus will be calculated for the 

1439 point farthest from the centroidal axis of the ellipse. 

1440 

1441 Returns 

1442 ======= 

1443 

1444 S_x, S_y: numbers or SymPy expressions 

1445 S_x is the section modulus with respect to the x-axis 

1446 S_y is the section modulus with respect to the y-axis 

1447 A negative sign indicates that the section modulus is 

1448 determined for a point below the centroidal axis. 

1449 

1450 Examples 

1451 ======== 

1452 

1453 >>> from sympy import Symbol, Ellipse, Circle, Point2D 

1454 >>> d = Symbol('d', positive=True) 

1455 >>> c = Circle((0, 0), d/2) 

1456 >>> c.section_modulus() 

1457 (pi*d**3/32, pi*d**3/32) 

1458 >>> e = Ellipse(Point2D(0, 0), 2, 4) 

1459 >>> e.section_modulus() 

1460 (8*pi, 4*pi) 

1461 >>> e.section_modulus((2, 2)) 

1462 (16*pi, 4*pi) 

1463 

1464 References 

1465 ========== 

1466 

1467 .. [1] https://en.wikipedia.org/wiki/Section_modulus 

1468 

1469 """ 

1470 x_c, y_c = self.center 

1471 if point is None: 

1472 # taking x and y as maximum distances from centroid 

1473 x_min, y_min, x_max, y_max = self.bounds 

1474 y = max(y_c - y_min, y_max - y_c) 

1475 x = max(x_c - x_min, x_max - x_c) 

1476 else: 

1477 # taking x and y as distances of the given point from the center 

1478 point = Point2D(point) 

1479 y = point.y - y_c 

1480 x = point.x - x_c 

1481 

1482 second_moment = self.second_moment_of_area() 

1483 S_x = second_moment[0]/y 

1484 S_y = second_moment[1]/x 

1485 

1486 return S_x, S_y 

1487 

1488 

1489class Circle(Ellipse): 

1490 """A circle in space. 

1491 

1492 Constructed simply from a center and a radius, from three 

1493 non-collinear points, or the equation of a circle. 

1494 

1495 Parameters 

1496 ========== 

1497 

1498 center : Point 

1499 radius : number or SymPy expression 

1500 points : sequence of three Points 

1501 equation : equation of a circle 

1502 

1503 Attributes 

1504 ========== 

1505 

1506 radius (synonymous with hradius, vradius, major and minor) 

1507 circumference 

1508 equation 

1509 

1510 Raises 

1511 ====== 

1512 

1513 GeometryError 

1514 When the given equation is not that of a circle. 

1515 When trying to construct circle from incorrect parameters. 

1516 

1517 See Also 

1518 ======== 

1519 

1520 Ellipse, sympy.geometry.point.Point 

1521 

1522 Examples 

1523 ======== 

1524 

1525 >>> from sympy import Point, Circle, Eq 

1526 >>> from sympy.abc import x, y, a, b 

1527 

1528 A circle constructed from a center and radius: 

1529 

1530 >>> c1 = Circle(Point(0, 0), 5) 

1531 >>> c1.hradius, c1.vradius, c1.radius 

1532 (5, 5, 5) 

1533 

1534 A circle constructed from three points: 

1535 

1536 >>> c2 = Circle(Point(0, 0), Point(1, 1), Point(1, 0)) 

1537 >>> c2.hradius, c2.vradius, c2.radius, c2.center 

1538 (sqrt(2)/2, sqrt(2)/2, sqrt(2)/2, Point2D(1/2, 1/2)) 

1539 

1540 A circle can be constructed from an equation in the form 

1541 `a*x**2 + by**2 + gx + hy + c = 0`, too: 

1542 

1543 >>> Circle(x**2 + y**2 - 25) 

1544 Circle(Point2D(0, 0), 5) 

1545 

1546 If the variables corresponding to x and y are named something 

1547 else, their name or symbol can be supplied: 

1548 

1549 >>> Circle(Eq(a**2 + b**2, 25), x='a', y=b) 

1550 Circle(Point2D(0, 0), 5) 

1551 """ 

1552 

1553 def __new__(cls, *args, **kwargs): 

1554 evaluate = kwargs.get('evaluate', global_parameters.evaluate) 

1555 if len(args) == 1 and isinstance(args[0], (Expr, Eq)): 

1556 x = kwargs.get('x', 'x') 

1557 y = kwargs.get('y', 'y') 

1558 equation = args[0].expand() 

1559 if isinstance(equation, Eq): 

1560 equation = equation.lhs - equation.rhs 

1561 x = find(x, equation) 

1562 y = find(y, equation) 

1563 

1564 try: 

1565 a, b, c, d, e = linear_coeffs(equation, x**2, y**2, x, y) 

1566 except ValueError: 

1567 raise GeometryError("The given equation is not that of a circle.") 

1568 

1569 if S.Zero in (a, b) or a != b: 

1570 raise GeometryError("The given equation is not that of a circle.") 

1571 

1572 center_x = -c/a/2 

1573 center_y = -d/b/2 

1574 r2 = (center_x**2) + (center_y**2) - e/a 

1575 

1576 return Circle((center_x, center_y), sqrt(r2), evaluate=evaluate) 

1577 

1578 else: 

1579 c, r = None, None 

1580 if len(args) == 3: 

1581 args = [Point(a, dim=2, evaluate=evaluate) for a in args] 

1582 t = Triangle(*args) 

1583 if not isinstance(t, Triangle): 

1584 return t 

1585 c = t.circumcenter 

1586 r = t.circumradius 

1587 elif len(args) == 2: 

1588 # Assume (center, radius) pair 

1589 c = Point(args[0], dim=2, evaluate=evaluate) 

1590 r = args[1] 

1591 # this will prohibit imaginary radius 

1592 try: 

1593 r = Point(r, 0, evaluate=evaluate).x 

1594 except ValueError: 

1595 raise GeometryError("Circle with imaginary radius is not permitted") 

1596 

1597 if not (c is None or r is None): 

1598 if r == 0: 

1599 return c 

1600 return GeometryEntity.__new__(cls, c, r, **kwargs) 

1601 

1602 raise GeometryError("Circle.__new__ received unknown arguments") 

1603 

1604 def _eval_evalf(self, prec=15, **options): 

1605 pt, r = self.args 

1606 dps = prec_to_dps(prec) 

1607 pt = pt.evalf(n=dps, **options) 

1608 r = r.evalf(n=dps, **options) 

1609 return self.func(pt, r, evaluate=False) 

1610 

1611 @property 

1612 def circumference(self): 

1613 """The circumference of the circle. 

1614 

1615 Returns 

1616 ======= 

1617 

1618 circumference : number or SymPy expression 

1619 

1620 Examples 

1621 ======== 

1622 

1623 >>> from sympy import Point, Circle 

1624 >>> c1 = Circle(Point(3, 4), 6) 

1625 >>> c1.circumference 

1626 12*pi 

1627 

1628 """ 

1629 return 2 * S.Pi * self.radius 

1630 

1631 def equation(self, x='x', y='y'): 

1632 """The equation of the circle. 

1633 

1634 Parameters 

1635 ========== 

1636 

1637 x : str or Symbol, optional 

1638 Default value is 'x'. 

1639 y : str or Symbol, optional 

1640 Default value is 'y'. 

1641 

1642 Returns 

1643 ======= 

1644 

1645 equation : SymPy expression 

1646 

1647 Examples 

1648 ======== 

1649 

1650 >>> from sympy import Point, Circle 

1651 >>> c1 = Circle(Point(0, 0), 5) 

1652 >>> c1.equation() 

1653 x**2 + y**2 - 25 

1654 

1655 """ 

1656 x = _symbol(x, real=True) 

1657 y = _symbol(y, real=True) 

1658 t1 = (x - self.center.x)**2 

1659 t2 = (y - self.center.y)**2 

1660 return t1 + t2 - self.major**2 

1661 

1662 def intersection(self, o): 

1663 """The intersection of this circle with another geometrical entity. 

1664 

1665 Parameters 

1666 ========== 

1667 

1668 o : GeometryEntity 

1669 

1670 Returns 

1671 ======= 

1672 

1673 intersection : list of GeometryEntities 

1674 

1675 Examples 

1676 ======== 

1677 

1678 >>> from sympy import Point, Circle, Line, Ray 

1679 >>> p1, p2, p3 = Point(0, 0), Point(5, 5), Point(6, 0) 

1680 >>> p4 = Point(5, 0) 

1681 >>> c1 = Circle(p1, 5) 

1682 >>> c1.intersection(p2) 

1683 [] 

1684 >>> c1.intersection(p4) 

1685 [Point2D(5, 0)] 

1686 >>> c1.intersection(Ray(p1, p2)) 

1687 [Point2D(5*sqrt(2)/2, 5*sqrt(2)/2)] 

1688 >>> c1.intersection(Line(p2, p3)) 

1689 [] 

1690 

1691 """ 

1692 return Ellipse.intersection(self, o) 

1693 

1694 @property 

1695 def radius(self): 

1696 """The radius of the circle. 

1697 

1698 Returns 

1699 ======= 

1700 

1701 radius : number or SymPy expression 

1702 

1703 See Also 

1704 ======== 

1705 

1706 Ellipse.major, Ellipse.minor, Ellipse.hradius, Ellipse.vradius 

1707 

1708 Examples 

1709 ======== 

1710 

1711 >>> from sympy import Point, Circle 

1712 >>> c1 = Circle(Point(3, 4), 6) 

1713 >>> c1.radius 

1714 6 

1715 

1716 """ 

1717 return self.args[1] 

1718 

1719 def reflect(self, line): 

1720 """Override GeometryEntity.reflect since the radius 

1721 is not a GeometryEntity. 

1722 

1723 Examples 

1724 ======== 

1725 

1726 >>> from sympy import Circle, Line 

1727 >>> Circle((0, 1), 1).reflect(Line((0, 0), (1, 1))) 

1728 Circle(Point2D(1, 0), -1) 

1729 """ 

1730 c = self.center 

1731 c = c.reflect(line) 

1732 return self.func(c, -self.radius) 

1733 

1734 def scale(self, x=1, y=1, pt=None): 

1735 """Override GeometryEntity.scale since the radius 

1736 is not a GeometryEntity. 

1737 

1738 Examples 

1739 ======== 

1740 

1741 >>> from sympy import Circle 

1742 >>> Circle((0, 0), 1).scale(2, 2) 

1743 Circle(Point2D(0, 0), 2) 

1744 >>> Circle((0, 0), 1).scale(2, 4) 

1745 Ellipse(Point2D(0, 0), 2, 4) 

1746 """ 

1747 c = self.center 

1748 if pt: 

1749 pt = Point(pt, dim=2) 

1750 return self.translate(*(-pt).args).scale(x, y).translate(*pt.args) 

1751 c = c.scale(x, y) 

1752 x, y = [abs(i) for i in (x, y)] 

1753 if x == y: 

1754 return self.func(c, x*self.radius) 

1755 h = v = self.radius 

1756 return Ellipse(c, hradius=h*x, vradius=v*y) 

1757 

1758 @property 

1759 def vradius(self): 

1760 """ 

1761 This Ellipse property is an alias for the Circle's radius. 

1762 

1763 Whereas hradius, major and minor can use Ellipse's conventions, 

1764 the vradius does not exist for a circle. It is always a positive 

1765 value in order that the Circle, like Polygons, will have an 

1766 area that can be positive or negative as determined by the sign 

1767 of the hradius. 

1768 

1769 Examples 

1770 ======== 

1771 

1772 >>> from sympy import Point, Circle 

1773 >>> c1 = Circle(Point(3, 4), 6) 

1774 >>> c1.vradius 

1775 6 

1776 """ 

1777 return abs(self.radius) 

1778 

1779 

1780from .polygon import Polygon, Triangle