12980 Reputation

19 Badges

12 years, 290 days

MaplePrimes Activity

These are replies submitted by tomleslie

model. There appears to be a transfer from the final state (ie resolved=recovered/dead) back to the susceptible population, which is presumably modelling the fact that infection only offers partial immunity


Im using Windows 7 with Excel 2010 (14.0.7268.5000) SP2 MSO (14.0.72682.5000).

Does anything different if you exclude the 'header' line in the matrix you are trying to output?

ie replace

ExcelTools:-Export( M1, "C:/Temp/M1.xlsx", 1, "B2");


ExcelTools:-Export( M1[2..,1..2], "C:/Temp/M1.xlsx", 1, "B3");

you may as well have the attached before I transfer it to the bitbucket for good

I have no idea whether it reflects what you want, but at least it functions

Aren't posted worksheets - nice?!!

May be an Excel or OS issue, because OP's code produces the attached in Maple 2022 on my machine


which is working just fine. Only thing I had a slight doubt about was the use (in Maple) of the subscripted variable V__b. However on export and within Excel this just renders as a literal name with two underscores in it - which I owuld regard as an acceptable compromise.

Since I can't make this fail, I can't fix it. Can only suggest that the OP changes the name in the header row of the matrix to be exported from V__b, to something like Vb, just to see if it makes any difference


the option in the textplot3d() command should have read font=..., not labelfont=....  That's what too much use of copy/paste will do to you!

Fixed in the attached


The OP's original optimisation commands were

`sol1≔Optimization:-Minimize`(new, t = -4.17 .. -4.10, x = 0.5 .. 0.8)
`sol2 ≔ Optimization:-Maximize`(new, t = -4.17 .. 0, x = 0.2 .. 2)

check the ranges on 'x' - this is what I used for my original response


that was why OP restricted the range of 'x' in the optimization problems to avoid the value x=0 (and also to avoid what would appear to be infinities at x~=+/-9 and x~=+/-19


still absolutely refuse to upload an executable worksheet, I am going to follow your example, and not upload the functional worksheet I created - after all who needs this - right?

So the problems in you code are

  1. You have an global variable 'cir' which is used in the draw() command. The variable 'cir' is defined nowhere
  2. For some values of 'i' in your final seq() command, the gradient of the line L3 becomes infinite and changes sign form +ve to -ve, this causes the 'order' of the results from some of your solve() commands to change, which means that between two frames of the animation, much of the your construction is reflected about the x-axis.
  3. There are similar animation glitches when the tangent lines are either horizontal or vertical
  4. For some some values of 'i' in your final seq() command the three points supplied to the triangle command are collinear, which means that the tranngle cannot be created, Since it cannot be created you will get an error from the draw() command, because you are trying to draw something which does not exist
  1. Execute your worksheet
  2. Save the worksheet (with the output, including the error) to a .mw file
  3. On the Mapleprimes site, locate the comment to whihc yu wish to reply, and click the 'Reply' at the lower right of the comment
  4. Enter a "Title"
  5. In the toolbar here, there is a big green up-arrow, third from the right in the top row. Click this icon. This will provide a pop-up where you can select the file saved at stage 2 above
  6. Click 'Upload' in this pop-up
  7. Click 'Insert Link' in this popup
  8. A link to your worksheet will now appear in the 'Reply' which you are composing


a complete worksheet using the big green up-arrow in the Mapleprimes toolbar


no idea what you want - maybe the attached,which does dra the relevant line on the relevant surface.

 f := transform((X, Y) -> [X, 0.5, Y]):

is used to transrom a 2D plot to a 3D plot - see the help at plottools/transform.

I replied to the question at


see the attached, whihc for some reason won't display inline here

You can browse books on Maple here

depends which documentation you are using!

I may be a bit old-fashioned, but I think that for learning Maple from scratch, the (pdf) manuals available here

are a better bet than simply using the built-in help.

These manuals are available in worksheet format in the built-in help, which I don't find as convenient as the pdf versions


now the analytical and numerical solutions agree and all graphs display- see the attached.


ode := diff(theta(y), y, y)+Br*((-3*We*(y+1+(1/2)*x^2)*((9*(1+(1/2)*x^2)^2+9*y^2)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^3/(4*(1+(1/2)*x^2)^7)+3*y*(k+1)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^2/(1+(1/2)*x^2)^5+3*(k+1)^2*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)/(2*(1+(1/2)*x^2)^3)-(-120*lambda*(1+(1/2)*x^2)^4+(63*k^2-36*k+63)*(k-1)*(1+(1/2)*x^2)^3+288*Q*(k^2-(11/8)*k+1)*(1+(1/2)*x^2)^2+486*Q^2*(k-1)*(1+(1/2)*x^2)+324*Q^3)/(15*(1+(1/2)*x^2)^5))/(16*(1+(1/2)*x^2)^2)-3*We*(y-1-(1/2)*x^2)*((9*(1+(1/2)*x^2)^2+9*y^2)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^3/(4*(1+(1/2)*x^2)^7)+3*y*(k+1)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^2/(1+(1/2)*x^2)^5+3*(k+1)^2*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)/(2*(1+(1/2)*x^2)^3)-(-120*lambda*(1+(1/2)*x^2)^4+(63*k^2-36*k+63)*(k-1)*(1+(1/2)*x^2)^3+288*Q*(k^2-(11/8)*k+1)*(1+(1/2)*x^2)^2+486*Q^2*(k-1)*(1+(1/2)*x^2)+324*Q^3)/(15*(1+(1/2)*x^2)^5))/(16*(1+(1/2)*x^2)^2)-3*We*(y-1-(1/2)*x^2)*(y+1+(1/2)*x^2)*(9*y*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^3/(2*(1+(1/2)*x^2)^7)+(3*(k+1))*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^2/(1+(1/2)*x^2)^5)/(16*(1+(1/2)*x^2)^2)+(-3*y*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)/(1+(1/2)*x^2)^2-k-1)/(x^2+2))^2-(1/2)*We*(-3*We*(y+1+(1/2)*x^2)*((9*(1+(1/2)*x^2)^2+9*y^2)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^3/(4*(1+(1/2)*x^2)^7)+3*y*(k+1)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^2/(1+(1/2)*x^2)^5+3*(k+1)^2*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)/(2*(1+(1/2)*x^2)^3)-(-120*lambda*(1+(1/2)*x^2)^4+(63*k^2-36*k+63)*(k-1)*(1+(1/2)*x^2)^3+288*Q*(k^2-(11/8)*k+1)*(1+(1/2)*x^2)^2+486*Q^2*(k-1)*(1+(1/2)*x^2)+324*Q^3)/(15*(1+(1/2)*x^2)^5))/(16*(1+(1/2)*x^2)^2)-3*We*(y-1-(1/2)*x^2)*((9*(1+(1/2)*x^2)^2+9*y^2)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^3/(4*(1+(1/2)*x^2)^7)+3*y*(k+1)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^2/(1+(1/2)*x^2)^5+3*(k+1)^2*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)/(2*(1+(1/2)*x^2)^3)-(-120*lambda*(1+(1/2)*x^2)^4+(63*k^2-36*k+63)*(k-1)*(1+(1/2)*x^2)^3+288*Q*(k^2-(11/8)*k+1)*(1+(1/2)*x^2)^2+486*Q^2*(k-1)*(1+(1/2)*x^2)+324*Q^3)/(15*(1+(1/2)*x^2)^5))/(16*(1+(1/2)*x^2)^2)-3*We*(y-1-(1/2)*x^2)*(y+1+(1/2)*x^2)*(9*y*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^3/(2*(1+(1/2)*x^2)^7)+(3*(k+1))*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^2/(1+(1/2)*x^2)^5)/(16*(1+(1/2)*x^2)^2)+(-3*y*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)/(1+(1/2)*x^2)^2-k-1)/(x^2+2))^4) = 0

bc := theta(-1-(1/2)*x^2) = 0, theta(1+(1/2)*x^2) = 1

Q__vals := [.5617, .4392, .2564, .1645, 0.659e-1]

`λ__vals` := [0.96e-2, 0.10e-2, 0.93e-2, 0.75e-2, 0.31e-2]; k__vals := [.1, .3, .5, .7, .9]

We__vals := [.1, .3, .5, .7, .9]

colors := [black, red, blue, green, orange, cyan, purple, yellow, navy]; for j to numelems(k__vals) do sol := dsolve(eval([ode, bc], [k = k__vals[j], We = .1, Br = .3, x = 0, Q = Q__vals[j], lambda = `λ__vals`[j]]), numeric); pltn[j] := plots:-odeplot(sol, [y, theta(y)], y = 0 .. 1, axes = boxed, color = colors[j]) end do; plots:-display([seq(pltn[j], j = 1 .. 5)], title = "Numerical Solution")



Analytical_sol := -(3*y+3+3*x^2*(1/2))*((-2/3+(-(5/12)*k^2+(1/6)*k-5/12)*Br)*(1+(1/2)*x^2)^5-Br*((-(1/12)*k^2+(1/6)*k-3/4)*y+Q*(k-1))*(1+(1/2)*x^2)^4-(((1/6)*k+7/6)*y+Q)*Br*((-(1/2)*k+1/2)*y+Q)*(1+(1/2)*x^2)^3+Br*((1/4)*(k-1)^2*y^2-(1/3)*Q*(k-5)*y+Q^2)*y*(1+(1/2)*x^2)^2-((1-k)*y+Q)*Br*y^2*Q*(1+(1/2)*x^2)+Q^2*y^3*Br)/(4*(1+(1/2)*x^2)^6)+We*(-15*lambda*(k-1)*(1+(1/2)*x^2)^9+(27*(-10*lambda*(k+1)*y*(1/27)+149*k^4*(1/432)-49*k^3*(1/108)+29*k^2*(1/72)-49*k*(1/108)+149/432-10*Q*lambda*(1/9)))*(1+(1/2)*x^2)^8-(1/2)*(9*(11*k-11))*(10*y^2*lambda*(1/33)-(1/36)*(5*(k^2-(14/11)*k+1))*(k+1)*y-Q*(k^2-(2/11)*k+1))*(1+(1/2)*x^2)^7+(27*(((1/8)*k^4-(1/2)*k^3+(3/4)*k^2-(1/2)*k+1/8-(10/9)*Q*lambda)*y^2+(1/108)*(145*k^2-250*k+145)*Q*(k+1)*y+23*Q^2*(k^2-(26/23)*k+1)*(1/6)))*(1+(1/2)*x^2)^6+(27*(4*k-4))*(-(1/32)*(k+1)*(k-1)^2*y^3+(1/4)*Q*(k-1)^2*y^2+5*Q^2*(k+1)*y*(1/8)+Q^3)*(1+(1/2)*x^2)^5+(27*(-(1/16)*(k-1)^4*y^4-3*Q*(k+1)*(k-1)^2*y^3*(1/4)+3*Q^2*(k-1)^2*y^2+5*Q^3*(k+1)*y*(1/3)+2*Q^4))*(1+(1/2)*x^2)^4+(27*(4*k-4))*Q*(-(1/8)*(k-1)^2*y^2-3*Q*(k+1)*y*(1/8)+Q^2)*y^2*(1+(1/2)*x^2)^3+(27*(-3*(k-1)^2*y^2*(1/2)-Q*(k+1)*y+2*Q^2))*Q^2*y^2*(1+(1/2)*x^2)^2-54*Q^3*y^4*(k-1)*(1+(1/2)*x^2)-27*Q^4*y^4)*Br*(y+1+(1/2)*x^2)*(y-1-(1/2)*x^2)/(20*(1+(1/2)*x^2)^12)


for j to numelems(k__vals) do sol[j] := eval(Analytical_sol, [k = k__vals[j], We = .1, Br = .3, x = 0, Q = Q__vals[j], lambda = `λ__vals`[j]]); plta[j] := plot(sol[j], y = 0 .. 1, axes = boxed, color = colors[j], linestyle = dash, thickness = 6) end do; plots:-display([seq(plta[j], j = 1 .. numelems(k__vals))], title = "Analytical Solution"); plots:-display([seq(plots:-display([pltn[j], plta[j]]), j = 1 .. numelems(k__vals))], title = "Numerical and Analytical Solutions")







The term "garbage" may have been a bit strong - but I was having a bad day.

However you did state

plot the arrow individually; it works. Plot the trail individually; it works

which for the code you posted,, is simply untrue. The "arrow" part works, but the "trail" part generates the error

Error, (in plots:-pointplot3d) points are not in the correct format

This tells me that you have not even tested the code which you posted. This fills me with great joy, because now I know I cannot trust any of your code and have to check it line by line.

Furthermore, the code overall is (I think) unnecessarily complicated. You don't need any "procedures", just correct use of the display() command - see the attached

When you post code here using the big green up-arrow in the Mapleprimes toolbar you have a choice between "Insert link" and "Insert contents". Using the latter means that code will be displayed (as in the attached) and animations will play. This is generally preferred behaviour - it gives potential responders a quick idea of what they might be getting into without downloading your worksheet and firing up Maple

  r:= rand(1 .. 10):
  testData:= Matrix( [ seq( [r(), r(), r()],
                            i = 1 .. 20
# Set the length of the "trail" and display the
# results as an animation
  lenTrail:= 5:
  ( [ seq
      ( display
        ( [ plottools:-sphere( [5,5,5],
            arrow( testData[j,..],
            pointplot3d( testData[max([1, j-lenTrail])..j,...],
        j=1..op([1,1], testData)





1 2 3 4 5 6 7 Last Page 1 of 199