kill(all)$

maxima :derivatives of the delta function , version 1.0a
peter.vlasschaert@gmail.com,18/02/2020

see : maxima implimentation @ moment only for , Laplace Transform.
'delta function'

load 'packages'
===========

load(lrats)$

defintion : delta function property 'integral'
==============================

p1:integrate(delta(x)*f(x),x,minf,inf) =f(0);

p1a:args(p1);

 p1 : only value when delta(0)  => f(0)

a) integration by parts
=================

p2a:integrate(diff(u(x)*v(x),x),x);

p2a1:args(p2a);

p2b:integrate('diff(u(x)*v(x)),x);

p2b1:args(p2b);

p3:p2b1[1] =p2a1[1];

p3a:lhs(p3) ='diff(u(x)*v(x),x);

p3b:part(rhs(p3),1);

p3c:subst([p3a],rhs(p3-p3b) =lhs(p3-p3b));

p3d:subst([%c1=0],integrate(lhs(p3c),x) =integrate(rhs(p3c)+p3b,x)-integrate(p3b,x));

definition :partial integration
====================

p4:lhs(p3d) =ev(rhs(p3d),integrate);

example :x^2*sin(x)
-----------------------

a)  example : partial integration ( use p4)
-------------------------------------------------

p5:[u(x)=sin(2*x),v(x)=x^2];
p5a:subst(p5,p4);

p5b:lhs(p5a) = ev(rhs(p5a),nouns);

p5b1:lhs(p5b) =expand(rhs(p5b));

p5b2:lhs(p5b)=collectterms(rhs(p5b1),sin(2*x),cos(2*x));

b)  example : direct integration ( use p4)
-------------------------------------------------

p5c:lhs(p5b) =ev(lhs(p5b),nouns);

p5c1:lhs(p5b) = expand(rhs(p5b));

p5c2:lhs(p5b)= collectterms(rhs(p5c1),sin(2*x),cos(2*x));

syntax: collectterms(expr, ...,....,...)

combination : two calculation methods
============================
1e) expand result
2e) collectterms(1e,..,...)
============================

check result : p5b2 ? p5c2
--------------------------------

p6:is(rhs(p5b2) = rhs(p5c2));

use :definition :p4
=============
to put in :integration limits: - ∞,∞

p6:p4;

p6a1:args(lhs(p6));
p6a2:args(rhs(p6));

p6a3:args(part(rhs(p6),2));

p6a4:args(p6a3[1]);

reverse function : args to integrate

p6b:[q(var1),var1,a1,a2] =integrate('q(var1),var1,a1,a2);

construct : reverse function.
====================

p7:[a,b]; /*integration limits*/
p7a:flatten(append(p6a1,p7)); /*function to integrate, args*/

p7b:[lhs(p6b)[1]=p7a[1],lhs(p6b)[2]=p7a[2],lhs(p6b)[3]=p7a[3],lhs(p6b)[4]=p7a[4]];

p7b1:makelist(lhs(p6b)[i]=p7a[i],i,length(p7a));

p7c:lratsubst(p7b1,p6b);

p7d:ev(lratsubst([q(x)=p7a[1]],p7c),nouns);

intargs(a1,a2,vari,ff):=block([],
    mi:integrate(ff,vari),
    m1:[a1,a2],
    m2:args(mi),
    m:flatten(append(m2,m1)),
    mp:[ff,vari,a1,a2] =integrate(ff,vari,a1,a2)
   );

p8a:intargs(a,b,x,p6a1[1]);
p8b:intargs(a,b,x,p6a4[1]);

what about :u(x)*v(x)
===============
limits ,[a1,a2]=[a,b]

p8c:[u(x)*v(x),a,b] =[u(b)*v(b)-u(a)*v(a)];

p8c1:rhs(p8c)[1];

definition :partial integration,with integration limits
====================================
a = lower limit
b = upper limit

p9:rhs(p8a) = p8c1-rhs(p8b);

a=minf = -∞
b=inf    =∞

p9a:[a=minf,b=inf]$
p9a1:subst(p9a,p9);

computation :
==========

p10:integrate(f(x)*diff(delta(x),x,n),x,minf,inf);

example : n=1

p10a:subst([n=1],args(lhs(p9a1))=args(lhs(p10)));

p10b:lhs(p10a)[1] =rhs(p10a)[1];

p10b1:args(lhs(p10b)) = args(rhs(p10b));

p10b2:map("=",lhs(p10b1),rhs(p10b1));

p10b2a:subst([%c1=0],integrate(p10b2[2],x));

use  :p9a1

p10b3:subst(p10b2,p9a1);

p10b4:subst([p10b2a],p10b3);

rem : only exit delta(x=0),see :definition 'p1'
rem : -u(-inf)*v(-inf)+u(inf)*v(inf) = 0,symmetric.

how to do substitution ,trough args

p10b51:args(lhs(p10b4));

p10b52a:args(part(rhs(p10b4),1));

p10b52b:args(p10b52a[1]);

p10b53c:p10b51[1]=-p10b52b[1];

p11:intargs(minf,inf,x,lhs(%))=intargs(minf,inf,x,rhs(%))$

p11a:rhs(lhs(p11)) = rhs(rhs(p11));

use :p1,only x=0 ,use p10b52a & p10b52b

final result :

p11b:lhs(p11a) =-subst([delta(x)=1,f(x)=f(0)],p10b52b[1]);

meaning :'diff(f(x),x,1) @ x=0 →  'diff(f(0),x,1)

how to do : maxima

p11c:limit(diff(f(x),x),x,0)=subst([f(x)=f(0)],diff(f(x),x));

combination result :

p11d:lhs(p11b) =-limit(subst([delta(x)=1],p10b52b[1]),x,0);

use  : p11c

p11d1:subst([p11c],p11d);

find nth power of dirac delta function.
===========================
use : n=2,use 'p10'

example : n=2

p12:subst([n=2],p10);

use :condition , symmetric ,see 'below p10b4'

how to find parts from : 'p10b4'
--------------------------------------

p12a:lhs(p10b4) = part(rhs(p10b4),1);

use:p9a1

p12a1:lhs(p9a1) = part(rhs(p9a1),1);

sign : rhs of eq: 'p12a'

first part :rhs
----------------

p12b1: part(rhs(p10b4),1,0);

second part:rhs
-------------------

p12b2:part(rhs(p10b4),1,1);

rhs
----

p12b3:part(rhs(p10b4),1);

for n=2
=====
use :p12a,p12a1

→ n=1

p12c:p11a;

→ n=2;

p12c1:p12=-integrate(diff(delta(x),x,1)*diff(f(x),x,1),x,minf,inf);

use :partial integration.

p14:p12a1;

p14a:[u(x)=diff(delta(x),x,1),v(x)=f(x)];

p14b:subst(p14a,p14);

p14c:[u(x)=delta(x),v(x)=diff(f(x),x)];
p14d:subst([p14c],p14);
p14e:(-1)*p14d;

combination:
=========

p15:lratsubst([p14e],p14b);

only for x=0

p16:lhs(p15) =limit(diff(f(x),x,2),x,0);

p16a:rhs(p16);

p16b:subst([f(x)=f(0)],args(p16a)[1]);

p17:lhs(p15) = p16b;

even: m=2  => coef  (-1)^m =(-1)^2 =1, 'rhs'
odd :  m=3  => coef  (-1)^m =(-1)^3 =-1, 'rhs'

first subst 
-----------:

p15a;p14a;

second subst
-----------------

p15b:p14c;

relation
=====  

p15c:args(lhs(p15)) = args(rhs(p15));

p15c:lhs(p15c)[1]=rhs(p15c)[1];