Advanced Search

dharr

Dr. David Harrington

8902 Reputation

22 Badges

21 years, 136 days
University of Victoria
Professor or university staff
Victoria, British Columbia, Canada

Social Networks and Content at Maplesoft.com

Maple Application Center
Contact Dr. David Harrington
I am a retired professor of chemistry at the University of Victoria, BC, Canada. My research areas are electrochemistry and surface science. I have been a user of Maple since about 1990.

MaplePrimes Activity


These are answers submitted by dharr

typos plus...

dharr 8902
January 06 2026
2 8

I don't know about (24) - perhaps the paper has a typo in the ode (missing m?).

For (25) and (26), I corrected typos in magenta and the tests are verified.

I don't know about (30) or (31).

There seem to be typos in (40) and (41) but that doesn't solve the problems.

ode-17.mw

You didn't provide explanation about the Weierstraa cases (what is f3?). Maple has the WeierstrassP and WeierstrassPPrime.

<> construction...

dharr 8902
December 29 2025
4 3

The <> constructors are useful here:

ga_zero := <IdentityMatrix(2), ZeroMatrix(2); ZeroMatrix(2), -IdentityMatrix(2) >;

or

ga_zero := <<IdentityMatrix(2) | ZeroMatrix(2)>, <ZeroMatrix(2) | -IdentityMatrix(2)>>;

 

One way...

dharr 8902
December 29 2025
2 1

This method only works since your envelope is always increasing in size. Then you can use the running minimum and maximum values as your envelope.

numsolve-1229.mw

If we know enough......

dharr 8902
December 28 2025
2 2

What are we allowed to know here? I don't think you could find dio1 from dio without knowing the form we are seeking. If we do know that form, then we can find the missing numbers using solve as follows:

> 

restart
> 

dio1 := ((x-1)*(y-2)*(z-3))^2+((x-4)*(y-5)*(z-6))^2+((x-7)*(y-8)*(z-9))^2

(x-1)^2*(y-2)^2*(z-3)^2+(x-4)^2*(y-5)^2*(z-6)^2+(x-7)^2*(y-8)^2*(z-9)^2

> 

dio := expand(dio1)

3*x^2*y^2*z^2-36*x^2*y^2*z-30*x^2*y*z^2-24*x*y^2*z^2+126*x^2*y^2+432*x^2*y*z+93*x^2*z^2+360*x*y^2*z+312*x*y*z^2+66*y^2*z^2-1692*x^2*y-1476*x^2*z-1440*x*y^2-5040*x*y*z-1104*x*z^2-1080*y^2*z-948*y*z^2+6120*x^2+21096*x*y+18576*x*z+4554*y^2+16056*y*z+3540*z^2-79848*x-69300*y-61272*z+268452

Suppose we are allowed to know that we want an expression of the following form

> 

dioguess := ((x-a)*(y-b)*(z-c))^2+((x-d)*(y-e)*(z-f))^2+((x-g)*(y-h)*(z-i))^2

(x-a)^2*(y-b)^2*(z-c)^2+(x-d)^2*(y-e)^2*(z-f)^2+(x-g)^2*(y-h)^2*(z-i)^2

Then we require the following to be zero

> 

q := collect(dioguess-dio, {x, y, z}, distributed)

-268452+(a^2+d^2+g^2-66)*y^2*z^2+(-2*a^2*c-2*d^2*f-2*g^2*i+1080)*y^2*z+(-2*a^2*b-2*d^2*e-2*g^2*h+948)*y*z^2+(4*a^2*b*c+4*d^2*e*f+4*g^2*h*i-16056)*y*z+d^2*e^2*f^2+g^2*h^2*i^2+a^2*b^2*c^2+(-2*a*c^2-2*d*f^2-2*g*i^2+1440)*x*y^2+(4*a*b*c^2+4*d*e*f^2+4*g*h*i^2-21096)*x*y+(-2*a*b^2-2*d*e^2-2*g*h^2+1104)*x*z^2+(4*a*b^2*c+4*d*e^2*f+4*g*h^2*i-18576)*x*z+(a^2*b^2+d^2*e^2+g^2*h^2-3540)*z^2+(-2*a^2*b^2*c-2*d^2*e^2*f-2*g^2*h^2*i+61272)*z+(-2*a*b^2*c^2-2*d*e^2*f^2-2*g*h^2*i^2+79848)*x+(a^2*c^2+d^2*f^2+g^2*i^2-4554)*y^2+(-2*a^2*b*c^2-2*d^2*e*f^2-2*g^2*h*i^2+69300)*y+(b^2*c^2+e^2*f^2+h^2*i^2-6120)*x^2+(c^2+f^2+i^2-126)*x^2*y^2+(-2*b*c^2-2*e*f^2-2*h*i^2+1692)*x^2*y+(b^2+e^2+h^2-93)*x^2*z^2+(-2*b^2*c-2*e^2*f-2*h^2*i+1476)*x^2*z+(-2*d-2*g-2*a+24)*x*y^2*z^2+(4*a*c+4*d*f+4*g*i-360)*x*y^2*z+(4*a*b+4*d*e+4*g*h-312)*x*y*z^2+(-8*a*b*c-8*d*e*f-8*g*h*i+5040)*x*y*z+(-2*f-2*i-2*c+36)*x^2*y^2*z+(-2*e-2*h-2*b+30)*x^2*y*z^2+(4*b*c+4*e*f+4*h*i-432)*x^2*y*z

So we want the coefficients of each term to be zero. There are 26 coefficients

> 

eqs := {coeffs(q, {x, y, z})}; nops(%)

{-2*d-2*g-2*a+24, -2*e-2*h-2*b+30, -2*f-2*i-2*c+36, a^2+d^2+g^2-66, b^2+e^2+h^2-93, c^2+f^2+i^2-126, 4*a*b+4*d*e+4*g*h-312, -2*a*b^2-2*d*e^2-2*g*h^2+1104, -2*a^2*b-2*d^2*e-2*g^2*h+948, a^2*b^2+d^2*e^2+g^2*h^2-3540, 4*a*c+4*d*f+4*g*i-360, -2*a*c^2-2*d*f^2-2*g*i^2+1440, -2*a^2*c-2*d^2*f-2*g^2*i+1080, a^2*c^2+d^2*f^2+g^2*i^2-4554, 4*b*c+4*e*f+4*h*i-432, -2*b*c^2-2*e*f^2-2*h*i^2+1692, -2*b^2*c-2*e^2*f-2*h^2*i+1476, b^2*c^2+e^2*f^2+h^2*i^2-6120, -8*a*b*c-8*d*e*f-8*g*h*i+5040, 4*a*b*c^2+4*d*e*f^2+4*g*h*i^2-21096, 4*a*b^2*c+4*d*e^2*f+4*g*h^2*i-18576, 4*a^2*b*c+4*d^2*e*f+4*g^2*h*i-16056, -2*a*b^2*c^2-2*d*e^2*f^2-2*g*h^2*i^2+79848, -2*a^2*b*c^2-2*d^2*e*f^2-2*g^2*h*i^2+69300, -2*a^2*b^2*c-2*d^2*e^2*f-2*g^2*h^2*i+61272, a^2*b^2*c^2+d^2*e^2*f^2+g^2*h^2*i^2-268452}

26

Solve these equations for the unknowns {a,b,c...i}. Of course there are multiple solutions corresponding to permutations of the three terms in dioguess

> 

sols := solve(eqs, {a, b, c, d, e, f, g, h, i})

{a = 1, b = 2, c = 3, d = 7, e = 8, f = 9, g = 4, h = 5, i = 6}, {a = 4, b = 5, c = 6, d = 7, e = 8, f = 9, g = 1, h = 2, i = 3}, {a = 1, b = 2, c = 3, d = 4, e = 5, f = 6, g = 7, h = 8, i = 9}, {a = 7, b = 8, c = 9, d = 4, e = 5, f = 6, g = 1, h = 2, i = 3}, {a = 4, b = 5, c = 6, d = 1, e = 2, f = 3, g = 7, h = 8, i = 9}, {a = 7, b = 8, c = 9, d = 1, e = 2, f = 3, g = 4, h = 5, i = 6}

> 

eval(dioguess, sols[1])

(x-1)^2*(y-2)^2*(z-3)^2+(x-4)^2*(y-5)^2*(z-6)^2+(x-7)^2*(y-8)^2*(z-9)^2

> 

NULL

Download dio.mw

Finding generators...

dharr 8902
December 25 2025
1 0

See

Download Generators.mw

AI garbage?...

dharr 8902
December 25 2025
1 9

No, I don't think you can get EllipticE(x,-1) from that integral. I think MMA is wrong, and the AI probably got the wrong result from MMA.

Edit - see below - MMA and Maple have different definitions.

Download EllipticE.mw

simple example...

dharr 8902
December 23 2025
3 6

Here's a simple example - remembering a mutable structure like a Matrix is problematic.

This is a serious bug. I submitted an SCR.

forget.mw

AIC...

dharr 8902
December 17 2025
2 9

I didn't know about symbolic regression, but it seems from your responses to @sand15 that you might also be interested in automatic selection from a class of models. He has provided a nice summary of the issues. Since he mentioned AIC and provided nice code for rational function fitting, I thought I would illustrate the principle of AIC for those functions.

Basically the idea is that more fitting parameters generally improve the fit but after a while the extra parameters are not meaningful (overfitting). AIC quantifies this tradeoff into one parameter, AIC. I recklessly ignored the full rank warnings; this is just to give you some idea of the method.

One can also use an F-ratio test to see whether a model with another parameter is statistically significant. Here you need to choose a confidence level, and it is not always clear what to use. I used to use F-ratio (at the 1% level), but more recently switched to AIC, though I have the impression that it still overfits more that I would like; that may just prove I'm too cautious. AIC can select between different classes of models, but you have to generate them in some way, which I guess is what symbolic regression does.

At the end of the day your data looks suspiciously noise-free, and the key issue may be more of finding what type of function would look like that (with less that the 16 parameters AIC decided on) than the regression itself. That still seems like a human endeavour.

@sand15's FitRationalFraction procedure and the dataset are in the startup code

> 

restart;
> 

kernelopts(version);

`Maple 2025.2, X86 64 WINDOWS, Nov 11 2025, Build ID 1971053`

> 

X := <map2(op, 1, dataset)>:
Y := <map2(op, 2, dataset)>:
n := nops(dataset);

101

Fit all rational fraction models with numerator/denominator degrees between 5 and 10

> 

AIC:=Array(5..10,5..10):
> 

for ij in indices(AIC) do
  sij := ij[];
  F := FitRationalFraction(ij, X, Y, AllResults=true):
  S := F[1]["residualsumofsquares"];
  p := numelems(F[1]["parametervalues"]);
  fit[sij] := F[2];
  AIC[sij]:=evalf(n*ln(2*Pi)-n*ln(n)+n+n*ln(S)+2*(p+1));
end do:

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

Warning, model is not of full rank

In all rows and columns, the AIC gets better (more negative) and then worse (as overfitting occurs).

> 

AIC;

`Array(5..10, 5..10, {})`

We want the fit with the most negative AIC

> 

min(AIC);
min[index](AIC);
bestfit:=fit[%];

HFloat(-1772.0130337707733)

7, 9

y = (HFloat(5.639522374600771e-5)-HFloat(0.01288317349838677)*x-HFloat(0.06604775314173278)*x^2+HFloat(0.15402604574862966)*x^3-HFloat(0.10061405755123216)*x^4+HFloat(0.023236449348067213)*x^5-HFloat(0.0019176140897419658)*x^6+HFloat(4.3721095224913685e-5)*x^7)/(-HFloat(2.4682505243289557e-5)*x^9-HFloat(0.029674002394674673)*x^7+HFloat(0.21654010768516482)*x^6+HFloat(0.0018072506978093039)*x^8-HFloat(2.232927596702025)*x+HFloat(5.086432887393217)*x^2-HFloat(4.445396816679889)*x^3+HFloat(2.781300306907042)*x^4-HFloat(0.9769206707691607)*x^5+1)

> 

plots:-display(
  plot(rhs(bestfit), x=(min..max)(X), color="SteelBlue")
  , Statistics:-ScatterPlot(X, Y, symbol=circle, symbolsize=10, color=black)
  , view=[(min..max)(X), (min..max)(Y)]
  , size=[1000, 400]
)

> 

NULL

Download RationalFractionFitAIC.mw

A way...

dharr 8902
December 12 2025
2 0

I didn't try to figure out which solutions you might want; most of them have many a[i] and b[i] zero, so you may want to solve for different variables.

> 

restart
> 

with(PDEtools)
> 

with(LinearAlgebra)
> 

NULL
> 

with(SolveTools)
> 

undeclare(prime)

`There is no more prime differentiation variable; all derivatives will be displayed as indexed functions`

 (1)
> 

ode := (diff(diff(diff(diff(U(xi), xi), xi), xi), xi))*`&epsilon;`-(-6*w^2*`&epsilon;`+alpha)*(diff(diff(U(xi), xi), xi))-U(xi)*(3*w^4*`&epsilon;`-alpha*w^2-4*mu*U(xi)^2+k)

(diff(diff(diff(diff(U(xi), xi), xi), xi), xi))*epsilon-(-6*epsilon*w^2+alpha)*(diff(diff(U(xi), xi), xi))-U(xi)*(3*w^4*epsilon-alpha*w^2-4*mu*U(xi)^2+k)

 (2)
> 

F := sum(a[i]*G(xi)^i, i = 0 .. 2)+sum(b[i]*G(xi)^(-i), i = 1 .. 2)

a[0]+a[1]*G(xi)+a[2]*G(xi)^2+b[1]/G(xi)+b[2]/G(xi)^2

 (3)
> 

Y := G(xi) = sech(h*xi)
> 

Y1 := (diff(G(xi), xi))^2 = h^2*G(xi)^2*(1-G(xi)^2)
> 

odetest(Y, Y1)

0

 (4)
> 

Y2 := diff(G(xi), xi) = sqrt(h^2*G(xi)^2*(1-G(xi)^2))

diff(G(xi), xi) = (h^2*G(xi)^2*(1-G(xi)^2))^(1/2)

 (5)

Look, ma, no square roots!

> 

UU := U(xi) = dsubs(Y2, F)

U(xi) = (a[2]*G(xi)^4+a[1]*G(xi)^3+a[0]*G(xi)^2+b[1]*G(xi)+b[2])/G(xi)^2

 (6)
> 

L := numer(normal(eval(ode, UU)))

12*G(xi)^6*mu*a[1]^2*b[2]+12*G(xi)^6*mu*a[2]*b[1]^2-3*G(xi)^4*epsilon*w^4*b[2]+G(xi)^5*alpha*w^2*b[1]+12*G(xi)^5*mu*a[0]^2*b[1]+12*G(xi)^5*mu*a[1]*b[1]^2+G(xi)^4*alpha*w^2*b[2]+12*G(xi)^4*mu*a[0]^2*b[2]+12*G(xi)^4*mu*a[0]*b[1]^2+12*G(xi)^4*mu*a[2]*b[2]^2+12*G(xi)^3*mu*a[1]*b[2]^2+12*G(xi)^2*mu*a[0]*b[2]^2+12*G(xi)^2*mu*b[1]^2*b[2]+12*G(xi)*mu*b[1]*b[2]^2+2*G(xi)^7*(diff(diff(diff(diff(G(xi), xi), xi), xi), xi))*epsilon*a[2]+6*G(xi)^6*(diff(diff(G(xi), xi), xi))^2*epsilon*a[2]+G(xi)^6*(diff(diff(diff(diff(G(xi), xi), xi), xi), xi))*epsilon*a[1]-G(xi)^4*(diff(diff(diff(diff(G(xi), xi), xi), xi), xi))*epsilon*b[1]+6*G(xi)^3*(diff(diff(G(xi), xi), xi))^2*epsilon*b[1]+24*G(xi)*(diff(G(xi), xi))^4*epsilon*b[1]-2*G(xi)^3*(diff(diff(diff(diff(G(xi), xi), xi), xi), xi))*epsilon*b[2]+18*G(xi)^2*(diff(diff(G(xi), xi), xi))^2*epsilon*b[2]-2*G(xi)^7*(diff(diff(G(xi), xi), xi))*alpha*a[2]-2*G(xi)^6*(diff(G(xi), xi))^2*alpha*a[2]-G(xi)^6*(diff(diff(G(xi), xi), xi))*alpha*a[1]+G(xi)^4*(diff(diff(G(xi), xi), xi))*alpha*b[1]-2*G(xi)^3*(diff(G(xi), xi))^2*alpha*b[1]+2*G(xi)^3*(diff(diff(G(xi), xi), xi))*alpha*b[2]-6*G(xi)^2*(diff(G(xi), xi))^2*alpha*b[2]+12*G(xi)^11*mu*a[1]*a[2]^2+12*G(xi)^10*mu*a[0]*a[2]^2+12*G(xi)^10*mu*a[1]^2*a[2]-3*G(xi)^8*epsilon*w^4*a[2]+12*G(xi)^9*mu*a[2]^2*b[1]-3*G(xi)^7*epsilon*w^4*a[1]+G(xi)^8*alpha*w^2*a[2]+12*G(xi)^8*mu*a[0]^2*a[2]+12*G(xi)^8*mu*a[0]*a[1]^2+12*G(xi)^8*mu*a[2]^2*b[2]-3*G(xi)^6*epsilon*w^4*a[0]+G(xi)^7*alpha*w^2*a[1]+12*G(xi)^7*mu*a[0]^2*a[1]+12*G(xi)^7*mu*a[1]^2*b[1]-3*G(xi)^5*epsilon*w^4*b[1]+G(xi)^6*alpha*w^2*a[0]+120*(diff(G(xi), xi))^4*epsilon*b[2]+4*G(xi)^12*mu*a[2]^3+4*G(xi)^9*mu*a[1]^3-G(xi)^8*k*a[2]+4*G(xi)^6*mu*a[0]^3-G(xi)^7*k*a[1]-G(xi)^6*k*a[0]-G(xi)^5*k*b[1]+4*G(xi)^3*mu*b[1]^3-G(xi)^4*k*b[2]+4*mu*b[2]^3+6*G(xi)^6*(diff(diff(G(xi), xi), xi))*epsilon*w^2*a[1]-6*G(xi)^4*(diff(diff(G(xi), xi), xi))*epsilon*w^2*b[1]+12*G(xi)^3*(diff(G(xi), xi))^2*epsilon*w^2*b[1]-12*G(xi)^3*(diff(diff(G(xi), xi), xi))*epsilon*w^2*b[2]+36*G(xi)^2*(diff(G(xi), xi))^2*epsilon*w^2*b[2]+24*G(xi)^9*mu*a[0]*a[1]*a[2]+24*G(xi)^8*mu*a[1]*a[2]*b[1]+24*G(xi)^7*mu*a[0]*a[2]*b[1]+24*G(xi)^7*mu*a[1]*a[2]*b[2]+24*G(xi)^6*mu*a[0]*a[1]*b[1]+24*G(xi)^6*mu*a[0]*a[2]*b[2]+24*G(xi)^5*mu*a[0]*a[1]*b[2]+24*G(xi)^5*mu*a[2]*b[1]*b[2]+24*G(xi)^4*mu*a[1]*b[1]*b[2]+24*G(xi)^3*mu*a[0]*b[1]*b[2]+8*G(xi)^6*(diff(G(xi), xi))*(diff(diff(diff(G(xi), xi), xi), xi))*epsilon*a[2]+8*G(xi)^3*(diff(G(xi), xi))*(diff(diff(diff(G(xi), xi), xi), xi))*epsilon*b[1]-36*G(xi)^2*(diff(G(xi), xi))^2*(diff(diff(G(xi), xi), xi))*epsilon*b[1]+24*G(xi)^2*(diff(G(xi), xi))*(diff(diff(diff(G(xi), xi), xi), xi))*epsilon*b[2]-144*G(xi)*(diff(G(xi), xi))^2*(diff(diff(G(xi), xi), xi))*epsilon*b[2]+12*G(xi)^7*(diff(diff(G(xi), xi), xi))*epsilon*w^2*a[2]+12*G(xi)^6*(diff(G(xi), xi))^2*epsilon*w^2*a[2]

 (7)

Still have derivatives of G, so

> 

L1 := dsubs(Y2, L)

120*G(xi)^12*epsilon*h^4*a[2]+24*G(xi)^11*epsilon*h^4*a[1]-120*G(xi)^10*epsilon*h^4*a[2]-20*G(xi)^9*epsilon*h^4*a[1]+6*G(xi)^10*alpha*h^2*a[2]+16*G(xi)^8*h^4*epsilon*a[2]+2*G(xi)^9*alpha*h^2*a[1]+G(xi)^7*epsilon*h^4*a[1]-4*G(xi)^8*h^2*alpha*a[2]-8*G(xi)^6*epsilon*h^4*b[2]-G(xi)^7*h^2*alpha*a[1]+G(xi)^5*h^4*epsilon*b[1]+2*G(xi)^6*alpha*h^2*b[2]+16*G(xi)^4*h^4*epsilon*b[2]-G(xi)^5*h^2*alpha*b[1]-4*G(xi)^4*h^2*alpha*b[2]-36*G(xi)^10*epsilon*h^2*w^2*a[2]-12*G(xi)^9*epsilon*h^2*w^2*a[1]+24*G(xi)^8*h^2*epsilon*w^2*a[2]+6*G(xi)^7*h^2*epsilon*w^2*a[1]-12*G(xi)^6*epsilon*h^2*w^2*b[2]+6*G(xi)^5*h^2*epsilon*w^2*b[1]+24*G(xi)^4*h^2*epsilon*w^2*b[2]+12*G(xi)^6*mu*a[1]^2*b[2]+12*G(xi)^6*mu*a[2]*b[1]^2-3*G(xi)^4*epsilon*w^4*b[2]+G(xi)^5*alpha*w^2*b[1]+12*G(xi)^5*mu*a[0]^2*b[1]+12*G(xi)^5*mu*a[1]*b[1]^2+G(xi)^4*alpha*w^2*b[2]+12*G(xi)^4*mu*a[0]^2*b[2]+12*G(xi)^4*mu*a[0]*b[1]^2+12*G(xi)^4*mu*a[2]*b[2]^2+12*G(xi)^3*mu*a[1]*b[2]^2+12*G(xi)^2*mu*a[0]*b[2]^2+12*G(xi)^2*mu*b[1]^2*b[2]+12*G(xi)*mu*b[1]*b[2]^2+12*G(xi)^11*mu*a[1]*a[2]^2+12*G(xi)^10*mu*a[0]*a[2]^2+12*G(xi)^10*mu*a[1]^2*a[2]-3*G(xi)^8*epsilon*w^4*a[2]+12*G(xi)^9*mu*a[2]^2*b[1]-3*G(xi)^7*epsilon*w^4*a[1]+G(xi)^8*alpha*w^2*a[2]+12*G(xi)^8*mu*a[0]^2*a[2]+12*G(xi)^8*mu*a[0]*a[1]^2+12*G(xi)^8*mu*a[2]^2*b[2]-3*G(xi)^6*epsilon*w^4*a[0]+G(xi)^7*alpha*w^2*a[1]+12*G(xi)^7*mu*a[0]^2*a[1]+12*G(xi)^7*mu*a[1]^2*b[1]-3*G(xi)^5*epsilon*w^4*b[1]+G(xi)^6*alpha*w^2*a[0]+4*G(xi)^12*mu*a[2]^3+4*G(xi)^9*mu*a[1]^3-G(xi)^8*k*a[2]+4*G(xi)^6*mu*a[0]^3-G(xi)^7*k*a[1]-G(xi)^6*k*a[0]-G(xi)^5*k*b[1]+4*G(xi)^3*mu*b[1]^3-G(xi)^4*k*b[2]+4*mu*b[2]^3+24*G(xi)^9*mu*a[0]*a[1]*a[2]+24*G(xi)^8*mu*a[1]*a[2]*b[1]+24*G(xi)^7*mu*a[0]*a[2]*b[1]+24*G(xi)^7*mu*a[1]*a[2]*b[2]+24*G(xi)^6*mu*a[0]*a[1]*b[1]+24*G(xi)^6*mu*a[0]*a[2]*b[2]+24*G(xi)^5*mu*a[0]*a[1]*b[2]+24*G(xi)^5*mu*a[2]*b[1]*b[2]+24*G(xi)^4*mu*a[1]*b[1]*b[2]+24*G(xi)^3*mu*a[0]*b[1]*b[2]

 (8)
> 

indets(L1)

{alpha, epsilon, h, k, mu, w, xi, a[0], a[1], a[2], b[1], b[2], G(xi)}

 (9)
> 

eqs := {coeffs(collect(L1, G(xi)), G(xi))}

{4*mu*b[2]^3, 12*mu*b[1]*b[2]^2, 120*epsilon*h^4*a[2]+4*mu*a[2]^3, 12*mu*a[0]*b[2]^2+12*mu*b[1]^2*b[2], 24*epsilon*h^4*a[1]+12*mu*a[1]*a[2]^2, 24*mu*a[0]*b[1]*b[2]+12*mu*a[1]*b[2]^2+4*mu*b[1]^3, -120*epsilon*h^4*a[2]-36*epsilon*h^2*w^2*a[2]+6*alpha*h^2*a[2]+12*mu*a[0]*a[2]^2+12*mu*a[1]^2*a[2], -20*epsilon*h^4*a[1]-12*epsilon*h^2*w^2*a[1]+2*alpha*h^2*a[1]+24*mu*a[0]*a[1]*a[2]+4*mu*a[1]^3+12*mu*a[2]^2*b[1], 16*epsilon*h^4*b[2]+24*epsilon*h^2*w^2*b[2]-3*epsilon*w^4*b[2]-4*alpha*h^2*b[2]+alpha*w^2*b[2]+12*mu*a[0]^2*b[2]+12*mu*a[0]*b[1]^2+24*mu*a[1]*b[1]*b[2]+12*mu*a[2]*b[2]^2-k*b[2], epsilon*h^4*b[1]+6*epsilon*h^2*w^2*b[1]-3*epsilon*w^4*b[1]-alpha*h^2*b[1]+alpha*w^2*b[1]+12*mu*a[0]^2*b[1]+24*mu*a[0]*a[1]*b[2]+12*mu*a[1]*b[1]^2+24*mu*a[2]*b[1]*b[2]-k*b[1], 16*epsilon*h^4*a[2]+24*epsilon*h^2*w^2*a[2]-3*epsilon*w^4*a[2]-4*alpha*h^2*a[2]+alpha*w^2*a[2]+12*mu*a[0]^2*a[2]+12*mu*a[0]*a[1]^2+24*mu*a[1]*a[2]*b[1]+12*mu*a[2]^2*b[2]-k*a[2], epsilon*h^4*a[1]+6*epsilon*h^2*w^2*a[1]-3*epsilon*w^4*a[1]-alpha*h^2*a[1]+alpha*w^2*a[1]+12*mu*a[0]^2*a[1]+24*mu*a[0]*a[2]*b[1]+12*mu*a[1]^2*b[1]+24*mu*a[1]*a[2]*b[2]-k*a[1], -8*epsilon*h^4*b[2]-12*epsilon*h^2*w^2*b[2]-3*epsilon*w^4*a[0]+2*alpha*h^2*b[2]+alpha*w^2*a[0]+4*mu*a[0]^3+24*mu*a[0]*a[1]*b[1]+24*mu*a[0]*a[2]*b[2]+12*mu*a[1]^2*b[2]+12*mu*a[2]*b[1]^2-k*a[0]}

 (10)
> 

indets(eqs); nops(eqs)

{alpha, epsilon, h, k, mu, w, a[0], a[1], a[2], b[1], b[2]}

  

13

 (11)
> 

solution1 := solve(eqs, [alpha, mu, w, a[0], a[1], a[2], b[1], b[2]], explicit); nops(%)

28

 (12)
> 

ans1 := solution1[1]

[alpha = (epsilon*h^4+6*epsilon*h^2*w^2-3*epsilon*w^4-k)/((h-w)*(h+w)), mu = 0, w = w, a[0] = 0, a[1] = 0, a[2] = 0, b[1] = b[1], b[2] = 0]

 (13)
> 

ode1 := eval(ode, ans1)

(diff(diff(diff(diff(U(xi), xi), xi), xi), xi))*epsilon-(-6*w^2*epsilon+(epsilon*h^4+6*epsilon*h^2*w^2-3*epsilon*w^4-k)/((h-w)*(h+w)))*(diff(diff(U(xi), xi), xi))-U(xi)*(3*w^4*epsilon-(epsilon*h^4+6*epsilon*h^2*w^2-3*epsilon*w^4-k)*w^2/((h-w)*(h+w))+k)

 (14)
> 

UU1 := eval(eval(UU, ans1), Y)

U(xi) = b[1]/sech(h*xi)

 (15)
> 

odetest(UU1, ode1)

0

 (16)

NULL

Download test-F-p2.mw

simpler...

dharr 8902
December 09 2025
2 3

@sand15 I think you meant rhs(sol) in unapply. Then the second step can be simpler:

Sol := unapply(rhs(sol), x);
eval(IC, y = Sol);

Collecting with symbolic exponents...

dharr 8902
December 07 2025
4 6

Edit: handling of terms like ln(r) improved.

https://www.mapleprimes.com/questions/242062-Collect-Using-Pattern

> 

restart;

We want to collect A with respect to r, including symbolic exponents, in order to get B. I've added a constant term to the OP's example, to further test the method.
It is not clear from the question if the order matters, or what it would be for more complicated exponents.

For simplicity, assume that the order of exponents is just the default sort order, which would work here.

> 

A:=r^(2*a)+r^2*(1+a+r^(2*a)) + r + a*r + a;
B:=a+(1+a)*r+(1+a)*r^2+r^(2*a)+r^(2+2*a);

r^(2*a)+r^2*(1+a+r^(2*a))+r+a*r+a

a+(1+a)*r+(1+a)*r^2+r^(2*a)+r^(2+2*a)

Notice that Maple's automatic simplification may frustrate any attempt to get a desired order, e.g., the following is already collected, but does not retain the input order, and is perhaps not in the expected order

> 

r^3 + a*r^2 + r;

a*r^2+r^3+r

The reason collect does not work is that it is for polynomials, which have integer exponents but not symbolic ones.
The following finds the exponent wrt r for any exponent not containing r

> 

ex := (y,r) -> r*diff(y,r)/y;

proc (y, r) options operator, arrow; r*(diff(y, r))/y end proc

Examples. The presence of r in the second case indicates we should probably just leave this alone (this is done in the procedure below).

> 

ex(7*r^(2*tan(a)), r);
ex(a*ln(r), r);

2*tan(a)

1/ln(r)

Expand, find the terms and their corresponding exponents and coefficients. (assumes it is type `+`)

> 

terms := [op(expand(A))];
exps := expand(ex~(terms, r));
cfs := zip((t,e)->expand(t/r^e), terms, exps);

[(r^a)^2, r^2, a*r^2, r^2*(r^a)^2, r, a*r, a]

[2*a, 2, 2, 2+2*a, 1, 1, 0]

[1, 1, a, 1, 1, a, a]

Sort them

> 

exps2, p := sort(exps, output = [sorted, permutation]);
cfs2 := cfs[p];

[0, 1, 1, 2, 2, 2*a, 2+2*a], [7, 5, 6, 2, 3, 1, 4]

[a, 1, a, 1, a, 1, 1]

Assemble the results

> 

s:=0: cf:=cfs2[1]: expt:=exps2[1]:
for i from 2 to nops(exps2) do
  if exps2[i] = expt then
    cf += cfs2[i]
  else
    s += cf*r^expt;
    expt := exps2[i];
    cf := cfs2[i]
  end if;
end do;
s += cf*r^expt;

a+(1+a)*r+(1+a)*r^2+r^(2*a)+r^(2+2*a)

As a procedure.

> 

Collect := proc(expr, x::name)
local zz, ex, term, cfs, exps, p, s, cf, i;
if expr::function then return expr end if;
for i,term in expand(expr)+zz do # zz forces type `+`
  ex := expand(x*diff(term,x)/term);
  if has(ex, x) then
     cfs[i] := term;
     exps[i] := 0;
  else
     exps[i] := ex;
     cfs[i] := expand(term/x^ex)
  end if;
end do;
exps := convert(exps, list);
cfs := convert(cfs, list);
exps, p := sort(exps, 'output' = ['sorted', 'permutation']);
cfs := cfs[p];
s := 0;
cf := cfs[1];
ex := exps[1]:
for i from 2 to nops(exps) do
  if exps[i] = ex then
    cf += cfs[i]
  else
    s += cf*x^ex;
    ex := exps[i];
    cf := cfs[i]
  end if;
end do;
s += cf*x^ex - zz;
end proc:

Test

> 

Collect(A + ln(r), r);

a+ln(r)+(1+a)*r+(1+a)*r^2+r^(2*a)+r^(2+2*a)

> 

Collect(r^3 + a*r^2 + r, r);

a*r^2+r^3+r

> 

Collect(ln(x)+3*x^2, x);

ln(x)+3*x^2

x^f(x) is left alone, i.e., is treated as a constant term.

> 

Collect(3*x^2 + 2*x^x, x);

3*x^2+2*x^x

> 

 

Download CollectSymbolicExponents.mw

simplify with side relations...

dharr 8902
December 04 2025
1 3 
sqrt(1-cos(x)^2);
simplify(%,{cos(x)^2=1-sin(x)^2});

expand...

dharr 8902
December 01 2025
1 0

For some reason, the original simplify gives a numerator/denominator, with sqrt(1-t^2) in both.  If expand is used first then we don't have numerator/denominator, Then simplify can deal with the diffs. As well the denominator sqrt(1-t^2) might cancel or be simpler. So simplify(expand(B)) works here. The collect serves the same purpose as expand - giving priority to the diffs means we get a polynomial in diffs and not numerator/denominator, which then simplifies. The collect can be just wrt diff, so simplify(collect(B,diff))works.

Error in paper...

dharr 8902
November 23 2025
0 1

As far as I can see, there is some error in the paper. Maybe one of the other case(s) works.

Download pde.mw

Numbers...

dharr 8902
November 19 2025
3 1

For the numbers.

> 

restart

On a 64 bit machine

> 

kernelopts(wordsize);

64

> 

evalhf(DBL_MIN);  #$MinMachineNumber

0.100000010000000008e-306

> 

evalhf(DBL_MAX);  #$MaxMachineNumber

0.99999999000000001e308

> 

Maple_floats(MAX_FLOAT); #$MaxNumber

0.1e9223372036854775807

> 

Maple_floats(MIN_FLOAT); #$MinNumber

0.1e-9223372036854775805

> 

NULL

Download Numbers.mw

For the calling of previous commands by referring to them by the equation labels on the right-hand side, use ctrl-L to insert the number in the appropriate place.

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