Advanced Search

vv

14122 Reputation

20 Badges

10 years, 203 days
Contact vv

MaplePrimes Activity


These are answers submitted by vv

sort...

vv 14122
April 27 2017
3 6

You need sort (not Groebner) for this.

f:= u[1]^2+u[2]+v[3]:
sort(f, order=plex(v[3], u[2],u[1]) );
        

I don't understand what you mean by u[i]=v[i]; u[i], v[j] are supposed to be the indeterminates.

 

Simp with evala...

vv 14122
April 27 2017
1 6

 

> 

restart;
> 

Simp:=proc(E)
convert(E,RootOf):
evala(%):
convert(%,radical):
simplify(%) assuming positive
end:
> 

f:= i -> -(2*(((A*sqrt(s)*a^(2*i+2)-a^(2*i)*B*L*dz*(sqrt(s)*c+s^(3/2)))*2^i*sqrt(s+c)*sqrt(L^2*s*(s+c)*dz^2+4*a^2)-(s+c)*(A*L*s*dz*2^i-B*2^(i+1))*a^(2*i+2)+(2^(i+1)*c*s+2^i*(c^2+s^2))*a^(2*i)*dz^2*L^2*B*s)*(1/(-dz*L*sqrt(s*(s+c)*(L^2*c*dz^2*s+L^2*dz^2*s^2+4*a^2))+L^2*s*(s+c)*dz^2+2*a^2))^i+(dz*L*sqrt(s*(s+c)*(L^2*c*dz^2*s+L^2*dz^2*s^2+4*a^2))+L^2*s*(s+c)*dz^2+2*a^2)^(-i)*((A*sqrt(s)*a^(2*i+2)-a^(2*i)*B*L*dz*(sqrt(s)*c+s^(3/2)))*2^i*sqrt(s+c)*sqrt(L^2*s*(s+c)*dz^2+4*a^2)+(s+c)*(A*L*s*dz*2^i-B*2^(i+1))*a^(2*i+2)-(2^(i+1)*c*s+2^i*(c^2+s^2))*a^(2*i)*dz^2*L^2*B*s)))*a^2/(sqrt(L^2*s*(s+c)*dz^2+4*a^2)*sqrt(s+c)*sqrt(s)*(dz*L*sqrt(s)*sqrt(s+c)*sqrt(L^2*s*(s+c)*dz^2+4*a^2)-L^2*s*(s+c)*dz^2-2*a^2)*(dz*L*sqrt(s)*sqrt(s+c)*sqrt(L^2*s*(s+c)*dz^2+4*a^2)+L^2*s*(s+c)*dz^2+2*a^2)):
> 

for i to 5 do f[i]=Simp(f(i)) od;

f[1] = A

  

f[2] = (A*L^2*c*dz^2*s+A*L^2*dz^2*s^2+B*L*c*dz+B*L*dz*s+A*a^2)/a^2

  

f[3] = (A*L^4*c^2*dz^4*s^2+2*A*L^4*c*dz^4*s^3+A*L^4*dz^4*s^4+B*L^3*c^2*dz^3*s+2*B*L^3*c*dz^3*s^2+B*L^3*dz^3*s^3+3*A*L^2*a^2*c*dz^2*s+3*A*L^2*a^2*dz^2*s^2+2*B*L*a^2*c*dz+2*B*L*a^2*dz*s+A*a^4)/a^4

  

f[4] = (A*L^6*c^3*dz^6*s^3+3*A*L^6*c^2*dz^6*s^4+3*A*L^6*c*dz^6*s^5+A*L^6*dz^6*s^6+B*L^5*c^3*dz^5*s^2+3*B*L^5*c^2*dz^5*s^3+3*B*L^5*c*dz^5*s^4+B*L^5*dz^5*s^5+5*A*L^4*a^2*c^2*dz^4*s^2+10*A*L^4*a^2*c*dz^4*s^3+5*A*L^4*a^2*dz^4*s^4+4*B*L^3*a^2*c^2*dz^3*s+8*B*L^3*a^2*c*dz^3*s^2+4*B*L^3*a^2*dz^3*s^3+6*A*L^2*a^4*c*dz^2*s+6*A*L^2*a^4*dz^2*s^2+3*B*L*a^4*c*dz+3*B*L*a^4*dz*s+A*a^6)/a^6

  

f[5] = (A*L^8*c^4*dz^8*s^4+4*A*L^8*c^3*dz^8*s^5+6*A*L^8*c^2*dz^8*s^6+4*A*L^8*c*dz^8*s^7+A*L^8*dz^8*s^8+B*L^7*c^4*dz^7*s^3+4*B*L^7*c^3*dz^7*s^4+6*B*L^7*c^2*dz^7*s^5+4*B*L^7*c*dz^7*s^6+B*L^7*dz^7*s^7+7*A*L^6*a^2*c^3*dz^6*s^3+21*A*L^6*a^2*c^2*dz^6*s^4+21*A*L^6*a^2*c*dz^6*s^5+7*A*L^6*a^2*dz^6*s^6+6*B*L^5*a^2*c^3*dz^5*s^2+18*B*L^5*a^2*c^2*dz^5*s^3+18*B*L^5*a^2*c*dz^5*s^4+6*B*L^5*a^2*dz^5*s^5+15*A*L^4*a^4*c^2*dz^4*s^2+30*A*L^4*a^4*c*dz^4*s^3+15*A*L^4*a^4*dz^4*s^4+10*B*L^3*a^4*c^2*dz^3*s+20*B*L^3*a^4*c*dz^3*s^2+10*B*L^3*a^4*dz^3*s^3+10*A*L^2*a^6*c*dz^2*s+10*A*L^2*a^6*dz^2*s^2+4*B*L*a^6*c*dz+4*B*L*a^6*dz*s+A*a^8)/a^8

 (1)
> 

 

Floats problem...

vv 14122
April 26 2017
1 0


 

Maple gives the correct result if the parameters are rational.
Let us examine what happens for floats.

The problem reduces to compute

> 

Int(f, t=0..6*Pi);

Int(f, t = 0 .. 6*Pi)

 (1)
> 

#where
> 

f:= sqrt((-.4000000000*sin(t)-.5333333334*sin(.6666666668*t))^2+(.4000000000*cos(t)-.5333333334*cos(.6666666668*t))^2);

((-.4000000000*sin(t)-.5333333334*sin(.6666666668*t))^2+(.4000000000*cos(t)-.5333333334*cos(.6666666668*t))^2)^(1/2)

 (2)
> 

int(f,t=0..6*Pi);  #due to wrong manipulation of floats

12.67823876+0.*I

 (3)
> 

int(f,t=0..6*Pi,numeric);  #ok

11.52567160

 (4)
> 

F:=int(f,t);  #discontinuous antiderivative and odd manipulation of floats!!!

-1.120000000*((0.2133333334e20*cos(.8333333334*t)^2-0.2177777778e20)*(cos(.8333333334*t)^2-1.))^(1/2)*(-1.*cos(.8333333334*t)^2+1.)^(1/2)*EllipticE(cos(.8333333334*t), .9897433186)/((0.2133333334e20*cos(.8333333334*t)^4-0.4311111112e20*cos(.8333333334*t)^2+0.2177777778e20)^(1/2)*sin(.8333333334*t))

 (5)
> 

plot(F,t=0..6*Pi);

 
> 

int(convert(f,rational),t=0..6*Pi);evalf(%); #ok, correct approach

(56/5)*EllipticE((4/7)*3^(1/2))

  

11.52567160

 (6)
> 

 

The roots...

vv 14122
April 26 2017
2 2
> 

restart;
> 

Det1 := (1/256)*(Aiso*(c+t)^2*(a^2+b^2)*(mu-1)*Pi^2-4*a^2*b^2*c*Gc)*(16*Aiso^2*Do^2*(c+t)^4*(a^2+b^2)^6*Pi^12+10*(a^2+b^2)^5*((c+t)^2*Aiso+4*Do)*Gc*(c+t)^2*c*a^2*Do*b^2*Aiso*Pi^10+(a^2+b^2)^4*((c+t)^2*Aiso+4*Do)^2*Gc^2*c^2*a^4*b^4*Pi^8-(1024/81)*a^6*Aiso^2*b^6*Tcr^2*(c+t)^4*(a^2+b^2)^2*Pi^4-(2560/81)*a^8*Aiso*b^8*c*Gc*Tcr^2*(c+t)^2*(a^2+b^2)*Pi^2-(1024/81)*a^10*b^10*c^2*Gc^2*Tcr^2)*(Aiso*(c+t)^2*(a^2+4*b^2)*(mu-1)*Pi^2-4*a^2*b^2*c*Gc)*(Aiso*(c+t)^2*(a^2+b^2)*(mu-1)*Pi^2-a^2*b^2*c*Gc)*(16*(c+t)^4*(a^2+(1/4)*b^2)^3*Do^2*(a^2+4*b^2)^3*Aiso^2*Pi^12+(10*(a^2+b^2))*((c+t)^2*Aiso+4*Do)*Gc*(c+t)^2*(a^2+(1/4)*b^2)^2*c*a^2*Do*(a^2+4*b^2)^2*b^2*Aiso*Pi^10+((c+t)^2*Aiso+4*Do)^2*Gc^2*(a^2+(1/4)*b^2)^2*c^2*a^4*(a^2+4*b^2)^2*b^4*Pi^8-(1024/81)*(c+t)^4*(a^2+(1/4)*b^2)*Tcr^2*a^6*(a^2+4*b^2)*b^6*Aiso^2*Pi^4-(2560/81)*a^8*Aiso*b^8*c*Gc*Tcr^2*(c+t)^2*(a^2+b^2)*Pi^2-(1024/81)*a^10*b^10*c^2*Gc^2*Tcr^2)*((mu-1)*(c+t)^2*(a^2+(1/4)*b^2)*Aiso*Pi^2-a^2*b^2*c*Gc)/(b^20*a^20*(c+t)^16):
> 

ec:=subs(Tcr^2=T,Det1):
> 

S:=simplify ~([solve(ec,T)]); #the roots squared

[(81/1024)*((a^2*b^2*c*Gc+8*Do*Pi^2*(a^2+b^2))*(c+t)^2*Aiso+4*a^2*b^2*c*Do*Gc)*Pi^8*(a^2+b^2)^4*((c+t)^2*(a^2*b^2*c*Gc+2*Do*Pi^2*(a^2+b^2))*Aiso+4*a^2*b^2*c*Do*Gc)/(((1/2)*Pi^2*(c+t)^2*(a^2+b^2)*Aiso+a^2*b^2*c*Gc)*(2*Pi^2*(c+t)^2*(a^2+b^2)*Aiso+a^2*b^2*c*Gc)*b^6*a^6), (81/1024)*(a^2+4*b^2)^2*((a^2*b^2*c*Gc+8*Do*Pi^2*(a^2+(1/4)*b^2))*(c+t)^2*Aiso+4*a^2*b^2*c*Do*Gc)*Pi^8*((a^2*b^2*c*Gc+2*Do*Pi^2*(a^2+4*b^2))*(c+t)^2*Aiso+4*a^2*b^2*c*Do*Gc)*(a^2+(1/4)*b^2)^2/(((1/2)*Pi^2*(c+t)^2*(a^2+4*b^2)*Aiso+a^2*b^2*c*Gc)*b^6*(2*Pi^2*(c+t)^2*(a^2+(1/4)*b^2)*Aiso+a^2*b^2*c*Gc)*a^6)]

 (1)
> 

AllSols:=sqrt~(S)[],-sqrt~(S)[]:
> 

 

Check...

vv 14122
April 26 2017
1 0


 

> 

restart;
> 

f:=x^(2-s)*(x^(-s)*GAMMA(3-s)*GAMMA(2-2*s)/(GAMMA(2-s)*GAMMA(3-2*s))+x^(s-2)*GAMMA(3-s)*GAMMA(2*s-2)/GAMMA(s))/(2-s)+(1/2)*(2*s*x-x+1)*(x+1)/(((x+1)^s)^2*(2*s^2-3*s+1))+x^(1-s)*(x^(-s)*GAMMA(2-s)*GAMMA(-2*s+1)/(GAMMA(1-s)*GAMMA(2-2*s))+x^(-1+s)*GAMMA(2-s)*GAMMA(2*s-1)/GAMMA(s))/(1-s)+(x+1)/((2*s-1)*((x+1)^s)^2)+x/((-1+s)*x^s)-(x+1)/((-1+s)*(x+1)^s);

x^(2-s)*(x^(-s)*GAMMA(3-s)*GAMMA(2-2*s)/(GAMMA(2-s)*GAMMA(3-2*s))+x^(s-2)*GAMMA(3-s)*GAMMA(2*s-2)/GAMMA(s))/(2-s)+(1/2)*(2*s*x-x+1)*(x+1)/(((x+1)^s)^2*(2*s^2-3*s+1))+x^(1-s)*(x^(-s)*GAMMA(2-s)*GAMMA(-2*s+1)/(GAMMA(1-s)*GAMMA(2-2*s))+x^(-1+s)*GAMMA(2-s)*GAMMA(2*s-1)/GAMMA(s))/(1-s)+(x+1)/((2*s-1)*((x+1)^s)^2)+x/((-1+s)*x^s)-(x+1)/((-1+s)*(x+1)^s)

 (1)
> 

f:=simplify(f) assuming s<1,s>1/2;

(1/8)*(8*(x+1)^(-2*s)*sin(Pi*s)*GAMMA(s)*s*x^2+16*(x+1)^(-2*s)*sin(Pi*s)*GAMMA(s)*s*x-4*(x+1)^(-2*s)*sin(Pi*s)*GAMMA(s)*x^2-16*(x+1)^(-s)*sin(Pi*s)*GAMMA(s)*s*x+2*Pi^(1/2)*GAMMA(s-1/2)*4^s*s^2+16*sin(Pi*s)*GAMMA(s)*x^(1-s)*s-8*x^(-2*s+1)*sin(Pi*s)*GAMMA(s)*s+8*(x+1)^(-2*s)*sin(Pi*s)*GAMMA(s)*s-8*(x+1)^(-2*s)*sin(Pi*s)*GAMMA(s)*x-16*(x+1)^(-s)*sin(Pi*s)*GAMMA(s)*s+8*(x+1)^(-s)*sin(Pi*s)*GAMMA(s)*x-8*x^(2-2*s)*sin(Pi*s)*GAMMA(s)*s-3*Pi^(1/2)*GAMMA(s-1/2)*4^s*s-8*sin(Pi*s)*GAMMA(s)*x^(1-s)+8*x^(-2*s+1)*sin(Pi*s)*GAMMA(s)-4*(x+1)^(-2*s)*sin(Pi*s)*GAMMA(s)+8*(x+1)^(-s)*sin(Pi*s)*GAMMA(s)+4*x^(2-2*s)*sin(Pi*s)*GAMMA(s)+Pi^(1/2)*4^s*GAMMA(s-1/2))/((2*s-1)*(-1+s)*sin(Pi*s)*GAMMA(s))

 (2)
> 

 

> 

F:=asympt(f,x,3) assuming s<1,s>1/2;

(1/8)*(4*sin(Pi*s)*GAMMA(s)-8*sin(Pi*s)*GAMMA(s)*s)*(1/x)^(2*s-2)/((2*s-1)*(-1+s)*sin(Pi*s)*GAMMA(s))+(1/8)*(16*sin(Pi*s)*GAMMA(s)*s-8*sin(Pi*s)*GAMMA(s))*(1/x)^(-1+s)/((2*s-1)*(-1+s)*sin(Pi*s)*GAMMA(s))+(1/8)*(2*Pi^(1/2)*GAMMA(s-1/2)*4^s*s^2+Pi^(1/2)*4^s*GAMMA(s-1/2)-3*Pi^(1/2)*GAMMA(s-1/2)*4^s*s)/((2*s-1)*(-1+s)*sin(Pi*s)*GAMMA(s))+O((1/x)^(2*s-1))

 (3)
> 

M:=MultiSeries:-asympt(f,x,3)  assuming s<1,s>1/2;

(1/8)*Pi^(1/2)*4^s*GAMMA(s-1/2)/(GAMMA(s)*sin(Pi*s))-2*s*(1/x)^(2*s-1)/(2*s-1)+(1/x)^s+O((1/x)^(2*s))

 (4)
> 

evalf(eval(F,[s=3/4,x=10000]));

162.6220576+O((1/10000)*10000^(1/2))

 (5)
> 

evalf(eval(M,[s=3/4,x=10000]));

2.593057555+O((1/100000000)*10000^(1/2))

 (6)
> 

evalf(eval(f,[s=3/4,x=10000]));

2.593057712

 (7)
> 

# So, it seems that asympt is wrong but MultiSeries:-asympt is OK.
# Note that it is easy to obtain the asymptotic by hand developing just (1+x)^(-s) and (1+x)^(-2*s).
> 

 

solve,identity...

vv 14122
April 24 2017
2 4

f:=-4*a*sin(2*x)+4*b*cos(2*x)=-8*sin(2*x)+20*cos(2*x):
solve( identity(f,x), {a,b} );
                         {a = 2, b = 5}


Ad hoc methods are also possible:
solve({ eval(f,x=0), eval(diff(f,x),x=0)});
or
solve({ eval(f,x=0), eval(f,x=Pi/4)});

convert,rational...

vv 14122
April 23 2017
2 1

restart;
x:= rand(10^1000..10^1200)()/rand(10^1000..10^1200)():
n:=100:
v:=convert(evalf[n](x), rational, n);

evalf[10](evalf[5*n](x-v));
     -2.407569284*10^(-99)

 

evalindets...

vv 14122
April 23 2017
2 0

For example:

evalindets(Tours, list(integer), u -> [1,u[]]);

color...

vv 14122
April 22 2017
0 0


 

> 

restart;
> 

with(plots):
with(plottools):

z := (m+I*n)/(p+I*q):
g := proc (z) options operator, arrow; (z-I)/(z+I) end proc:
bz := simplify(evalc(Im(z))):
a := simplify(evalc(Re(g(z)))):
b := simplify(evalc(Im(g(z)))):

r := 15;
Lista := Vector():
Listb := Vector():
C := Array(datatype = float[8]):
j := 1:
for m from -r to r do for n from -r to r do for p from -r to r do for q from -r to r do
if p <> 0 and q <> 0 and m^2-2*m*q+n^2+2*n*p+p^2+q^2 <> 0 and bz >= 0 then
   Lista(j) := a; Listb(j) := b;
   C(j) := a^2+b^2;  # color (HUE)
   j := j+1
end if
end do end do end do end do; j;

pointplot(Lista, Listb, symbol = point, symbolsize = 1, size = [1200, 1200], color = COLOR(HUE, C));
> 

 


 

Download C.mw

strange bug...

vv 14122
April 15 2017
1 5

It seems bo be a strange bug which crashes Maple and freezes Windows.
It seems to be caused by the fact that piecewise has floats and the rest of the expression does not.
A workaround:

fe11_1 := evalf(abs(-piecewise(1 <= p and p <= 9, 9.310043871-1.372270968*p+0.6222709675e-1*p^2, 0)+exp(-(1/40)*p*ln(3)-(1/4)*p*ln(2)+(1/40)*ln(3)+(1/10)*ln(p)+(13/4)*ln(2)))):
evalf(Int(fe11_1, p = 1 .. 11/3)+Int(fe11_1, p = 11/3 .. 19/3)+Int(fe11_1, p = 19/3 .. 9));

      1.342268852

If evalf in fe11_1 is removed ==> crash/freeze in Maple&Windows.

Combine symbolic+numeric...

vv 14122
April 13 2017
1 0

You final double integral is improper and it will be hard to compute it with desired accuracy.
A better idea is to use some approximations and symbolics.

ee:= t -> add(t^k/k!,k=0..12):
G1:=-0.9445379894:
f:= (x) -> 0.9/abs(x-0.4)^(1/3)+0.1/abs(x-0.6)^(1/2):
U1 := unapply(-ee(-x)*((int(f(t)*ee(t), t = 0 .. x))+G1)/2-ee(x)*((int(f(t)*ee(-t), t = 0 .. x))+G1)/2, x):
U:= unapply(-ee(x)/2*((int(f(t)*ee(-t),t=0..x))+G1)+ee(-x)/2*((int(f(t)*ee(t),t=0..x))+G1), x):
evalf(Int(U1(x)^2+U(x)^2-2*f(x)*U(x), x=0..1));

         -0.3753314046

Edit. For ee:= t -> add(t^k/k!,k=0..25):  and Digits:=25 ==>

-0.3753314045949344185803002

LaTeX...

vv 14122
April 11 2017
0 4

In my opinion Maple should be used for mathematics, not for typesetting. For this, use LaTeX.

IsIsomorphic from scratch...

vv 14122
April 10 2017
0 0

Here is a simple IsI. Not optimized; it will be slow for large graphs but could be used to test IsIsomorphic.

IsI:=proc(G1::GRAPHLN,G2::GRAPHLN,phi::name)
local A,A1,B,B1,f,k,n:=nops(op(3,G1));
A:=convert(op(4,G1),list);
B:=convert(op(4,G2),list);
if sort(nops~(A)) <> sort(nops~(B)) then return false fi;
for f in combinat:-permute(n) do
  A1:=subsindets(A, posint, k->f[k]);
  B1:=[seq(B[f[k]],k=1..n)];
  if A1=B1 then if nargs=3 then phi:=[seq(k=f[k],k=1..n)] fi;
     return true fi
od;
false
end:

with(GraphTheory):
G1:=Graph(8,Trail(1,2,3,4,5,6,7,8)):
G2:=Graph(8,Trail(2,3,4,5,6,7,8,1)):
IsI(G1,G2,'ff');ff;
                              true
    [1 = 1, 2 = 8, 3 = 7, 4 = 6, 5 = 5, 6 = 4, 7 = 3, 8 = 2]
IsIsomorphic(G1,G2,'ff');ff;
                              true
    [1 = 1, 2 = 8, 3 = 7, 4 = 6, 5 = 5, 6 = 4, 7 = 3, 8 = 2]

 

capped...

vv 14122
April 09 2017
1 0

A capped version.

ST:=proc(c::list(realcons), ur::range, r::realcons:=1, R::realcons:= 2)
local p1,p2,u,v,t,S:=[c[1]+(R+r*cos(u))*cos(v),c[2]+(R+r*cos(u))*sin(v),c[3]+r*sin(u)];
p1:=plot3d(S,u=ur,v=0..2*Pi);
p2:=plot3d(t*~eval(S,u=lhs(ur))+(1-t)*~eval(S,u=rhs(ur)), t=0..1,v=0..2*Pi);
plots:-display(p1,p2,_rest);
end proc:

plots:-display(ST([0, 0, 0], 0 .. Pi, 1, 2), lightmodel = light4, orientation = [-140, 60], scaling = constrained, style = patchnogrid);

 

plot(u)...

vv 14122
April 04 2017
1 2

Actually if u is undefined,
plot(u)  <==>  plot(u, u = -10 .. 10)

because -10 .. 10  is the default range.

First 91 92 93 94 95 96 97 Last Page 93 of 121