Maple 2022 Questions and Posts

These are Posts and Questions associated with the product, Maple 2022

I'm not sure if this is fixed in newer Maples or if there's a work around. 

Whenever I size a 3d plot, that I'm trying to stretch out the width while keeping the height, never seems to work.  For example using the size=[3000,800] produces a plot area that's bigger but NOT actually a plot stretched in the x axis.  Going to size=[3000,3000] of course then makes the plot and the area bigger and so scales both x and y bigger.  However I don't want the y axis scaled up - I'm trying to scale the plot up - not the area.  And what happened to the window zoom, icon - we've lost zoom control to just magnify + and magnify - (at least in 2022) this seems like a regression. 

Is this plot size a bug or just a plot command that fails to function like it should? 

Hi everyone, I am trying to plot graphs for dp/dx versus x from my ordinary differential equation numerically. My file is working, but the outcome is straight lines, which means I am doing something wrong. Could anyone  please have a look on my file.

Help-dpdx.mw

the expected  results should be  look like this

I am very confused by the y-value of the rightmost point on the plot below.

restart

I'd like to find the values of x for which x^2/(10^(-8)-x) = 5*10^(-3).

So I ask Maple to solve this equation.

evalf(solve(x^2/(10^(-8)-x) = 5*10^(-3)))

-0.5000010000e-2, 0.10000e-7

(1)

Do these solutions work?

eval(x^2/(10^(-8)-x), x = 1.0000*10^(-8)) = Float(infinity) 

eval(x^2/(10^(-8)-x), x = -0.5000010000e-2) = 0.5000000000e-2 NULL

Suppose I define the function

f := proc (x) options operator, arrow; x^2/(1/100000000-x) end proc = proc (x) options operator, arrow; x^2/(1/100000000-x) end proc 

f(10^(-8))

Error, (in f) numeric exception: division by zero

 

f(.999999*10^(-8)) = 0.9999980000e-2NULL

f(.99999*10^(-8)) = 0.9999800001e-3NULL

Now, the function seems continuous between these two points

plot(f, .99999*10^(-8) .. .999999*10^(-8))

 

It is late, and perhaps I am just tired and not seeing things clearly. I expected the topmost point on the right to have a y-value of 0.00999998, ie almost 0.01.

I expected that the bottom leftmost point to be 0.0009999800001, ie almost 0.001, and it is.

And I thus expected to show that there must be some x for which we have f(x)=0.005, which if I am not mistaken is between the two numbers. After all, 0.999998e-2-0.5e-2 = 0.499998e-2NULL

0.5e-2-0.9999800001e-3 = 0.4000020000e-2NULL

what am i missing here?


Download plotq.mw

I have a question about a calculation I was just trying to do. For some context, I am trying to calculate the pH of a solution of a weak acid in water. When the concentration of the weak acid is low enough, we need to consider the effect of ionization of the water itself (ie, autoprotolysis of water), since the order of magnitude of this ionization is the same as the order of magnitude of the ions due to the weak acid in this case of very low concentration.

There is a specific formula that can be used for this, which is shown below.

K_a is a constant that depends on the weak acid being considered, and K_w is a constant as well (for the autoprotolysis of water). I am interested in solving for the concentration of hydronium, given an initial concentration of the weak acid [HA]_initial.

I'd like to solve the equation for different values of [HA]_initial.

In the calculations below, I am trying to solve for [H__3*LinearAlgebra:-Transpose(O)] in the expression

K__a = ([H__3*LinearAlgebra:-Transpose(O)])([H__3*LinearAlgebra:-Transpose(O)]-K__w/[H__3*LinearAlgebra:-Transpose(O)])/(HA__initial+[-H__3*LinearAlgebra:-Transpose(O)]+K__w/[H__3*LinearAlgebra:-Transpose(O)])

Below, I use x as the concentration of hydronium [H__3*LinearAlgebra:-Transpose(O)] and K__a = 5*10^(-3), K__w = 10^(-14), HA__initial = 10^(-3)

evalf(solve((x^2-10^(-14))/(10^(-3)-x+10^(-14)/x) = 5*10^(-3), x))

0.854101976e-3-0.4e-11*I, -0.5854101969e-2-0.465063510e-12*I, -0.9e-11+0.3865063510e-11*I

(1)

If I use a manual simplification (noting that K__w is extremely small compared to reasonable values of [H__3*LinearAlgebra:-Transpose(O)], then I am able to get real solutions.

evalf(solve(x^2/(10^(-3)-x) = 5*10^(-3), x))

-0.5854101966e-2, 0.854101966e-3

(2)

Now, looking at the complex solutions, the imaginary parts are very small and the real parts of two of the three solutions match the real solutions above, so I guess perhaps I could have just ignored those complex parts.

On the other hand, if I try to do similar calculations but with HA__initial = 10^(-6), I get

evalf(solve((x^2-10^(-14))/(10^(-6)-x+10^(-14)/x) = 5*10^(-3), x))

0.1009702e-5-0.3e-11*I, -0.5000999802e-2-0.332050808e-12*I, -0.9902e-8+0.3132050808e-11*I

(3)

evalf(solve(x^2/(10^(-6)-x) = 5*10^(-3), x))

-0.5000999800e-2, 0.999800e-6

(4)

So here, I don't see a complex solution that looks like the positive real solution above. Which means that I am not sure if the equation with the complex solutions (which is not simplified) is useful to me (I am doing various such calculations with different values of "`HA__initial`)".

Download concentrations.mw

Today I was working on some plots involving pH, which is defined as -log_10 ([hydronium]), that is, the negative of the base 10 logarithm of the concentration of hydronium in a solution.

I reached an expression for a variable x that is a function of an initial concentration C_i.

I wanted to plot the pH given by -log_10 (0.0001+x).

Note that x(0)=0, and so for this latter plot we should have the point (0, 4).

I am not able to see any part of the plot near (0,4), as can be seen below.

plot(-log[10](0.1e-3+x))

 

x := -0.2550000000e-2+0.5000000000e-4*sqrt(2601.+(2.000000000*10^6)*C__i)

-0.2550000000e-2+0.5000000000e-4*(2601.+2000000.000*C__i)^(1/2)

(1)

plot(-log[10](0.1e-3+x))

 

I want to see the plot being 4 at C__i = 0.

 

Note that subs({C__i = 0}, x) = 0. and evalf(subs({C__i = 0}, -log[10](0.1e-3+x))) = 4.000000000 

NULL

plot(-log[10](0.1e-3+x), C__i = 0 .. 1)

 

plot(-log[10](0.1e-3+x), C__i = 0 .. 1, view = [0 .. 1, 1 .. 4])

 

plot(-log[10](0.1e-3+x), C__i = 0 .. 1, view = [0 .. .1, 1 .. 4])

 

NULL

Download plotatzero.mw

I aim to conduct a numerical frequency sweep analysis on a nonlinear, coupled two-degree-of-freedom vibration model and compare the results with the analytical solution. However, I am currently facing two main difficulties:

  1. I am unsure how to determine whether the computed response has reached a steady state. If I simply use the maximum value to represent the steady-state amplitude, it can be misleading—since transient responses prior to reaching steady state may yield a higher peak, as seen in the uploaded code.

  2. I do not know how to properly select the steady-state result from the previous frequency as the initial condition for the next frequency step.
    2.mw

Hi!

I am studying Burger's equation, and I would like to see the steps that Maple takes to solve this.  "ShowSteps" doesn't seem to work.

Unfortunately, I am unable to share the worksheet I made.

Server Error - MaplePrimes

 
 

MaplePrimes
 
 
 
Connect with Maplesoft:
Questions  |  Posts  |  Tags  |  Users  |  Unanswered  |  Maplesoft Blog  |  Badges  |  Recent
© Maplesoft, a division of Waterloo Maple Inc. .  |  maplesoft.com  |  Terms of Use  |  Privacy  |  Consent Preferences  | Trademarks

 

loading

Error occurred during PDF generation. Please refresh the page and try again

How do we simplify the arguments of the exponential in (1)? Further how to express (1) into hyperbolic/trig functions? 
 

restart

with(LinearAlgebra)

Bans := -8*delta2*delta4*delta1*(delta3^2+delta4^2)*exp((-(2*I)*y*a*(delta3^2+delta4^2)*delta1^3+((2*I)*y*delta3^3*a+(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+((2*I)*y*delta4^2*a-I*t)*delta3+delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta1^2-(2*I)*(y*delta2^2*a-(1/2)*t)*(delta3^2+delta4^2)*delta1+(2*(I*y*delta3^3*a+(1/2)*(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-I*t*(1/2)+I*y*delta4^2*a)*delta3+(1/2)*delta4^2*(x+y+t)^2*(B2+B1+B3+B4)))*delta2^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))-8*delta2*delta3*delta4*(delta1^2+delta2^2)*exp(((2*I)*y*a*(delta3^2+delta4^2)*delta1^3+(-(2*I)*y*delta3^3*a+(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-(2*I)*y*delta4^2*a+I*t)*delta3+delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta1^2+(2*I)*(y*delta2^2*a-(1/2)*t)*(delta3^2+delta4^2)*delta1-(2*(I*y*delta3^3*a-(1/2)*(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-I*t*(1/2)+I*y*delta4^2*a)*delta3-(1/2)*delta4^2*(x+y+t)^2*(B2+B1+B3+B4)))*delta2^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))+(delta2^2-2*delta4*delta2+delta4^2+(delta3-delta1)^2)*((delta3+I*delta4)*delta2^2+(I*delta3^2+I*delta4^2)*delta2+delta1*(delta4^2+I*delta4*delta1+delta3*(delta3+delta1)))*exp((-2*y*a*(delta3^2+delta4^2)*delta2^3+(-2*y*delta4^3*a+2*(x+y+t)^2*(B1+B3)*delta4^2+(-2*a*delta3^2*y-t)*delta4+2*delta3^2*(x+y+t)^2*(B1+B3))*delta2^2-(2*(delta3^2+delta4^2))*(y*a*delta1^2+(1/2)*t)*delta2+2*delta1^2*(-y*delta4^3*a+(x+y+t)^2*(B1+B3)*delta4^2+(-(1/2)*t-y*a*delta3^2)*delta4+delta3^2*(x+y+t)^2*(B1+B3)))/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))+(delta2^2+2*delta4*delta2+delta4^2+(delta3-delta1)^2)*((delta3+I*delta4)*delta2^2+(-I*delta3^2-I*delta4^2)*delta2+delta1*(delta4^2+I*delta4*delta1+delta3*(delta3+delta1)))*exp((2*y*a*(delta3^2+delta4^2)*delta2^3+(-2*y*delta4^3*a+2*(x+y+t)^2*(B2+B3)*delta4^2+(-2*a*delta3^2*y-t)*delta4+2*delta3^2*(x+y+t)^2*(B2+B3))*delta2^2+(2*(delta3^2+delta4^2))*(y*a*delta1^2+(1/2)*t)*delta2+2*delta1^2*(-y*delta4^3*a+(x+y+t)^2*(B2+B3)*delta4^2+(-(1/2)*t-y*a*delta3^2)*delta4+delta3^2*(x+y+t)^2*(B2+B3)))/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))-(delta2^2+2*delta4*delta2+delta4^2+(delta3-delta1)^2)*((I*delta4-delta3)*delta2^2+(-I*delta3^2-I*delta4^2)*delta2+delta1*(-delta4^2+I*delta4*delta1-delta3*(delta3+delta1)))*exp((-2*y*a*(delta3^2+delta4^2)*delta2^3+(2*y*delta4^3*a+2*(x+y+t)^2*(B1+B4)*delta4^2+(2*a*delta3^2*y+t)*delta4+2*delta3^2*(x+y+t)^2*(B1+B4))*delta2^2-(2*(delta3^2+delta4^2))*(y*a*delta1^2+(1/2)*t)*delta2+(2*(y*delta4^3*a+(x+y+t)^2*(B1+B4)*delta4^2+(y*a*delta3^2+(1/2)*t)*delta4+delta3^2*(x+y+t)^2*(B1+B4)))*delta1^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))-(delta2^2-2*delta4*delta2+delta4^2+(delta3-delta1)^2)*((I*delta4-delta3)*delta2^2+(I*delta3^2+I*delta4^2)*delta2+delta1*(-delta4^2+I*delta4*delta1-delta3*(delta3+delta1)))*exp((2*y*a*(delta3^2+delta4^2)*delta2^3+(2*y*delta4^3*a+2*(x+y+t)^2*(B2+B4)*delta4^2+(2*a*delta3^2*y+t)*delta4+2*delta3^2*(x+y+t)^2*(B2+B4))*delta2^2+(2*(delta3^2+delta4^2))*(y*a*delta1^2+(1/2)*t)*delta2+(2*(y*delta4^3*a+(x+y+t)^2*(B2+B4)*delta4^2+(y*a*delta3^2+(1/2)*t)*delta4+delta3^2*(x+y+t)^2*(B2+B4)))*delta1^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))

-8*delta2*delta4*delta1*(delta3^2+delta4^2)*exp((-(2*I)*y*a*(delta3^2+delta4^2)*delta1^3+((2*I)*y*delta3^3*a+(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+((2*I)*y*delta4^2*a-I*t)*delta3+delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta1^2-(2*I)*(y*delta2^2*a-(1/2)*t)*(delta3^2+delta4^2)*delta1+2*(I*y*delta3^3*a+(1/2)*(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-((1/2)*I)*t+I*y*delta4^2*a)*delta3+(1/2)*delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta2^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))-8*delta2*delta3*delta4*(delta1^2+delta2^2)*exp(((2*I)*y*a*(delta3^2+delta4^2)*delta1^3+(-(2*I)*y*delta3^3*a+(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-(2*I)*y*delta4^2*a+I*t)*delta3+delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta1^2+(2*I)*(y*delta2^2*a-(1/2)*t)*(delta3^2+delta4^2)*delta1-2*(I*y*delta3^3*a-(1/2)*(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-((1/2)*I)*t+I*y*delta4^2*a)*delta3-(1/2)*delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta2^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))+(delta2^2-2*delta4*delta2+delta4^2+(delta3-delta1)^2)*((delta3+I*delta4)*delta2^2+(I*delta3^2+I*delta4^2)*delta2+delta1*(delta4^2+I*delta4*delta1+delta3*(delta3+delta1)))*exp((-2*y*a*(delta3^2+delta4^2)*delta2^3+(-2*y*delta4^3*a+2*(x+y+t)^2*(B1+B3)*delta4^2+(-2*a*delta3^2*y-t)*delta4+2*delta3^2*(x+y+t)^2*(B1+B3))*delta2^2-2*(delta3^2+delta4^2)*(y*a*delta1^2+(1/2)*t)*delta2+2*delta1^2*(-y*delta4^3*a+(x+y+t)^2*(B1+B3)*delta4^2+(-(1/2)*t-y*a*delta3^2)*delta4+delta3^2*(x+y+t)^2*(B1+B3)))/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))+(delta2^2+2*delta4*delta2+delta4^2+(delta3-delta1)^2)*((delta3+I*delta4)*delta2^2+(-I*delta3^2-I*delta4^2)*delta2+delta1*(delta4^2+I*delta4*delta1+delta3*(delta3+delta1)))*exp((2*y*a*(delta3^2+delta4^2)*delta2^3+(-2*y*delta4^3*a+2*(x+y+t)^2*(B2+B3)*delta4^2+(-2*a*delta3^2*y-t)*delta4+2*delta3^2*(x+y+t)^2*(B2+B3))*delta2^2+2*(delta3^2+delta4^2)*(y*a*delta1^2+(1/2)*t)*delta2+2*delta1^2*(-y*delta4^3*a+(x+y+t)^2*(B2+B3)*delta4^2+(-(1/2)*t-y*a*delta3^2)*delta4+delta3^2*(x+y+t)^2*(B2+B3)))/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))-(delta2^2+2*delta4*delta2+delta4^2+(delta3-delta1)^2)*((I*delta4-delta3)*delta2^2+(-I*delta3^2-I*delta4^2)*delta2+delta1*(-delta4^2+I*delta4*delta1-delta3*(delta3+delta1)))*exp((-2*y*a*(delta3^2+delta4^2)*delta2^3+(2*y*delta4^3*a+2*(x+y+t)^2*(B1+B4)*delta4^2+(2*a*delta3^2*y+t)*delta4+2*delta3^2*(x+y+t)^2*(B1+B4))*delta2^2-2*(delta3^2+delta4^2)*(y*a*delta1^2+(1/2)*t)*delta2+2*(y*delta4^3*a+(x+y+t)^2*(B1+B4)*delta4^2+(y*a*delta3^2+(1/2)*t)*delta4+delta3^2*(x+y+t)^2*(B1+B4))*delta1^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))-(delta2^2-2*delta4*delta2+delta4^2+(delta3-delta1)^2)*((I*delta4-delta3)*delta2^2+(I*delta3^2+I*delta4^2)*delta2+delta1*(-delta4^2+I*delta4*delta1-delta3*(delta3+delta1)))*exp((2*y*a*(delta3^2+delta4^2)*delta2^3+(2*y*delta4^3*a+2*(x+y+t)^2*(B2+B4)*delta4^2+(2*a*delta3^2*y+t)*delta4+2*delta3^2*(x+y+t)^2*(B2+B4))*delta2^2+2*(delta3^2+delta4^2)*(y*a*delta1^2+(1/2)*t)*delta2+2*(y*delta4^3*a+(x+y+t)^2*(B2+B4)*delta4^2+(y*a*delta3^2+(1/2)*t)*delta4+delta3^2*(x+y+t)^2*(B2+B4))*delta1^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))

(1)

argexp1 := simplify((-(2*I)*y*a*(delta3^2+delta4^2)*delta1^3+((2*I)*y*delta3^3*a+(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+((2*I)*y*delta4^2*a-I*t)*delta3+delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta1^2-(2*I)*(y*delta2^2*a-(1/2)*t)*(delta3^2+delta4^2)*delta1+(2*(I*y*delta3^3*a+(1/2)*(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-I*t*(1/2)+I*y*delta4^2*a)*delta3+(1/2)*delta4^2*(x+y+t)^2*(B2+B1+B3+B4)))*delta2^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)), size)

(-(2*I)*y*a*(delta3^2+delta4^2)*delta1^3+((2*I)*y*delta3^3*a+(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+((2*I)*y*delta4^2*a-I*t)*delta3+delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta1^2-(2*I)*(y*delta2^2*a-(1/2)*t)*(delta3^2+delta4^2)*delta1+2*(y*delta3^3*a*I+(1/2)*(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-((1/2)*I)*t+y*delta4^2*a*I)*delta3+(1/2)*delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta2^2)/((delta1+delta2*I)*(delta4*I-delta3)*(delta3+delta4*I)*(delta2*I-delta1))

(2)

 

argexp2 := ((2*I)*y*a*(delta3^2+delta4^2)*delta1^3+(-(2*I)*y*delta3^3*a+(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-(2*I)*y*delta4^2*a+I*t)*delta3+delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta1^2+(2*I)*(y*delta2^2*a-(1/2)*t)*(delta3^2+delta4^2)*delta1-(2*(I*y*delta3^3*a-(1/2)*(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-I*t*(1/2)+I*y*delta4^2*a)*delta3-(1/2)*delta4^2*(x+y+t)^2*(B2+B1+B3+B4)))*delta2^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1))

((2*I)*y*a*(delta3^2+delta4^2)*delta1^3+(-(2*I)*y*delta3^3*a+(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-(2*I)*y*delta4^2*a+t*I)*delta3+delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta1^2+(2*I)*(y*delta2^2*a-(1/2)*t)*(delta3^2+delta4^2)*delta1-2*(y*delta3^3*a*I-(1/2)*(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-((1/2)*I)*t+y*delta4^2*a*I)*delta3-(1/2)*delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta2^2)/((delta1+delta2*I)*(delta4*I-delta3)*(delta3+delta4*I)*(delta2*I-delta1))

(3)

simplify(argexp1+argexp2)

2*(x+y+t)^2*(B2+B1+B3+B4)

(4)

terms := op(Bans)

-8*delta2*delta4*delta1*(delta3^2+delta4^2)*exp((-(2*I)*y*a*(delta3^2+delta4^2)*delta1^3+((2*I)*y*delta3^3*a+(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+((2*I)*y*delta4^2*a-I*t)*delta3+delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta1^2-(2*I)*(y*delta2^2*a-(1/2)*t)*(delta3^2+delta4^2)*delta1+2*(I*y*delta3^3*a+(1/2)*(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-((1/2)*I)*t+I*y*delta4^2*a)*delta3+(1/2)*delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta2^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1))), -8*delta2*delta3*delta4*(delta1^2+delta2^2)*exp(((2*I)*y*a*(delta3^2+delta4^2)*delta1^3+(-(2*I)*y*delta3^3*a+(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-(2*I)*y*delta4^2*a+I*t)*delta3+delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta1^2+(2*I)*(y*delta2^2*a-(1/2)*t)*(delta3^2+delta4^2)*delta1-2*(I*y*delta3^3*a-(1/2)*(x+y+t)^2*(B2+B1+B3+B4)*delta3^2+(-((1/2)*I)*t+I*y*delta4^2*a)*delta3-(1/2)*delta4^2*(x+y+t)^2*(B2+B1+B3+B4))*delta2^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1))), (delta2^2-2*delta4*delta2+delta4^2+(delta3-delta1)^2)*((delta3+I*delta4)*delta2^2+(I*delta3^2+I*delta4^2)*delta2+delta1*(delta4^2+I*delta4*delta1+delta3*(delta3+delta1)))*exp((-2*y*a*(delta3^2+delta4^2)*delta2^3+(-2*y*delta4^3*a+2*(x+y+t)^2*(B1+B3)*delta4^2+(-2*a*delta3^2*y-t)*delta4+2*delta3^2*(x+y+t)^2*(B1+B3))*delta2^2-2*(delta3^2+delta4^2)*(y*a*delta1^2+(1/2)*t)*delta2+2*delta1^2*(-y*delta4^3*a+(x+y+t)^2*(B1+B3)*delta4^2+(-(1/2)*t-y*a*delta3^2)*delta4+delta3^2*(x+y+t)^2*(B1+B3)))/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1))), (delta2^2+2*delta4*delta2+delta4^2+(delta3-delta1)^2)*((delta3+I*delta4)*delta2^2+(-I*delta3^2-I*delta4^2)*delta2+delta1*(delta4^2+I*delta4*delta1+delta3*(delta3+delta1)))*exp((2*y*a*(delta3^2+delta4^2)*delta2^3+(-2*y*delta4^3*a+2*(x+y+t)^2*(B2+B3)*delta4^2+(-2*a*delta3^2*y-t)*delta4+2*delta3^2*(x+y+t)^2*(B2+B3))*delta2^2+2*(delta3^2+delta4^2)*(y*a*delta1^2+(1/2)*t)*delta2+2*delta1^2*(-y*delta4^3*a+(x+y+t)^2*(B2+B3)*delta4^2+(-(1/2)*t-y*a*delta3^2)*delta4+delta3^2*(x+y+t)^2*(B2+B3)))/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1))), -(delta2^2+2*delta4*delta2+delta4^2+(delta3-delta1)^2)*((I*delta4-delta3)*delta2^2+(-I*delta3^2-I*delta4^2)*delta2+delta1*(-delta4^2+I*delta4*delta1-delta3*(delta3+delta1)))*exp((-2*y*a*(delta3^2+delta4^2)*delta2^3+(2*y*delta4^3*a+2*(x+y+t)^2*(B1+B4)*delta4^2+(2*a*delta3^2*y+t)*delta4+2*delta3^2*(x+y+t)^2*(B1+B4))*delta2^2-2*(delta3^2+delta4^2)*(y*a*delta1^2+(1/2)*t)*delta2+2*(y*delta4^3*a+(x+y+t)^2*(B1+B4)*delta4^2+(y*a*delta3^2+(1/2)*t)*delta4+delta3^2*(x+y+t)^2*(B1+B4))*delta1^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1))), -(delta2^2-2*delta4*delta2+delta4^2+(delta3-delta1)^2)*((I*delta4-delta3)*delta2^2+(I*delta3^2+I*delta4^2)*delta2+delta1*(-delta4^2+I*delta4*delta1-delta3*(delta3+delta1)))*exp((2*y*a*(delta3^2+delta4^2)*delta2^3+(2*y*delta4^3*a+2*(x+y+t)^2*(B2+B4)*delta4^2+(2*a*delta3^2*y+t)*delta4+2*delta3^2*(x+y+t)^2*(B2+B4))*delta2^2+2*(delta3^2+delta4^2)*(y*a*delta1^2+(1/2)*t)*delta2+2*(y*delta4^3*a+(x+y+t)^2*(B2+B4)*delta4^2+(y*a*delta3^2+(1/2)*t)*delta4+delta3^2*(x+y+t)^2*(B2+B4))*delta1^2)/((delta1+I*delta2)*(I*delta4-delta3)*(delta3+I*delta4)*(I*delta2-delta1)))

(5)

NULL


 

Download argument.mw

I am taking a course in quantum mechanics and trying to do the calculations in Maple.

In the worksheet below I try to solve a system of equations to find some constants. 

It isn't very relevant that the domain is physics. 

The issues I have are with manipulations of expressions (specifically, simplifications).

restart

Levine - Ch. 2.4 - Particle in a Rectangular Well

NULL

Define

s__1 := sqrt(2*m*(V__0-E))/`ℏ` = 2^(1/2)*(m*(V__0-E))^(1/2)/`ℏ`NULL

s__2 := -s__1 = -2^(1/2)*(m*(V__0-E))^(1/2)/`ℏ`NULL

NULLs := sqrt(2*m*E)/`ℏ` 

2^(1/2)*(m*E)^(1/2)/`ℏ`

(1)

`ψ__1`, `ψ__2` and `ψ__3` are wavefunctions.

`ψ__1` := proc (x) options operator, arrow; C*exp('s__1'*x) end proc = proc (x) options operator, arrow; C*exp('s__1'*x) end procNULL

`ψ__2` := proc (x) options operator, arrow; A*cos('s'*x)+B*sin('s'*x) end proc = proc (x) options operator, arrow; A*cos('s'*x)+B*sin('s'*x) end procNULL

`ψ__3` := proc (x) options operator, arrow; G*exp('s__2'*x) end proc = proc (x) options operator, arrow; G*exp('s__2'*x) end procNULL

NULL

assume(G <> 0, C <> 0)

NULL

I would also like to assume that A and B cannot both be zero simultaneously.


As can be seen above, there are four unknown constants, "A,B,C,"and G.

 

We can obtain values for these unknowns by imposing boundary conditions.

 

Some of the boundary conditions involve the first derivatives of the wavefunctions.

`&psi;__1,d` := D(`&psi;__1`) = proc (z) options operator, arrow; C*s__1*exp(s__1*z) end procNULL

`&psi;__2,d` := D(`&psi;__2`) = proc (z) options operator, arrow; -A*s*sin(s*z)+B*s*cos(s*z) end procNULL

`&psi;__3,d` := D(`&psi;__3`) = proc (z) options operator, arrow; G*s__2*exp(s__2*z) end procNULL

NULL

The boundary conditions are

NULL

`&psi;__1`(0) = `&psi;__2`(0)*`&psi;__2`(l) and `&psi;__2`(0)*`&psi;__2`(l) = `&psi;__3`(l)*(D(`&psi;__1`))(0) and `&psi;__3`(l)*(D(`&psi;__1`))(0) = (D(`&psi;__2`))(0)*(D(`&psi;__2`))(l) and (D(`&psi;__2`))(0)*(D(`&psi;__2`))(l) = (D(`&psi;__3`))(l)

NULL

My question is how to find the constants in Maple.

NULLNULL

In this problem, E < V__0 and m > 0.

 

assume(E < V__0, m > 0)

 

Solving the first boundary condition is easy.

expr1 := solve(`&psi;__1`(0) = `&psi;__2`(0), {C}) = {C = A}

NULL

Next, I solve the third boundary condition and sub in the result from solving the first boundary condition. We get B in terms of A.

expr3 := eval(solve(`&psi;__1,d`(0) = `&psi;__2,d`(0), {B}), expr1) = {B = A*(m*(V__0-E))^(1/2)/(m*E)^(1/2)}NULL

NULL

Already here it is not clear to me why Maple does not cancel the m's.

 

Next, consider the second boundary condition

expr2 := solve(`&psi;__2`(l) = `&psi;__3`(l), {G}) = {G = (A*cos(s*l)+B*sin(s*l))/exp(s__2*l)}

 

We can obtain G in terms of just A by subbing in previously found relationships, as follows

expr2 := simplify(subs(expr3, expr2))

{G = exp(2^(1/2)*m^(1/2)*(V__0-E)^(1/2)*l/`&hbar;`)*A*(sin(2^(1/2)*m^(1/2)*E^(1/2)*l/`&hbar;`)*(V__0-E)^(1/2)+cos(2^(1/2)*m^(1/2)*E^(1/2)*l/`&hbar;`)*E^(1/2))/E^(1/2)}

(4)

NULL

Finally, consider the fourth boundary condition.

solve(`&psi;__2,d`(l) = `&psi;__3,d`(l), {G})

{G = (m*E)^(1/2)*(A*sin(2^(1/2)*(m*E)^(1/2)*l/`&hbar;`)-B*cos(2^(1/2)*(m*E)^(1/2)*l/`&hbar;`))/((m*(V__0-E))^(1/2)*exp(-2^(1/2)*(m*(V__0-E))^(1/2)*l/`&hbar;`))}

(5)

expr4 := simplify(subs(expr3, solve(`&psi;__2,d`(l) = `&psi;__3,d`(l), {G})))

{G = (sin(2^(1/2)*m^(1/2)*E^(1/2)*l/`&hbar;`)*E^(1/2)-cos(2^(1/2)*m^(1/2)*E^(1/2)*l/`&hbar;`)*(V__0-E)^(1/2))*exp(2^(1/2)*m^(1/2)*(V__0-E)^(1/2)*l/`&hbar;`)*A/(V__0-E)^(1/2)}

(6)

result := simplify(subs(expr2, expr4))

{exp(2^(1/2)*m^(1/2)*(V__0-E)^(1/2)*l/`&hbar;`)*A*(sin(2^(1/2)*m^(1/2)*E^(1/2)*l/`&hbar;`)*(V__0-E)^(1/2)+cos(2^(1/2)*m^(1/2)*E^(1/2)*l/`&hbar;`)*E^(1/2))/E^(1/2) = (sin(2^(1/2)*m^(1/2)*E^(1/2)*l/`&hbar;`)*E^(1/2)-cos(2^(1/2)*m^(1/2)*E^(1/2)*l/`&hbar;`)*(V__0-E)^(1/2))*exp(2^(1/2)*m^(1/2)*(V__0-E)^(1/2)*l/`&hbar;`)*A/(V__0-E)^(1/2)}

(7)

At this point, we have this complicated expression.

 

The calculations above come from the book "Quantum Chemistry" by Levine. The variable "result" above should, after a few more steps of manipulating the expression, give us

 

(2*E-V__0)*sin(sqrt(2*mE)*l/`&hbar;`) = 2*sqrt(-E^2+E*V__0)*cos(sqrt(2*mE)*l/`&hbar;`)

(2*E-V__0)*sin(2^(1/2)*mE^(1/2)*l/`&hbar;`) = 2*(-E^2+E*V__0)^(1/2)*cos(2^(1/2)*mE^(1/2)*l/`&hbar;`)

(8)

 

I've done the calculation by hand. It is not difficult. The exponentials and the A cancel on each side and we are left with an equation involving sin(sl) and cos(sl) which easily comes out to the desired equation above.

I would like to be able to do it here in Maple.

The first thing I would do at this point is to start eliminating some of the clutter. For example, I would like to replace sqrt(2*mE)/`&hbar;` with s.

 

This is a first task (that I have asked about in a separate question).

 

Then, I would want Maple to do the cancellations for me, but I have already used "simplify".

 

Finally, I think I would like to collect terms involving cos and sin.

 

Here is what happens when I try to get Maple to solve the equations for me in one go

 

solve({`&psi;__1`(0) = `&psi;__2`(0), `&psi;__1,d`(0) = `&psi;__2,d`(0), `&psi;__2`(l) = `&psi;__3`(l), `&psi;__2,d`(l) = `&psi;__3,d`(l)}, {A, B, C, G})

{A = 0, B = 0, C = 0, G = 0}

(9)

Download solve_constants.mw

In a calculation, I obtained the following expression from Maple.

expr1:=A*sqrt(-m*(-V__0 + E))/sqrt(m*E)

A*(-m*(-V__0+E))^(1/2)/(m*E)^(1/2)

(1)

The first question I have is why the "m's don't cancel."

Next, consider the expression

expr2:=exp(sqrt(2)*sqrt(m)*sqrt(-V__0 + E)*l*I/`&hbar;`)

exp(I*2^(1/2)*m^(1/2)*(-V__0+E)^(1/2)*l/`&hbar;`)

(2)

simplify(expr2,[E-V__0=y])

exp(I*2^(1/2)*m^(1/2)*y^(1/2)*l/`&hbar;`)

(3)

simplify(expr2,[2*m=y])

exp(I*y^(1/2)*(-V__0+E)^(1/2)*l/`&hbar;`)

(4)

Why doesn't the simplification below work?

simplify(expr2,[2*m*(E-V__0)=y])

exp(I*2^(1/2)*m^(1/2)*(-V__0+E)^(1/2)*l/`&hbar;`)

(5)

I also tried with other commands

eval(expr2,2*m=y)

exp(I*2^(1/2)*m^(1/2)*(-V__0+E)^(1/2)*l/`&hbar;`)

(6)

eval(expr2,E-V__0=y)

exp(I*2^(1/2)*m^(1/2)*y^(1/2)*l/`&hbar;`)

(7)

eval(expr2,2*m*(E-V__0)=y)

exp(I*2^(1/2)*m^(1/2)*(-V__0+E)^(1/2)*l/`&hbar;`)

(8)

algsubs(2*m=y,expr2)

exp(I*y^(1/2)*(-V__0+E)^(1/2)*l/`&hbar;`)

(9)

algsubs(E-V__0=y,expr2)

exp(I*2^(1/2)*m^(1/2)*y^(1/2)*l/`&hbar;`)

(10)

algsubs(2*m*(E-V__0)=y,expr2)

exp(I*2^(1/2)*m^(1/2)*(-V__0+E)^(1/2)*l/`&hbar;`)

(11)
 

NULL

In reality, I would like to the entirety of sqrt(2m(E-V__0)/hbar) a variable by the name of s.

Download simplify_roots.mw

Below is a simple demonstration of something that I am very confused about in another worksheet I am working on.

Basically, I have a function of x called psi_I and I want to define another function of x called psi_I,d which is just the derivative of psi_I.

In, the snippet below, this works if I use the variable names psi_1 and psi_1,d, but not if I replace the "1" with a "I" as I want to do.

In addition, I don't know why the output (6) is different from (7) or (8). All three are from the same command in Maple.

restart

`&psi;__1` := proc (x) options operator, arrow; C*exp(a*x) end proc

proc (x) options operator, arrow; C*exp(a*x) end proc

(1)

`&psi;__1,d` := unapply(diff(`&psi;__1`(x), x), x)

proc (x) options operator, arrow; C*a*exp(a*x) end proc

(2)

`&psi;__1,d`(x)

C*a*exp(a*x)

(3)

NULL

restart

`&psi;__I` := proc (x) options operator, arrow; C*exp(a*x) end proc

proc (x) options operator, arrow; C*exp(a*x) end proc

(4)

`&psi;__I,d` := unapply(diff(`#msub(mi("psi",fontstyle = "normal"),mi("1"))`(x), x), x)

D(psi__1)

(5)

`#msub(mi("psi",fontstyle = "normal"),mi("1,d"))`(x)

`psi__1,d`(x)

(6)

`#msub(mi("psi",fontstyle = "normal"),mi("I,d"))`(x)

(D(psi__1))(x)

(7)

`&psi;__I,d`(x)

(D(psi__1))(x)

(8)

NULL

Download psi_subscript.mw

restart

with(LinearAlgebra)

A := `<,>`(`<|>`(-2*`&omega;__0`^2, `&omega;__0`^2), `<|>`(`&omega;__0`^2, -`&omega;__0`^2))

Matrix(%id = 36893488151921005020)

(1)

AA := subs(`&omega;__0` = 3, A)

Matrix(%id = 36893488151921002244)

(2)

evalf([Eigenvectors(AA)])

[Vector[column](%id = 36893488151920990196), Matrix(%id = 36893488151920990316)]

(3)

NULL

Below is an eigenvector equation for AA for the eigenvalue -3.4376...

 

S := AA+3.43769410*IdentityMatrix(2)

Matrix(%id = 36893488151920987180)

(4)

When I solve the system below, I expect to get the vector `<,>`(.6180, 1)

Vector[column](%id = 36893488151920981284)

(5)

, just like the result of the Eigenvectors command used above.

 

LinearSolve(S, `<,>`(0, 0))

Vector[column](%id = 36893488151920975988)

(6)

Instead I get the zero vector.

 

In addition, when I compute the reduced row echelon form, I expect to have one of the rows be zero since matrix S has rank 1.

 

Rank(S)

1

(7)

ReducedRowEchelonForm(S)

Matrix(%id = 36893488151959135044)

(8)

NULL

Therefore, I guess my question is maybe about calculations and not about Maple. That being said, I did the calculations manually, then I asked chatgpt, which agreed with me.

 

Why are the calculations not coinciding with Maple?

 

I can do the calculations "manually" with Maple:

 

solve(-14.5623059000000*x+9*y = 0)

{x = x, y = 1.618033989*x}

(9)

solve(9*x-5.56230590000000*y = 0)

{x = x, y = 1.618033988*x}

(10)

 

Now, if I try to solve them simultaneously, it seems I have the same issue as with LinearSolve

NULL

solve({9*x-5.56230590000000*y = 0, -14.5623059000000*x+9*y = 0})

{x = 0., y = 0.}

(11)

I tried something else that worked, but the questions above still remain.

 

We start with A again

 

A = Matrix(%id = 36893488151921005020)NULL

In my problem I have the variable

 

`&omega;__a` := (1/2)*`&omega;__0`*(sqrt(5)-1) = (1/2)*omega__0*(5^(1/2)-1)NULL

NULL

and it turns out that -`#msub(mi("omega",fontstyle = "normal"),mi("a"))`^2 is an eigenvalue of A.

In particular, for `&omega;__0` = 3 we get the values I showed previously.

 

For example, here is the eigenvalue we were looking at beforeevalf(subs(`&omega;__0` = 3, -`&omega;__a`^2))

-3.437694099

(12)

Okay so the eigenvector equation without subbing in a value for `&omega;__0` is (A+I*`&omega;__a`^2)*v = 0.

 

Here is what the matrix on the left looks like

NULL

A+`&omega;__a`^2*IdentityMatrix(2)

Matrix(%id = 36893488151964015780)

(13)

 

And here we solve for the eigenvectors

sol := LinearSolve(A+`&omega;__a`^2*IdentityMatrix(2), `<,>`(0, 0))

Vector[column](%id = 36893488151964001572)

(14)

evalf(subs(_t2[1] = 1, sol))

Vector[column](%id = 36893488151964010244)

(15)

Which is the expected result.

 

So, Maple can solve the eigenvector equation in this more general case but not if I specify `&omega;__0` before asking it to solve.NULL

Download evcalc.mw

The document below is long because I go through the context of my question step-by-step. 

The context is solving a first-order system of differential equations where we have complex eigenvalues.

But essentially, my question has to do with the end of this reasoning. I find four solutions x1, x2, x3, and x4, and I would like to check whether they are indeed solutions by subbing them into the original system x'=Ax.

When we make such a substitution, we get an equation of the form yl=yr, where yl and yr are 4x1 vectors. If the equality is true then we can say we have a solution.

My question is about how to check for this equality. I amusing Maple's "Equals" but this doesn't seem to work for all four solutions.

with(LinearAlgebra)

NULL

I explain the entire reasoning that leads to the final vectors that I compare

NULL

The following matrix has complex eigenvalues 3+`&+-`(4*I), each with multiplicity 2.

 

A := `<,>`(`<|>`(3, -4, 1, 0), `<|>`(4, 3, 0, 1), `<|>`(0, 0, 3, -4), `<|>`(0, 0, 4, 3))

Matrix(%id = 36893488151974184588)

(1)

NULL

Consider the eigenvalue `&lambda;__1` = 3-4*I. The associated eigenvectors are

evec := LinearSolve(A-(3-4*I)*IdentityMatrix(4), `<,>`(0, 0, 0, 0), free = c)

Vector[column](%id = 36893488151974164828)

(2)

We see this gives us only one linearly independent eigenvector.

 

We can find another one by forming a rank 2 generalized eigenvector from one of the eigenvectors above.

 

The pair of such eigenvectors is a chain of generalized eigenvectors of length 2 and they satisfy

 

(A-I*`&lambda;__1`)*w__2 = w__1*(A-I*`&lambda;__1`)*w__1 and w__1*(A-I*`&lambda;__1`)*w__1 = 0

NULL

These equations come from trying the solution w__1*e^(`&lambda;__1`*t) and (t*w__1+w__2)*exp(`&lambda;__1`*t), where w__1 is an eigenvector for `&lambda;__1`.

The equations above lead to

 

(A-I*`&lambda;__1`)^2*w__2 = 0

 

which we can solve for the generalized eigenvector w__2 as follows

NULL

`w__2,gen` := LinearSolve((A-(3-4*I)*IdentityMatrix(4))^2, `<,>`(0, 0, 0, 0), free = c)

Vector[column](%id = 36893488151974138212)

(3)

One such generalized eigenvector is

 

w__2 := subs({c[2] = I, c[4] = I}, `w__2,gen`)

Vector[column](%id = 36893488151982985812)

(4)

which gives us

 

w__1 := (A-(3-4*I)*IdentityMatrix(4)).w__2

Vector[column](%id = 36893488151982979180)

(5)

which is indeed an eigenvector for `&lambda;__1`.

Thus, we have the two complex solutions

 

s__1 := w__1*exp((3-4*I)*t)

Vector[column](%id = 36893488151982971708)

(6)

s__2 := (t*w__1+w__2)*exp((3-4*I)*t)

Vector[column](%id = 36893488151982963516)

(7)

But we want real solutions.

NULL

x__1 := `assuming`([Re(s__1)], [t > 0])

Vector[column](%id = 36893488151974145788)

(8)

x__2 := `assuming`([Im(s__1)], [t > 0])

Vector[column](%id = 36893488151982991948)

(9)

x__3 := `assuming`([Re(s__2)], [t > 0])

Vector[column](%id = 36893488151982977980)

(10)

x__4 := `assuming`([Im(s__2)], [t > 0])

Vector[column](%id = 36893488151982967252)

(11)

Let's check that each of the vectors above (the real and imaginary parts of our complex solutions) are real solutions.

 

I would like to substitute, say, x__1 into the first order system of differential equations diff(x, x) = Ax.

 

At first I tried simply

NULL

subs({x = x__1}, diff(x, t) = A.x)

0 = Matrix(%id = 36893488151974184588).Vector[column](%id = 36893488151982948468)

(12)

This doesn't work because the diff(x, t) command evaluates to zero before we do the desired substitution.

 

Then I tried

 

subs({x = x__1}, 'diff(x, t)' = A.x)

diff(Vector[column](%id = 36893488151982948468), t) = Matrix(%id = 36893488151974184588).Vector[column](%id = 36893488151982948468)

(13)

 

   

result := eval(subs({x = x__1}, 'diff(x, t)' = A.x))

Vector[column](%id = 36893488151991655348) = Vector[column](%id = 36893488151982967612)

(14)

NULL

NULL

By inspection, it seems the two sides are equal and so x__1 indeed satisfies the system.

 

But how do I get maple to tell me this?

 

evalb(result)

false

(15)

Equal(rhs(result), lhs(result))

true

(16)

Okay, but does this use of Equals always work?

 

sols := [x__1, x__2, x__3, x__4]; for i to 4 do result := eval(subs({x = sols[i]}, 'diff(x, t)' = A.x)); print(simplify(result))*print(Equal(rhs(result), lhs(result))) end do

true

(17)

 

Apparently not.``

 

What is going on here?

It sure looks like the righthand and lefthand sides of each expression are the same.

Download evs_equal.mw

Dear Maple Community,

I come to you with a question about the reduced involutive form (rif) package. Namely, I decided to try the classic example from the "LONG GUIDE TO THE STANDARD FORM PACKAGE", which dates back to 1993. Here is the link to the complete documentation:

https://wayback.cecm.sfu.ca/~wittkopf/files/standard_manual.txt

So, the example is the following:

2.1 SIMPLE EXAMPLES

EXAMPLE A

Consider the system of nonlinear PDEs:       

y Zxxx - x Zxyy  =  Zyy - y Zy

                        2     2    2
2 y x Zxxx Zxyy + x Zxxx + x y Zxyy  =  0

                  2    2
y Zxyy - x W + 2 x  y Z  =  0

                 2    2
Zyy - y Zy  + 2 x  y W  =  x W

where the dependent variables W and Z are functions of the
independent variables x and y, and Zxxx denotes the third partial
derivative of Z with respect to x etc.

After making computations back in 1993 with Maple V, they obtain the following involutive form:

In our original notation the (considerably) simplified system is:
                                            2
  Zxxx = 0, Zxy = 0, Zyy = y Zy, W = 2 x y Z

So, I tried the same system of PDEs in the modern Maple and the modern rifsimp() command. You can find the result below:

demo_question.mw

The problem is that nowadays [Maple 2022.1] , I get only the trivial solution: z = 0 and w = 0.

Could someone clarify, please, where the truth is and what am I doing wrong?

Thanks a lot in advance for any help and clarification!

Best regards,

Dr. Denys D.
 

restart:

with(DETools):

PDE1 := y*diff(z(x,y), x$3) - x*diff(z(x,y),x,y$2) = diff(z(x,y),y$2) - y*diff(z(x,y), y);

y*(diff(diff(diff(z(x, y), x), x), x))-x*(diff(diff(diff(z(x, y), x), y), y)) = diff(diff(z(x, y), y), y)-y*(diff(z(x, y), y))

(1)

PDE2 := 2*x*y*diff(z(x,y),x$3)*diff(z(x,y),x,y$2) + x*(diff(z(x,y),x$3))^2 + x*y^2*(diff(z(x,y),x,y$2))^2 = 0;

2*x*y*(diff(diff(diff(z(x, y), x), x), x))*(diff(diff(diff(z(x, y), x), y), y))+x*(diff(diff(diff(z(x, y), x), x), x))^2+x*y^2*(diff(diff(diff(z(x, y), x), y), y))^2 = 0

(2)

PDE3 := y*diff(z(x,y),x,y$2) - x*w(x,y) + 2*x^2*y*z(x,y)^2 = 0;

y*(diff(diff(diff(z(x, y), x), y), y))-x*w(x, y)+2*x^2*y*z(x, y)^2 = 0

(3)

PDE4 := diff(z(x,y), y$2) - y*diff(z(x,y),y) + 2*x^2*y*w(x,y)^2 = x*w(x,y);

diff(diff(z(x, y), y), y)-y*(diff(z(x, y), y))+2*x^2*y*w(x, y)^2 = x*w(x, y)

(4)

sys := [PDE1, PDE2, PDE3, PDE4]:

rif := rifsimp(sys, [[w], [z]], indep = [x,y]);

table( [( Case ) = [[z(x, y)*(8*z(x, y)^2*y^2*x^2-1) = 0, diff(z(x, y), x), "false split"]], ( Solved ) = [w(x, y) = 0, z(x, y) = 0] ] )

(5)
 

 

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