Advanced Search

vv

14077 Reputation

20 Badges

10 years, 69 days
Contact vv

MaplePrimes Activity


These are answers submitted by vv

Int...

vv 14077
September 07 2021
2 0

The computation is slow on my (old) computer. It needs several hours.
To  interrupt and then continue later, just save the partial results. 

restart;
read "d:/temp/Jfile.mpl";
max(indices(J,nolist));
smax := 20;  R := 1680;  k := 1/4;  gam := 1/2;    vmax:=R/2;
n:=0;
Digits:=15:
f:=proc(v,s) 
  evalf(Re(  sinh(s)*coth(1/2*s)^(2*I/k)/(gam-I*k*cosh(s))^5*exp(-(k^2*sinh(s)^2/(gam-I*
  k*cosh(s))+I*k*(1-cosh(s)))*v)*Hypergeom([I/k],[1],I*k*v)*v^n  ))
end:
V0:=0;  dV:=20;  # <-----
for V from V0 to vmax-dV by dV do
  J[V]:=evalf(1/R * Int(f, [V..V+dV, 0..smax], method=_CubaCuhre, epsilon=1e-5));
od;
save J, "d:/temp/Jfile.mpl";
add(entries(J,nolist));
max(indices(J,nolist));

Similar for the imaginary part if needed.

simplify/JacobiP...

vv 14077
September 07 2021
3 4

 

> 

restart;
> 

`simplify/JacobiP`:=proc(ex::algebraic)
  evalindets(ex, specfunc({int,Int}),
    proc(J)
      local J11,J12, n,a,b,x, m;
      if not type(op(1,J), `*`) then return J fi;
      J11:=select(type, [op(op(1,J))], specfunc(JacobiP)^2 );
      J12:=select(type, [op(op(1,J))], specfunc(JacobiP));
      if not(J11<>[] xor J12<>[]) or nops(J11)>1 then return J fi;
      if J11<>[] then J12:=[op([1,1],J11)$2] fi;
      if nops(J12)<>2 or op(J12[1])[2..] <> op(J12[2])[2..] then return J fi;
      n,a,b,x := op(J12[1]); m:=op(1,J12[2]);
      if  not mul(J12)*(1-x)^a*(1+x)^b = op(1,J) or not (op(-1,J) = (x=-1..1))   then return  J fi;
      piecewise(m<>n, 0, 2^(a+b+1)*GAMMA(n+a+1)*GAMMA(n+b+1)/((2*n+a+b+1)*GAMMA(n+a+b+1)*n!))
    end proc
  )
end proc:
> 

 

> 

# Examples
> 

J:=Int(JacobiP(n, a, b, x)*JacobiP(m, a, b, x)*(1 - x)^a*(1 + x)^b, x=-1..1):
> 

simplify(J, JacobiP);

piecewise(m <> n, 0, 2^(a+b+1)*GAMMA(n+a+1)*GAMMA(n+b+1)/((2*n+a+b+1)*GAMMA(n+a+b+1)*factorial(n)))

 (1)
> 

J:=Int(JacobiP(n, a, b, x)*JacobiP(m, a, b, x)*(1 - x)^a*(1 + x)^b, x=-1..100):
> 

simplify(J, JacobiP);

Int(JacobiP(n, a, b, x)*JacobiP(m, a, b, x)*(1-x)^a*(1+x)^b, x = -1 .. 100)

 (2)
> 

J:=Int(JacobiP(n, a, b, x)^2*(1 - x)^a*(1 + x)^b, x=-1..1);

Int(JacobiP(n, a, b, x)^2*(1-x)^a*(1+x)^b, x = -1 .. 1)

 (3)
> 

simplify(J, JacobiP);

2^(a+b+1)*GAMMA(n+a+1)*GAMMA(n+b+1)/((2*n+a+b+1)*GAMMA(n+a+b+1)*factorial(n))

 (4)
> 

J:=int(JacobiP(3, a, b, x)*JacobiP(n, a, b, x)*(1 - x)^a*(1 + x)^b, x=-1..1);

int(JacobiP(3, a, b, x)*JacobiP(n, a, b, x)*(1-x)^a*(1+x)^b, x = -1 .. 1)

 (5)
> 

simplify(J, JacobiP);

piecewise(n <> 3, 0, 2^(a+b+1)*GAMMA(4+a)*GAMMA(4+b)/((6*(7+a+b))*GAMMA(4+a+b)))

 (6)
> 

simplify((J+1)^2+exp(7*x), JacobiP);

(piecewise(n <> 3, 0, 2^(a+b+1)*GAMMA(4+a)*GAMMA(4+b)/((6*(7+a+b))*GAMMA(4+a+b)))+1)^2+exp(7*x)

 (7)
> 

simplify((J+1)^2+exp(7*x), JacobiP) assuming n>3;

1+exp(7*x)

 (8)
> 

 


Download SimpJP-vv.mw

By hand...

vv 14077
September 06 2021
2 0

                      

Solution...

vv 14077
September 02 2021
3 2
> 

 

> 

restart;
> 

d2:=(A,B)-> (A[1]-B[1])^2 + (A[2]-B[2])^2:
area:=(A,B,C) -> LinearAlgebra:-Determinant( <Matrix([A,B,C])|<1,1,1>> )/2:
> 

A:=[0,0];
M:=[0,-7/sqrt(2)];
N:=[0,7/sqrt(2)];
S:=[7/sqrt(2),0];
P:=M+11/7*(S-M);

[0, 0]

  

[0, -(7/2)*2^(1/2)]

  

[0, (7/2)*2^(1/2)]

  

[(7/2)*2^(1/2), 0]

  

[(11/2)*2^(1/2), 2*2^(1/2)]

 (1)
> 

T:=P+k*~(N-P);
k:=solve( T[1]^2+T[2]^2-49/2 )[2];

[-(11/2)*k*2^(1/2)+(11/2)*2^(1/2), (3/2)*k*2^(1/2)+2*2^(1/2)]

  

44/65

 (2)
> 

B:=[x,y]: solve([d2(B,T)=d2(B,P),d2(B,P)=d2(B,S)],[x,y]):
B:=eval([x,y],%[]);

[(7/2)*2^(1/2), 2*2^(1/2)]

 (3)
> 

Q:=2*B-P;

[(3/2)*2^(1/2), 2*2^(1/2)]

 (4)
> 

#plot([M,N,P,Q,M]);
> 

Area__MNPQ = area(Q,N,M)+area(Q,P,N);

Area__MNPQ = 33/2

 (5)
> 

 

 

0...

vv 14077
August 31 2021
1 0

F = 0, G = 0 is obviously a solution (gamma = anything).
Depending on your coefficients, this could be the only one.

Answer...

vv 14077
August 27 2021
2 0 
A:=<-2,3,-5>: B:=<-6,1,-1>: C:=<2,-3,7>:
with(LinearAlgebra): local D:
k:=Norm(B-A,2) / Norm(C-A,2):
D:=(B + k*C)/(1+k);
eqAD:=A + t*(D-A); # parametric equation

Definite int...

vv 14077
August 26 2021
3 2

The first problem is that you cannot have t=1:

fN(0.5, 1);
Error, (in solnproc) unable to compute solution for t>HFloat(0.2):
Newton iteration is not converging

So, let's take  t := 0.1;

It is not possible to compute the indefinite integral numerically. In your case, what you want is a definite integral.

t := 0.1;
                            t := 0.1

A1:=x*R(z)*R(z)*(fN)(x, t);
                A1 := 0.8692871388 x fN(x, 0.1)

A2:= X -> int(A1, x=0..X):
A2(0); # 0, as you want
                               0.

plot(A2, 0..1);

 

Answer...

vv 14077
August 25 2021
1 2

Actually the equation  chi(i) = a (a>=0) may have up to 4 solutions. (The first one for i<=0.)

> 

restart
> 

p := proc (S) options operator, arrow; .32*exp(S)/(1.2+1.901885*exp(S)-S^3) end proc

proc (S) options operator, arrow; .32*exp(S)/(1.2+1.901885*exp(S)-S^3) end proc

 (1)
> 

eta := proc (S) options operator, arrow; (D(p))(S)*S/p(S) end proc

proc (S) options operator, arrow; (eval(diff(p(x), x), x = S))*S/p(S) end proc

 (2)
> 

f := proc (k) options operator, arrow; A*k^alpha end proc

proc (k) options operator, arrow; A*k^alpha end proc

 (3)
> 

w := proc (k) options operator, arrow; f(k)-(D(f))(k)*k end proc

proc (k) options operator, arrow; f(k)-(eval(diff(f(x), x), x = k))*k end proc

 (4)
> 

beta := 1

1

 (5)
> 

A := 28

28

 (6)
> 

alpha := .33

.33

 (7)
> 

chi := proc (i) options operator, arrow; i*((1+beta)/beta+eta(i))/(1+eta(i)) end proc

proc (i) options operator, arrow; i*((1+beta)/beta+eta(i))/(1+eta(i)) end proc

 (8)
> 

#ss:=simplify(chi(i));
> 

#plot(ss, i=-2..18)
> 

#sort(select(u -> Im(u)=0,[solve(ss,i),solve(1/ss,i)]));
> 

#select(u -> (u>6), [solve(diff(ss,i),i)]);
> 

ii:=[-1.059280978, -0.7553094684, 0, 3.612421353, 6.399747910, 9.270348450, 9.270348450, 50]:
> 

chin:=proc(a)
  local b:= evalf(RootOf(200000*_Z^5 - 200000*_Z^4*a - 200000*_Z^4 + 400000*_Z^3*a - 760754*exp(_Z)*_Z + 380377*exp(_Z)*a - 240000*_Z^2 + 240000*_Z*a - 480000*_Z + 240000*a, _Z, ii[2*n-1] .. ii[2*n]));
`if`(b::numeric,b,undefined)
end:
> 

NULL
> 

Hn:=k -> chin(w(k)):
Vn := k -> ln((w(k) - Hn(k))*(A*alpha)^beta*(p(Hn(k))*Hn(k))^(alpha*beta)):
> 

for n from 1 to 4 do
  plot(Vn, 1 .. 2, title=Plot||n);
od;

Warning, unable to evaluate the function to numeric values in the region; see the plotting command's help page to ensure the calling sequence is correct

  

 

 

 

 


Download Worksheet_1_vv.mw

evaln...

vv 14077
August 22 2021
2 4

If you really want to create the (global) m||i variables even when they exist, just replace

M := [seq(m || i, i = 2 .. l)];

with

M := [seq(evaln(m || i), i = 2 .. l)];


 

quadratic...

vv 14077
August 21 2021
2 2

It's a quadratic equation. Just use solve.

solve((a*y-x*(a+1))/(x^(2)-b*a^(2)*(y-x)^(2))=A, y);

 

Answer...

vv 14077
August 20 2021
1 0

With a little help, Maple can compute the coefficients.

1. Expand to a series around z=0 the function 1/(1 + z)^2.

2. Substitute here z = exp(x-5*t) and expand around t=0.

3. Compute the sum for the coefficient of t^k.

 

One obtains the coefficient

 

> 

c := k -> (-1)^k*(polylog(-1 - k, -exp(x)) + polylog(-k, -exp(x)))*5^k/k!;

proc (k) options operator, arrow; (-1)^k*(polylog(-1-k, -exp(x))+polylog(-k, -exp(x)))*5^k/factorial(k) end proc

 (1)
> 

add(map(combine,simplify(c(k)))*t^k,k=0..6); # check

-exp(x)*(2+exp(x))/(1+exp(x))^2+10*exp(x)*t/(1+exp(x))^3+(50*exp(2*x)-25*exp(x))*t^2/(1+exp(x))^4+(125/3)*exp(x)*(4*exp(2*x)-7*exp(x)+1)*t^3/(1+exp(x))^5+(625/12)*exp(x)*(8*exp(3*x)-33*exp(2*x)+18*exp(x)-1)*t^4/(1+exp(x))^6+(625/12)*exp(x)*(16*exp(4*x)-131*exp(3*x)+171*exp(2*x)-41*exp(x)+1)*t^5/(1+exp(x))^7+(3125/72)*exp(x)*(32*exp(5*x)-473*exp(4*x)+1208*exp(3*x)-718*exp(2*x)+88*exp(x)-1)*t^6/(1+exp(x))^8

 (2)

 

asympt...

vv 14077
August 20 2021
1 0

@Michael_Watson It seems that you want an asymptotic expansion. Note that you missed each time a paranthesis.

> 

restart;
> 

L/(z*L*sqrt((z/L)^2 + 1)): % = asympt(%,L,2);

1/(z*(z^2/L^2+1)^(1/2)) = 1/z+O(1/L^2)

 (1)
> 

L/(z*z*sqrt(1+(L/z)^2)): % = asympt(%,L,2) assuming z>0;

L/(z^2*(1+L^2/z^2)^(1/2)) = 1/z+O(1/L^2)

 (2)
> 

L/(z*sqrt(z^2 + L^2)): % = asympt(%,L,2);

L/(z*(L^2+z^2)^(1/2)) = 1/z+O(1/L^2)

 (3)


Download Lz.mw

Answer...

vv 14077
August 19 2021
2 3 
EQSYS:=proc(S1::set(`=`), S2::set(`=`), {vars:={x,y,z}, param:={t}}) 
  local A:=eliminate(S1,param)[2], B:=eliminate(S2,param)[2];
  evalb( eval(B,solve(A,vars)) union eval(A,solve(B,vars)) = {0});
end:

sys1:={16*x-2*y-11*z=0, 14*x-y-10*z=3}:
sys2:={(x-2)/3=(y-5)/2, (y-5)/2=(z-2)/4}:
sys2a:={ (x-2)/3=t, (y-5)/2=t, (z-2)/4=t }:
sys3:={x+y-z=2, x-z=77}:

EQSYS(sys1, sys2);
                              true

EQSYS(sys1, sys2a);
                              true

EQSYS(sys1, sys3);
                             false

 

inert...

vv 14077
August 17 2021
3 4

You may use inert operators

ex := L/(z*sqrt(z^2 + L^2)):
subs(L=(L%/z), numer(ex)) %* (1/denom(ex*z));

But, why do you need this?

Answer...

vv 14077
August 17 2021
2 2

It makes sense to use Maple only for numeric purposes here. 
For the two examples, a few numeric computations will suggest the general result, but the proof must be given by hand.

a) The limit  (M⊗N)(a,b) = sqrt(a*b). The proof is simple, as mmcdara showed.

M := (a,b) -> (a+b)/2:
N := (a,b) -> 2*a*b/(a+b):
MN:=proc(M::procedure, N::procedure, a, b, n::nonnegint, nmax::nonnegint:=n)
  local an:=a, bn:=b, ab:=Vector(1), k; 
  ab[1]:=[a,b]; 
  for k to nmax do  an, bn := evalf(M(an,bn)), evalf(N(an,bn));
    if k>=n then ab(k+1-n):=[an,bn] fi 
  od;
 seq(ab);
end:
MN(M,N,1,2, 2,5);             # the sequence [a1,b1],[a2,b2], ...

    [1.416666666, 1.411764706], [1.414215686, 1.414211438], [1.414213562, 1.414213562], [1.414213562, 1.414213562]

MN(M,N, 1,2, 10) = sqrt(2.);  # [a10,b10]

            [1.414213562, 1.414213562] = 1.414213562

b)  The limit M⊗N is the logarithmic mean; this can be proven similarly to the Gauss' proof for the AGM mean.

L := (a, b) -> (b-a)/ln(b/a):  # the logarithmic mean
a:=2.:  b:=8.:
MN( (a,b) -> (a+sqrt(a*b))/2, (a,b) -> (b+sqrt(a*b))/2,  a, b, 50 ) = L(a,b);

        [4.328085124, 4.328085124] = 4.328085123

First 24 25 26 27 28 29 30 Last Page 26 of 121