Advanced Search

vv

12453 Reputation

19 Badges

9 years, 286 days
Contact vv

MaplePrimes Activity


These are replies submitted by vv

Digits...

vv 12453
August 16 2019

@mmcdara 

It works in Maple 2019. But not properly for Digits>15 (try 25 and 20). Probably Carl will fix this.

Unfortunately you have not included a corresponding Your_NormInv to compare the accuracy.

z1=z...

vv 12453
August 16 2019

@mmcdara 

It works in newer versions.
For Maple 2015 just replace z1=z  with  abs(z1-z)<10^(-Digits+1)   (or something similar).

Risch...

vv 12453
August 11 2019
0 0

@DJJerome1976 

AFAIK the Risch algorithm is not fully implemented in Maple.
You can see the methods which are tried using
infolevel[int]:=5;

BTW, rewriting just a bit the integrand, Maple integrates easily:

int(sin(x)^(1/3)*(1-sin(x)^2)*cos(x), x);
   

simpler example...

vv 12453
August 09 2019
0 0

@nm 

Take

ode:=y(x)*diff(y(x),x) - y(x) = 0:
dsolve(ode);

   y(x) = 0, y(x) = x + _C1

This is obviously the correct answer (two solutions; actually there are other solutions by combining these on intervals).

If we solve ode wrt y'  ==>  diff(y(x),x) = 1
[obtained  dividing by y (supposed to be <> 0; that is how solve works)].

Of course, dsolve( diff(y(x),x) = 1 ) ==>  y(x) = x + _C1,
so, y=0 disappears. But y(x)=0 satisfies (of course) the ode.

asympt...

vv 12453
July 28 2019
0 0

@Rouben Rostamian  

OP wants an asymptotic expansion for G0(s). He should use asympt(G0(s),s)

The question does not seem to have much ...

vv 12453
July 28 2019
0 0

The question does not seem to have much sense.

ok, last change...

vv 12453
July 28 2019
0 0

@9009134 

# K[r, s]
for p from 1 to 3 do  
for  s from 1 to 3  do
ff := L[p, s](r, theta, phi)*F[p, 1](r, theta):
Aij:=[indets([F[p,s](r,theta),ff],indexed)[]];
# print(indets=Aij);
g1:=expand(ff, Aij, distributed);
g2:=coeffs(g1, Aij, 'T'):
C1:=int~([g2], theta = 0 .. 2*Pi, r = 0.5 .. 1, epsilon=1e-8, numeric):
kk:=add(C1 .~ 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;

 

@Carl Love Unfortunately, as I see,...

vv 12453
July 28 2019
0 0

@Carl Love 

Unfortunately, as I see, sw_zgeevx_  also fails sometimes.

The bug seems to occur only for multiple (repeated) eigenvalues in the non-symmetric case, so, in practice should be extremely rare.
The examples were found using e.g.

LinearAlgebra:-CompanionMatrix((x-15)^3, x)^+;

 

 

 

@9009134  1. f[1,1] is not define...

vv 12453
July 27 2019
0 0

@9009134 

1.

f[1,1] is not defined.

2. You must rename my f  to ff  (because you have introduced a f)

ff := L[p, s](r, theta, phi)*F[p, 1](r, theta):
Aij:=[indets(F[p,s](r,theta),indexed)[]];
g1:=expand(ff, Aij, distributed);

3. Now you have B[i,j] etc

Sorry, but I cannot continue with so many changes.


 

asympt...

vv 12453
July 26 2019
0 0

@michaelvio 

With series(...,r=r0) you may approximate only near r0.
In your example you have a huge r. Use asympt instead.
[Note that you have not actually used the power of slode; only a few terms, just like the standard series]

> 

restart;
R:=370; Order:=6:
ode:= diff(g(r),r$2)- r/R*g(r)=0:
ic:=g(2*R)=0, D[1](g)(0)=R:
G:=rhs(dsolve([ode,ic],g(r)));

370

  

-740*370^(1/3)*Pi*AiryBi(2*370^(2/3))*AiryAi((1/370)*370^(2/3)*r)/(GAMMA(2/3)*(AiryAi(2*370^(2/3))*3^(2/3)+3^(1/6)*AiryBi(2*370^(2/3))))+740*370^(1/3)*AiryAi(2*370^(2/3))*Pi*AiryBi((1/370)*370^(2/3)*r)/(GAMMA(2/3)*(AiryAi(2*370^(2/3))*3^(2/3)+3^(1/6)*AiryBi(2*370^(2/3))))

 (1)
> 

asympt(G, r):
> 

Gapp:=convert(evalf[50](%),polynom);

(0.47817456179569551515317081621201242285756903386552e-602*(1/r)^(1/4)+0.95811159840225010761653791669736639828136520861439e-602*(1/r)^(7/4)+0.14782110748219797338014145566276295493936785779979e-600*(1/r)^(13/4)+0.43638226822222969691007895829878914061306256867382e-599*(1/r)^(19/4))*(exp(19.235384061671344751855362921213373873781867091120/(1/r)^(3/2)))^(1/555)

 (2)
> 

evalf[50](eval([G,Gapp], r=400));

[-0.40702431977939449083274686710027843869975497791315e-117, 0.27836877438958079551079281963767137932752590618428e-482]

 (3)
> 

evalf[50](eval([G,Gapp], r=1000));

[0.81928332104566902988587362191002268450537068088687e-127, 0.81928332104563623792440214614327162747404927396426e-127]

 (4)
> 

evalf[50](eval([G,Gapp], r=10^10));

[0.16056840058550626226900357319617922990909434242673e15051930008898, 0.16056840058550626226900357319617922981745765695502e15051930008898]

 (5)

details...

vv 12453
July 26 2019
0 0
> 

restart;

The integral

> 

J:=Int( exp(I*m*omega + I*b*cos(omega) ), omega=0..2*Pi) ;#          (m integer,  b positive constant)

Int(exp(I*m*omega+I*b*cos(omega)), omega = 0 .. 2*Pi)

 (1)
> 

value(J);  # wrong!

0

 (2)

 

> 

f:=expand(op(1,J)) assuming m::integer;

exp(I*m*omega)*exp(I*b*cos(omega))

 (3)

Prepare the change of variables

> 

g := expand(z^m * exp(I*b*(z+1/z)/2));  z = exp(I*omega);

z^m*exp(((1/2)*I)*b*z)*exp(((1/2)*I)*b/z)

  

z = exp(I*omega)

 (4)
> 

G:=simplify(g, [exp(I*omega)=z]) assuming m::integer

z^m*exp(((1/2)*I)*b*z)*exp(((1/2)*I)*b/z)

 (5)
> 

answer = 2*Pi*residue(G/z, z=0); # Residue theorem, z=0 essential singularity

answer = 2*Pi*residue(z^m*exp(((1/2)*I)*b*z)*exp(((1/2)*I)*b/z)/z, z = 0)

 (6)
> 

fii:=convert(G, Sum, dummy=i);  # we need the constant (coeff of z^0):

z^m*(Sum((((1/2)*I)*b*z)^i/factorial(i), i = 0 .. infinity))*(Sum((((1/2)*I)*b/z)^i/factorial(i), i = 0 .. infinity))

 (7)
> 

ai:=eval(op([2,1], fii), z=1); ai_bis:=eval(op([3,1], fii), z=1);

(((1/2)*I)*b)^i/factorial(i)

  

(((1/2)*I)*b)^i/factorial(i)

 (8)
> 

Sum(ai * subs(i=i+m, ai), i=0..infinity);  # coeff of z^(-m) in fii

Sum((((1/2)*I)*b)^i*(((1/2)*I)*b)^(i+m)/(factorial(i)*factorial(i+m)), i = 0 .. infinity)

 (9)
> 

answer := simplify(2*Pi*value(%)) assuming m::nonnegint;

2*Pi*BesselJ(m, b)*exp(((1/2)*I)*m*Pi)

 (10)
> 

evalf(eval([J,answer],[m=3,b=7/3])); #check1

[-0.2869900475e-12-1.168527755*I, -1.168527755*I]

 (11)
> 

evalf(eval([J,answer],[m=2,b=5]));   #check2

[-.2925772544-0.3974126502e-12*I, -.2925772544]

 (12)
> 

evalf(eval([J,answer],[m=0,b=1]));   #check3

[4.807878861+0.3488081236e-12*I, 4.807878862]

 (13)
> 

evalf(eval([J,answer],[m=7,b=4]));   #check4

[-0.2695364704e-12-0.9535405641e-1*I, -0.9535405642e-1*I]

 (14)

assuming......

vv 12453
July 23 2019
0 0

@Carl Love 

I assumed OP knew what he was doing. The original problem was not stated.

constants...

vv 12453
July 23 2019
0 0

@9009134 

Your constants are not assigned. After giving them numeric values, you may try to compute the integrals, but the integrands being very "odd", you may need to fine-tune int as above.

@9009134  Try Lps:=expand(L[p, s...

vv 12453
July 22 2019
0 0

@9009134 

Try

Lps:=expand(L[p, s](r, theta, phi)*F[p,1](r, theta));
kk:=int(Lps, [theta=0..2*Pi,r=0.5..1], numeric, epsilon=1e-7);
# kk:=int( L[p, s](r, theta, phi)*F[p,1](r, theta), [theta=0..2*Pi,r=0.5..1], numeric);

 

?...

vv 12453
July 21 2019
0 0

@9009134 

You keep modifying the worksheet and have eliminated my code.

Now you have functions of 3 variables and want to integrate numerically wrt 2 of them.
It wil not work of course.

First 48 49 50 51 52 53 54 Last Page 50 of 166