14 years, 0 days

## rule for cosine...

Just expand the exponentials and use the rule for the cosine:

exp2cos:=A::algebraic*exp(I*B::algebraic)+A::algebraic*1/exp(I*B::algebraic)=2*A*cos(B):applyrule(exp2cos,expand(r));

2      2
1/2   1/2       k  sigma
sigma 2    Pi    exp(- ---------) cos(k theta)
2



Or, without expanding:

exp2cosh:=A::algebraic*exp(a::algebraic)+A::algebraic*exp(b::algebraic)=2*A*cosh((a-b)/2)*exp((a+b)/2):
applyrule(exp2cosh,r);

2      2
1/2   1/2       k  sigma
sigma 2    Pi    exp(- ---------) cos(k theta)
2


If the last one is E(x)*E(y) = -1/n if x<>y and E(x)^2 = 1/n, you can apply such rules by applyrule as in these examples:

r:=[E(x::name)*E(y::name) = -1/n,E(x::name)^2 = 1/n]:
applyrule(r,E(u1)*E(u2));
1
- -
n

map2(applyrule,r,[E(u1)*E(u2)*E(u3)*E(u4),E(u1)^2*E(u2)*E(u3),E(u1)^2*E(u3)^2]);
[1     1   1 ]
[--, - --, --]
[ 2     2   2]
[n     n   n ]



## get maple to give the Frechet derivative...

For expressions similar to u(x)*u(y), the following procedure will do it, by expanding F(u+v)-F(u), and using rules to remove nonlinear terms in v:

fd:=proc(F,u::name,v::name)
local r:
r:=[A::algebraic*v(a::name)*v(b::name)=0,conditional(A::algebraic*v(b::name)^n::posint=0,n>1)]:
subs(u=u+v,F)-F;
simplify(%);
applyrule(r,%);
end proc:

fd(u(x)*u(y),u,v);
u(x) v(y) + v(x) u(y)

fd(u(x)*u(y)^2,u,v);
2
2 u(x) u(y) v(y) + v(x) u(y)

fd(u(x)*u(y)*u(z),u,v);
u(x) u(y) v(z) + u(x) v(y) u(z) + v(x) u(y) u(z)


## the problem here...

Tracing shows that int(sqrt(1 + sec(x)^2), x) is transformed by the change of variable u=tan(x/2), producing an integrand with a fourth power:

[...]
{--> enter int/algebraic2/algebraic, args = 1/(_X^2-1)*(2*_X^4+2)^(1/2)/(1+_X^2)
[...]


and probably, because of it, it is dispatched to the elliptic integration routines. On the other hand, int(sqrt(2 + tan(x)^2), x) produces a quadratic term:

[...]
{--> enter int/algebraic2/algebraic, args = (2+_X^2)^(1/2)/(1+_X^2)
[...]


and it is integrated in terms of trigonometrics.

## inert and applyrule...

First, for such symbolic computation it is advisable to use the inert form of the differential functions like %diff and %Laplacian, instead of the active forms diff and Laplacian that you would use for explicit calculations. Then, the required transformation rules, like permuting the order of the derivatives, can be performed by applyrule:

e2 := %diff(u(x, y), x)-%diff(v(x, y), y);
G2:=map(%Laplacian,e2):
Sxx := map(%diff,G2, x):
Syy := map(%diff,G2, y):
r1:=%Laplacian(%diff(y::function,x::name))=%diff(%Laplacian(y),x):
r2:=%Laplacian(-1*A::algebraic)=-%Laplacian(A):
r3:=%diff(-1*A::algebraic,x::name)=-%diff(A,x):
s:=%Laplacian(u(x,y))=P,%Laplacian(v(x,y))=Q:applyrule([r1,r2,r3,s],Sxx);

/ 2   \   /  2    \
|d    |   | d     |
|--- P| - |----- Q|
|  2  |   \dx dy  /
\dx   /

applyrule([r1,r2,r3,s],Syy);

/  2    \   / 2   \
| d     |   |d    |
|----- P| - |--- Q|
\dy dx  /   |  2  |
\dy   /



Note that linearity of the inert differential functions has to be imposed explicitly, e.g. by using map. And no package is needed for these symbolic computations.

## how can get the result...

Several steps are involved for arriving to the answer, the same as you would do by hand. So, the main issue is implementing every of these steps as a rule. Also, as these symbolic computations do not make use the inner computational routines associated with function calls like sum(...) or int(...), it is more efficient avoiding them altogether by using the inert forms like Sum(...) or Int(...):

eq := I*hbar*(Sum((diff(c[n](t), t))*f[n](r)*exp(-I*omega[n]*t), n = l .. k)) = (1/2)*E[0](e_.r_)*e*(Sum(c[n](t)*f[n](r)*omega[n]*(exp(I*t*(-omega[n]+Omega))+exp(-I*t*(omega[n]+Omega))), n = l .. k)):
req := Int(conjugate(f[m](r))*rhs(eq), r):
leq := Int(conjugate(f[m](r))*lhs(eq), r):

#rules
sumint:=Int(A::algebraic*Sum(v::algebraic,r::equation),u::name)=Sum(Int(A*v,u),r):
intsum:=Sum(v::algebraic*Int(A::algebraic,u::name),r::equation)=Int(Sum(A*v,r),u):
sup:= Int(A::algebraic*conjugate(f[m](r))*f[n](r), r) = A*delta[m, n]:
kron:=[Sum(F::Not(anyindex(identical(n)))*a::anyindex(identical(n))*delta[m,n],r::equation)=F*a,Sum(F::Not(anyindex(identical(n)))*delta[m,n],r::equation)=F]:

#lhs
applyrule([sumint,sup],Int(conjugate(f[m](r))*leq,r)):
applyrule(intsum,%):
subs(n=m,applyrule(kron,%)):
LHS:=IntegrationTools:-Expand(%):

#rhs
applyrule([sumint,sup],Int(conjugate(f[m](r))*req,r)):
applyrule(intsum,%):
subs(n=m,applyrule(kron,%)):
IntegrationTools:-Expand(%):
RHS:=collect(%,[Int,E[0],c[m]],simplify):

(LHS=RHS)/op(indets(LHS,specfunc(anything,Int)));
/d         \
hbar |-- c[m](t)| I
\dt        /
------------------- = 1/2 e omega[m]
exp(omega[m] t I)

(exp(t (-omega[m] + Omega) I) + exp(-I t (omega[m] + Omega)))

c[m](t) E[0](e_ . r_)


## Error and forbid an unwanted assignment...

Actually, an error can be raised for trying to assign 0, by a typed assignment with assertion failure:

restart:
kernelopts(assertlevel=2):
a::Non(0):=0;Error, assertion failed in assignment, expected Non(0), got 0

a;
a


## anyfunc...

I would do it in two stages, one for the function call using anyfunc, and then for its arguments, like:

g(u+v*w):
patmatch(%,f::anyfunc(+(algebraic)),'p'):p;
patmatch(op(subs(%,f)),a::name*b::name +c::name,'p'):p;

[f = g(u + v w)]
[a = v, b = w, c = u]


Note that the type +(algebraic) is simple but a bit more general than the pattern in the next line. A more strict (but longer) type version is:

g(u+v*w): patmatch(%,f::anyfunc(&+(name,&*(name,name))),'p'):p;
patmatch(op(subs(%,f)),a::name*b::name +c::name,'p'):p;

[f = g(u + v w)]
[a = v, b = w, c = u]


Variants may be needed, depending on details. Also, if you do know the name of a variable, v say, and you want linearity in it, you could do instead something like:

g(u+v*w):
patmatch(%,f::anyfunc(+(linear(v))),'p'):p;
patmatch(op(subs(%,f)),a::name*v +c::name,'p'):p;

[f = g(u + v w)]
[a = w, c = u]


## simple transformation...

I think that such a simple algebraic transformation should be possible to be dealt with also simply, in uniform way, without much need of creativity or fine tuning, as is the case of using applyrule for the above examples:

csq:=conditional(a::algebraic+2*b::algebraic,a=b^2)=(b+1)^2-1:
applyrule(csq,x^2+2*x);
2
(x + 1)  - 1

applyrule(csq,x^4+2*x^2);
2     2
(x  + 1)  - 1

applyrule(csq,x^4/sin(x)^4+2*x^2/sin(x)^2);

/   2       \2
|  x        |
|------- + 1|  - 1
|      2    |
\sin(x)     /


## assuming real...

Just assuming real:

int((r^2+r*cos(theta))*r*(a^2-r^2)^(1/2),[r = 0 .. a,theta = 0 .. 2*Pi]) assuming real;
5
4/15 a  Pi signum(a)



## problems with a radiation boundary condi...

See this thread, where problems were also reported when a radiation boundary condition is present.

## shortest form...

You can easily make the transformations from exp to cosh and sinh with applyrule:

applyrule~(1/2*exp(a::algebraic)+1/2*exp(b::algebraic)=cosh((a-b)/2)*exp((a+b)/2),p):
applyrule~(1/2*exp(a::algebraic)-1/2*exp(b::algebraic)=sinh((a-b)/2)*exp((a+b)/2),%);


## by applyrule...

I would use:

ee:=itail*(exp(vd/Vth)-1)/(1+exp(vd/Vth)):
applyrule((exp(z::algebraic)-1)/(exp(z::algebraic)+1)=tanh(z/2),ee);

vd
tanh(-----) itail
2 Vth

 First 27 28 29 Page 29 of 29
﻿