Preben Alsholm

13728 Reputation

22 Badges

20 years, 249 days

MaplePrimes Activity


These are replies submitted by Preben Alsholm

@Markiyan Hirnyk

limit(J,a=R/sqrt(2),left) assuming R>0;

Answer: 5/3*sqrt(2)*Pi*R^3-4/3*R^3*Pi

(Now done as it happens in Maple 16).

@Christopher2222 Here it is:

restart;
L := .22:
lambda := 3.4:
rho := 2950:
cp := 1050:
a := lambda/(rho*cp):
h := 9:
T0 := 298:
sigma := 5.6697*10^(-8): epsilon:= .85:
PDE := diff(T(x, t), t) = lambda*(diff(T(x, t), x, x))/(rho*cp);
IBC := {-lambda*(D[1](T))(L, t) = h*(T(L, t)-T0)+sigma*epsilon*(T(L, t)^4-T0^4), T(0, t) = 1573, T(x, 0) = 298};
pde:=diff(u(x,t),t)=(1+u(x,t)^2)/(1+t^2);
ibc:={u(x,0)=0};
SOL:=pdsolve({PDE,pde}, IBC union ibc, numeric,timestep=100);

#SOL := pdsolve(PDE, IBC, numeric, timestep = 100);
SOL:-plot(T,t = 80*3600, thickness = 3);
SOL:-animate(T,t = 0 .. 80*3600, frames = 200, thickness = 3);
p1 := SOL:-value(T);
p1(.22, 80*3600);

# In the uploaded worksheet I have also given another method: Add the term T(x,t)^2*piecewise(x<0,1,0) to the right hand side of the PDE to make the equation nonlinear. But since x>0 in our problem this doesn't hurt.

pde_workaround.mw

@Christopher2222 Here it is:

restart;
L := .22:
lambda := 3.4:
rho := 2950:
cp := 1050:
a := lambda/(rho*cp):
h := 9:
T0 := 298:
sigma := 5.6697*10^(-8): epsilon:= .85:
PDE := diff(T(x, t), t) = lambda*(diff(T(x, t), x, x))/(rho*cp);
IBC := {-lambda*(D[1](T))(L, t) = h*(T(L, t)-T0)+sigma*epsilon*(T(L, t)^4-T0^4), T(0, t) = 1573, T(x, 0) = 298};
pde:=diff(u(x,t),t)=(1+u(x,t)^2)/(1+t^2);
ibc:={u(x,0)=0};
SOL:=pdsolve({PDE,pde}, IBC union ibc, numeric,timestep=100);

#SOL := pdsolve(PDE, IBC, numeric, timestep = 100);
SOL:-plot(T,t = 80*3600, thickness = 3);
SOL:-animate(T,t = 0 .. 80*3600, frames = 200, thickness = 3);
p1 := SOL:-value(T);
p1(.22, 80*3600);

# In the uploaded worksheet I have also given another method: Add the term T(x,t)^2*piecewise(x<0,1,0) to the right hand side of the PDE to make the equation nonlinear. But since x>0 in our problem this doesn't hurt.

pde_workaround.mw

In the procedure `pdsolve/default/construct_method` (lines  176 - 181 in Maple 12) there is a check on the value of a local table INFO produced by `pdsolve/numeric` and given as argument to `pdsolve/default/construct_method`.
If INFO["linear"] = true then .... else .... end if;
The problem seems to be that somehow INFO["linear"] has got the value 'true' (in Maple 12) at that point, whereas at the similar point in Maple 16 it has the value 'false'.

If you introduce a peaceful nonlinear pde quite independent of the other pde it appears that you are able to work your way around the problem

restart;
PDE := diff(T(x, t), t) = diff(T(x, t), x, x);
IBC_4:={D[1](T)(1, t) = -T(1,t)^4, T(0, t) = 1, T(x, 0) = 2};
pde:=diff(u(x,t),t)=(1+u(x,t)^2)/(1+t^2);
pdsolve(pde);
ibc:={u(x,0)=0};
SOL:=pdsolve({PDE,pde}, IBC_4 union ibc, numeric,spacestep=.005);
SOL:-plot(T,t=2);
SOL:-plot(u,t=2);
SOL:-value(T)(1,1);

#This workaround also works on the original problem from 2010:

http://www.mapleprimes.com/questions/97391-Whats-Wrong-With-Maple-Solution

PS. Instead you could add the term T(x,t)^2*piecewise(x<0,1,0) to the right hand side of the PDE to make the equation nonlinear. Since x>0 in our problem this doesn't hurt.

@Alejandro Jakubi Here is one:

restart;
PDE := diff(T(x, t), t) = diff(T(x, t), x, x);
IBC_4:={D[1](T)(1, t) = -T(1,t)^4, T(0, t) = 1, T(x, 0) = 2};
IBC_0:={D[1](T)(1, t) = 0, T(0, t) = 1, T(x, 0) = 2};
IBC_1:={D[1](T)(1, t) = -T(1,t), T(0, t) = 1, T(x, 0) = 2};
SOL_4 := pdsolve(PDE, IBC_4, numeric,spacestep=.005);
SOL_0 := pdsolve(PDE, IBC_0, numeric,spacestep=.005);
SOL_1 := pdsolve(PDE, IBC_1, numeric,spacestep=.005);
p4:=SOL_4:-plot(t = 1, color = red):
p0:=SOL_0:-plot(t = 1, color = blue):
p1:=SOL_1:-plot(t = 1, color = black):
plots:-display(p0,p4,p1);
SOL_4:-value()(1,1);
SOL_0:-value()(1,1);
SOL_1:-value()(1,1);

#There is no difference between SOL_4 and SOL_0 in Maple 12. However, SOL_1 is clearly different from the other two.   In Maple 16 all 3 are different and SOL_0 and SOL_1 are as in Maple 12.
This issue came up in 2010:

http://www.mapleprimes.com/questions/97391-Whats-Wrong-With-Maple-Solution

What I said was that there was no convergence in the original problem after reformulation neither in Maple 12 nor in Maple 16.

The original code stripped to essentials:

restart;
L := .22:
lambda := 3.4:
rho := 2950:
cp := 1050:
a := lambda/(rho*cp):
h := 9:
T0 := 298:
sigma := 5.6697*10^(-8): epsilon:= .85:
PDE := diff(T(x, t), t) = lambda*(diff(T(x, t), x, x))/(rho*cp);
IBC := {-lambda*(D[1](T))(L, t) = h*(T(L, t)-T0)+sigma*epsilon*(T(L, t)^4-T0^4), T(0, t) = 1573, T(x, 0) = 298};
SOL := pdsolve(PDE, IBC, numeric, timestep = 100);
SOL:-plot(t = 80*3600, thickness = 3);
SOL:-animate(t = 0 .. 80*3600, frames = 200, thickness = 3);
p1 := SOL:-value();
p1(.22, 80*3600);
##That "worked", but the result is wrong in Maple 12, but fine in Maple 16.
####################################################
##And the reformulated version:
restart;
L := .22:
lambda := 3.4:
rho := 2950:
cp := 1050:
a := lambda/(rho*cp):
h := 9:
T0 := 298:
sigma := 5.6697*10^(-8): epsilon:= .85:
sys:= {diff(T(x, t), t) = lambda*(diff(T(x, t), x, x))/(rho*cp),T4(x,t)=T(x,t)^4};
IBC := {-lambda*(D[1](T))(L, t) = h*(T(L, t)-T0)+sigma*epsilon*(T4(L, t)-T0^4), T(0, t) = 1573, T(x, 0) = 298};
SOL := pdsolve(sys, IBC, numeric, timestep = 100);
SOL:-plot(t = 80*3600, thickness = 3);
##Here there is an error message in Maple 12 and 16:

Error, (in pdsolve/numeric/plot) unable to compute solution for t>HFloat(0.0):
Newton iteration is not converging
#The only difference between Maple  12 and 16 is that Maple 12 says 't>0.' instead of 't>HFloat(0.0)'.



What I said was that there was no convergence in the original problem after reformulation neither in Maple 12 nor in Maple 16.

The original code stripped to essentials:

restart;
L := .22:
lambda := 3.4:
rho := 2950:
cp := 1050:
a := lambda/(rho*cp):
h := 9:
T0 := 298:
sigma := 5.6697*10^(-8): epsilon:= .85:
PDE := diff(T(x, t), t) = lambda*(diff(T(x, t), x, x))/(rho*cp);
IBC := {-lambda*(D[1](T))(L, t) = h*(T(L, t)-T0)+sigma*epsilon*(T(L, t)^4-T0^4), T(0, t) = 1573, T(x, 0) = 298};
SOL := pdsolve(PDE, IBC, numeric, timestep = 100);
SOL:-plot(t = 80*3600, thickness = 3);
SOL:-animate(t = 0 .. 80*3600, frames = 200, thickness = 3);
p1 := SOL:-value();
p1(.22, 80*3600);
##That "worked", but the result is wrong in Maple 12, but fine in Maple 16.
####################################################
##And the reformulated version:
restart;
L := .22:
lambda := 3.4:
rho := 2950:
cp := 1050:
a := lambda/(rho*cp):
h := 9:
T0 := 298:
sigma := 5.6697*10^(-8): epsilon:= .85:
sys:= {diff(T(x, t), t) = lambda*(diff(T(x, t), x, x))/(rho*cp),T4(x,t)=T(x,t)^4};
IBC := {-lambda*(D[1](T))(L, t) = h*(T(L, t)-T0)+sigma*epsilon*(T4(L, t)-T0^4), T(0, t) = 1573, T(x, 0) = 298};
SOL := pdsolve(sys, IBC, numeric, timestep = 100);
SOL:-plot(t = 80*3600, thickness = 3);
##Here there is an error message in Maple 12 and 16:

Error, (in pdsolve/numeric/plot) unable to compute solution for t>HFloat(0.0):
Newton iteration is not converging
#The only difference between Maple  12 and 16 is that Maple 12 says 't>0.' instead of 't>HFloat(0.0)'.



Nonlinear terms in the IBCs seem to get ignored (or complained about) in Maple 12.

I tried the following simplified version. First the obvious way, which compared to Maple 16 gives the wrong result:

restart;
PDE := diff(T(x, t), t) = diff(T(x, t), x, x);
IBC:={D[1](T)(1, t) = -T(1,t)^4, T(0, t) = 1, T(x, 0) = 2};
SOL := pdsolve(PDE, IBC, numeric,spacestep=.005);
SOL:-plot(t = 1, thickness = 3, colour = red);
SOL:-animate(t=0..1);
SOL:-value()(1,1);
#Now introduce T4(x,t) = T(x,t)^4:
restart;
PDE := diff(T(x, t), t) = diff(T(x, t), x, x);
IBC := {D[1](T)(1, t) = -T4(1, t), T(0, t) = 1, T(x, 0) = 2};
sys:={PDE,T4(x,t)=T(x,t)^4};
SOL := pdsolve(sys, IBC, numeric,spacestep=.005);
SOL:-plot(T,t = 1, thickness = 3, colour = red);
SOL:-animate(T,t=0..1);
SOL:-value(T)(1,1);
#Now the results agree with those from Maple 16.

#However, unfortunately when trying this last version on the original problem Maple 12 as well as Maple 16 complains:
Error, (in pdsolve/numeric/plot) unable to compute solution for t>HFloat(0.0):
Newton iteration is not converging


Nonlinear terms in the IBCs seem to get ignored (or complained about) in Maple 12.

I tried the following simplified version. First the obvious way, which compared to Maple 16 gives the wrong result:

restart;
PDE := diff(T(x, t), t) = diff(T(x, t), x, x);
IBC:={D[1](T)(1, t) = -T(1,t)^4, T(0, t) = 1, T(x, 0) = 2};
SOL := pdsolve(PDE, IBC, numeric,spacestep=.005);
SOL:-plot(t = 1, thickness = 3, colour = red);
SOL:-animate(t=0..1);
SOL:-value()(1,1);
#Now introduce T4(x,t) = T(x,t)^4:
restart;
PDE := diff(T(x, t), t) = diff(T(x, t), x, x);
IBC := {D[1](T)(1, t) = -T4(1, t), T(0, t) = 1, T(x, 0) = 2};
sys:={PDE,T4(x,t)=T(x,t)^4};
SOL := pdsolve(sys, IBC, numeric,spacestep=.005);
SOL:-plot(T,t = 1, thickness = 3, colour = red);
SOL:-animate(T,t=0..1);
SOL:-value(T)(1,1);
#Now the results agree with those from Maple 16.

#However, unfortunately when trying this last version on the original problem Maple 12 as well as Maple 16 complains:
Error, (in pdsolve/numeric/plot) unable to compute solution for t>HFloat(0.0):
Newton iteration is not converging


I find it desirable that the printed output distinguishes between different objects.

I find it unfortunate (but can easily live with it) that 'Pi' and 'pi' prints the same.
'e' and exp(1) are clearly distinguishable in print.

With typesetting set to 'Extended' as you describe the output from these two commands look the same, but they are not:
evalf(-a);
-a;
#Both are products, but the factors are not the same:

evalf(-a);
op(%);
eval(%%,a=5/6);
-a;
op(%);
eval(%%,a=5/6);

I find it desirable that the printed output distinguishes between different objects.

I find it unfortunate (but can easily live with it) that 'Pi' and 'pi' prints the same.
'e' and exp(1) are clearly distinguishable in print.

With typesetting set to 'Extended' as you describe the output from these two commands look the same, but they are not:
evalf(-a);
-a;
#Both are products, but the factors are not the same:

evalf(-a);
op(%);
eval(%%,a=5/6);
-a;
op(%);
eval(%%,a=5/6);

@AliKhan

I am not using 'assume' but 'assuming' as seen below.

restart;
W := m*Pi/b; K := 2*n*Pi/t;
A := simplify(2*int(sin(W*(x+1/2*b))*cos(K*x), x = -1/2*b .. 1/2*b)/t) assuming m::integer;
Aodd:=simplify(A) assuming m::odd;
simplify(A) assuming m::even;
#Supposed answer:
ANS:=4*m*b*t*cos(Pi*n*b/t)*sin((1/2)*m*Pi)^2/(Pi*(-m^2*t^2+4*n^2*b^2));
simplify(ANS) assuming m::even;
simplify(ANS) assuming m::odd;
#We notice that ANS has the wrong sign. To make a quick check:
#Numerical integration with the special values t=1,b=1,n=1,m=7:
evalf(eval(2*Int(sin(W*(x+1/2*b))*cos(K*x), x = -1/2*b .. 1/2*b)/t,{t=1,b=1,n=1,m=7}));
#Compare with the next two outputs:
evalf(eval(Aodd,{t=1,b=1,n=1,m=7}));
evalf(eval(ANS,{t=1,b=1,n=1,m=7}));
#So Maple is right: ANS has the wrong sign.
#The special case m*t = 2*n*b (notice 'Int' instead of 'int'):
Aspec := eval(2*Int(sin(W*(x+(1/2)*b))*cos(K*x), x = -1/2*b .. 1/2*b)/t, n=m*t/(2*b));
A1:=value(Aspec);
#Answers don't simplify easily
A1even:=simplify(A1) assuming m::even;
eval(A1even,m=2*p);
simplify(%) assuming p::integer;
simplify(A1) assuming m::odd;
#Proposed answer:
ANS1:=b*sin((1/2)*m*Pi)/t;
simplify(ANS1) assuming m::odd;
#At least the last two 'simplies' to the same thing! But somewhat nicer:
eval(A1,m=2*p+1);
simplify(%) assuming p::integer;
eval(ANS1,m=2*p+1);
simplify(%) assuming p::integer;


@AliKhan

I am not using 'assume' but 'assuming' as seen below.

restart;
W := m*Pi/b; K := 2*n*Pi/t;
A := simplify(2*int(sin(W*(x+1/2*b))*cos(K*x), x = -1/2*b .. 1/2*b)/t) assuming m::integer;
Aodd:=simplify(A) assuming m::odd;
simplify(A) assuming m::even;
#Supposed answer:
ANS:=4*m*b*t*cos(Pi*n*b/t)*sin((1/2)*m*Pi)^2/(Pi*(-m^2*t^2+4*n^2*b^2));
simplify(ANS) assuming m::even;
simplify(ANS) assuming m::odd;
#We notice that ANS has the wrong sign. To make a quick check:
#Numerical integration with the special values t=1,b=1,n=1,m=7:
evalf(eval(2*Int(sin(W*(x+1/2*b))*cos(K*x), x = -1/2*b .. 1/2*b)/t,{t=1,b=1,n=1,m=7}));
#Compare with the next two outputs:
evalf(eval(Aodd,{t=1,b=1,n=1,m=7}));
evalf(eval(ANS,{t=1,b=1,n=1,m=7}));
#So Maple is right: ANS has the wrong sign.
#The special case m*t = 2*n*b (notice 'Int' instead of 'int'):
Aspec := eval(2*Int(sin(W*(x+(1/2)*b))*cos(K*x), x = -1/2*b .. 1/2*b)/t, n=m*t/(2*b));
A1:=value(Aspec);
#Answers don't simplify easily
A1even:=simplify(A1) assuming m::even;
eval(A1even,m=2*p);
simplify(%) assuming p::integer;
simplify(A1) assuming m::odd;
#Proposed answer:
ANS1:=b*sin((1/2)*m*Pi)/t;
simplify(ANS1) assuming m::odd;
#At least the last two 'simplies' to the same thing! But somewhat nicer:
eval(A1,m=2*p+1);
simplify(%) assuming p::integer;
eval(ANS1,m=2*p+1);
simplify(%) assuming p::integer;


int(simplify(f),r); doesn't work for me in Maple 16, 15, and 12.

Even when specifying n success doesn't come often with or without simplify:

[seq(int(f,r),n=1..10)];
has~(%,int);

int(simplify(f),r); doesn't work for me in Maple 16, 15, and 12.

Even when specifying n success doesn't come often with or without simplify:

[seq(int(f,r),n=1..10)];
has~(%,int);

First 198 199 200 201 202 203 204 Last Page 200 of 230