Advanced Search

dharr

Dr. David Harrington

8637 Reputation

22 Badges

21 years, 73 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 replies submitted by dharr

solve for w...

dharr 8637
6 hours ago
0 0

@Andiguys Not following what tou are doing here, but you didn't include w in the variables to solve for, so maybe that is it?

native isolve...

dharr 8637
Yesterday at 2:54 PM

@Alfred_F Maple's isolve does this directly and correctly.

Ellipse.mw

simplified cubic equation...

dharr 8637
Yesterday at 9:16 AM
0 0

@Andiguys Solve for w separately, then Pn is found by solving a cubic polynomial, which can be given in a simple way.

Inside_root.mw

Inversion in Maple....

dharr 8637
November 05 2025
0 0

@dharr I realized that for the cases where Algfield gives one number as a polynomial of the other, inversion can be done by simply solving that polynomial. 

> 

restart;
> 

local gamma:

Example 1.

> 

a := RootOf(4*_Z^3 - 12*_Z^2 - 24*_Z - 25, index=1);
b := RootOf(2*_Z^3 - 3, index=1);
alias(alpha = a, beta = b);

RootOf(4*_Z^3-12*_Z^2-24*_Z-25, index = 1)

RootOf(2*_Z^3-3, index = 1)

alpha, beta

Algfield finds a in terms of b as a polynomial, so they a is in the same field defined by n

> 

evala(Algfield({a,b}));
> 

eqa := %[1][];

[[alpha = beta^2+2*beta+1], [], {beta}, true, {}]

alpha = beta^2+2*beta+1

So to invert, we just treat beta as an unknown x and solve
roots and evala(Roots()) work in the same way

> 

eq := alpha-subs(beta = x, rhs(eqa));
roots(eq);

alpha-x^2-2*x-1

[[(2/13)*alpha^2-(12/13)*alpha-28/13, 1], [-(2/13)*alpha^2+(12/13)*alpha+2/13, 1]]

Example 2

> 

c := RootOf(_Z^4 + 1, index=1);
d := RootOf(_Z^4 + 6*_Z^2 + 1, index=1);
alias(gamma = c, delta = d);

RootOf(_Z^4+1, index = 1)

RootOf(_Z^4+6*_Z^2+1, index = 1)

alpha, beta, gamma, delta

> 

evala(Algfield({c,d}));
> 

eqd := %[1][];
eq := delta-subs(gamma = x, rhs(eqd));
roots(eq);

[[delta = gamma^3-gamma^2+gamma], [], {gamma}, true, {}]

delta = gamma^3-gamma^2+gamma

delta-x^3+x^2-x

[[(1/4)*delta^3+(1/4)*delta^2+(7/4)*delta+3/4, 1]]

> 

NULL

Download evala2.mw

General procedure...

dharr 8637
November 04 2025

@Alfred_F Here is a general procedure, tested on your two examples. If it isn't clear or needs more features, please let me know. Choose more _Z1 values for more solutions, but they get very large very quickly.

Procedure QuadraticIsolve is in the startup code edit region.

It accepts a quadratic in two variables and returns integer solutions.
If there is no second argument it returns all solutions if there are a finite number of solutions

or general solutions parameterized by variable _Z1 if there are an infinite number of solutions.

If a set or range of integers is given as a second argument, then solutions corresponding only

to these values of _Z1 are given.
Use infolevel[QuadraticIsolve]:=2; to generate qualitative information about the solving process

and infolevel[QuadraticIsolve]:=3; for more quantitative information.

> 

restart
> 

infolevel[QuadraticIsolve] := 3

3

First case  -  there are 128 general solutions involving an arbitrary integer parameter _Z1

> 

P := 14*x^2-5*x*y-17*y^2+2*x+11*y-12; sols := QuadraticIsolve(P); nops(%)

14*x^2-5*x*y-17*y^2+2*x+11*y-12

QuadraticIsolve: disc = 977, transformed equation: X^2-977*Y^2 = -559328

QuadraticIsolve: hyperbola, solved by transformation

QuadraticIsolve: returning general solutions containing parameter _Z1

QuadraticIsolve: solutions may be non-integer for some values of the parameter

128

Look at one of the general solutions

> 

sols[92]

{x = (151471/1954)*977^(1/2)*((108832847723078562849+3481871275306470280*977^(1/2))^_Z1-(108832847723078562849-3481871275306470280*977^(1/2))^_Z1)+(4734529/1954)*(108832847723078562849+3481871275306470280*977^(1/2))^_Z1+(4734529/1954)*(108832847723078562849-3481871275306470280*977^(1/2))^_Z1-13/977, y = 318/977-(2524409/977)*(108832847723078562849+3481871275306470280*977^(1/2))^_Z1-(2524409/977)*(108832847723078562849-3481871275306470280*977^(1/2))^_Z1-(80763/977)*977^(1/2)*((108832847723078562849+3481871275306470280*977^(1/2))^_Z1-(108832847723078562849-3481871275306470280*977^(1/2))^_Z1)}

Choose a specific _Z1. (Some solutions may not be integer)

> 

evala(eval(sols[92], _Z1 = 1))

{x = 1054805055873880264681484, y = -1124825473059920639559012}

Specify chosen _Z1 values as the second argument (set or range); here 0 and 1. Only integer solutions are returned.

> 

sols := QuadraticIsolve(P, {0, 1}); nops(%)

QuadraticIsolve: disc = 977, transformed equation: X^2-977*Y^2 = -559328

QuadraticIsolve: hyperbola, solved by transformation

36

> 

sols

{{x = -4354274635632925848630, y = 4643321530916357188489}, {x = -842369127639866026312480, y = -650531927181484066810521}, {x = -86052394265244921080883408, y = -66455224963899021648732681}, {x = -1834717862750735136472687617, y = -1416887808357334226606681712}, {x = -1935541499837998630, y = -1494750342474165031}, {x = -10004971671472080, y = 10669125001482999}, {x = -97938928148608, y = 104440342385799}, {x = -4593556006657, y = 4898486956848}, {x = -710258662, y = -548507680}, {x = -225037758, y = 239976289}, {x = -10554777, y = 11255428}, {x = -103321, y = 110180}, {x = -4846, y = 5168}, {x = -1632, y = -1260}, {x = -16, y = -12}, {x = -1, y = 0}, {x = 20, y = -21}, {x = 3870, y = 2989}, {x = 395342, y = 305309}, {x = 8429063, y = 6509468}, {x = 54514467908, y = -58133265300}, {x = 172057228892, y = 132873721299}, {x = 3668421709343, y = 2832992527848}, {x = 374748386313759, y = 289404943697160}, {x = 7989987545117724, y = 6170385197335548}, {x = 23725301048449838, y = -25300241809328680}, {x = 2423663086517669214, y = -2584551489063081992}, {x = 51674773213364371719, y = -55105065055676535380}, {x = 163094654023171420679, y = 125952241260259919180}, {x = 3477331195842375008654, y = 2685420073047531974168}, {x = 355227494974617759355998, y = 274329648738604371091960}, {x = 1054805055873880264681484, y = -1124825473059920639559012}, {x = 22489434422549442376595999, y = -23982335477370792469719000}, {x = 2297412873087832090660316703, y = -2449920492316976033389932792}, {x = 48983000093817740901123542892, y = -52234605764925810264202912221}, {x = 154598945607990784875891494868, y = 119391305695605386194795009740}}

Check

> 

map2(eval, P, sols)

{0}

Second case

> 

QuadraticIsolve(3*x^2-5*x*y+7*y^2+2*x-11*y+19)

QuadraticIsolve: disc = -59, (a+c)*delta = 2100

QuadraticIsolve: imaginary ellipse; has no solutions

{}

 

``

Download DiophantineConics3.mw

syntax fixed...

dharr 8637
November 04 2025
0 0

@Andiguys Data needs to be in a set, not list if you use union. 

new_Q_solve_sand15.mw

remark...

dharr 8637
November 04 2025
0 0

@sursumCorda I took this idea from a much longer algorithm for finding when one field is a subfield of another, so I think it is much harder to properly distinguish the two cases. Glad you like it.

Using IntegerRelations:-LinearDependency...

dharr 8637
November 03 2025
0 0

@sursumCorda Here's an algorithm using IntegerRelations:-LinearDependency, based on the LLL algorithm. It uses floats as an intermediate step but checks the obtained result is exact. It does the original examples, but not their reverse (expressing degree 6 roots in terms of degree 2 or 3 roots).

> 

restart;
> 

local gamma:
> 

algconvert := proc(a::algnumext, b::algnumext)
description "writes algebraic number a as a polynomial in algebraic number b";
local Bs, i, Banum, cfs;
Bs := seq(b^i, i = 0..degree(op(1,b), ':-_Z')-1);
Banum:=evalf([Bs,a]);
cfs := IntegerRelations:-LinearDependency(Banum);
if evala(add(cfs*~[Bs,a])) <> 0 then return FAIL end if;
-evala(add(cfs[1..-2]*~[Bs]))/cfs[-1];
end proc:

Example 1.

> 

a := RootOf(4*_Z^3 - 12*_Z^2 - 24*_Z - 25, index=1);
b := RootOf(2*_Z^3 - 3, index=1);
alias(alpha = a, beta = b);

RootOf(4*_Z^3-12*_Z^2-24*_Z-25, index = 1)

RootOf(2*_Z^3-3, index = 1)

alpha, beta

> 

algconvert(a, b);
evala(a - %);
algconvert(b, a);
evala(b - %);

beta^2+2*beta+1

0

-(2/13)*alpha^2+(12/13)*alpha+2/13

0

Example 2

> 

c := RootOf(_Z^4 + 1, index=1);
d := RootOf(_Z^4 + 6*_Z^2 + 1, index=1);
alias(gamma = c, delta = d);

RootOf(_Z^4+1, index = 1)

RootOf(_Z^4+6*_Z^2+1, index = 1)

alpha, beta, gamma, delta

> 

algconvert(c, d);
evala(c - %);
algconvert(d, c);
evala(d - %);

(1/4)*delta^3+(1/4)*delta^2+(7/4)*delta+3/4

0

gamma^3-gamma^2+gamma

0

Example 3

> 

a1 := convert(3^(1/3), RootOf);
a2 := convert(4^(1/4), RootOf);
3^(1/3) + 2^(1/2) + 1^(1/1);
B := RootOf(evala(Minpoly(%, x)), index=2);
evala(%% - B);
alias(alpha__1 = a1, alpha__2 = a2, beta = B);
 

RootOf(_Z^3-3, index = 1)

RootOf(_Z^2-2, index = 1)

3^(1/3)+2^(1/2)+1

RootOf(_Z^6-6*_Z^5+9*_Z^4-2*_Z^3+9*_Z^2-60*_Z+50, index = 2)

0

alpha, beta, gamma, delta, alpha__1, alpha__2

> 

algconvert(a1, B);

FAIL

Succeeds with higher digits

> 

Digits := 50;

50

> 

algconvert(a1, B);
evala(a1 - %);
algconvert(a2, B);
evala(a2 - %);

277/151-(232/755)*beta-(174/755)*beta^2-(52/755)*beta^3+(213/755)*beta^4-(48/755)*beta^5

0

-428/151+(987/755)*beta+(174/755)*beta^2+(52/755)*beta^3-(213/755)*beta^4+(48/755)*beta^5

0

These fail

> 

algconvert(B, a1);
algconvert(B, a2);

FAIL

FAIL

> 

NULL

Download evala.mw

remarks...

dharr 8637
November 03 2025
0 0

@sand15 That's a nice analysis. Both fits basically treat the first two points the same, so really many functions can fit equally "well" through two points. Condition numbers this high just confirm meaningless data. As you probably deduced, a+b*nu_c*nu^2 might make some sense since log(1+x) is approx x for small x = A*c*nu^2.

But really, better data more spread out is needed before one spends more time.

plot shape...

dharr 8637
November 03 2025
0 0

@MichaelVio Here's some plots to emphasize @sand15's point that your data do not fit to a parabola shape centered at positive nu, if that is what your are expecting.

MLfit.mw

nice...

dharr 8637
November 02 2025
0 0

@sand15 Nice analysis. The wildly different orders of magnitudes of the various quantities makes me nervous, aside from the likely poorly posed problem.

isolve...

dharr 8637
November 02 2025

@Alfred_F @janhardo Maple's isolve finds the general formula for the Pell's equation, so no extra procedure is required. (It expresses it in terms of an arbitrary integer parameter _Z1 rather than the equivalent recursion.) My point to @janhardo was that his procedure testing if integers on a grid satisfied this equation is not a realistic way since very large numbers can be needed even for equations with relatively small coefficients.

large numbers...

dharr 8637
November 02 2025

@janhardo The test-on-a-grid method fails fairly quickly - try to find a solution for x^2 - 151*y^2 = 1 aside from (1,0) and (-1,0).

example...

dharr 8637
October 28 2025
0 0

It's not clear what you want to do "on line". If you just want help with your integration, please upload a worksheet showing the example and what you mean by unstable. Use the green up-arrow, choose your file, click upload and then one of insert link or insert contents.

another workaround...

dharr 8637
October 26 2025
0 0

@GunnerMunk If you put simplify in the definition of the velocity function, then it works.

This is a bug and I have submitted an SCR ("bug report") to Maple.

 

restart  

Opgave b

 

Udbredelsesfarten er givet ved

 

"v(d):=simplify(sqrt(g*d)):"

 

med

 

g := 9.82*Unit('m'/'s'^2)

 

Vanddybden er 2,0 km. Dermed er farten af tsunamien

 

v__s := v(2.0*Unit('km'))"(=)"140.1427843*Units:-Unit(m/s)

 

Bølgelængden er

lambda := 120*Unit('km')

 

Dermed er perioden``

T := lambda/v__s

856.2695582*Units:-Unit(s)

T := lambda/v__s"(=)"856.2695582*Units:-Unit(s)

 

Dermed er tsunamiens fart v = 141 m/s og perioden

 

Download Error_2.mw

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