vv

12453 Reputation

19 Badges

9 years, 285 days

MaplePrimes Activity


These are answers submitted by vv

Seems to be an old bug.

with(GraphTheory):
G:=Graph( Trail(1,2,3,4,5) ):
DrawGraph(G): # ok
DrawPlanar(G); # error

Probably the graph G being "collinear", it's "too planar"  :-)

 

Actually you do not have integral equations, but expressions containing integrals depending on the unknowns.
You should proceed as in the next example:

restart;
Digits:=15:
eq1 := proc(x,y) local t; int(sin(x*t)/(t^2+7), t=0..(x+y), numeric)  - 1/10. end:
eq2 := proc(x,y) local t; int(exp(x*t)/(t^2+3), t=0..(x-y), numeric) + 0.13 end:
solxy:=fsolve([eq1,eq2], [1 .. 2.5,  1.5 .. 4]);

          solxy := [1.474881121279997081543, 2.084760793500574599373]
eq1(solxy[]), eq2(solxy[]); # check
        
-.4e-15, .1e-14

restart;

alias(p = p(x, t), q = q(x, t), u = u(x, t));

p, q, u

(1)

p_x := diff(p,x) = a*p + u*q; #   (1)

diff(p, x) = a*p+u*q

(2)

q_x := diff(q,x) =  -conjugate(u)*p - a*q; # (2)     # where a is complex parameter, p_{x} mean derivative of p w.r.t x

diff(q, x) = -conjugate(u)*p-a*q

(3)

`diff/conjugate` := (z, x) ->  conjugate(diff(z,x));  # needed, x :: real

proc (z, x) options operator, arrow; conjugate(diff(z, x)) end proc

(4)

#Define, u=p^2+conjugate(q)^2,  (3)  
#Now take the derivative of (3) w.r.t x and by using (1) and (2), we get

u := p^2+conjugate(q)^2;

p^2+conjugate(q)^2

(5)

ux:=diff(u,x);

2*p*(diff(p, x))+2*conjugate(q)*conjugate(diff(q, x))

(6)

eval(ux, [p_x, q_x]);

2*p*(a*p+(p^2+conjugate(q)^2)*q)+2*conjugate(q)*conjugate(-(q^2+conjugate(p)^2)*p-a*q)

(7)

ux:=expand(%);

2*a*p^2+2*q*p^3+2*p*q*conjugate(q)^2-2*conjugate(q)^3*conjugate(p)-2*conjugate(q)*p^2*conjugate(p)-2*conjugate(q)^2*conjugate(a)

(8)

# check wrt your result
u__x := 2*(a*p^2 - conjugate(a)*conjugate(q)^2) + 2*(p^2 + conjugate(q)^2)*(p*q - conjugate(p)*conjugate(q));
simplify(ux - u__x);

2*a*p^2-2*conjugate(q)^2*conjugate(a)+2*(p^2+conjugate(q)^2)*(p*q-conjugate(p)*conjugate(q))

 

0

(9)


Download diff_conj.mw

plots:-display(seq(seq(plottools:-polygon([[m,n],[m+1/2,n+1/2],[m+1,n],[m+1/2,n-1/2]],
               color=`if`(m::odd,red,blue)), m=0..3),n=0..3), axes=none);

Note that ignoring the color, another possibility is

plots:-inequal(abs(x-floor(x)-1/2) + abs(y-floor(y)-1/2)<=1/2, x=0..4, y=0..4);

Yes, ChainTools  and also SolveTools:-SemiAlgebraic  are too slow for this problem.
But after inspection and some numeric experiments it seems that the system has solution iff:

( a=1  and  1 <= b <= sqrt(4/3) )    or    (a=b>1).

In general, such a brute-force solution is used only if a direct one is not available. This is not the case.
For t=0 the number of the roots of the equation is either <= 4 or infinity (so 6 is out of the question).

It remains t=1 which was solved in your previous post [or, just form the 4 quadratic equations (obtained from abs(u) = +- u) and find the conditions to have integer roots]. The case t=2 is similar.

I see two problems in your worksheet.

1. The distribution corresponding to the RV MIXTURE is not correctly defined. It should be:

pdf := (p,x) -> p*PDF(U1, x)*Heaviside(x+3) + (1-p)*PDF(U2, x)*Heaviside(x-1);
cdf := (p,x) -> int(pdf(p,y), y=-infinity..x);

MIXTURE := p -> RandomVariable( Distribution( PDF = (x -> pdf(p,x)), CDF = (x -> cdf(p,x)), Conditions=[p >=0, p <=1] )):

2. It seems that Sample works only for distributions for which the support is an interval and the PDF is differentiable here.
The support of  MIXTURE is (-3,-1) union (1,3).
[I don't know how to declare such a support, or even if it's possible].
So,  Sample(MIXTURE(1/3), 10);  will fail, but

Sample(MIXTURE(1/3), 10, method=[envelope, range=1..3]);  # works.

Probably it's possible to use somehow method=custom.

I have set interface(typesetting=standard), because typesetting=extended  exhibits (here too) its weakness, not related to our problem.

1. Let's take the following version of your example:

restart;
interface(typesetting=standard):
#interface(typesetting=extended):
j0:=1; n0:=0;
V:=Vector(n0):
for j from 1 to 5 do
  try
    print('j'=j);
    V(j):=j/(j-j0)
  catch:
    print("Ante",lastexception); 
    V(j):=17; # clears lastexception apparently
    print("Post",lastexception);
  end try
end do:
V, lastexception;

It is equivalent to the original one.
V(j):=17  will clear lastexception iff j0 > n0. So, it seems that the extension of a rtable is inside of some try ... end try and lastexception is reset, but only locally.
Outside try .. end try,  lastexception has always the correct value.


2. For the second example, a workaround would be:

restart;
interface(typesetting=standard);
lastexception;
5/0;
lastexception;

5/5;
lastexception;
x:=47;
lastexception;
table();
lastexception;
[1,2,3];
lastexception;
# rtable();
try  rtable(); catch : end try ;
lastexception; # Now NOT cleared!

 

The arithmetic with infinity is very delicate.
The automatic simplification must be lite, so there are compromises.
When an expression is suspected to contain infinity,  some kind of simplification is recommended.
It is strange that even is must be sometimes preceded by simplify:

eq := Pi*infinity + I*infinity = infinity + I*infinity:   #  sqrt(2) too
is(eq), is(simplify(eq));

                          false, true

(Maybe this is a consequence of the fact that a constant for complex infinity is missing in Maple.)

Edit. I add some other stange results about arithmetic with oo:

limit( (x^2+1) + x*I, x=infinity );
                            infinity
limit( (x^2+1) + x^2*I, x=infinity );
                        (1 + I) infinity
limit( exp(x) + exp(x+1)*I, x=infinity );
                   -infinity I (-exp(1) + I)
simplify(%);
                   -infinity I (-exp(1) + I)
is(%=infinity+I*infinity);
                             false

The function f is strictly convex (because x -> x^2 is strictly convex and the other 4 terms are convex).
So f has at most 2 real solutions (and you want 6 integer ones).

msolve(x^3+y = 123, 124);
       
{x = x, y = 123*x^3 + 123}

The equation is trivial because x is obviously arbitrary.
Take a nontrivial one:

msolve(x^3+y^2 = 123, 124);
   {x = 7, y = 20}, {x = 7, y = 42}, {x = 7, y = 82},  {x = 7, y = 104}, {x = 11, y = 30}, {x = 11, y = 32},
    {x = 11, y = 92}, {x = 11, y = 94}, {x = 27, y = 30}, {x = 27, y = 32}, {x = 27, y = 92}, {x = 27, y = 94},
    {x = 35, y = 20}, {x = 35, y = 42}, {x = 35, y = 82}, {x = 35, y = 104}, {x = 39, y = 18}, {x = 39, y = 44},
    {x = 39, y = 80}, {x = 39, y = 106}, {x = 43, y = 10}, {x = 43, y = 52}, {x = 43, y = 72}, {x = 43, y = 114},
    {x = 51, y = 20}, {x = 51, y = 42}, {x = 51, y = 82}, {x = 51, y = 104}, {x = 55, y = 30}, {x = 55, y = 32},
    {x = 55, y = 92}, {x = 55, y = 94}, {x = 71, y = 18}, {x = 71, y = 44}, {x = 71, y = 80}, {x = 71, y = 106},
    {x = 83, y = 10}, {x = 83, y = 52}, {x = 83, y = 72}, {x = 83, y = 114}, {x = 91, y = 10}, {x = 91, y = 52},
    {x = 91, y = 72}, {x = 91, y = 114}, {x = 99, y = 0}, {x = 99, y = 62}, {x = 107, y = 18}, {x = 107, y = 44},
    {x = 107, y = 80}, {x = 107, y = 106}, {x = 119, y = 0}, {x = 119, y = 62}, {x = 123, y = 0}, {x = 123, y = 62}

Inside a loop, the two statement separators (semicolon and colon) are equivalent. Ony the terminator after end counts (semicolon in your case). See:
?for
?;
?printlevel

I don't know what you mean by "not well posed" in this context. I think that it is more important to state clearly the class of the function.

So, if u is supposed to be differentiable in (0, oo) x (0, oo) and continuous in [0, oo) x [0, oo) then the problem has a unique solution:

u(x, t) =  piecewise(x>t, exp(-(x-t)^2), 1);

Maple does not find it of course, but does half of the job, by solving the pde alone.
(bc is interpreted as D[1](u)(0+, t) = 0 for t>0).

In Maple 2019 does
(1/2)*(int(cosh(cos(2*t)), t = 0 .. 2*Pi) + int(cosh(-cos(2*t)), t = 0 .. 2*Pi));
produce 2*Pi*BesselI(0, 1)?  Normally, it should simplify cosh(-cos(2*t)) to cosh(cos(2*t)).

In Maple 2018 I have to manipulate like this:

 

J:=Int(cosh(cos(2*t)), t = 0 .. 2*Pi);

Int(cosh(cos(2*t)), t = 0 .. 2*Pi)

(1)

J1:=IntegrationTools:-Change(J, 2*t=u);

(1/2)*(Int(cosh(cos(u)), u = 0 .. 4*Pi))

(2)

J2:=simplify(convert(J1,exp));

(1/4)*(Int(exp(cos(u))+exp(-cos(u)), u = 0 .. 4*Pi))

(3)

IntegrationTools:-Expand(J2);

(1/4)*(Int(exp(cos(u)), u = 0 .. 4*Pi))+(1/4)*(Int(1/exp(cos(u)), u = 0 .. 4*Pi))

(4)

value(simplify(%));

2*Pi*BesselI(0, 1)

(5)

 

 


 

restart

f := ((1 - a)^2 + a^2*((1 - exp(-y))*(1 - exp(-x)) - 2 + exp(-x) + exp(-y)) + a*(2 - exp(-x) - exp(-y) + (1 - exp(-y))*(1 - exp(-x))))/(1 - a*exp(-x)*exp(-y))^3;

((1-a)^2+a^2*((1-exp(-y))*(1-exp(-x))-2+exp(-x)+exp(-y))+a*(2-exp(-x)-exp(-y)+(1-exp(-y))*(1-exp(-x))))/(1-a*exp(-x)*exp(-y))^3

(1)

F:=simplify(eval(f, a=3/10));

(-390*exp(-x-y)+600*exp(-x)+600*exp(-y)-1300)/(-10+3*exp(-x-y))^3

(2)

#s := 2*evalf(int((int(f*exp(-x)*exp(-y), x = 0 .. y + t,AllSolutions)), y = 0 .. infinity,AllSolutions)) assuming real ;

S:=2*Int(F*exp(-x)*exp(-y), x = 0 .. y + t,y = 0 .. infinity);

2*(Int((-390*exp(-x-y)+600*exp(-x)+600*exp(-y)-1300)*exp(-x)*exp(-y)/(-10+3*exp(-x-y))^3, x = 0 .. y+t, y = 0 .. infinity))

(3)

plot(S, t=-1..4)

 

#solve(-10 + 3*exp(-t))

ans:=value(S) assuming t>-ln(10/3);

(1/3)*(-9*exp(-(1/2)*t)*30^(1/2)*arccoth((1/3)*exp((1/2)*t)*30^(1/2))-160*exp((3/2)*t)*30^(1/2)*arccoth((1/3)*exp((1/2)*t)*30^(1/2))+100*exp((5/2)*t)*30^(1/2)*arccoth((1/3)*exp((1/2)*t)*30^(1/2))+69*exp((1/2)*t)*30^(1/2)*arccoth((1/3)*exp((1/2)*t)*30^(1/2))+300*exp(2*t)-180*exp(t)+27)/(100*exp(2*t)-60*exp(t)+9)

(4)

# Check

t0:=1;
eval(ans,t=t0):
evalf(%) = evalf(eval(S, t=t0));

1

 

1.657047690 = 1.657047690

(5)

t0:=-ln(10/3)+1e-2;
eval(ans,t=t0):
evalf(%) = evalf[15](eval(S, t=t0));

-ln(10/3)+0.1e-1

 

-5.925318866 = -5.92522395230400

(6)

 

First 35 36 37 38 39 40 41 Last Page 37 of 111