Coverage for /usr/lib/python3/dist-packages/sympy/simplify/sqrtdenest.py: 6%

1from sympy.core import Add, Expr, Mul, S, sympify

2from sympy.core.function import _mexpand, count_ops, expand_mul

3from sympy.core.sorting import default_sort_key

4from sympy.core.symbol import Dummy

5from sympy.functions import root, sign, sqrt

6from sympy.polys import Poly, PolynomialError

9def is_sqrt(expr):

10 """Return True if expr is a sqrt, otherwise False."""

12 return expr.is_Pow and expr.exp.is_Rational and abs(expr.exp) is S.Half

15def sqrt_depth(p):

16 """Return the maximum depth of any square root argument of p.

18 >>> from sympy.functions.elementary.miscellaneous import sqrt

19 >>> from sympy.simplify.sqrtdenest import sqrt_depth

21 Neither of these square roots contains any other square roots

22 so the depth is 1:

24 >>> sqrt_depth(1 + sqrt(2)*(1 + sqrt(3)))

25 1

27 The sqrt(3) is contained within a square root so the depth is

28 2:

30 >>> sqrt_depth(1 + sqrt(2)*sqrt(1 + sqrt(3)))

31 2

32 """

33 if p is S.ImaginaryUnit:

34 return 1

35 if p.is_Atom:

36 return 0

37 elif p.is_Add or p.is_Mul:

38 return max([sqrt_depth(x) for x in p.args], key=default_sort_key)

39 elif is_sqrt(p):

40 return sqrt_depth(p.base) + 1

41 else:

42 return 0

45def is_algebraic(p):

46 """Return True if p is comprised of only Rationals or square roots

47 of Rationals and algebraic operations.

49 Examples

50 ========

52 >>> from sympy.functions.elementary.miscellaneous import sqrt

53 >>> from sympy.simplify.sqrtdenest import is_algebraic

54 >>> from sympy import cos

55 >>> is_algebraic(sqrt(2)*(3/(sqrt(7) + sqrt(5)*sqrt(2))))

56 True

57 >>> is_algebraic(sqrt(2)*(3/(sqrt(7) + sqrt(5)*cos(2))))

58 False

59 """

61 if p.is_Rational:

62 return True

63 elif p.is_Atom:

64 return False

65 elif is_sqrt(p) or p.is_Pow and p.exp.is_Integer:

66 return is_algebraic(p.base)

67 elif p.is_Add or p.is_Mul:

68 return all(is_algebraic(x) for x in p.args)

69 else:

70 return False

73def _subsets(n):

74 """

75 Returns all possible subsets of the set (0, 1, ..., n-1) except the

76 empty set, listed in reversed lexicographical order according to binary

77 representation, so that the case of the fourth root is treated last.

79 Examples

80 ========

82 >>> from sympy.simplify.sqrtdenest import _subsets

83 >>> _subsets(2)

84 [[1, 0], [0, 1], [1, 1]]

86 """

87 if n == 1:

88 a = [[1]]

89 elif n == 2:

90 a = [[1, 0], [0, 1], [1, 1]]

91 elif n == 3:

92 a = [[1, 0, 0], [0, 1, 0], [1, 1, 0],

93 [0, 0, 1], [1, 0, 1], [0, 1, 1], [1, 1, 1]]

94 else:

95 b = _subsets(n - 1)

96 a0 = [x + [0] for x in b]

97 a1 = [x + [1] for x in b]

98 a = a0 + [[0]*(n - 1) + [1]] + a1

99 return a

100

101

102def sqrtdenest(expr, max_iter=3):

103 """Denests sqrts in an expression that contain other square roots

104 if possible, otherwise returns the expr unchanged. This is based on the

105 algorithms of [1].

106

107 Examples

108 ========

109

110 >>> from sympy.simplify.sqrtdenest import sqrtdenest

111 >>> from sympy import sqrt

112 >>> sqrtdenest(sqrt(5 + 2 * sqrt(6)))

113 sqrt(2) + sqrt(3)

114

115 See Also

116 ========

117

118 sympy.solvers.solvers.unrad

119

120 References

121 ==========

122

123 .. [1] https://web.archive.org/web/20210806201615/https://researcher.watson.ibm.com/researcher/files/us-fagin/symb85.pdf

124

125 .. [2] D. J. Jeffrey and A. D. Rich, 'Symplifying Square Roots of Square Roots

126 by Denesting' (available at https://www.cybertester.com/data/denest.pdf)

127

128 """

129 expr = expand_mul(expr)

130 for i in range(max_iter):

131 z = _sqrtdenest0(expr)

132 if expr == z:

133 return expr

134 expr = z

135 return expr

136

137

138def _sqrt_match(p):

139 """Return [a, b, r] for p.match(a + b*sqrt(r)) where, in addition to

140 matching, sqrt(r) also has then maximal sqrt_depth among addends of p.

141

142 Examples

143 ========

144

145 >>> from sympy.functions.elementary.miscellaneous import sqrt

146 >>> from sympy.simplify.sqrtdenest import _sqrt_match

147 >>> _sqrt_match(1 + sqrt(2) + sqrt(2)*sqrt(3) + 2*sqrt(1+sqrt(5)))

148 [1 + sqrt(2) + sqrt(6), 2, 1 + sqrt(5)]

149 """

150 from sympy.simplify.radsimp import split_surds

151

152 p = _mexpand(p)

153 if p.is_Number:

154 res = (p, S.Zero, S.Zero)

155 elif p.is_Add:

156 pargs = sorted(p.args, key=default_sort_key)

157 sqargs = [x**2 for x in pargs]

158 if all(sq.is_Rational and sq.is_positive for sq in sqargs):

159 r, b, a = split_surds(p)

160 res = a, b, r

161 return list(res)

162 # to make the process canonical, the argument is included in the tuple

163 # so when the max is selected, it will be the largest arg having a

164 # given depth

165 v = [(sqrt_depth(x), x, i) for i, x in enumerate(pargs)]

166 nmax = max(v, key=default_sort_key)

167 if nmax[0] == 0:

168 res = []

169 else:

170 # select r

171 depth, _, i = nmax

172 r = pargs.pop(i)

173 v.pop(i)

174 b = S.One

175 if r.is_Mul:

176 bv = []

177 rv = []

178 for x in r.args:

179 if sqrt_depth(x) < depth:

180 bv.append(x)

181 else:

182 rv.append(x)

183 b = Mul._from_args(bv)

184 r = Mul._from_args(rv)

185 # collect terms comtaining r

186 a1 = []

187 b1 = [b]

188 for x in v:

189 if x[0] < depth:

190 a1.append(x[1])

191 else:

192 x1 = x[1]

193 if x1 == r:

194 b1.append(1)

195 else:

196 if x1.is_Mul:

197 x1args = list(x1.args)

198 if r in x1args:

199 x1args.remove(r)

200 b1.append(Mul(*x1args))

201 else:

202 a1.append(x[1])

203 else:

204 a1.append(x[1])

205 a = Add(*a1)

206 b = Add(*b1)

207 res = (a, b, r**2)

208 else:

209 b, r = p.as_coeff_Mul()

210 if is_sqrt(r):

211 res = (S.Zero, b, r**2)

212 else:

213 res = []

214 return list(res)

215

216

217class SqrtdenestStopIteration(StopIteration):

218 pass

219

220

221def _sqrtdenest0(expr):

222 """Returns expr after denesting its arguments."""

223

224 if is_sqrt(expr):

225 n, d = expr.as_numer_denom()

226 if d is S.One: # n is a square root

227 if n.base.is_Add:

228 args = sorted(n.base.args, key=default_sort_key)

229 if len(args) > 2 and all((x**2).is_Integer for x in args):

230 try:

231 return _sqrtdenest_rec(n)

232 except SqrtdenestStopIteration:

233 pass

234 expr = sqrt(_mexpand(Add(*[_sqrtdenest0(x) for x in args])))

235 return _sqrtdenest1(expr)

236 else:

237 n, d = [_sqrtdenest0(i) for i in (n, d)]

238 return n/d

239

240 if isinstance(expr, Add):

241 cs = []

242 args = []

243 for arg in expr.args:

244 c, a = arg.as_coeff_Mul()

245 cs.append(c)

246 args.append(a)

247

248 if all(c.is_Rational for c in cs) and all(is_sqrt(arg) for arg in args):

249 return _sqrt_ratcomb(cs, args)

250

251 if isinstance(expr, Expr):

252 args = expr.args

253 if args:

254 return expr.func(*[_sqrtdenest0(a) for a in args])

255 return expr

256

257

258def _sqrtdenest_rec(expr):

259 """Helper that denests the square root of three or more surds.

260

261 Explanation

262 ===========

263

264 It returns the denested expression; if it cannot be denested it

265 throws SqrtdenestStopIteration

266

267 Algorithm: expr.base is in the extension Q_m = Q(sqrt(r_1),..,sqrt(r_k));

268 split expr.base = a + b*sqrt(r_k), where `a` and `b` are on

269 Q_(m-1) = Q(sqrt(r_1),..,sqrt(r_(k-1))); then a**2 - b**2*r_k is

270 on Q_(m-1); denest sqrt(a**2 - b**2*r_k) and so on.

271 See [1], section 6.

272

273 Examples

274 ========

275

276 >>> from sympy import sqrt

277 >>> from sympy.simplify.sqrtdenest import _sqrtdenest_rec

278 >>> _sqrtdenest_rec(sqrt(-72*sqrt(2) + 158*sqrt(5) + 498))

279 -sqrt(10) + sqrt(2) + 9 + 9*sqrt(5)

280 >>> w=-6*sqrt(55)-6*sqrt(35)-2*sqrt(22)-2*sqrt(14)+2*sqrt(77)+6*sqrt(10)+65

281 >>> _sqrtdenest_rec(sqrt(w))

282 -sqrt(11) - sqrt(7) + sqrt(2) + 3*sqrt(5)

283 """

284 from sympy.simplify.radsimp import radsimp, rad_rationalize, split_surds

285 if not expr.is_Pow:

286 return sqrtdenest(expr)

287 if expr.base < 0:

288 return sqrt(-1)*_sqrtdenest_rec(sqrt(-expr.base))

289 g, a, b = split_surds(expr.base)

290 a = a*sqrt(g)

291 if a < b:

292 a, b = b, a

293 c2 = _mexpand(a**2 - b**2)

294 if len(c2.args) > 2:

295 g, a1, b1 = split_surds(c2)

296 a1 = a1*sqrt(g)

297 if a1 < b1:

298 a1, b1 = b1, a1

299 c2_1 = _mexpand(a1**2 - b1**2)

300 c_1 = _sqrtdenest_rec(sqrt(c2_1))

301 d_1 = _sqrtdenest_rec(sqrt(a1 + c_1))

302 num, den = rad_rationalize(b1, d_1)

303 c = _mexpand(d_1/sqrt(2) + num/(den*sqrt(2)))

304 else:

305 c = _sqrtdenest1(sqrt(c2))

306

307 if sqrt_depth(c) > 1:

308 raise SqrtdenestStopIteration

309 ac = a + c

310 if len(ac.args) >= len(expr.args):

311 if count_ops(ac) >= count_ops(expr.base):

312 raise SqrtdenestStopIteration

313 d = sqrtdenest(sqrt(ac))

314 if sqrt_depth(d) > 1:

315 raise SqrtdenestStopIteration

316 num, den = rad_rationalize(b, d)

317 r = d/sqrt(2) + num/(den*sqrt(2))

318 r = radsimp(r)

319 return _mexpand(r)

320

321

322def _sqrtdenest1(expr, denester=True):

323 """Return denested expr after denesting with simpler methods or, that

324 failing, using the denester."""

325

326 from sympy.simplify.simplify import radsimp

327

328 if not is_sqrt(expr):

329 return expr

330

331 a = expr.base

332 if a.is_Atom:

333 return expr

334 val = _sqrt_match(a)

335 if not val:

336 return expr

337

338 a, b, r = val

339 # try a quick numeric denesting

340 d2 = _mexpand(a**2 - b**2*r)

341 if d2.is_Rational:

342 if d2.is_positive:

343 z = _sqrt_numeric_denest(a, b, r, d2)

344 if z is not None:

345 return z

346 else:

347 # fourth root case

348 # sqrtdenest(sqrt(3 + 2*sqrt(3))) =

349 # sqrt(2)*3**(1/4)/2 + sqrt(2)*3**(3/4)/2

350 dr2 = _mexpand(-d2*r)

351 dr = sqrt(dr2)

352 if dr.is_Rational:

353 z = _sqrt_numeric_denest(_mexpand(b*r), a, r, dr2)

354 if z is not None:

355 return z/root(r, 4)

356

357 else:

358 z = _sqrt_symbolic_denest(a, b, r)

359 if z is not None:

360 return z

361

362 if not denester or not is_algebraic(expr):

363 return expr

364

365 res = sqrt_biquadratic_denest(expr, a, b, r, d2)

366 if res:

367 return res

368

369 # now call to the denester

370 av0 = [a, b, r, d2]

371 z = _denester([radsimp(expr**2)], av0, 0, sqrt_depth(expr))[0]

372 if av0[1] is None:

373 return expr

374 if z is not None:

375 if sqrt_depth(z) == sqrt_depth(expr) and count_ops(z) > count_ops(expr):

376 return expr

377 return z

378 return expr

379

380

381def _sqrt_symbolic_denest(a, b, r):

382 """Given an expression, sqrt(a + b*sqrt(b)), return the denested

383 expression or None.

384

385 Explanation

386 ===========

387

388 If r = ra + rb*sqrt(rr), try replacing sqrt(rr) in ``a`` with

389 (y**2 - ra)/rb, and if the result is a quadratic, ca*y**2 + cb*y + cc, and

390 (cb + b)**2 - 4*ca*cc is 0, then sqrt(a + b*sqrt(r)) can be rewritten as

391 sqrt(ca*(sqrt(r) + (cb + b)/(2*ca))**2).

392

393 Examples

394 ========

395

396 >>> from sympy.simplify.sqrtdenest import _sqrt_symbolic_denest, sqrtdenest

397 >>> from sympy import sqrt, Symbol

398 >>> from sympy.abc import x

399

400 >>> a, b, r = 16 - 2*sqrt(29), 2, -10*sqrt(29) + 55

401 >>> _sqrt_symbolic_denest(a, b, r)

402 sqrt(11 - 2*sqrt(29)) + sqrt(5)

403

404 If the expression is numeric, it will be simplified:

405

406 >>> w = sqrt(sqrt(sqrt(3) + 1) + 1) + 1 + sqrt(2)

407 >>> sqrtdenest(sqrt((w**2).expand()))

408 1 + sqrt(2) + sqrt(1 + sqrt(1 + sqrt(3)))

409

410 Otherwise, it will only be simplified if assumptions allow:

411

412 >>> w = w.subs(sqrt(3), sqrt(x + 3))

413 >>> sqrtdenest(sqrt((w**2).expand()))

414 sqrt((sqrt(sqrt(sqrt(x + 3) + 1) + 1) + 1 + sqrt(2))**2)

415

416 Notice that the argument of the sqrt is a square. If x is made positive

417 then the sqrt of the square is resolved:

418

419 >>> _.subs(x, Symbol('x', positive=True))

420 sqrt(sqrt(sqrt(x + 3) + 1) + 1) + 1 + sqrt(2)

421 """

422

423 a, b, r = map(sympify, (a, b, r))

424 rval = _sqrt_match(r)

425 if not rval:

426 return None

427 ra, rb, rr = rval

428 if rb:

429 y = Dummy('y', positive=True)

430 try:

431 newa = Poly(a.subs(sqrt(rr), (y**2 - ra)/rb), y)

432 except PolynomialError:

433 return None

434 if newa.degree() == 2:

435 ca, cb, cc = newa.all_coeffs()

436 cb += b

437 if _mexpand(cb**2 - 4*ca*cc).equals(0):

438 z = sqrt(ca*(sqrt(r) + cb/(2*ca))**2)

439 if z.is_number:

440 z = _mexpand(Mul._from_args(z.as_content_primitive()))

441 return z

442

443

444def _sqrt_numeric_denest(a, b, r, d2):

445 r"""Helper that denest

446 $\sqrt{a + b \sqrt{r}}, d^2 = a^2 - b^2 r > 0$

447

448 If it cannot be denested, it returns ``None``.

449 """

450 d = sqrt(d2)

451 s = a + d

452 # sqrt_depth(res) <= sqrt_depth(s) + 1

453 # sqrt_depth(expr) = sqrt_depth(r) + 2

454 # there is denesting if sqrt_depth(s) + 1 < sqrt_depth(r) + 2

455 # if s**2 is Number there is a fourth root

456 if sqrt_depth(s) < sqrt_depth(r) + 1 or (s**2).is_Rational:

457 s1, s2 = sign(s), sign(b)

458 if s1 == s2 == -1:

459 s1 = s2 = 1

460 res = (s1 * sqrt(a + d) + s2 * sqrt(a - d)) * sqrt(2) / 2

461 return res.expand()

462

463

464def sqrt_biquadratic_denest(expr, a, b, r, d2):

465 """denest expr = sqrt(a + b*sqrt(r))

466 where a, b, r are linear combinations of square roots of

467 positive rationals on the rationals (SQRR) and r > 0, b != 0,

468 d2 = a**2 - b**2*r > 0

469

470 If it cannot denest it returns None.

471

472 Explanation

473 ===========

474

475 Search for a solution A of type SQRR of the biquadratic equation

476 4*A**4 - 4*a*A**2 + b**2*r = 0 (1)

477 sqd = sqrt(a**2 - b**2*r)

478 Choosing the sqrt to be positive, the possible solutions are

479 A = sqrt(a/2 +/- sqd/2)

480 Since a, b, r are SQRR, then a**2 - b**2*r is a SQRR,

481 so if sqd can be denested, it is done by

482 _sqrtdenest_rec, and the result is a SQRR.

483 Similarly for A.

484 Examples of solutions (in both cases a and sqd are positive):

485

486 Example of expr with solution sqrt(a/2 + sqd/2) but not

487 solution sqrt(a/2 - sqd/2):

488 expr = sqrt(-sqrt(15) - sqrt(2)*sqrt(-sqrt(5) + 5) - sqrt(3) + 8)

489 a = -sqrt(15) - sqrt(3) + 8; sqd = -2*sqrt(5) - 2 + 4*sqrt(3)

490

491 Example of expr with solution sqrt(a/2 - sqd/2) but not

492 solution sqrt(a/2 + sqd/2):

493 w = 2 + r2 + r3 + (1 + r3)*sqrt(2 + r2 + 5*r3)

494 expr = sqrt((w**2).expand())

495 a = 4*sqrt(6) + 8*sqrt(2) + 47 + 28*sqrt(3)

496 sqd = 29 + 20*sqrt(3)

497

498 Define B = b/2*A; eq.(1) implies a = A**2 + B**2*r; then

499 expr**2 = a + b*sqrt(r) = (A + B*sqrt(r))**2

500

501 Examples

502 ========

503

504 >>> from sympy import sqrt

505 >>> from sympy.simplify.sqrtdenest import _sqrt_match, sqrt_biquadratic_denest

506 >>> z = sqrt((2*sqrt(2) + 4)*sqrt(2 + sqrt(2)) + 5*sqrt(2) + 8)

507 >>> a, b, r = _sqrt_match(z**2)

508 >>> d2 = a**2 - b**2*r

509 >>> sqrt_biquadratic_denest(z, a, b, r, d2)

510 sqrt(2) + sqrt(sqrt(2) + 2) + 2

511 """

512 from sympy.simplify.radsimp import radsimp, rad_rationalize

513 if r <= 0 or d2 < 0 or not b or sqrt_depth(expr.base) < 2:

514 return None

515 for x in (a, b, r):

516 for y in x.args:

517 y2 = y**2

518 if not y2.is_Integer or not y2.is_positive:

519 return None

520 sqd = _mexpand(sqrtdenest(sqrt(radsimp(d2))))

521 if sqrt_depth(sqd) > 1:

522 return None

523 x1, x2 = [a/2 + sqd/2, a/2 - sqd/2]

524 # look for a solution A with depth 1

525 for x in (x1, x2):

526 A = sqrtdenest(sqrt(x))

527 if sqrt_depth(A) > 1:

528 continue

529 Bn, Bd = rad_rationalize(b, _mexpand(2*A))

530 B = Bn/Bd

531 z = A + B*sqrt(r)

532 if z < 0:

533 z = -z

534 return _mexpand(z)

535 return None

536

537

538def _denester(nested, av0, h, max_depth_level):

539 """Denests a list of expressions that contain nested square roots.

540

541 Explanation

542 ===========

543

544 Algorithm based on <http://www.almaden.ibm.com/cs/people/fagin/symb85.pdf>.

545

546 It is assumed that all of the elements of 'nested' share the same

547 bottom-level radicand. (This is stated in the paper, on page 177, in

548 the paragraph immediately preceding the algorithm.)

549

550 When evaluating all of the arguments in parallel, the bottom-level

551 radicand only needs to be denested once. This means that calling

552 _denester with x arguments results in a recursive invocation with x+1

553 arguments; hence _denester has polynomial complexity.

554

555 However, if the arguments were evaluated separately, each call would

556 result in two recursive invocations, and the algorithm would have

557 exponential complexity.

558

559 This is discussed in the paper in the middle paragraph of page 179.

560 """

561 from sympy.simplify.simplify import radsimp

562 if h > max_depth_level:

563 return None, None

564 if av0[1] is None:

565 return None, None

566 if (av0[0] is None and

567 all(n.is_Number for n in nested)): # no arguments are nested

568 for f in _subsets(len(nested)): # test subset 'f' of nested

569 p = _mexpand(Mul(*[nested[i] for i in range(len(f)) if f[i]]))

570 if f.count(1) > 1 and f[-1]:

571 p = -p

572 sqp = sqrt(p)

573 if sqp.is_Rational:

574 return sqp, f # got a perfect square so return its square root.

575 # Otherwise, return the radicand from the previous invocation.

576 return sqrt(nested[-1]), [0]*len(nested)

577 else:

578 R = None

579 if av0[0] is not None:

580 values = [av0[:2]]

581 R = av0[2]

582 nested2 = [av0[3], R]

583 av0[0] = None

584 else:

585 values = list(filter(None, [_sqrt_match(expr) for expr in nested]))

586 for v in values:

587 if v[2]: # Since if b=0, r is not defined

588 if R is not None:

589 if R != v[2]:

590 av0[1] = None

591 return None, None

592 else:

593 R = v[2]

594 if R is None:

595 # return the radicand from the previous invocation

596 return sqrt(nested[-1]), [0]*len(nested)

597 nested2 = [_mexpand(v[0]**2) -

598 _mexpand(R*v[1]**2) for v in values] + [R]

599 d, f = _denester(nested2, av0, h + 1, max_depth_level)

600 if not f:

601 return None, None

602 if not any(f[i] for i in range(len(nested))):

603 v = values[-1]

604 return sqrt(v[0] + _mexpand(v[1]*d)), f

605 else:

606 p = Mul(*[nested[i] for i in range(len(nested)) if f[i]])

607 v = _sqrt_match(p)

608 if 1 in f and f.index(1) < len(nested) - 1 and f[len(nested) - 1]:

609 v[0] = -v[0]

610 v[1] = -v[1]

611 if not f[len(nested)]: # Solution denests with square roots

612 vad = _mexpand(v[0] + d)

613 if vad <= 0:

614 # return the radicand from the previous invocation.

615 return sqrt(nested[-1]), [0]*len(nested)

616 if not(sqrt_depth(vad) <= sqrt_depth(R) + 1 or

617 (vad**2).is_Number):

618 av0[1] = None

619 return None, None

620

621 sqvad = _sqrtdenest1(sqrt(vad), denester=False)

622 if not (sqrt_depth(sqvad) <= sqrt_depth(R) + 1):

623 av0[1] = None

624 return None, None

625 sqvad1 = radsimp(1/sqvad)

626 res = _mexpand(sqvad/sqrt(2) + (v[1]*sqrt(R)*sqvad1/sqrt(2)))

627 return res, f

628

629 # sign(v[1])*sqrt(_mexpand(v[1]**2*R*vad1/2))), f

630 else: # Solution requires a fourth root

631 s2 = _mexpand(v[1]*R) + d

632 if s2 <= 0:

633 return sqrt(nested[-1]), [0]*len(nested)

634 FR, s = root(_mexpand(R), 4), sqrt(s2)

635 return _mexpand(s/(sqrt(2)*FR) + v[0]*FR/(sqrt(2)*s)), f

636

637

638def _sqrt_ratcomb(cs, args):

639 """Denest rational combinations of radicals.

640

641 Based on section 5 of [1].

642

643 Examples

644 ========

645

646 >>> from sympy import sqrt

647 >>> from sympy.simplify.sqrtdenest import sqrtdenest

648 >>> z = sqrt(1+sqrt(3)) + sqrt(3+3*sqrt(3)) - sqrt(10+6*sqrt(3))

649 >>> sqrtdenest(z)

650 0

651 """

652 from sympy.simplify.radsimp import radsimp

653

654 # check if there exists a pair of sqrt that can be denested

655 def find(a):

656 n = len(a)

657 for i in range(n - 1):

658 for j in range(i + 1, n):

659 s1 = a[i].base

660 s2 = a[j].base

661 p = _mexpand(s1 * s2)

662 s = sqrtdenest(sqrt(p))

663 if s != sqrt(p):

664 return s, i, j

665

666 indices = find(args)

667 if indices is None:

668 return Add(*[c * arg for c, arg in zip(cs, args)])

669

670 s, i1, i2 = indices

671

672 c2 = cs.pop(i2)

673 args.pop(i2)

674 a1 = args[i1]

675

676 # replace a2 by s/a1

677 cs[i1] += radsimp(c2 * s / a1.base)

678

679 return _sqrt_ratcomb(cs, args)