Coverage for /usr/lib/python3/dist-packages/sympy/functions/special/delta

1from sympy.core import S, diff

2from sympy.core.function import Function, ArgumentIndexError

3from sympy.core.logic import fuzzy_not

4from sympy.core.relational import Eq, Ne

5from sympy.functions.elementary.complexes import im, sign

6from sympy.functions.elementary.piecewise import Piecewise

7from sympy.polys.polyerrors import PolynomialError

8from sympy.polys.polyroots import roots

9from sympy.utilities.misc import filldedent

12###############################################################################

13################################ DELTA FUNCTION ###############################

14###############################################################################

17class DiracDelta(Function):

18 r"""

19 The DiracDelta function and its derivatives.

21 Explanation

22 ===========

24 DiracDelta is not an ordinary function. It can be rigorously defined either

25 as a distribution or as a measure.

27 DiracDelta only makes sense in definite integrals, and in particular,

28 integrals of the form ``Integral(f(x)*DiracDelta(x - x0), (x, a, b))``,

29 where it equals ``f(x0)`` if ``a <= x0 <= b`` and ``0`` otherwise. Formally,

30 DiracDelta acts in some ways like a function that is ``0`` everywhere except

31 at ``0``, but in many ways it also does not. It can often be useful to treat

32 DiracDelta in formal ways, building up and manipulating expressions with

33 delta functions (which may eventually be integrated), but care must be taken

34 to not treat it as a real function. SymPy's ``oo`` is similar. It only

35 truly makes sense formally in certain contexts (such as integration limits),

36 but SymPy allows its use everywhere, and it tries to be consistent with

37 operations on it (like ``1/oo``), but it is easy to get into trouble and get

38 wrong results if ``oo`` is treated too much like a number. Similarly, if

39 DiracDelta is treated too much like a function, it is easy to get wrong or

40 nonsensical results.

42 DiracDelta function has the following properties:

44 1) $\frac{d}{d x} \theta(x) = \delta(x)$

45 2) $\int_{-\infty}^\infty \delta(x - a)f(x)\, dx = f(a)$ and $\int_{a-

46 \epsilon}^{a+\epsilon} \delta(x - a)f(x)\, dx = f(a)$

47 3) $\delta(x) = 0$ for all $x \neq 0$

48 4) $\delta(g(x)) = \sum_i \frac{\delta(x - x_i)}{\|g'(x_i)\|}$ where $x_i$

49 are the roots of $g$

50 5) $\delta(-x) = \delta(x)$

52 Derivatives of ``k``-th order of DiracDelta have the following properties:

54 6) $\delta(x, k) = 0$ for all $x \neq 0$

55 7) $\delta(-x, k) = -\delta(x, k)$ for odd $k$

56 8) $\delta(-x, k) = \delta(x, k)$ for even $k$

58 Examples

59 ========

61 >>> from sympy import DiracDelta, diff, pi

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

64 >>> DiracDelta(x)

65 DiracDelta(x)

66 >>> DiracDelta(1)

67 0

68 >>> DiracDelta(-1)

69 0

70 >>> DiracDelta(pi)

71 0

72 >>> DiracDelta(x - 4).subs(x, 4)

73 DiracDelta(0)

74 >>> diff(DiracDelta(x))

75 DiracDelta(x, 1)

76 >>> diff(DiracDelta(x - 1), x, 2)

77 DiracDelta(x - 1, 2)

78 >>> diff(DiracDelta(x**2 - 1), x, 2)

79 2*(2*x**2*DiracDelta(x**2 - 1, 2) + DiracDelta(x**2 - 1, 1))

80 >>> DiracDelta(3*x).is_simple(x)

81 True

82 >>> DiracDelta(x**2).is_simple(x)

83 False

84 >>> DiracDelta((x**2 - 1)*y).expand(diracdelta=True, wrt=x)

85 DiracDelta(x - 1)/(2*Abs(y)) + DiracDelta(x + 1)/(2*Abs(y))

87 See Also

88 ========

90 Heaviside

91 sympy.simplify.simplify.simplify, is_simple

92 sympy.functions.special.tensor_functions.KroneckerDelta

94 References

95 ==========

97 .. [1] https://mathworld.wolfram.com/DeltaFunction.html

99 """

100

101 is_real = True

102

103 def fdiff(self, argindex=1):

104 """

105 Returns the first derivative of a DiracDelta Function.

106

107 Explanation

108 ===========

109

110 The difference between ``diff()`` and ``fdiff()`` is: ``diff()`` is the

111 user-level function and ``fdiff()`` is an object method. ``fdiff()`` is

112 a convenience method available in the ``Function`` class. It returns

113 the derivative of the function without considering the chain rule.

114 ``diff(function, x)`` calls ``Function._eval_derivative`` which in turn

115 calls ``fdiff()`` internally to compute the derivative of the function.

116

117 Examples

118 ========

119

120 >>> from sympy import DiracDelta, diff

121 >>> from sympy.abc import x

122

123 >>> DiracDelta(x).fdiff()

124 DiracDelta(x, 1)

125

126 >>> DiracDelta(x, 1).fdiff()

127 DiracDelta(x, 2)

128

129 >>> DiracDelta(x**2 - 1).fdiff()

130 DiracDelta(x**2 - 1, 1)

131

132 >>> diff(DiracDelta(x, 1)).fdiff()

133 DiracDelta(x, 3)

134

135 Parameters

136 ==========

137

138 argindex : integer

139 degree of derivative

140

141 """

142 if argindex == 1:

143 #I didn't know if there is a better way to handle default arguments

144 k = 0

145 if len(self.args) > 1:

146 k = self.args[1]

147 return self.func(self.args[0], k + 1)

148 else:

149 raise ArgumentIndexError(self, argindex)

150

151 @classmethod

152 def eval(cls, arg, k=S.Zero):

153 """

154 Returns a simplified form or a value of DiracDelta depending on the

155 argument passed by the DiracDelta object.

156

157 Explanation

158 ===========

159

160 The ``eval()`` method is automatically called when the ``DiracDelta``

161 class is about to be instantiated and it returns either some simplified

162 instance or the unevaluated instance depending on the argument passed.

163 In other words, ``eval()`` method is not needed to be called explicitly,

164 it is being called and evaluated once the object is called.

165

166 Examples

167 ========

168

169 >>> from sympy import DiracDelta, S

170 >>> from sympy.abc import x

171

172 >>> DiracDelta(x)

173 DiracDelta(x)

174

175 >>> DiracDelta(-x, 1)

176 -DiracDelta(x, 1)

177

178 >>> DiracDelta(1)

179 0

180

181 >>> DiracDelta(5, 1)

182 0

183

184 >>> DiracDelta(0)

185 DiracDelta(0)

186

187 >>> DiracDelta(-1)

188 0

189

190 >>> DiracDelta(S.NaN)

191 nan

192

193 >>> DiracDelta(x - 100).subs(x, 5)

194 0

195

196 >>> DiracDelta(x - 100).subs(x, 100)

197 DiracDelta(0)

198

199 Parameters

200 ==========

201

202 k : integer

203 order of derivative

204

205 arg : argument passed to DiracDelta

206

207 """

208 if not k.is_Integer or k.is_negative:

209 raise ValueError("Error: the second argument of DiracDelta must be \

210 a non-negative integer, %s given instead." % (k,))

211 if arg is S.NaN:

212 return S.NaN

213 if arg.is_nonzero:

214 return S.Zero

215 if fuzzy_not(im(arg).is_zero):

216 raise ValueError(filldedent('''

217 Function defined only for Real Values.

218 Complex part: %s found in %s .''' % (

219 repr(im(arg)), repr(arg))))

220 c, nc = arg.args_cnc()

221 if c and c[0] is S.NegativeOne:

222 # keep this fast and simple instead of using

223 # could_extract_minus_sign

224 if k.is_odd:

225 return -cls(-arg, k)

226 elif k.is_even:

227 return cls(-arg, k) if k else cls(-arg)

228 elif k.is_zero:

229 return cls(arg, evaluate=False)

230

231 def _eval_expand_diracdelta(self, **hints):

232 """

233 Compute a simplified representation of the function using

234 property number 4. Pass ``wrt`` as a hint to expand the expression

235 with respect to a particular variable.

236

237 Explanation

238 ===========

239

240 ``wrt`` is:

241

242 - a variable with respect to which a DiracDelta expression will

243 get expanded.

244

245 Examples

246 ========

247

248 >>> from sympy import DiracDelta

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

250

251 >>> DiracDelta(x*y).expand(diracdelta=True, wrt=x)

252 DiracDelta(x)/Abs(y)

253 >>> DiracDelta(x*y).expand(diracdelta=True, wrt=y)

254 DiracDelta(y)/Abs(x)

255

256 >>> DiracDelta(x**2 + x - 2).expand(diracdelta=True, wrt=x)

257 DiracDelta(x - 1)/3 + DiracDelta(x + 2)/3

258

259 See Also

260 ========

261

262 is_simple, Diracdelta

263

264 """

265 wrt = hints.get('wrt', None)

266 if wrt is None:

267 free = self.free_symbols

268 if len(free) == 1:

269 wrt = free.pop()

270 else:

271 raise TypeError(filldedent('''

272 When there is more than 1 free symbol or variable in the expression,

273 the 'wrt' keyword is required as a hint to expand when using the

274 DiracDelta hint.'''))

275

276 if not self.args[0].has(wrt) or (len(self.args) > 1 and self.args[1] != 0 ):

277 return self

278 try:

279 argroots = roots(self.args[0], wrt)

280 result = 0

281 valid = True

282 darg = abs(diff(self.args[0], wrt))

283 for r, m in argroots.items():

284 if r.is_real is not False and m == 1:

285 result += self.func(wrt - r)/darg.subs(wrt, r)

286 else:

287 # don't handle non-real and if m != 1 then

288 # a polynomial will have a zero in the derivative (darg)

289 # at r

290 valid = False

291 break

292 if valid:

293 return result

294 except PolynomialError:

295 pass

296 return self

297

298 def is_simple(self, x):

299 """

300 Tells whether the argument(args[0]) of DiracDelta is a linear

301 expression in *x*.

302

303 Examples

304 ========

305

306 >>> from sympy import DiracDelta, cos

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

308

309 >>> DiracDelta(x*y).is_simple(x)

310 True

311 >>> DiracDelta(x*y).is_simple(y)

312 True

313

314 >>> DiracDelta(x**2 + x - 2).is_simple(x)

315 False

316

317 >>> DiracDelta(cos(x)).is_simple(x)

318 False

319

320 Parameters

321 ==========

322

323 x : can be a symbol

324

325 See Also

326 ========

327

328 sympy.simplify.simplify.simplify, DiracDelta

329

330 """

331 p = self.args[0].as_poly(x)

332 if p:

333 return p.degree() == 1

334 return False

335

336 def _eval_rewrite_as_Piecewise(self, *args, **kwargs):

337 """

338 Represents DiracDelta in a piecewise form.

339

340 Examples

341 ========

342

343 >>> from sympy import DiracDelta, Piecewise, Symbol

344 >>> x = Symbol('x')

345

346 >>> DiracDelta(x).rewrite(Piecewise)

347 Piecewise((DiracDelta(0), Eq(x, 0)), (0, True))

348

349 >>> DiracDelta(x - 5).rewrite(Piecewise)

350 Piecewise((DiracDelta(0), Eq(x, 5)), (0, True))

351

352 >>> DiracDelta(x**2 - 5).rewrite(Piecewise)

353 Piecewise((DiracDelta(0), Eq(x**2, 5)), (0, True))

354

355 >>> DiracDelta(x - 5, 4).rewrite(Piecewise)

356 DiracDelta(x - 5, 4)

357

358 """

359 if len(args) == 1:

360 return Piecewise((DiracDelta(0), Eq(args[0], 0)), (0, True))

361

362 def _eval_rewrite_as_SingularityFunction(self, *args, **kwargs):

363 """

364 Returns the DiracDelta expression written in the form of Singularity

365 Functions.

366

367 """

368 from sympy.solvers import solve

369 from sympy.functions.special.singularity_functions import SingularityFunction

370 if self == DiracDelta(0):

371 return SingularityFunction(0, 0, -1)

372 if self == DiracDelta(0, 1):

373 return SingularityFunction(0, 0, -2)

374 free = self.free_symbols

375 if len(free) == 1:

376 x = (free.pop())

377 if len(args) == 1:

378 return SingularityFunction(x, solve(args[0], x)[0], -1)

379 return SingularityFunction(x, solve(args[0], x)[0], -args[1] - 1)

380 else:

381 # I don't know how to handle the case for DiracDelta expressions

382 # having arguments with more than one variable.

383 raise TypeError(filldedent('''

384 rewrite(SingularityFunction) does not support

385 arguments with more that one variable.'''))

386

387

388###############################################################################

389############################## HEAVISIDE FUNCTION #############################

390###############################################################################

391

392

393class Heaviside(Function):

394 r"""

395 Heaviside step function.

396

397 Explanation

398 ===========

399

400 The Heaviside step function has the following properties:

401

402 1) $\frac{d}{d x} \theta(x) = \delta(x)$

403 2) $\theta(x) = \begin{cases} 0 & \text{for}\: x < 0 \\ \frac{1}{2} &

404 \text{for}\: x = 0 \\1 & \text{for}\: x > 0 \end{cases}$

405 3) $\frac{d}{d x} \max(x, 0) = \theta(x)$

406

407 Heaviside(x) is printed as $\theta(x)$ with the SymPy LaTeX printer.

408

409 The value at 0 is set differently in different fields. SymPy uses 1/2,

410 which is a convention from electronics and signal processing, and is

411 consistent with solving improper integrals by Fourier transform and

412 convolution.

413

414 To specify a different value of Heaviside at ``x=0``, a second argument

415 can be given. Using ``Heaviside(x, nan)`` gives an expression that will

416 evaluate to nan for x=0.

417

418 .. versionchanged:: 1.9 ``Heaviside(0)`` now returns 1/2 (before: undefined)

419

420 Examples

421 ========

422

423 >>> from sympy import Heaviside, nan

424 >>> from sympy.abc import x

425 >>> Heaviside(9)

426 1

427 >>> Heaviside(-9)

428 0

429 >>> Heaviside(0)

430 1/2

431 >>> Heaviside(0, nan)

432 nan

433 >>> (Heaviside(x) + 1).replace(Heaviside(x), Heaviside(x, 1))

434 Heaviside(x, 1) + 1

435

436 See Also

437 ========

438

439 DiracDelta

440

441 References

442 ==========

443

444 .. [1] https://mathworld.wolfram.com/HeavisideStepFunction.html

445 .. [2] https://dlmf.nist.gov/1.16#iv

446

447 """

448

449 is_real = True

450

451 def fdiff(self, argindex=1):

452 """

453 Returns the first derivative of a Heaviside Function.

454

455 Examples

456 ========

457

458 >>> from sympy import Heaviside, diff

459 >>> from sympy.abc import x

460

461 >>> Heaviside(x).fdiff()

462 DiracDelta(x)

463

464 >>> Heaviside(x**2 - 1).fdiff()

465 DiracDelta(x**2 - 1)

466

467 >>> diff(Heaviside(x)).fdiff()

468 DiracDelta(x, 1)

469

470 Parameters

471 ==========

472

473 argindex : integer

474 order of derivative

475

476 """

477 if argindex == 1:

478 return DiracDelta(self.args[0])

479 else:

480 raise ArgumentIndexError(self, argindex)

481

482 def __new__(cls, arg, H0=S.Half, **options):

483 if isinstance(H0, Heaviside) and len(H0.args) == 1:

484 H0 = S.Half

485 return super(cls, cls).__new__(cls, arg, H0, **options)

486

487 @property

488 def pargs(self):

489 """Args without default S.Half"""

490 args = self.args

491 if args[1] is S.Half:

492 args = args[:1]

493 return args

494

495 @classmethod

496 def eval(cls, arg, H0=S.Half):

497 """

498 Returns a simplified form or a value of Heaviside depending on the

499 argument passed by the Heaviside object.

500

501 Explanation

502 ===========

503

504 The ``eval()`` method is automatically called when the ``Heaviside``

505 class is about to be instantiated and it returns either some simplified

506 instance or the unevaluated instance depending on the argument passed.

507 In other words, ``eval()`` method is not needed to be called explicitly,

508 it is being called and evaluated once the object is called.

509

510 Examples

511 ========

512

513 >>> from sympy import Heaviside, S

514 >>> from sympy.abc import x

515

516 >>> Heaviside(x)

517 Heaviside(x)

518

519 >>> Heaviside(19)

520 1

521

522 >>> Heaviside(0)

523 1/2

524

525 >>> Heaviside(0, 1)

526 1

527

528 >>> Heaviside(-5)

529 0

530

531 >>> Heaviside(S.NaN)

532 nan

533

534 >>> Heaviside(x - 100).subs(x, 5)

535 0

536

537 >>> Heaviside(x - 100).subs(x, 105)

538 1

539

540 Parameters

541 ==========

542

543 arg : argument passed by Heaviside object

544

545 H0 : value of Heaviside(0)

546

547 """

548 if arg.is_extended_negative:

549 return S.Zero

550 elif arg.is_extended_positive:

551 return S.One

552 elif arg.is_zero:

553 return H0

554 elif arg is S.NaN:

555 return S.NaN

556 elif fuzzy_not(im(arg).is_zero):

557 raise ValueError("Function defined only for Real Values. Complex part: %s found in %s ." % (repr(im(arg)), repr(arg)) )

558

559 def _eval_rewrite_as_Piecewise(self, arg, H0=None, **kwargs):

560 """

561 Represents Heaviside in a Piecewise form.

562

563 Examples

564 ========

565

566 >>> from sympy import Heaviside, Piecewise, Symbol, nan

567 >>> x = Symbol('x')

568

569 >>> Heaviside(x).rewrite(Piecewise)

570 Piecewise((0, x < 0), (1/2, Eq(x, 0)), (1, True))

571

572 >>> Heaviside(x,nan).rewrite(Piecewise)

573 Piecewise((0, x < 0), (nan, Eq(x, 0)), (1, True))

574

575 >>> Heaviside(x - 5).rewrite(Piecewise)

576 Piecewise((0, x < 5), (1/2, Eq(x, 5)), (1, True))

577

578 >>> Heaviside(x**2 - 1).rewrite(Piecewise)

579 Piecewise((0, x**2 < 1), (1/2, Eq(x**2, 1)), (1, True))

580

581 """

582 if H0 == 0:

583 return Piecewise((0, arg <= 0), (1, True))

584 if H0 == 1:

585 return Piecewise((0, arg < 0), (1, True))

586 return Piecewise((0, arg < 0), (H0, Eq(arg, 0)), (1, True))

587

588 def _eval_rewrite_as_sign(self, arg, H0=S.Half, **kwargs):

589 """

590 Represents the Heaviside function in the form of sign function.

591

592 Explanation

593 ===========

594

595 The value of Heaviside(0) must be 1/2 for rewriting as sign to be

596 strictly equivalent. For easier usage, we also allow this rewriting

597 when Heaviside(0) is undefined.

598

599 Examples

600 ========

601

602 >>> from sympy import Heaviside, Symbol, sign, nan

603 >>> x = Symbol('x', real=True)

604 >>> y = Symbol('y')

605

606 >>> Heaviside(x).rewrite(sign)

607 sign(x)/2 + 1/2

608

609 >>> Heaviside(x, 0).rewrite(sign)

610 Piecewise((sign(x)/2 + 1/2, Ne(x, 0)), (0, True))

611

612 >>> Heaviside(x, nan).rewrite(sign)

613 Piecewise((sign(x)/2 + 1/2, Ne(x, 0)), (nan, True))

614

615 >>> Heaviside(x - 2).rewrite(sign)

616 sign(x - 2)/2 + 1/2

617

618 >>> Heaviside(x**2 - 2*x + 1).rewrite(sign)

619 sign(x**2 - 2*x + 1)/2 + 1/2

620

621 >>> Heaviside(y).rewrite(sign)

622 Heaviside(y)

623

624 >>> Heaviside(y**2 - 2*y + 1).rewrite(sign)

625 Heaviside(y**2 - 2*y + 1)

626

627 See Also

628 ========

629

630 sign

631

632 """

633 if arg.is_extended_real:

634 pw1 = Piecewise(

635 ((sign(arg) + 1)/2, Ne(arg, 0)),

636 (Heaviside(0, H0=H0), True))

637 pw2 = Piecewise(

638 ((sign(arg) + 1)/2, Eq(Heaviside(0, H0=H0), S.Half)),

639 (pw1, True))

640 return pw2

641

642 def _eval_rewrite_as_SingularityFunction(self, args, H0=S.Half, **kwargs):

643 """

644 Returns the Heaviside expression written in the form of Singularity

645 Functions.

646

647 """

648 from sympy.solvers import solve

649 from sympy.functions.special.singularity_functions import SingularityFunction

650 if self == Heaviside(0):

651 return SingularityFunction(0, 0, 0)

652 free = self.free_symbols

653 if len(free) == 1:

654 x = (free.pop())

655 return SingularityFunction(x, solve(args, x)[0], 0)

656 # TODO

657 # ((x - 5)**3*Heaviside(x - 5)).rewrite(SingularityFunction) should output

658 # SingularityFunction(x, 5, 0) instead of (x - 5)**3*SingularityFunction(x, 5, 0)

659 else:

660 # I don't know how to handle the case for Heaviside expressions

661 # having arguments with more than one variable.

662 raise TypeError(filldedent('''

663 rewrite(SingularityFunction) does not

664 support arguments with more that one variable.'''))

Coverage for /usr/lib/python3/dist-packages/sympy/functions/special/delta_functions.py: 22%

131 statements