MaplePrimes Questions

This ode can be solved by just looking at it


We see the solution is y=-x and y=c__1 because either (x+y)=0 or y'=0

But for some reason ODESteps(ode) says it cannot compute integral.

Any idea why?




`Standard Worksheet Interface, Maple 2024.0, Windows 10, March 01 2024 Build ID 1794891`


ode := (x+y(x))*(diff(y(x), x)) = 0

"[[,,"Let's solve"],[,,(x+y(x)) ((ⅆ)/(ⅆx) y(x))=0],["•",,"Highest derivative means the order of the ODE is" 1],[,,(ⅆ)/(ⅆx) y(x)],["•",,"Integrate both sides with respect to" x],[,,∫(x+y(x)) ((ⅆ)/(ⅆx) y(x)) ⅆx=∫0 ⅆx+`c__1`],["•",,"Cannot compute integral"],[,,∫(x+y(x)) ((ⅆ)/(ⅆx) y(x)) ⅆx=`c__1`]]"





Here is another simpler example that also confused it



`Standard Worksheet Interface, Maple 2024.0, Windows 10, March 01 2024 Build ID 1794891`


ode := x*(diff(y(x), x)) = 0

"[[,,"Let's solve"],[,,x ((ⅆ)/(ⅆx) y(x))=0],["•",,"Highest derivative means the order of the ODE is" 1],[,,(ⅆ)/(ⅆx) y(x)],["•",,"Integrate both sides with respect to" x],[,,∫x ((ⅆ)/(ⅆx) y(x)) ⅆx=∫0 ⅆx+`c__1`],["•",,"Cannot compute integral"],[,,∫x ((ⅆ)/(ⅆx) y(x)) ⅆx=`c__1`]]"






Here is another one which it gets wrong. 



`Standard Worksheet Interface, Maple 2024.0, Windows 10, March 01 2024 Build ID 1794891`


ode := y(x)*(diff(y(x), x)) = 0

"[[,,"Let's solve"],[,,y(x) ((ⅆ)/(ⅆx) y(x))=0],["•",,"Highest derivative means the order of the ODE is" 1],[,,(ⅆ)/(ⅆx) y(x)],["•",,"Integrate both sides with respect to" x],[,,∫y(x) ((ⅆ)/(ⅆx) y(x)) ⅆx=∫0 ⅆx+`c__1`],["•",,"Evaluate integral"],[,,((y(x))^2)/2=`c__1`],["•",,"Solve for" y(x)],[,,{y(x)=sqrt(2) sqrt(`c__1`),y(x)=-sqrt(2) sqrt(`c__1`)}]]"


y(x) = 0, y(x) = -c__1


The correct solution is given by dsolve, which is y=0 and y=constant (I do not know why dsolve put minus sign in front of the constant, but it is still correct).




I am stuck this command works seemlessly in Maple:

ThermophysicalData:-CoolProp:-Property(D, T = 20*Unit('degC'), P = 760*Unit('mmHg'), water)

but it does not work in Maple Flow. Does anyone knows why? Thank you so much for your help in the matter.

Why does 'simplify' not work when calculating Eigenvectors? Further, how can we express (2) in a more simplified form by using 'simplify'?

Is there a shortcut in Maple that will comment out a block of highlighted code? I tried searching on Mapleprimes but everytime I search I get a page generation error.

Hello :) 

I have a math problem, where I first need to use Linear regression to find the equation based on a set of data. I did that, no problem. 

However, in the next part of the problem I need to check if the residuals are under "normal distribution". Usually, I check if a dataset is normally distributed via "QQ-plot", and there will be no problems. But this time, because I need to check the residuals, I need to use the "residualQQplot(data,LinReg)" command to make it happen. But when I read the mean-value, mu, it says "-0," and nothing else? I know it should be "-3,2752*10^-15. 

The standard deviation is correct.

How do I fix this, so the residualQQplot shows me the right result? 

I have attached the worksheet here.

Thank you! 

what is the mathematics behind isolve? How can one selct special solutions?

Hello :) 

I don't know how to further simplify a complex equation as a solution to a differential equation. I have already tried "simplify(%)" etc. but it only repeats the same equation. I also know that it can be simplified to a much more "nice" equation in WordMat and TI-Nspire. I have attached a screenshot of my work. 

Thank you! 





I need to find the solution to the differential equation: y'=0.0768*f(x)^2/3-0.0102*f(x) when f(1)=59. 


dsolve({diff(f(x), x) = 0.768e-1*f(x)^(2/3)-0.102e-1*f(x), f(1) = 59}, f(x))

f(x) = 2097152/4913+(384/4913)*(16384-4352*59^(1/3)+289*59^(2/3))*exp(17/2500-(17/2500)*x)+(1/4913)*(835584*59^(1/3)-110976*59^(2/3)-1807285)*exp(51/5000-(51/5000)*x)+(49152/4913)*(17*59^(1/3)-128)*exp(17/5000-(17/5000)*x)



f(x) = 2097152/4913+(384/4913)*(16384-4352*59^(1/3)+289*59^(2/3))*exp(17/2500-(17/2500)*x)+(1/4913)*(835584*59^(1/3)-110976*59^(2/3)-1807285)*exp(51/5000-(51/5000)*x)+(49152/4913)*(17*59^(1/3)-128)*exp(17/5000-(17/5000)*x)




I can't understand this behavior. Any idea why it happens?

Solve is able to solve equation   f(y)=x+A for y, but can't solve   f(y)=x for y.

This is unexpected for me. I do not see why it can solve it when RHS is x+A but not when RHS is just x.




`Standard Worksheet Interface, Maple 2024.0, Windows 10, March 01 2024 Build ID 1794891`


`The "Physics Updates" version in the MapleCloud is 1745. The version installed in this computer is 1744 created 2024, April 17, 19:33 hours Pacific Time, found in the directory C:\Users\Owner\maple\toolbox\2024\Physics Updates\lib\`




(sin(y)+1)^(1/2)*(-2*sin(y)+2)^(1/2)*(-sin(y))^(1/2)*EllipticF((sin(y)+1)^(1/2), (1/2)*2^(1/2))/(cos(y)*sin(y)^(1/2))

Warning, solutions may have been lost

{%}; #to eliminate duplicates

(sin(y)+1)^(1/2)*(-2*sin(y)+2)^(1/2)*(-sin(y))^(1/2)*EllipticF((sin(y)+1)^(1/2), (1/2)*2^(1/2))/(cos(y)*sin(y)^(1/2))

{arctan(JacobiSN(((1/2)*I)*2^(1/2)*(x+b), (1/2)*2^(1/2))^2-1, -(1/2)*JacobiSN(((1/2)*I)*2^(1/2)*(x+b), (1/2)*2^(1/2))*(4-2*JacobiSN(((1/2)*I)*2^(1/2)*(x+b), (1/2)*2^(1/2))^2)^(1/2)*2^(1/2)), arctan(JacobiSN(((1/2)*I)*2^(1/2)*(x+b), (1/2)*2^(1/2))^2-1, (1/2)*JacobiSN(((1/2)*I)*2^(1/2)*(x+b), (1/2)*2^(1/2))*(4-2*JacobiSN(((1/2)*I)*2^(1/2)*(x+b), (1/2)*2^(1/2))^2)^(1/2)*2^(1/2))}


I can trick it to solve  f(y)=x for y  by asking it to solve f(y)=x+A for y and then set A=0 in the solution. But one should not have to do this. Is this a bug or Am I missing something?


Dear Colleague. 

I am trying to improve the results of abs(res[i] - exy) in the following codes.

Digits := 30:

# Define the function
f := proc(n)
    -0.5*y[n] + 0.5*sin(x[n] - Pi)
end proc:

# Define equations
e1 := y[n+2] = 2*h*delta[n] + y[n] - h^2*(-2*sin(u)*f(n)*u^2 - 2*sin(u)*f(n+2)*u^2 + 2*sin(2*u)*f(n+1)*u^2 + 2*cos(u)*f(n)*u - 2*cos(u)*f(n+2)*u + 2*cos(2*u)*f(n+1)*u - 2*cos(2*u)*f(n)*u - 2*sin(u)*f(n) + 2*sin(u)*f(n+2) + sin(2*u)*f(n) - sin(2*u)*f(n+2) - 2*f(n+1)*u + 2*f(n+2)*u)/((2*sin(u) - sin(2*u))*u^2):
e2 := y[n+1] = h*delta[n] + y[n] - (1/2)*h^2*(-sin(u)*f(n)*u^2 - sin(u)*f(n+2)*u^2 + sin(2*u)*f(n+1)*u^2 + 2*cos(u)*f(n)*u - 2*cos(u)*f(n+2)*u + 2*cos(2*u)*f(n+1)*u - 2*cos(2*u)*f(n)*u + 4*sin(u)*f(n+1) - 4*sin(u)*f(n) - 2*sin(2*u)*f(n+1) + 2*sin(2*u)*f(n) - 2*f(n+1)*u + 2*f(n+2)*u)/((2*sin(u) - sin(2*u))*u^2):
e3 := h*delta[n+2] = h*delta[n] + h^2*(2*sin(u)*f(n)*u + 2*sin(u)*f(n+2)*u - 2*sin(2*u)*f(n+1)*u - 2*cos(2*u)*f(n+1) + cos(2*u)*f(n) + cos(2*u)*f(n+2) + 2*f(n+1) - f(n) - f(n+2))/(u*(2*sin(u) - sin(2*u))):

epsilon := 10^(-10):
inx := 0:
ind := 1:
iny := 0:
h := 0.01:
n := 0:
omega := 1:
u := omega * h:
tol := 1e-4:
N := solve(h * p = 8 * Pi, p):

err := Vector(round(N)):
exy_lst := Vector(round(N)):

c := 1:
for j from 0 to 2 do
    t[j] := inx + j * h:
end do:

vars := y[n+1], y[n+2], delta[n+2]:

step := [seq(eval(x, x = c * h), c = 1 .. N)]:
printf("%6s%15s%15s%16s%15s%15s%15s\n", "h", "Num.y", "Num.z", "Ex.y", "Ex.z", "Error y", "Error z");

st := time():
for k from 1 to N / 2 do
    par1 := x[0] = t[0], x[1] = t[1], x[2] = t[2]:
    par2 := y[n] = iny, delta[n] = ind:    
    res := eval(<vars>, fsolve(eval({e1, e2, e3}, [par1, par2]), {vars}));

    for i from 1 to 2 do
        exy := eval(sin(c * h)):
        exz := eval(cos(c * h)):
        printf("%6.5f%17.9f%15.9f%15.9f%15.9f%13.5g%15.5g\n", h * c, res[i], res[i+1], exy, exz, abs(res[i] - exy), abs(res[i+1] - exz));
        err[c] := abs(evalf(res[i] - exy));
        if Norm(err) <= tol then 
            h := 0.1 * h * (c + 1) * (tol/Norm(err))^(0.2);
        end if;
        exy_lst[c] := exy;
        numerical_y1[c] := res[i];
        c := c + 1;
    end do;
    iny := res[2];
    ind := res[3];
    inx := t[2];
    for j from 0 to 2 do
        t[j] := inx + j * h;
    end do;
end do:
v := time() - st;
v / 4;
printf("Maximum error is %.13g\n", max(err));
NFE = evalf((N / 4 * 3) + 1);

# Get array of numerical and exact solutions for y1
numerical_array_y1 := [seq(numerical_y1[i], i = 1 .. N)]:
exact_array_y1 := [seq(exy_lst[i], i = 1 .. N)]:

# Get array of time steps
time_t := [seq(step[i], i = 1 .. N)]:

# Display graphs for y1
numerical_plot_y1 := plot(time_t, numerical_array_y1, style = point, symbol = asterisk, color = blue, symbolsize = 20, legend = ["TFIBF"]);
exact_plot_y1 := plot(time_t, exact_array_y1, style = point, symbol = box, color = red, symbolsize = 20, legend = ["EXACT"]);

display({numerical_plot_y1, exact_plot_y1});
Error_plot_y1 := plot(time_t, err, style = line, symbol = box, tickmarks = [piticks, decimalticks], color = navy, labels = [`h=Pi/8`, typeset(`Absolute Errors`)]);

I am suspecting that I didnt update the new h properly (I may be wrong, though). Please kindly help modify the code to allow the values of abs(res[i] - exy) to about 10^(-11). Thank you and best regards.

I am trying to solve the following recursion for any n, given a constant c. Here is my code for it:

c := 2:

A[i] := rsolve({a(0) = 1/n, a(i) = ((n - i + 1)/(n - i) + 1/(c*(n - i)))*a(i - 1)}, a);

total := evala(Simplify(sum(eval(A[i], i=k), k=0..n-1)));
evalf(eval(total, n = 6));

For c = 1, I get a valid (and correct) output, however for c = 2 for example, rsolve is returning A[i] = -GAMMA(-n)*GAMMA(-n + i - 1/2)/(GAMMA(-n - 1/2)*GAMMA(-n + i + 1)), which does not make sense when n is an integer. Is there something I am doing wrong here? Not sure why this is happening. Thanks!

Include print level in procedure
In the procedure code printlevel is not accepted
-enviroment variable
-interface variable 
error message : Error, (in interface) unknown interface variable, printlevel

More convenient in my opinion is to include on the printlevel depth in the procedure call?


fac := proc(n::integer)
    local previous_printlevel, result;
    previous_printlevel := interface('printlevel');  # Correct way to get the current printlevel
    interface('printlevel' = 3);  # Correct way to set the printlevel

    # De recursieve berekening
    if n = 0 then
        result := 1;  # Basisgeval
        result := n * fac(n - 1);  # Recursieve aanroep
    end if;

    interface('printlevel' = previous_printlevel);  # Restore the original printlevel
    return result;
end proc;


proc (n::integer) local previous_printlevel, result; previous_printlevel := interface('printlevel'); interface('printlevel' = 3); if n = 0 then result := 1 else result := n*fac(n-1) end if; interface('printlevel' = previous_printlevel); return result end proc


Error, (in interface) unknown interface variable, printlevel





I am very new to Maple. Many this are great, but I do not understand how maple deals with, especially, radiological units. In particular regarding joules [J] wich Maple seems to have alt least three types:

1) J - in relation to work

2) J(radiation) in relation to Gy (Gray) J/kg

3) J(dose_equivalent_index) in relation to Sv also J/kg

Why does Maple distinguish between these "joules"? As a phycisist they are all (well maby not entirely for Sv) equal to me. How can I make Maple treat them all at the "same joule"?

I have tried the following first:

with(Units[Standard]) and




Let S1: x^2 + (y - 2)^2 + (z + 1)^2 = 29 be a sphere and two points A(0, 0, 4), B(6, -2, 6); the line d passing through point C(4, -8, 4) and have direction v=(1, -1, 2). Find the point M such that M lies on the sphere S1, the angle AMB equals to 90 degree and distance from M to the line d is minimum

If I understand correctly, both of 

int(RETURN(is(y::positive)), y = 0 .. x) assuming 0 <= x, x < 1;
int(RETURN(coulditbe(y = 1)), y = 0 .. x) assuming 0 < x, x < 1;

should output `not`(true). However, Maple simply returns true for the second one. 
Isn't this result incorrect? Or am I missing something?

I would like to remove isomorphs from some graphs. That is to filter out non-isomorphic graphs.

graph_list := [GraphTheory:-CompleteGraph(3), GraphTheory:-PathGraph(3),Graph({{a,b},{b,c},{c,a}})]:

# Create a table to store non-isomorphic graphs
non_isomorphic_graphs := table():

# Counter for indexing the table
counter := 1:

# Iterate over each graph and check if it is isomorphic to any of the stored graphs
for g in graph_list do
    is_isomorphic := false:
    for key in indices(non_isomorphic_graphs,'nolist') do
        if GraphTheory:-IsIsomorphic(g, non_isomorphic_graphs[key]) then
            is_isomorphic := true:
        end if:
    end do:
    if not is_isomorphic then
        non_isomorphic_graphs[counter] := g:
        counter := counter + 1:
    end if:
end do:
DrawGraph~(non_isomorphic_graphs,  layoutoptions = [neutral_color = "pink", initial = spring])


A canonical form is a labeled graph Canon(G) that is isomorphic to G, such that every graph that is isomorphic to G has the same canonical form as G. I noticed that Maple has a function called CanonicalGraph. Can this function achieve the effect I want? I can easily achieve this by combining the  canonical form and property of sets  in  Sage.

graph_list = [Graph([(0, "a"), ("a", 2), (2, 0)]),graphs.PathGraph(3), graphs.CompleteGraph(3)]
non_isomorphic_graphs_labels = {g.canonical_label().copy(immutable=True) for g in graph_list}



An underlying motivation:My collaborators and I designed generation rules (algorithms) for 1-planar 4-trees;see Since the generating process is based on 1-planar embeddings, it will ultimately require filtering non-isomorphic graphs among a list of embeddings. I would be especially delighted to see that someone implement our algorithm in the future. Currently, I am stuck on handling some labeling details. It is somewhat similar to generating Apollonian networks (planar 3-trees). However, since its simplicial vertices are only two, the growth rate will not be too fast as the number of vertices increases.

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