Preben Alsholm

MaplePrimes Activity


These are replies submitted by Preben Alsholm

@tomleslie You wrote in your reply to Carl Love:
" I thought I was being somewhat more realistic than other "solutions" here because I didn't use the exact analytic solution to set up the shooting method. "

Well, in my answer, which predate yours, I only used the exact solution for comparison. I used a shooting method for the numerical solution. That I handled another ode (Problem 6) is irrelevant. Please have a look.

@tomleslie If you only have lines, points, and one filled area there doesn't seem to be covering in Maple 8 nor in Maple 2018.
But if you have two filled plots the one on top will cover the bottom one:
Define a smaller pink rectangle and display in two orders:
 

squPINK:=polygon([[1,1], [1,side], [side,side],[side,1]], color=pink):
plots[display]([squPINK,squ,poinC,tetA,lin[1],lin[2],lin[3],lin[4]], scaling=constrained,axes=normal);
plots[display]([squ,squPINK,poinC,tetA,lin[1],lin[2],lin[3],lin[4]], scaling=constrained,axes=normal);

In the first you see the pink rectangle, in the second you don't.
In Maple 8 you do see the black sides of the pink rectangle in the white plot though. Not in Maple 2018.

@9009134 I don't know if your system has a solution or not with the homogeneous BCS as given in my reply "Another look ...".
With BCS modified to BCS1 as below and using an approximate solution I made up a result turns up:
 

Digits:=15:
BCS1 := {theta(0) = 0, theta(1/1000) = 0, u(0) = 0, u(1/1000) = 0, w(0) = 0, w(1/1000) = 0, D(w)(0) = -0.1, D(w)(1/1000) = 0.1};
res:=dsolve(SYS union BCS1,numeric,approxsoln=[w(x)=100*x*(x-1/1000),u(x)=-3*10^(-7)*sin(2*Pi*1000*x),theta(x)=-200*(x-1/2000)],initmesh=2048,abserr=1e-7);

Plotting:
 

plots:-odeplot(res,[x,w(x)]);
plots:-odeplot(res,[x,u(x)]);
plots:-odeplot(res,[x,theta(x)]);

The first of these plots:

For the pure plot version I don't get red axes by the command:
 

plot(sin(x),x=-Pi..Pi,axis=[tickmarks=[color=red]], color=blue);
## Easier to see here:
plot(sin(x),x=-Pi..Pi,axis=[tickmarks=[color=red],thickness=3], color=blue);

I'm using Maple 2018.1 64bit with Windows 10.
Note: I can force the axes to be red (if I wanted to):

plot(sin(x),x=-Pi..Pi,axis=[tickmarks=[color=red],thickness=3,color=red], color=blue);


 

@9009134 I had another look at your system. I used the first ode to bring the total apparent differential order down from 10 to 8. That doesn't mean that I have a solution, but at least it clarifies the situation a bit:
 

restart;
dsys3 := {4.320502750*(diff(u(x), x, x))+4.320502750*(diff(w(x), x))*(diff(w(x), x, x)), 3.600418958*10^(-13)*(diff(theta(x), x, x))-1.206140351*theta(x)-1.206140351*(diff(w(x), x)), 1.206140351*(diff(theta(x), x))+1.206140351*(diff(w(x), x, x))+(4.320502750*(diff(u(x), x)+(1/2)*(diff(w(x), x))^2))*(diff(w(x), x, x)-8.379501384*10^20*(diff(w(x), x, x, x, x)))+(diff(w(x), x))*(4.320502750*(diff(u(x), x, x))+(diff(w(x), x))*(diff(w(x), x, x))-3.620365877*10^21*(diff(u(x), x, x, x, x))-1.086109763*10^22*(diff(w(x), x, x))*(diff(w(x), x, x, x))-3.620365877*10^21*(diff(w(x), x))*(diff(w(x), x, x, x, x)))+8.854000000*10^(-12)*(1.032500000-13000.00*w(x))/(2.500000000*10^(-6)-w(x))^2+6.503431892*10^(-32)/(2.500000000*10^(-6)-w(x))^4-4.451526315*10^10*(1.032500000-13000.00*w(x))*(diff(w(x), x))^2/(2.500000000*10^(-6)-w(x))^4+3.857989473*10^14*(diff(w(x), x))^2/(2.500000000*10^(-6)-w(x))^3-1.483842105*10^10*(1.032500000-13000.00*w(x))*(diff(w(x), x, x))/(2.500000000*10^(-6)-w(x))^3+9.644973683*10^13*(diff(w(x), x, x))/(2.500000000*10^(-6)-w(x))^2-1.089910330*10^(-9)*(diff(w(x), x))^2/(2.500000000*10^(-6)-w(x))^6-2.179820662*10^(-10)*(diff(w(x), x, x))/(2.500000000*10^(-6)-w(x))^5, theta(0) = 0, theta(0.1e-2) = 0, u(0) = 0, u(0.1e-2) = 0, w(0) = 0, w(0.1e-2) = 0, (D(u))(0) = 0, (D(u))(0.1e-2) = 0, (D(w))(0) = 0, (D(w))(0.1e-2) = 0}:
dsys3:=convert(dsys3,rational):
sys,bcs:=selectremove(has,dsys3,diff);
indets(sys,specfunc(diff));
sys[1];
eq:=solve(%,{diff(u(x),x,x)});
sys[3];
ode3:=eval(sys[3],eq);
indets(ode3,specfunc(diff));
sys[1..2];
SYS:={ode3,sys[1],sys[2]};
indets(SYS,specfunc(diff));
# Total order 2+2+4 = 8. So 8 boundary conditions, not 10.
bcs;
BCS:=remove(has,bcs,D(u)); # A natural choice 
dsolve(SYS union BCS,numeric); #No success
dsolve(SYS union BCS,numeric,method=bvp[midrich]); #No success

 

Instead of an image could you upload a worksheet? Use the big green arrow in the MaplePrimes editor.

@samen You may want to compare your results to what you get by applying Maple's pdsolve/numeric.
So what is the pde?
Example which may be close to yours.
 

restart;
pde:=diff(v(x,t),t)=a*diff(v(x,t),x,x);
res:=pdsolve(eval(pde,a=10^(-6)),{v(x,0)=1+x*0.9,v(0,t)=0.1,v(1,t)=0.5},numeric,timestep=.0001,spacestep=0.1);
res:-plot3d(t=1);

You can also give a method (including CrankNicholson). See the help page
?pdsolve,numeric,method

Note: In your code you should change the definition of the matrix nn to
 

nn := Matrix(N+2, N+2,(i, j)-> u[i-1, j-1]):

otherwise you won't see the boundary at x=1.
 

@9009134 I think that the immediate problem of not being able to convert to an explicit first order system can be overcome if you are willing to change the boundary conditions somewhat.

But you also have a possibly serious problem in the interior of your interval 0 .. 1/1000.
The denominators (1/400000 - w(x) )^n, n= 2..6 may easily become zero since w(0)=0 and w(1/1000) =0.

@schapplm There is a procedure sqrt in Maple. You can see the code by doing
 

showstat(sqrt);

In fact it appears that at least since Maple 15 sqrt is actually a module that has ModuleApply as a local procedure. This means that the module sqrt can be used as a procedure. That procedure is the procedure ModuleApply.
Details:
 

restart;
showstat(sqrt);
kernelopts(opaquemodules=false);
showstat(sqrt:-ModuleApply); # same procedure as the one above
op(3,eval(sqrt:-ModuleApply)); # option builtin, but can be seen
op(2,eval(sqrt)); # locals include 3 procedures and one large integer.
pr:=sqrt:-Primes; # The product of the first 35 primes:
ifactor(pr);
nops(%);
ithprime(35);
showstat(sqrt:-Product); 
showstat(sqrt:-Nested); 

 

@tomleslie Yes, just by replacing fsolve by solve (and using sigma1=1/2) makes sol = NULL.

About the boundary conditions: You have
D(f)(0)=1, D[1,1](f)(0), D[1,1,1](f)(6)=2
in your call to dsolve, but in the discrete version we find
 

bcEqs:=[subs(k=0,dfdef(h))=rhs(bc1),subs(k=0,d2fde2f(h))=rhs(bc2),
subs(k=n-1,d2fde2b(h))=rhs(bc3)];

which are discretized versions of D(f)(0) = rhs(bc1), (D@@2)(f)(0)=rhs(bc2), (D@@2)(f)(6) = rhs(bc3).
The notations d2fde2f and d2fde2b refer to forward and backward discretized second derivatives, respectively.

On my Windows 10 computer I get that sound when using the arrow keys to move from one input line to the next but only when there are no more input lines in the direction of the arrow.
That doesn't annoy me and is nothing new. I have the sound as early as Maple 12 (Standard Worksheet Interface), but not in Maple 8 (Classical Worksheet).
Could you give an example?

## I can produce the same sound by using the Delete key or the back delete key (or whatever it's called) when there is nothing to delete.
 

@tomleslie I'm using 64-bit Windows 10.

@farah adanan In the paper you link to the odes (13), (14), and (15) together with the boundary conditions (16) really will be a boundary value problem, thus you will not be using Runge-Kutta-Fehlberg.
Infinity in (16) must be replaced by a finite number of an appropriate size (not huge).
Mapleprimes has had quite a lot of problems in that ball park.
You should try to write the equations using Maple syntax. Then we could take it from there.
Boundary value problems are often difficult to solve.

@tomleslie For the same commands I get:

CUDA:-ComputeLevel(id=0);
                               5.
CUDA:-HasDoubleSupport(id=0);
                              true
CUDA:-Properties();
                              [table(["Total Global Memory" = 2147483648, "Minor" = 0, "Name" = "GeForce GTX 860M", "Texture Alignment" = 512, "MultiProcessor Count" = 5, "Device Overlap" = 1, "Clock Rate" = 1019500, "Max Threads Per Block" = 1024, "Max Grid Size" = [2147483647, 65535, 65535], "Max Threads Dimensions" = [1024, 1024, 64], "Major" = 5, "Memory Pitch" = 2147483647, "Resisters Per Block" = 65536, "Total Constant Memory" = 65536, "Kernel Exec Timeout Enabled" = true, "Shared Memory Per Block" = 49152, "Warp Size" = 32, "ID" = 0])]

 

@nm
solve evaluates its arguments before proceeding as is the general rule in Maple.
So try these:

restart;
x>a and x<b and x<>c;
x>a and x<>c;
restart;
x>-infinity and x<infinity and x<>1/2;

The output discards the '<>' part in both cases.  Thus solve never sees x<>1/2.

In the help page for solve/details we find under Parameters that the first argument can be an
"expression, equation, inequality, or procedure; or set or list of expressions, equations, or inequalities".
Admittedly the term 'expression' suggests something vague.

## Whether the 'generic' handling of the case 'x>a and x<>c' is a bug or not I shall find out when I report it as such by submitting anSCR.

First 8 9 10 11 12 13 14 Last Page 10 of 196