Preben Alsholm

13471 Reputation

22 Badges

20 years, 210 days

MaplePrimes Activity


These are answers submitted by Preben Alsholm

This is a bug.
 

restart;
A := abs(x*Heaviside(x) + x*Heaviside(-x));
plot(A,x=-1..1); #OK
int(A,x=-1..1); # Wrong: should be 1
int(A,x=-1..1,numeric); # 1.000000000
int(A,x=-1..1,method=_RETURNVERBOSE);
int(A,x=-1..1,method=_VERIFYFLOAT);
## Workaround:
B:= x*Heaviside(x) + x*Heaviside(-x);
BP:=convert(B,piecewise);
plot(BP,x=-1..1);
plot(abs(BP),x=-1..1);
int(abs(BP),x=-1..1); # 1

I shall submit a bug report.

In Maple 2024 use adaptive=true, or simply adaptive.
 

res := t^3 - 3*t^2*sqrt(t^2*(12*sqrt(2)*ln(t) + 9*ln(t)^2 + 8)^(1/3))*ln(t)/(12*sqrt(2)*ln(t) + 9*ln(t)^2 + 8)^(2/3) - 2*t^2*sqrt(t^2*(12*sqrt(2)*ln(t) + 9*ln(t)^2 + 8)^(1/3))*sqrt(2)/(12*sqrt(2)*ln(t) + 9*ln(t)^2 + 8)^(2/3);

evalf(eval(res,t=-5));
plot(res,t=-5..1,adaptive,thickness=3);

 

Inspired by Rouben and using his example:
 

restart;
plot([10 + sin(x),0], x=3*Pi..4*Pi);
plot([10 + sin(x),0], x=3*Pi..4*Pi,color=[blue,white],thickness=[3,1]);
plot([-2 + sin(x),0], x=3*Pi..4*Pi);
plot([-2 + sin(x),0], x=3*Pi..4*Pi,color=[blue,white],thickness=[3,1]);
plot([sin(x),0], x=Pi..4*Pi,color=[blue,white],thickness=[3,1]);

 

I'm not sure what you mean by "Kronecker not recognized".
The issue seems to be the following:
 

restart;
A:=Int(
         Sum(p__n*exp(-I*theta)^n, n = 0 .. infinity)*Sum(p__m*exp(theta*I)^m, m = 0 .. infinity), 
       theta = 0 .. 2*Pi);
eval(A,Int=int); # 0 Wrong!
# Counter example:
B:=subs({p__n=exp(-n),p__m=exp(-m)},A); # Don't use eval instead of subs
value(B); #  2*Pi*exp(2)/(exp(2) - 1)

 

A very useful feature of display is overrideoptions:
 

restart;
pt:=plottools:-line([1,2], [3,4],color=red,thickness=2,linestyle=dash); 
plots:-display(pt);
plots:-display(pt,overrideoptions,linestyle=1,thickness=5,color=blue);

restart;
with(IntegrationTools):
J1 := Int(h(x), x=-infinity..a) - Int(h(x), x=-infinity..b);
Ja:=Int(h(x), x=-infinity..a);
Jb:=Int(h(x), x=-infinity..b);
# Case 1: a>b.
Sab:=Split(Ja,b) assuming a>b;
Combine(Sab-Jb);
# Case 2: a<b. 
Sba:=Split(Jb,a) assuming a<b;
Combine(Ja-(Sba)) assuming a<b;
simplify(%);

MaplePrimes24-02-27_IntegrationTools_Combine.mw

NOTE: You don't need the assumption at the end, but it doesn't hurt.
You could add Flip(%); after simplify(%);
Then the results from both cases look the same.

restart;
###
with(LinearAlgebra):
###
v := <x, y>;
w:= <a,b>;
v^%T.v;
v^%T.w;
BilinearForm(v, v) assuming real;
##
## So you could e.g. do this:
M:=module() option package; export BilinearForm;
      BilinearForm:=proc(v::Vector,w::Vector) v^%T.w  end  proc;
end module;  
##             
with(M);
BilinearForm(<x, y>,<a,b>);
BilinearForm(<x, y,z>,<a,b,c>);
##
LinearAlgebra:-BilinearForm(<x, y>,<a,b>);

 

In the first procedure I would do this:
 

restart;
J1 := proc()
  local z:
  z := proc(u) if not u::realcons then return 'procname'(_passed) end if;
  fsolve(sqrt(x)=u, x) end proc:
  evalf(Int(z(u), u=0..1))
end proc:
##############
J1(); # 0.3333333333

In the next I would do likewise:

J2 := proc()
  local z:
   z := proc(q1, q2) if not [q1,q2]::list(realcons) then return 'procname'(_passed) end if;
     exp(
       2*(
         fsolve(1/2+(1/2)*erf((1/2)*x*sqrt(2)) = q1, x)
         *
         fsolve(1/2+(1/2)*erf((1/2)*x*sqrt(2)) = q2, x)
       )
       -
       fsolve(1/2+(1/2)*erf((1/2)*x*sqrt(2)) = sqrt(q1), x)^2
       -
       fsolve(1/2+(1/2)*erf((1/2)*x*sqrt(2)) = sqrt(q2), x)^2
    )
  end proc:
  evalf(Int(z(q1, q2), q1=0..1, q2=0..1))
end proc;

The execution of J2(), however takes an enormous amount of time (infinity?),
but fsolve has to work 4 times inside evalf/int for each value of q1=0..1, q2=0..1.

## Note:  If you comment out evalf(Int(z(q1, q2), q1=0..1, q2=0..1))  and replace it with
z(1/2,1/5) the execution works and returns 0.7300930564.
##
Note 2: I tried using narrower intervals:
evalf(Int(z(q1, q2), q1=0.1..0.9, q2=0.1..0.9))
J2() returned 0.4412675858 in a short time.

You are defining the procedure f to be the same as the procedure F so why not simply start with f:=F  ?

restart;
f := F;
((f@@(-1))@f)(x);
# Let us assume that y is defined this way

y := ((f@@(-1))@g)(x);
# When g is identical to f I would like to get y=x

'y' = eval(y, g=f); # y=x

A concrete example where F=ln:
 

restart;
f := ln;
((f@@(-1))@f)(x);
# Let us assume that y is defined this way

y := ((f@@(-1))@g)(x);
# When g is identical to f I would like to get y=x

'y' = eval(y, g=f); #y=x

 

I corrected a few errors, and here is what I have:
 

restart;
with(plots):
ode1 := diff(f(x), x$3) + 
        f(x)*diff(f(x), x, x) - lambda*diff(f(x), x) - diff(f(x), x)^2 + 2*beta*g(x) - Fr*diff(f(x), x)^2 
        = 0;
ode2 := diff(g(x), x$2) - diff(f(x), x)*g(x) + f(x)*diff(g(x), x) - 2*beta*diff(f(x), x) - lambda*g(x) - Fr*g(x)^2 
        = 0;
# You have 6 boundary conditions, but you can only have 5, unless one of the constants e.g lambda has to be determined.
# Your boundary conditions all have values 0, 
# which means that f(x) = 0, g(x)= 0 for all x on the interval is a (trivial) solution.
# To get a possible nontrivial solution you need to give an approximate solution:
# See ?dsolve,bvp
bcs := f(0) = 0, D(f)(0) = 0, f(5) = 0, D(f)(5) = 0, g(0) = 0, g(5) = 0;
a1 := [lambda = 0.2, beta = 0.2, Fr = 0.1];
####
sys1:=eval({ode1,ode2,bcs[1..5]},a1);# I arbitrarily removed the sixth bcs (g(5)=0).
res1:=dsolve(sys1,numeric);
odeplot(res1,[x,f(x)],thickness=5,view=-0.01..0.01);# As expected the zero solution.

To find a suitable approximate solution is not easy if you have no idea how f and g are expected to look.
You might have an idea?

Here is a proof that y(t) -> infinity as t -> infinity:
MaplePrimes24-02-07_odesys_infinity.mw

Here is the code directly:
 

restart;

sys_ode := diff(x(t), t) = -x(t)^2/(4*Pi*y(t)*(x(t)^2 + y(t)^2)), diff(y(t), t) = y(t)^2/(4*Pi*x(t)*(x(t)^2 + y(t)^2));
ics := x(0) = 1, y(0) = 1:

ode1,ode2:=sys_ode;

ode2/ode1;

ode:=diff(y(x),x)=-(y(x)/x)^3;

sol:=dsolve({ode,y(1)=1});

#From ode1 it follows that x(t) is decreasing from its initial value 1.
# sol shows that x cannot get lower than 1/sqrt(2), which is lower than 1. 
# Since x(t) > 1/sqrt(2), but decreasing we get by replacing x(t) by 1/sqrt(2):
ineq:=rhs(ode2) > subs(x(t)=1/sqrt(2),rhs(ode2));
ode3:=diff(y(t),t)=lhs(ineq);
sol3:=dsolve({ode3,y(0)=1});
limit(rhs(sol3),t=infinity); #infinity
# Thus y(t) from sys_ode will go to infinity too.

####### Graphical illustrations:
plot(rhs(sol),x=1/sqrt(2)..5,color=red,view=0..3); p1:=%: #Saving the plot

Digits:=15;
Sol := dsolve({sys_ode, ics}, {x(t),y(t)}, numeric,abserr=1e-16,relerr=1e-13,events=[[x(t)=5,halt]]);

Sol(100000);
1/sqrt(2.);
plots:-odeplot(Sol,[t,x(t)-1/sqrt(2)],99000..100000,caption=typeset("The difference between x(t) and ",1/sqrt(2)));
#Since x(t) is decreasing from 1 at t=0, you have to go backwards to get near to 1/sqrt(2)
plots:-odeplot(Sol,[x(t),y(t)], t=-43.1..100, color=blue,view=0..3,numpoints=10000);  p2:=%: # Saving the plot

plots:-display(p1,p2);# The last one is on top, so the color is blue
plots:-display(Array([p1,p2]));

 

Here is a simpler version. Let me know, if the result is what you expected:
 

restart;
with(LinearAlgebra):

xA := 1;
yA := 0;
xB := 0;
yB := 0;
xC := 0;
yC := 1;
Mat := Matrix(3, 3, [xA, xB, xC, yA, yB, yC, 1, 1, 1]);

Miv:=Mat^(-1); #Simpler
#Miv := MatrixInverse(Mat);

Miv^%T; #Simpler
#Transpose(Miv);
Miv^%T. <x, y, 1>;

phi := unapply(%,x,y);

for i to 6 do
    B || i := phi(xA || i, yA || i);
end do;

MaplePrimes24-02-05_LinearAlgebra.mw

Try this:

restart;

foo1:=overload([                   
                  proc(P1::list,P2::list,P3::list,$)
                      option overload(callseq_only);
                    [P1,P2,P3]
                  end proc,
                  proc(P1::list,P2::list,a::algebraic:=4,$)
                      option overload;
                    [P1,P2,a]
                  end proc
                       ]):
foo1([1,2],[3,4]); 
foo1([1,2],[3,4],x);
foo1([1,2],[3,4],[4,7]);

See ?overload where we read:
"If p1 raises an exception during the parameter type-checking phase only (not inside p1), and option overload(callseq_only) has been specified, then execution will proceed to p2(x,y)."

Note: See below for my comment on Ronan's original foo1.
Also try running my example above with option overload instead of option overload(callseq_only):
The results are the same. Thus the answer to the (revised) title is it didn't matter here.
Option overload(callseq_only) was introduced in Maple 16, but no example where it matters was given.

restart;
f :=piecewise(0 <= x and x <= 1.588125352, (-1)*0.39*x^2 + 1.459*x - y, 1.588125352 < x and x < 4, (x - 1.81)^2 + (y - 0.42)^2 + (-1)*0.94^2);

plots:-implicitplot(f, x = 0 .. 3, y = 0 .. 1.5, scaling = constrained,gridrefine=4,signchange=false);

In order to see that a sign change takes place you could try this:

plot3d([f,0],x = 0 .. 3, y = 0 .. 1.5);

 

Try this:

restart;
interface(typesetting=extended);
exp(-a); # No visible minus sign
interface(typesetting=standard);
exp(-a); #Fine also in Maple 2023.2
1 2 3 4 5 6 7 Last Page 1 of 158