Advanced Search

tomleslie

13876 Reputation

20 Badges

15 years, 163 days
Contact tomleslie

MaplePrimes Activity


These are answers submitted by tomleslie

Can't replicate...

tomleslie 13876
April 27 2022
1 0

wrapping your final execution group with

t1:=time[real]():
.
.
time[real]()-t1

I get

8.462e-008 in 0.205 secs in Maple 2015
8.462e-008 in 0.253 secs in Maple 2020
9.493e-008 in  0.197 secs in Maple 2022

So, no significant speed variation: slightly interesting that Maple 2015 and Maple 2020 produce the same answer, whereas Maple 2022 is *different*? .

All of the above performed on a machine running Windows 7.

The" brain-dead" way...

tomleslie 13876
April 26 2022
1 0

two equations in two unknowns - maybe a simple solve() command would work?

See the attached

> 

  restart;
#
# Define equations
#
  eq1:= m= 10^100 + 598*10^50 + 67497:
  eq2:= phi__m=10^100 + 596*10^50 + 66900:
#
# Obviously 'p' and 'q' are interchangeable in these
# equations - so we will get two solutions, either
#
# p=x, q=y, or
# p-y, q=x
#
  sol:= solve( [ eval( eq1, m=p*q),
                 eval( eq2, phi__m=(p-1)*(q-1))
               ],
               [p, q]
             );
#
# Check that these solutions are prime (we only need to check
# one of them)
#
  isprime~(rhs~(sol[1]));

[[p = 100000000000000000000000000000000000000000000000447, q = 100000000000000000000000000000000000000000000000151], [p = 100000000000000000000000000000000000000000000000151, q = 100000000000000000000000000000000000000000000000447]]

  

[true, true]

 (1)
> 

 

Download pq.mw

 

I agree with acer...

tomleslie 13876
April 25 2022
1 0

The attached worksheet shows that " big-O" expressions are returned for your examples.

Which begs the obvious question - what have you fiddled with (that you really ought not to fiddle with, because you don't understand the consequences)?

> 

  restart;
  interface(version);
  series(x/(1-x-x^2), x=0);
  taylor(sin(x),x=Pi);

`Standard Worksheet Interface, Maple 2022.0, Windows 7, March 8 2022 Build ID 1599809`

  

series(x+x^2+2*x^3+3*x^4+5*x^5+O(x^6),x,6)

  

series(-(x-Pi)+(1/6)*(x-Pi)^3-(1/120)*(x-Pi)^5+O((x-Pi)^7),x = Pi,7)

 (1)
> 

 

Download series.mw

 

The "short" way...

tomleslie 13876
April 25 2022
1 1

all you really need to do is apply the function

f:= L-> is(L=sort(L)):

as in the attached,

> 

  restart;
#
# Some lists for test purposes
#
  L__1:= [1, 1, 1, 1, 1]:  
  L__2:= [1,2,3,4,5,4,3,2,1]:
  L__3:= [1, 2, 2, 3, 4, 5, 5, 5]:
  L__4:= [5, 5, 4, 4, 3, 3, 1]:
  L__5:= [1, 0, 1, 0]:
  L__6:= [1, 1, 2, 2, 3, 4, 5, 2, 2, 6, 4, 2, 2, 1, 0]:
#
# Define the function whhc does the "work"
#
  f:= L-> is(L=sort(L)):
#
# Check unimodality on individual lists
#
  f(L__1);
  f(L__2);
  f(L__3);
  f(L__4);
  f(L__5);
  f(L__6);
#
# Check all of the above lists individually
#
  f~([L__1, L__2, L__3, L__4, L__5, L__6]);
#
# Or check them all the "smart" way
#
  f~([seq( eval(cat(L__, j)), j=1..6)]);

true

  

false

  

true

  

false

  

false

  

false

  

[true, false, true, false, false, false]

  

[true, false, true, false, false, false]

 (1)
> 

 

> 

 

> 

 

Download sortList.mw

 

A workaround(?)...

tomleslie 13876
April 25 2022
1 0

and, (before you ask) I have no idea why this appears to "work". My thought process was just

  1. Probably best to avoid using commands which are deprecated (in this case '||' )
  2. Probably best to avoid using indexed names (ie c[k]) for things entities which are not *really* indexable

See the attached,

> 

  restart;
  interface(version);
  alias(seq(cat(c__, k) = cat(_C, k), k = 1..10));
  ode:=x*diff(y(x),x$2)-cos(x)*diff(y(x),x)+sin(x)*y(x)=2;

  sol:=dsolve(ode);


  lprint(sol);
  latex(sol)

`Standard Worksheet Interface, Maple 2022.0, Windows 7, March 8 2022 Build ID 1599809`

  

c__1, c__2, c__3, c__4, c__5, c__6, c__7, c__8, c__9, c__10

  

x*(diff(diff(y(x), x), x))-cos(x)*(diff(y(x), x))+sin(x)*y(x) = 2

  

y(x) = (c__2+Int((c__1+2*x)/(exp(Ci(x))*x^2), x))*exp(Ci(x))*x

  

y(x) = (c__2+Int((c__1+2*x)/exp(Ci(x))/x^2,x))*exp(Ci(x))*x
y \! \left(x \right) =
\left(c_{2} +\textcolor{gray}{\int}\frac{c_{1} +2 x}{{\mathrm e}^{\mathrm{Ci}\left(x \right)} x^{2}}\textcolor{gray}{d}x \right) {\mathrm e}^{\mathrm{Ci}\left(x \right)} x

  
> 

 

Download oddDisp.mw

The code...

tomleslie 13876
April 23 2022
1 1 
  restart;
  n:=3:
  L:=[ x$n, y$n, z$n ]:
  [seq( combinat:-choose(L, j)[], j=1..n)];

will produce

[[x], [y], [z], [x, x], [x, y], [x, z], [y, y], [y, z], [z, z], [x, x, x], [x, x, y], [x, x, z], [x, y, y], [x, y, z], [x, z, z], [y, y, y], [y, y, z], [y, z, z], [z, z, z]]

 

Much easier with matrices...

tomleslie 13876
April 22 2022
2 3

Rather than use a lot of "for loops", it is generally much simpler to organise the calculation in matrix form as shown in the attached.

I have chosen to use 21 points in the time range t=0..1. You can adjust this just by changing the value of 'npts': with npts=51, the agreement between the rk2 method and Maple's default IVP solver (rkf45) gets significantly better.

> 

  restart;
  with(LinearAlgebra):
  with(plots):
#
# Define equations and ICS
#
  sys := [ diff(A[2](t), t) = -4/45*(-1323*A[7](t)*Pi^5 - 945*A[6](t)*Pi^4 - 630*A[5](t)*Pi^3 - 513*A[4](t)*Pi^2 - 189*A[3](t)*Pi + 3626*A[2](t))/Pi^2,
           diff(A[3](t), t) = 4/45*(-9744*A[7](t)*Pi^5 - 6960*A[6](t)*Pi^4 - 4415*A[5](t)*Pi^3 - 2784*A[4](t)*Pi^2 - 1392*A[3](t)*Pi + 14544*A[2](t))/Pi^3,
           diff(A[4](t), t) = -2/9*(-14112*A[7](t)*Pi^5 - 10215*A[6](t)*Pi^4 - 6720*A[5](t)*Pi^3 - 4032*A[4](t)*Pi^2 - 2016*A[3](t)*Pi + 12992*A[2](t))/Pi^4,
           diff(A[5](t), t) = 2/45*(-133455*A[7](t)*Pi^5 - 96000*A[6](t)*Pi^4 - 64000*A[5](t)*Pi^3 - 38400*A[4](t)*Pi^2 - 19200*A[3](t)*Pi + 81408*A[2](t))/Pi^5,
           diff(A[6](t), t) = -4096/45*(-63*A[7](t)*Pi^5 - 45*A[6](t)*Pi^4 - 30*A[5](t)*Pi^3 - 18*A[4](t)*Pi^2 - 9*A[3](t)*Pi + 26*A[2](t))/Pi^6,
           diff(A[7](t), t) = 32768/315*(-21*A[7](t)*Pi^5 - 15*A[6](t)*Pi^4 - 10*A[5](t)*Pi^3 - 6*A[4](t)*Pi^2 - 3*A[3](t)*Pi + 6*A[2](t))/Pi^7]:
  ICS := [A[2](0) = -0.001210685373, A[3](0) = -0.1636380032, A[4](0) = -0.003838287851, A[5](0) = 0.01100619795, A[6](0) = -0.001005637627, A[7](0) = -0.00002775982849]:
> 

#
# Convert to Matrix form and run the RK2 method for npts
# in the interval t=0..1
#
#  Rows in the results matrix are
#
#  resMat[1,..] <--> A[2](t)
#  resMat[2,..] <--> A[3](t)
#  resMat[3,..] <--> A[4](t)
#  resMat[4,..] <--> A[5](t)
#  resMat[5,..] <--> A[6](t)
#  resMat[6,..] <--> A[7](t)
#
  f, diffs:= GenerateMatrix
             ( rhs~(sys)-~lhs~(sys),
               [A[2](t), A[3](t), A[4](t), A[5](t), A[6](t), A[7](t)]
             ):
  npts:= 21:
  resMat:= Matrix(6,npts):
  resMat[1..6, 1]:= Vector(rhs~(ICS)):
  h:=1/(npts-1):
#
# The RK2 method
#
  for j from 2 to npts do
      k1:= f.resMat[..,j-1];
      k2:= f.(resMat[..,j-1]+h*k1);
      resMat[..,j]:= resMat[..,j-1]+h/2*(k1+k2);
  od:
> 

#
# Get the solution for the system using Maple's default solvers,
# just for comparison purposes
#
  sol:=dsolve([sys[], ICS[]], numeric):
> 

#
# Plot the solutions from the RK2 method above and Maple's
# default solvers for a couple of the independent variables
#
# I randoml picked A[2](t) and A[7](t).
#
# Solutions from Maple's default solvers are solid red lines
# Solutions from "crude" RK2 method above
#
  display( [ odeplot(sol, [t, A[2](t)], t=0..1),
             pointplot([seq( [h*(j-1),  resMat[1,j] ], j=1..npts)], symbol=solidbox)
           ],
           color=[ red, blue],
           legend=[ "MapleDefault", "RK2"],
           title= "Solutions for A[2](t)",
           titlefont=[times, bold, 16]
         );
  display( [ odeplot(sol, [t, A[7](t)], t=0..1),
             pointplot([seq( [h*(j-1),  resMat[6,j] ], j=1..npts)], symbol=solidbox)
           ],
           color=[ red, blue ],
           legend=[ "MapleDefault", "RK2"],
           title= "Solutions for A[7](t)",
           titlefont=[times, bold, 16]
         );

 

 
> 

 

> 

 

Download rk2.mw

 

Maybe...

tomleslie 13876
April 21 2022
0 0

as shown in the attached

> 

  restart;

  U[0](x):= lambda;
  n := 7;

  for k from 0 to n do
      U[k+1](x) := solve((k+1)*U[k+1](x)-(diff(U[k](x), x, x))+2*U[k](x)+3*(add(U[r](x)*U[k-r](x), r = 0 .. k)),U[k+1](x) )
  end do;

lambda

  

7

  

-3*lambda^2-2*lambda

  

9*lambda^3+9*lambda^2+2*lambda

  

-27*lambda^4-36*lambda^3-14*lambda^2-(4/3)*lambda

  

135*lambda^4+75*lambda^3+15*lambda^2+(2/3)*lambda+81*lambda^5

  

-351*lambda^4-108*lambda^3-(62/5)*lambda^2-(4/15)*lambda-486*lambda^5-243*lambda^6

  

630*lambda^4+(602/5)*lambda^3+(42/5)*lambda^2+(4/45)*lambda+1512*lambda^5+1701*lambda^6+729*lambda^7

  

-(8/315)*lambda-(508/105)*lambda^2-2187*lambda^8-5832*lambda^7-6156*lambda^6-3240*lambda^5-(4374/5)*lambda^4-(552/5)*lambda^3

  

(2/315)*lambda+(17/7)*lambda^2+6561*lambda^9+19683*lambda^8+24057*lambda^7+15309*lambda^6+(26811/5)*lambda^5+999*lambda^4+(605/7)*lambda^3

 (1)
> 

 


 

Download serSol2.mw

Difficult...

tomleslie 13876
April 21 2022
2 0

to determine what you want from your attempt at "cut-and-paste". in future please use the big green -up-arrow in the Mapleprimes toolbar to upload an actual worksheet. That way responders will not have to "guess" what you are after. My best guess(?!?) is shown in the attached

> 

  restart;
  U[0](x):=exp(x);
  n:=7:
  for k from 0 to n do
      U[k+1](x):=solve((k+1)*U[k+1](x)-diff(U[k](x),x$2)-2*U[k](x), U[k+1](x));
  end do;

exp(x)

  

3*exp(x)

  

(9/2)*exp(x)

  

(9/2)*exp(x)

  

(27/8)*exp(x)

  

(81/40)*exp(x)

  

(81/80)*exp(x)

  

(243/560)*exp(x)

  

(729/4480)*exp(x)

 (1)
> 

 

Download serSol.mw

See the attached...

tomleslie 13876
April 21 2022
1 1

II guessed a value for 'delaymax' - Note also that the form of the equations/ics means that x(t)=y(t), is that intentional?

> 

restart;

k[1] := 1; k[2] := 1; p := 1/2; q := 1/2; tau[1] := 1/2; tau[2] := 1/2;

eq1 := diff(x(t), t) = x(k[1]-x(t))-p*x(t-tau[1])*y(t);

eq2 := diff(y(t), t) = y(k[2]-y(t))-q*y(t-tau[2])*x(t);

ics := {x(0) = 2/3, y(0) = 2/3};

sys := {eq1, eq2};
sol:=dsolve(`union`(ics, sys), numeric, delaymax=1):
plots:-odeplot( sol, [t,y(t)], t=0..10);
plots:-odeplot( sol, [t,y(t)], t=0..10);

1

  

1

  

1/2

  

1/2

  

1/2

  

1/2

  

diff(x(t), t) = x(1-x(t))-(1/2)*x(t-1/2)*y(t)

  

diff(y(t), t) = y(1-y(t))-(1/2)*y(t-1/2)*x(t)

  

{x(0) = 2/3, y(0) = 2/3}

  

{diff(x(t), t) = x(1-x(t))-(1/2)*x(t-1/2)*y(t), diff(y(t), t) = y(1-y(t))-(1/2)*y(t-1/2)*x(t)}

  

 

 
> 

 


 

Download deldiff.mw

How I'd...

tomleslie 13876
April 20 2022
4 1

do it is shown in the attached, The precise details will vary a bit depending on whether all of your strings are entries a single array (which I assume in the attached), or you actually have a 1*1 Array for each string.

> 

  restart:
#
# Both the Units() package and the StringTools() package
# have a "Split()" command, so the following code use
# longform commands for the StringTools{} package, with
# an alias (just to save typing)
#
  alias(STR=StringTools):
  with(Units):
#
# Test Array - I only changed the last term in the three
# entries because I wanted to check that 'units' were being
# handled correctly
#
  arr:=Array(["0,00244140625;8,07751097960625;-3,555908203125E-05;0,002166748046875;1;1;0,0056941700604248;-0,0101856454842773;100 µA",
              "0,00244140625;8,07751097960625;-3,555908203125E-05;0,002166748046875;1;1;0,0056941700604248;-0,0101856454842773;200 mA",
              "0,00244140625;8,07751097960625;-3,555908203125E-05;0,002166748046875;1;1;0,0056941700604248;-0,0101856454842773;2 A"  ]);
  ArrayDims(arr);

Automatically loading the Units[Simple] subpackage
 

  

Vector[row](3, {(1) = "0,00244140625;8,07751097960625;-3,555908203125E-05;0,002166748046875;1;1;0,0056941700604248;-0,0101856454842773;100 µA", (2) = "0,00244140625;8,07751097960625;-3,555908203125E-05;0,002166748046875;1;1;0,0056941700604248;-0,0101856454842773;200 mA", (3) = "0,00244140625;8,07751097960625;-3,555908203125E-05;0,002166748046875;1;1;0,0056941700604248;-0,0101856454842773;2 A"})

  

1 .. 3

 (1)
> 

#
# Define a function which performs all necessary substitutions
# Note that these are done in the order specified by the list
# supplied to STR:-Subs()command, and this order is significant
#
  f:= z-> parse~
          ( STR:-Split
            ( STR:-Subs
              ( [ ","=".",
                  ";"=", ",
                  " µA"="*Unit(uA)",
                  " mA"="*Unit(mA)",
                  " A"="*Unit(A)"  
                ],
                z
              ),
              ","
            )
          ):
#
# Apply the above function elementwise to the test array
# and build the output array
#
  ans:=Array(convert(f~(arr), listlist));
  ArrayDims(ans);

Matrix(3, 9, {(1, 1) = 0.244140625e-2, (1, 2) = 8.07751097960625, (1, 3) = -0.3555908203125e-4, (1, 4) = 0.2166748046875e-2, (1, 5) = 1, (1, 6) = 1, (1, 7) = 0.56941700604248e-2, (1, 8) = -0.101856454842773e-1, (1, 9) = Units:-Simple:-`*`(100, Units:-Unit(uA)), (2, 1) = 0.244140625e-2, (2, 2) = 8.07751097960625, (2, 3) = -0.3555908203125e-4, (2, 4) = 0.2166748046875e-2, (2, 5) = 1, (2, 6) = 1, (2, 7) = 0.56941700604248e-2, (2, 8) = -0.101856454842773e-1, (2, 9) = Units:-Simple:-`*`(200, Units:-Unit(mA)), (3, 1) = 0.244140625e-2, (3, 2) = 8.07751097960625, (3, 3) = -0.3555908203125e-4, (3, 4) = 0.2166748046875e-2, (3, 5) = 1, (3, 6) = 1, (3, 7) = 0.56941700604248e-2, (3, 8) = -0.101856454842773e-1, (3, 9) = Units:-Simple:-`*`(2, Units:-Unit(A))})

  

1 .. 3, 1 .. 9

 (2)
> 

 

> 

 

Download arrStr.mw

You cannot...

tomleslie 13876
April 19 2022
1 0

use the Secant() command to solve systems of equations. It is intended for finding roots of a single, univariate expression.

However you have a number of other possibilities, see the attached. (The last couple of commands in the attached probably won't work for you, since the DirectSearch() package is an "add-on" , available from the Maple Applicarion Center)

> 

  restart;
> 

  f1 := x^2 - x*y + 2;
  f2 := x^2*y^2 - y - 43;
  F := [f1, f2];
  init_pt := {x = 2.0, y = 3.7};

  fsolve(F, init_pt);
  Optimization:-NLPSolve(f1^2 + f2^2, initialpoint = init_pt);
#
# OP (probably) won't have the DirectSearch() package
# installed, so the following is just for interest
#
  DirectSearch:-SolveEquations(F, initialpoint = init_pt);
  DirectSearch:-SolveEquations(F, AllSolutions=true);

f1 := x^2-x*y+2

  

f2 := x^2*y^2-y-43

  

F := [x^2-x*y+2, x^2*y^2-y-43]

  

init_pt := {x = 2.0, y = 3.7}

  

{x = 2.188578493, y = 3.102413663}

  

[2.11141742846706446*10^(-13), [x = 2.18857856657139, y = 3.10241354555813]]

  

[4.395192586061875*10^(-19), Vector(2, {(1) = -0.4384581587e-9, (2) = 0.4972662282e-9}), [x = 2.1885784931152497, y = 3.102413663628889], 104]

  

Matrix(%id = 36893488148118587860)

 (1)
> 

``

NULL

Download roots.mw

Running Syrup...

tomleslie 13876
April 14 2022
1 7

The attached shows Syrup running in Maple 2022.

(I'm seriously not impressed by the help page you reference)

> 

  restart;
  interface(version);
  with(Syrup);
  divider := "
              V 1 0
              R1 1 2
              R2 2 0
             .end":
  Solve(divider, 'returnall');

`Standard Worksheet Interface, Maple 2022.0, Windows 7, March 8 2022 Build ID 1599809`

  

[Draw, Library, Print, Solve, ToModelica, ToSpice, Version]

  

{v[1] = V, v[2] = V*R2/(R1+R2)}, {i[R1] = V/(R1+R2), i[R2] = V/(R1+R2), i[V] = -V/(R1+R2), v[R1] = V*R1/(R1+R2), v[R2] = V*R2/(R1+R2), v[V] = V}

 (1)
> 

  bridge := "
            V a 0
            R1 a b
            R2 a c
            R3 b 0
            R4 c 0
            R5 b c
         .end":
  Solve(bridge, 'returnall');

{v[a] = V, v[b] = V*(R1*R4+R2*R4+R2*R5+R4*R5)*R3/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5), v[c] = V*(R1*R3+R1*R5+R2*R3+R3*R5)*R4/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5)}, {i[R1] = V*(R2*R3+R2*R4+R2*R5+R4*R5)/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5), i[R2] = V*(R1*R3+R1*R4+R1*R5+R3*R5)/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5), i[R3] = V*(R1*R4+R2*R4+R2*R5+R4*R5)/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5), i[R4] = V*(R1*R3+R1*R5+R2*R3+R3*R5)/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5), i[R5] = -V*(R1*R4-R2*R3)/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5), i[V] = -V*(R1*R3+R1*R4+R1*R5+R2*R3+R2*R4+R2*R5+R3*R5+R4*R5)/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5), v[R1] = V*R1*(R2*R3+R2*R4+R2*R5+R4*R5)/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5), v[R2] = V*R2*(R1*R3+R1*R4+R1*R5+R3*R5)/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5), v[R3] = V*(R1*R4+R2*R4+R2*R5+R4*R5)*R3/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5), v[R4] = V*(R1*R3+R1*R5+R2*R3+R3*R5)*R4/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5), v[R5] = -V*(R1*R4-R2*R3)*R5/(R1*R2*R3+R1*R2*R4+R1*R2*R5+R1*R3*R4+R1*R4*R5+R2*R3*R4+R2*R3*R5+R3*R4*R5), v[V] = V}

 (2)
> 

 

Download syrup.mw

With...

tomleslie 13876
April 13 2022
1 1

a couple of fxies in the 'pellsolve' proc (which I guessed), the attached seems to function


 

> 

pellsolve := proc(D::posint)
                  local P, Q, a, A, B, i;
                  if   type(sqrt(D), integer)
                  then error "D must be a nonsquare integer";
                  end if;
                  P := 0;
                  Q := 1;
                  a := floor(sqrt(D));
                  A := 1, a;
                  B := 0, 1;
                  for i do
                      P := a*Q - P;
                      Q := (D - P^2)/Q;
                      a := floor((P + sqrt(D))/Q);
                      A := A[2], a*A[2] + A[1];
                      B := B[2], a*B[2] + B[1];
                      if Q = 1 and i mod 2 = 0
                      then break;
                      fi:
                  end do;
                  return A[1], B[1];
              end proc:
genpellsolve := proc(D::posint, N::integer)
                     local t, u, L1, L2, sols, x, y;
                     if   type(sqrt(D), integer)
                     then error "D must be a nonsquare integer";
                     end if:
                     t, u := pellsolve(D);
                     if   0 < N
                     then L1 := 0;
                          L2 := floor(sqrt(1/2*N*(t - 1)/D));
                     elif N < 0
                     then L1 := ceil(sqrt(-N/D));
                          L2 := floor(sqrt(-1/2*N*(t + 1)/D));
                     else return {[0, 0]};
                     end if;
                     sols := {};
                     for y from L1 to L2 do
                         x := sqrt(N + D*y^2);
                         if   type(x, integer)
                         then sols := sols union {[x, y]};
                               if   (x^2 + D*y^2) mod N <> 0 or (2*x*y) mod N <> 0
                               then sols := sols union {[-x, y]};
                               end if;
                         end if;
                     end do;
                     return sols;
                 end proc:
> 

seq( genpellsolve(10, j), j=1..10);

{[1, 0]}, {}, {}, {[2, 0]}, {}, {[-4, 1], [4, 1]}, {}, {}, {[-7, 2], [3, 0], [7, 2]}, {[10, 3]}

 (1)
> 

 

Download psolve.mw

Some observations...

tomleslie 13876
April 11 2022
2 0
  1. I can reproduce the difference in behaviour which you are seeing - so it's not you, something has changed between Maple 2021 and Maple 2022.
  2. The option 'legend" in any Maple plotting command is "normally" used when one wants to distinguish between multiple curves/datasets occurring in the same graph. When only a single curve/dataset occurs in a graph, then the "normal" approach for annotation would be to use either the "title" or "caption" options. So using 'legend' with a single curve/dataset is "unusual". The trouble with using "unusual" code is that when software is updated, "unusual" code is more likely to break
  3. The Histogram() command has two cateogries of options;those which are applicable to the generation of the histogram itself, and those which are applicable to the display of the the histogram. According to the help, the latter (including 'legend') are just passed to the plots:-display() command. However both the Histogram() command and the display() command were updated betwee Maple 2021 and Maple 2022 - so somewhere in these updates your "unusual" code broke
  4. The fact that the option 'legend' would normally be asociated with multiple curves/datasets is a reasonable indication of why a list is expected
  5. Can only suggest that rather than using the option 'legend' you go with either 'title' or 'caption'
First 27 28 29 30 31 32 33 Last Page 29 of 207