Dirac.mw
Functions with discontinuities. Dirac and continuous approximations.
Clara MASSE, John MASSE, François OLLIVIER
The derivative of the Heaviside function is the Dirac function
| > |
 |
 |
(1) |
The numerical evaluation of this function using "piecewise" is the same, but its derivative different.
| > |
 |
 |
(2) |
Using Heaviside, one may invert "int" and "diff" where the function is differentiable.
| > |
 |
 |
(3) |
| > |
 |
 |
(4) |
This no longer works with to successive integrations and derivations.
| > |
 |
| > |
 |
 |
(6) |
In a numerical environment, one may use some approximation of Heaviside, e.g. with "erf" function.
| > |
 |
 |
(7) |
We see that one needs to play with the value of parameter "b" in order to obtained the best precision.
| > |
![Digits := 10; 1; resu2 := dsolve([diff(Z(y), y) = diff(`*`(`+`(erf(`*`(`+`(y, `-`(a)), `*`(b))), 1), `/`(1, 2)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [Z(y), Z1(y)], numeric, param...](images/Dirac_28.gif)
![Digits := 10; 1; resu2 := dsolve([diff(Z(y), y) = diff(`*`(`+`(erf(`*`(`+`(y, `-`(a)), `*`(b))), 1), `/`(1, 2)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [Z(y), Z1(y)], numeric, param...](images/Dirac_29.gif) |
 |
 |
(9) |
The situation remains unchanged, even with the stiff option and a great number of digits.
We may observe the same behaviour with "arctan".
![Digits := 20; 1; evalf(subs([a = 0, b = 1, x = 1], resu)); 1; evalf(subs([a = 0, b = 5, x = 1], resu)); 1; evalf(subs([a = 0, b = 10, x = 1], resu)); 1; evalf(subs([a = 0, b = 15, x = 1], resu)); 1; e...](images/Dirac_18.gif){kind=small}
![Digits := 20; 1; evalf(subs([a = 0, b = 1, x = 1], resu)); 1; evalf(subs([a = 0, b = 5, x = 1], resu)); 1; evalf(subs([a = 0, b = 10, x = 1], resu)); 1; evalf(subs([a = 0, b = 15, x = 1], resu)); 1; e...](images/Dirac_19.gif){kind=small}
![Digits := 20; 1; evalf(subs([a = 0, b = 1, x = 1], resu)); 1; evalf(subs([a = 0, b = 5, x = 1], resu)); 1; evalf(subs([a = 0, b = 10, x = 1], resu)); 1; evalf(subs([a = 0, b = 15, x = 1], resu)); 1; e...](images/Dirac_20.gif){kind=small}
![resu2(parameters = [b = 1., a = 0.]); 1; resu2(1.); 1; resu2(parameters = [b = 5., a = 0.]); 1; resu2(1.); 1; resu2(parameters = [b = 10., a = 0.]); 1; resu2(1.); 1; resu2(parameters = [b = 15., a = 0...](images/Dirac_32.gif){kind=small}
![resu2(parameters = [b = 1., a = 0.]); 1; resu2(1.); 1; resu2(parameters = [b = 5., a = 0.]); 1; resu2(1.); 1; resu2(parameters = [b = 10., a = 0.]); 1; resu2(1.); 1; resu2(parameters = [b = 15., a = 0...](images/Dirac_33.gif){kind=small}
![resu2(parameters = [b = 1., a = 0.]); 1; resu2(1.); 1; resu2(parameters = [b = 5., a = 0.]); 1; resu2(1.); 1; resu2(parameters = [b = 10., a = 0.]); 1; resu2(1.); 1; resu2(parameters = [b = 15., a = 0...](images/Dirac_34.gif){kind=small}
![[a = 0., b = 1.]](images/Dirac_35.gif){kind=small}
![[y = 1., Z(y) = HFloat(7.707120112899687e-8), Z1(y) = HFloat(0.4275932988977108)]](images/Dirac_36.gif){kind=small}
![[a = 0., b = 5.]](images/Dirac_37.gif){kind=small}
![[y = 1., Z(y) = HFloat(-3.063579152425269e-10), Z1(y) = HFloat(1.0000000026453448)]](images/Dirac_38.gif){kind=small}
![[a = 0., b = 10.]](images/Dirac_39.gif){kind=small}
![[y = 1., Z(y) = HFloat(6.118405317357082e-9), Z1(y) = HFloat(1.000000006523538)]](images/Dirac_40.gif){kind=small}
![[a = 0., b = 15.]](images/Dirac_41.gif){kind=small}
![[y = 1., Z(y) = HFloat(6.8122416690601316e-9), Z1(y) = HFloat(1.0000000070661805)]](images/Dirac_42.gif){kind=small}
![[a = 0., b = 20.]](images/Dirac_43.gif){kind=small}
![[y = 1., Z(y) = HFloat(6.542206472095037e-9), Z1(y) = HFloat(1.0000000066936452)]](images/Dirac_44.gif){kind=small}
![[a = 0., b = 25.]](images/Dirac_45.gif){kind=small}
![[y = 1., Z(y) = HFloat(7.0008734049088085e-9), Z1(y) = HFloat(1.000000006903479)]](images/Dirac_46.gif){kind=small}
![Digits := 30; 1; resu3 := dsolve([diff(Z(y), y) = diff(`*`(`+`(erf(`*`(`+`(y, `-`(a)), `*`(b))), 1), `/`(1, 2)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [Z(y), Z1(y)], numeric, param...](images/Dirac_47.gif){kind=small}
![Digits := 30; 1; resu3 := dsolve([diff(Z(y), y) = diff(`*`(`+`(erf(`*`(`+`(y, `-`(a)), `*`(b))), 1), `/`(1, 2)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [Z(y), Z1(y)], numeric, param...](images/Dirac_48.gif){kind=small}
![Digits := 30; 1; resu3 := dsolve([diff(Z(y), y) = diff(`*`(`+`(erf(`*`(`+`(y, `-`(a)), `*`(b))), 1), `/`(1, 2)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [Z(y), Z1(y)], numeric, param...](images/Dirac_49.gif){kind=small}
![Digits := 30; 1; resu3 := dsolve([diff(Z(y), y) = diff(`*`(`+`(erf(`*`(`+`(y, `-`(a)), `*`(b))), 1), `/`(1, 2)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [Z(y), Z1(y)], numeric, param...](images/Dirac_50.gif){kind=small}
![Digits := 30; 1; resu3 := dsolve([diff(Z(y), y) = diff(`*`(`+`(erf(`*`(`+`(y, `-`(a)), `*`(b))), 1), `/`(1, 2)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [Z(y), Z1(y)], numeric, param...](images/Dirac_51.gif){kind=small}
![[a = 0., b = 1.]](images/Dirac_54.gif){kind=small}
![[y = 1., Z(y) = 0.77071199756160432299961122e-7, Z1(y) = .427593298897711602227789841677]](images/Dirac_55.gif){kind=small}
![[a = 0., b = 5.]](images/Dirac_56.gif){kind=small}
![[y = 1., Z(y) = -0.306357317549821454599000978430e-9, Z1(y) = 1.00000000264534810273850931227]](images/Dirac_57.gif){kind=small}
![[a = 0., b = 10.]](images/Dirac_58.gif){kind=small}
![[y = 1., Z(y) = 0.611840460259303704725888191778e-8, Z1(y) = 1.00000000652353924881966008655]](images/Dirac_59.gif){kind=small}
![[a = 0., b = 15.]](images/Dirac_60.gif){kind=small}
![[y = 1., Z(y) = 0.681223300593295428677411939526e-8, Z1(y) = 1.00000000706617521022333591176]](images/Dirac_61.gif){kind=small}
![[a = 0., b = 20.]](images/Dirac_62.gif){kind=small}
![[y = 1., Z(y) = 0.654220273959467641614373984506e-8, Z1(y) = 1.00000000669364210921004599322]](images/Dirac_63.gif){kind=small}
![[a = 0., b = 25.]](images/Dirac_64.gif){kind=small}
![[y = 1., Z(y) = 0.700088728875248589962619058857e-8, Z1(y) = 1.00000000690349332906039390659]](images/Dirac_65.gif){kind=small}
![Digits := 10; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_66.gif)
![Digits := 10; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_67.gif)
![Digits := 10; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_68.gif){kind=small}
![Digits := 10; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_69.gif){kind=small}
![Digits := 10; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_70.gif){kind=small}
![Digits := 10; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_71.gif){kind=small}
![Digits := 10; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_72.gif){kind=small}
![[a = 0., b = 1.]](images/Dirac_75.gif){kind=small}
![[y = 1., Z(y) = HFloat(-1.5245301486120025e-8), Z1(y) = HFloat(0.1816901153099545)]](images/Dirac_76.gif){kind=small}
![[a = 0., b = 10.]](images/Dirac_77.gif){kind=small}
![[y = 1., Z(y) = HFloat(-2.460041122392137e-8), Z1(y) = HFloat(0.8735173022884068)]](images/Dirac_78.gif){kind=small}
![[a = 0., b = 100.]](images/Dirac_79.gif){kind=small}
![[y = 1., Z(y) = HFloat(-3.085020925671462e-8), Z1(y) = HFloat(0.9872684529312092)]](images/Dirac_80.gif){kind=small}
![[a = 0., b = 1000.]](images/Dirac_81.gif){kind=small}
![[y = 1., Z(y) = HFloat(2.8364723288904784e-8), Z1(y) = HFloat(0.998726779043347)]](images/Dirac_82.gif){kind=small}
![[a = 0., b = 10000.]](images/Dirac_83.gif){kind=small}
![[y = 1., Z(y) = HFloat(-2.4139345964587916e-7), Z1(y) = HFloat(0.9998724244008227)]](images/Dirac_84.gif){kind=small}
![[a = 0., b = 100000.]](images/Dirac_85.gif){kind=small}
![[y = 1., Z(y) = HFloat(-3.6387096236949705e-6), Z1(y) = HFloat(0.9999836423651798)]](images/Dirac_86.gif){kind=small}
![[a = 0., b = 1000000.]](images/Dirac_87.gif){kind=small}
![[y = 1., Z(y) = HFloat(4.808290045290797e-6), Z1(y) = HFloat(1.0000035443669102)]](images/Dirac_88.gif){kind=small}
![[a = 0., b = 10000000.]](images/Dirac_89.gif){kind=small}
![[y = 1., Z(y) = HFloat(-3.279999071460797e-4), Z1(y) = HFloat(0.9996718756692161)]](images/Dirac_90.gif){kind=small}
![[a = 0., b = 100000000.]](images/Dirac_91.gif){kind=small}
![[y = 1., Z(y) = HFloat(-0.0024588343981494626), Z1(y) = HFloat(0.9975411540917314)]](images/Dirac_92.gif){kind=small}
![[a = 0., b = 1000000000.]](images/Dirac_93.gif){kind=small}
![[y = 1., Z(y) = HFloat(-0.04425362101111301), Z1(y) = HFloat(0.95574637761905)]](images/Dirac_94.gif){kind=small}
![Digits := 30; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_95.gif)
![Digits := 30; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_96.gif)
![Digits := 30; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_97.gif){kind=small}
![Digits := 30; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_98.gif){kind=small}
![Digits := 30; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_99.gif){kind=small}
![Digits := 30; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_100.gif){kind=small}
![Digits := 30; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], [...](images/Dirac_101.gif){kind=small}
![[a = 0., b = 1.]](images/Dirac_104.gif){kind=small}
![[y = 1., Z(y) = -0.15245301448503933195556427e-7, Z1(y) = .181690115309954389644926674467]](images/Dirac_105.gif){kind=small}
![[a = 0., b = 10.]](images/Dirac_106.gif){kind=small}
![[y = 1., Z(y) = -0.24600408807241527176729426334e-7, Z1(y) = .873517302288407692584220154609]](images/Dirac_107.gif){kind=small}
![[a = 0., b = 100.]](images/Dirac_108.gif){kind=small}
![[y = 1., Z(y) = -0.308502133124344907419083812e-7, Z1(y) = .987268452931203875115777812640]](images/Dirac_109.gif){kind=small}
![[a = 0., b = 1000.]](images/Dirac_110.gif){kind=small}
![[y = 1., Z(y) = 0.2836475094709602303435591685e-7, Z1(y) = .998726779043371895864235525693]](images/Dirac_111.gif){kind=small}
![[a = 0., b = 10000.]](images/Dirac_112.gif){kind=small}
![[y = 1., Z(y) = -0.241394928954764151618836507080e-6, Z1(y) = .999872424399353124789816947009]](images/Dirac_113.gif){kind=small}
![[a = 0., b = 100000.]](images/Dirac_114.gif){kind=small}
![[y = 1., Z(y) = -0.363868198771966774819400989911e-5, Z1(y) = .999983642392815517759802800082]](images/Dirac_115.gif){kind=small}
![[a = 0., b = 1000000.]](images/Dirac_116.gif){kind=small}
![[y = 1., Z(y) = 0.480825113717448439195606867696e-5, Z1(y) = 1.00000354432799753592827610560]](images/Dirac_117.gif){kind=small}
![[a = 0., b = 10000000.]](images/Dirac_118.gif){kind=small}
![[y = 1., Z(y) = -0.328003133546162168967901529826e-3, Z1(y) = .999671872442812120176265886944]](images/Dirac_119.gif){kind=small}
![[a = 0., b = 100000000.]](images/Dirac_120.gif){kind=small}
![[y = 1., Z(y) = -0.245881039555435603059004806481e-2, Z1(y) = .997541178094320977793014094960]](images/Dirac_121.gif){kind=small}
![[a = 0., b = 1000000000.]](images/Dirac_122.gif){kind=small}
![[y = 1., Z(y) = -0.442536838492834757756285778904e-1, Z1(y) = .955746314780879979378933291837]](images/Dirac_123.gif){kind=small}
![Digits := 100; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], ...](images/Dirac_124.gif)
![Digits := 100; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], ...](images/Dirac_125.gif)
![Digits := 100; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], ...](images/Dirac_126.gif){kind=small}
![Digits := 100; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], ...](images/Dirac_127.gif){kind=small}
![Digits := 100; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], ...](images/Dirac_128.gif){kind=small}
![Digits := 100; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], ...](images/Dirac_129.gif){kind=small}
![Digits := 100; 1; resu4 := dsolve([diff(Z(y), y) = diff(`/`(`*`(`+`(arctan(`*`(`+`(y, `-`(a)), `*`(b))), `*`(`/`(1, 2), `*`(Pi)))), `*`(Pi)), a, a), diff(Z1(y), y) = Z(y), Z(-1.) = 0., Z1(-1.) = 0.], ...](images/Dirac_130.gif){kind=small}
![[a = 0., b = 1.]](images/Dirac_133.gif){kind=small}
![[a = 0., b = 10.]](images/Dirac_139.gif){kind=small}
![[a = 0., b = 100.]](images/Dirac_145.gif){kind=small}
![[a = 0., b = 1000.]](images/Dirac_151.gif){kind=small}
![[a = 0., b = 10000.]](images/Dirac_157.gif){kind=small}
![[a = 0., b = 100000.]](images/Dirac_163.gif){kind=small}
![[a = 0., b = 1000000.]](images/Dirac_169.gif){kind=small}
![[a = 0., b = 10000000.]](images/Dirac_175.gif){kind=small}
![[a = 0., b = 100000000.]](images/Dirac_181.gif){kind=small}
![[a = 0., b = 1000000000.]](images/Dirac_187.gif){kind=small}
