vv

6117 Reputation

18 Badges

3 years, 328 days

MaplePrimes Activity


These are answers submitted by vv

If  F, G are distributions, F*G cannot be defined in general.
It is possible to define it when F is a test function and G is a distribution (in D'(R) or S'(R), i.e. tempered).
So, Dirac(x)^2  or Heaviside(x)*Dirac(x) do not exist as distributions.
[However, in some formal contexts, the latter could appear, I think]

Note that Maple acts symbolically (formally), so the user must check that the desired computations make sense.
For example,
int(x^2*Dirac(x)^2,x=-infinity..infinity);  #==>  not evaluated
int(x*Dirac(x)^2,x=-infinity..infinity);     # ==> 0
even if both are nonsense.

To obtain a simple and general solution it's better to change the strategy and also help Maple a little.

 

de := x^4*diff(y(x), x $ 2) + omega^2*y(x) = 0;

x^4*(diff(diff(y(x), x), x))+omega^2*y(x) = 0

(1)

# bc := y(a) = 0, y(b) = 0, D(y)(a) = 1;

dsolve({de, y(a) = 0, D(y)(a) = 1},y(x)):

Y:=simplify(rhs(%));

a*x*(sin(omega/a)*cos(omega/x)-cos(omega/a)*sin(omega/x))/omega

(2)

omega/a=omega/b+k*Pi;  #y(b)=0

omega/a = omega/b+k*Pi

(3)

omega=solve(%, omega);  # k in Z \ {0}

omega = -k*Pi*b*a/(a-b)

(4)

y(x)=(factor@combine)(eval(Y,%));

y(x) = -x*sin(k*Pi*b*(a-x)/((a-b)*x))*(a-b)/(k*Pi*b)

(5)

 

To compute the second integral:

int(sin(x)^(1/3)*cos(x)^3, x):
eval(IntegrationTools:-Change(%, cos(x)=u), u=cos(x));

Note that it is valid only for sin(x)>0. Check:

simplify(diff(%-%%,x)) assuming x>0,x<Pi;  # or sin(x)>0
     0

The integral given by Maple for your 1st example is valid in R.

After isolating y', the new ode is not exactly equivalent with the original, because the solutions given by the algebraic equation

2*x^(5/2) - 3*y(x)^(5/3) = 0      (*)

are lost [the ode was divided by the lhs of (*), which is seen "generically" <> 0].

The first 5 solutions given by dsolve are from (*).

 

Each solution is local, i.e. is valid in some intervals.

For example, sol[2] is valid for x > -ln(_C1)   if  _C1 > 0.

Note that the Cauchy problem  
dsolve({ode, y(0)=1});

gives  y(x) = 1, y(x) = 2*exp(x)-1
but the first one is nowhere valid.

Edit.

The correct (standard) general solution is:

y(x) = piecewise(0 <= x+c, exp(x+c)-1, 1-exp(-x-c))

 

f:=a*x+b;
pf:=unapply(f, a,b,x);
g:=(a,b)->curry(pf, a,b);

 

A:-B:-foo();           

In ?type,protected it is recommended:

select(type, {unames(), anames(anything)}, protected);

Executing this command several times, the number of names increases.
Probably the simple evaluation of some of these names generates other names which are not "seen" initially.

 

labels = ["x", ""], title=cat("y", " "$220)

restart;
a := 1/(i*sqrt(i+1)+(i+1)*sqrt(i)):
evalf(Sum(a, i=1..infinity));  # approx

                          1.000000000
 

b := (expand@simplify@expand@rationalize)(a):
sum(b, i=1..infinity); #exact
                               1

 

I posted an answer to your deleted question.
I think that the moderator who deleted it was wrong because:
1. It was answered.
2. The content was distinct from the previous question. Not to mention that the old thread is difficult to follow (and loads slowly) being very long.

So, here is the answer again.

restart;

Digits:=15;

15

(1)

F[1, 1] := (r, theta) -> A[1, 1] + A[1, 2]*cos(theta) + A[2, 1]*(2*r - 1) + A[2, 2]*(2*r - 1)*cos(theta);
L[1, 1] := (r, theta, phi) -> ((2.803059644*10^12*cos(theta)^2*r^2 + 7.474825716*10^11*r^3*cos(theta)^3 + 7.474825716*10^10*r^4*cos(theta)^4 + 4.671766072*10^12*r*cos(theta) + 2.919853795*10^12)*diff(F[1, 1](r, theta), theta, theta) + (1.584823993*10^13*cos(theta)*r^3 + 9.508943959*10^12*cos(theta)^2*r^4 + 2.535718389*10^12*r^5*cos(theta)^3 + 2.535718389*10^11*r^6*cos(theta)^4 + 9.905149957*10^12*r^2)*diff(F[1, 1](r, theta), r, r) + (1.981029991*10^13*r^2*cos(theta) + 1.426341594*10^13*r^3*cos(theta)^2 + 4.437507181*10^12*r^4*cos(theta)^3 + 5.071436778*10^11*r^5*cos(theta)^4 + 9.905149957*10^12*r)*diff(F[1, 1](r, theta), r) + ((-1)*1.167941518*10^12*sin(theta)*r + (-1)*5.606119287*10^11*sin(theta)*r^3*cos(theta)^2 + (-1)*7.474825716*10^10*sin(theta)*r^4*cos(theta)^3 + (-1)*1.401529822*10^12*sin(theta)*r^2*cos(theta))*diff(F[1, 1](r, theta), theta) + ((-1)*5.071436778*10^11*r^4*cos(theta)^4 + (-1)*3.803577584*10^12*r^3*cos(theta)^3 + ((-1)*1.109376795*10^13*r^2 + (-1)*7.474825716*10^10*r^4)*cos(theta)^2 + ((-1)*1.584823993*10^13*r + (-1)*3.737412858*10^11*r^3)*cos(theta) + (-1)*9.905149957*10^12 + (-1)*4.671766072*10^11*r^2)*F[1, 1](r, theta))/((39.0625 + 62.5*r*cos(theta) + 37.5*cos(theta)^2*r^2 + 10.*r^3*cos(theta)^3 + r^4*cos(theta)^4)*r^2):

proc (r, theta) options operator, arrow; A[1, 1]+A[1, 2]*cos(theta)+A[2, 1]*(2*r-1)+A[2, 2]*(2*r-1)*cos(theta) end proc

(2)

for p to 1 do  for s to 1  do
f := L[p, s](r, theta, phi)*F[p, 1](r, theta):
Aij:=[indets(F[p,s](r,theta),indexed)[]];
g1:=expand(f, Aij, distributed);
g2:=coeffs(g1, Aij, 'T'):
C:=int~([g2], theta = 0 .. 2*Pi, r = 0.5 .. 1, epsilon=1e-8, numeric):
kk:=add(C .~ T);
kkk := evalf( int(F[p, 1](r, theta)^2, theta = 0 .. 2*Pi, r = 0.5 .. 1) );
k[p, s] := kk/kkk;
print([p, s] = %);
od od;

[1, 1] = (-1716301349235.42*A[1, 1]^2+745939588799.615*A[1, 1]*A[2, 1]+513374393378.777*A[2, 1]^2+161631455946.711*A[1, 1]*A[1, 2]+427274927473.220*A[2, 1]*A[1, 2]+135421805163.233*A[2, 2]*A[1, 2]+427274927473.220*A[1, 1]*A[2, 2]+406304482715.359*A[2, 1]*A[2, 2]-1120768744270.57*A[1, 2]^2+178063137194.727*A[2, 2]^2)/(1.04719755119660*A[2, 1]^2+.523598775598299*A[2, 2]^2+3.14159265358979*A[1, 1]*A[2, 1]+1.57079632679490*A[2, 2]*A[1, 2]+3.14159265358979*A[1, 1]^2+1.57079632679490*A[1, 2]^2)

(3)
f:=exp(-(sqrt(4*x^2+4*y^2+4*z^2)^3));
F:=eval(f, [x=r*sin(v)*cos(u), y=r*sin(v)*sin(u), z=r*cos(v)]);
int(F*r^2*sin(v), r=0..2, u=0..2*Pi, v=0..Pi);
evalf(%);

4*(1/24 - exp(-64)/24)*Pi

0.5235987756

To obtain a numerical answer you must change the variables to integrate in a cube.
Of course the spherical one is the best here because the new integrand is continuous.
Or, use piecewise:

JJ:=Int(exp(-(sqrt(4*x^2+4*y^2+4*z^2)^3))*piecewise(x^2+y^2+z^2<4,1,0), x=-2..2,y=-2..2,z=-2..2):
evalf(JJ);

0.5235987756

restart;

Digits:=15;

15

(1)

F[1, 1] := (r, theta) -> A[1, 1] + A[1, 2]*cos(theta) + A[2, 1]*(2*r - 1) + A[2, 2]*(2*r - 1)*cos(theta);
L[1, 1] := (r, theta, phi) -> ((2.803059644*10^12*cos(theta)^2*r^2 + 7.474825716*10^11*r^3*cos(theta)^3 + 7.474825716*10^10*r^4*cos(theta)^4 + 4.671766072*10^12*r*cos(theta) + 2.919853795*10^12)*diff(F[1, 1](r, theta), theta, theta) + (1.584823993*10^13*cos(theta)*r^3 + 9.508943959*10^12*cos(theta)^2*r^4 + 2.535718389*10^12*r^5*cos(theta)^3 + 2.535718389*10^11*r^6*cos(theta)^4 + 9.905149957*10^12*r^2)*diff(F[1, 1](r, theta), r, r) + (1.981029991*10^13*r^2*cos(theta) + 1.426341594*10^13*r^3*cos(theta)^2 + 4.437507181*10^12*r^4*cos(theta)^3 + 5.071436778*10^11*r^5*cos(theta)^4 + 9.905149957*10^12*r)*diff(F[1, 1](r, theta), r) + ((-1)*1.167941518*10^12*sin(theta)*r + (-1)*5.606119287*10^11*sin(theta)*r^3*cos(theta)^2 + (-1)*7.474825716*10^10*sin(theta)*r^4*cos(theta)^3 + (-1)*1.401529822*10^12*sin(theta)*r^2*cos(theta))*diff(F[1, 1](r, theta), theta) + ((-1)*5.071436778*10^11*r^4*cos(theta)^4 + (-1)*3.803577584*10^12*r^3*cos(theta)^3 + ((-1)*1.109376795*10^13*r^2 + (-1)*7.474825716*10^10*r^4)*cos(theta)^2 + ((-1)*1.584823993*10^13*r + (-1)*3.737412858*10^11*r^3)*cos(theta) + (-1)*9.905149957*10^12 + (-1)*4.671766072*10^11*r^2)*F[1, 1](r, theta))/((39.0625 + 62.5*r*cos(theta) + 37.5*cos(theta)^2*r^2 + 10.*r^3*cos(theta)^3 + r^4*cos(theta)^4)*r^2):

proc (r, theta) options operator, arrow; A[1, 1]+A[1, 2]*cos(theta)+A[2, 1]*(2*r-1)+A[2, 2]*(2*r-1)*cos(theta) end proc

(2)

for p to 1 do  for s to 1  do
f := L[p, s](r, theta, phi)*F[p, 1](r, theta):
Aij:=[indets(F[p,s](r,theta),indexed)[]];
g1:=expand(f, Aij, distributed);
g2:=coeffs(g1, Aij, 'T'):
C:=int~([g2], theta = 0 .. 2*Pi, r = 0.5 .. 1, epsilon=1e-8, numeric):
kk:=add(C .~ T);
kkk := evalf( int(F[p, 1](r, theta)^2, theta = 0 .. 2*Pi, r = 0.5 .. 1) );
k[p, s] := kk/kkk;
print([p, s] = %);
od od;

[1, 1] = (-1716301349235.42*A[1, 1]^2+745939588799.615*A[1, 1]*A[2, 1]+513374393378.777*A[2, 1]^2+161631455946.711*A[1, 1]*A[1, 2]+427274927473.220*A[2, 1]*A[1, 2]+135421805163.233*A[2, 2]*A[1, 2]+427274927473.220*A[1, 1]*A[2, 2]+406304482715.359*A[2, 1]*A[2, 2]-1120768744270.57*A[1, 2]^2+178063137194.727*A[2, 2]^2)/(1.04719755119660*A[2, 1]^2+.523598775598299*A[2, 2]^2+3.14159265358979*A[1, 1]*A[2, 1]+1.57079632679490*A[2, 2]*A[1, 2]+3.14159265358979*A[1, 1]^2+1.57079632679490*A[1, 2]^2)

(3)

For a solution given implicitely, as in this case, odetest may fail.

Actually the implicit solution is

Int( 1/(y^3 + 1)^(2/3), y) + Int( 1/(x^3 + 1)^(2/3), x) + c = 0;

Int(1/(y^3+1)^(2/3), y)+Int(1/(x^3+1)^(2/3), x)+c = 0

(1)

The problem is that the integrals are computed in terms of hypergeom, and Maple is not able to further manipulate.

f:=1/(x^3 + 1)^(2/3);

1/(x^3+1)^(2/3)

(2)

int(f, x)

x*hypergeom([1/3, 2/3], [4/3], -x^3)

(3)

Z:=diff(% ,x) - f; # Z should be 0

hypergeom([1/3, 2/3], [4/3], -x^3)-(1/2)*x^3*hypergeom([4/3, 5/3], [7/3], -x^3)-1/(x^3+1)^(2/3)

(4)

is(Z=0) assuming x>0, x<1;

FAIL

(5)

# It seems that any usual convertion fails here.

Salut Viorel, bine ai venit!

If you are interested in an approximate solution based on series, here it is:

restart;
ode:= diff(g(r),r$2)- r/R*g(r)=0;
ic:=g(2*R)=0, D[1](g)(0)=R;
dsolve([ode,ic],g(r)) assuming R>0;
G:=rhs(%);
Order:=6;
simplify(series(G, r=0));
convert(%, polynom):
g1:=convert(series(%, R=0), polynom);  
R:=1/3;  #############
dsol:=dsolve([ode,ic],g(r), numeric);
plots:-odeplot(dsol, r=0..2*R);
plot(g1, r=0..2*R);

So, for small R and r, an approximation is

1 2 3 4 5 6 7 Last Page 1 of 70