tomleslie

13876 Reputation

20 Badges

15 years, 163 days

MaplePrimes Activity


These are answers submitted by tomleslie

which part of my original response did you not understabd?

If you are incapable of rational thought and need Maple to demonstrate it for you, see the attached

  restart;
  w:=1/(4*(exp(1/2)+exp(-1/2))+2*(exp(1/2)+exp(-1/2))^2):
  g1:=x->sum(a[i]*x^i,i=0..50):
  for j from 0 to 50 do
      r[j]:=(1/2+w)*(g1(j+1)+g1(j))-w*(g1(j-1)+g1(j+2))-g1(2*j+1):
  end do:

#
# Performming the solution which the OP does not want to believe
#
# Step 1 - convert the equations to a list
#
  equns:=convert(r, list):
#
# Step 2 - generate s set of all of the unknowns
#
  vars:=`union`(indets~(equns)[]):
#
# Step 3 - Genrerate the augmented matrix from the equations, and
# convert all coefficients to floating point numbers. The latter
# conversion is a matter of convenience whihc greatly speeds up
# all subsequent calculations
#
  M:=evalf~(LinearAlgebra:-GenerateMatrix( equns, vars, augmented=true)):
#
# Solve the equations implied by the above augmented matrix
#
  oldRT:=interface(rtablesize=50):
  LinearAlgebra:-LinearSolve(M);
  interface(rtablesize=oldRT):

_rtable[36893488148090695540]

(1)

 


 

Download linSolve.mw

there are a number of issues here

  1. Produce an "efficient" Maple code. The version provided by the OP "works", and is fairlyy "efficient" although it is hard-wired to only work for complete sub-graphs of order 4. The attached shows a simpler(?) , aka "more understandable" version which has the benefit of working for any clique size (provided as a parameter). In my test tests the two versions run in pretty much the same time.
  2. Can Maple "beat" optimised, compiled code on the same problem? No! But, if this is a complaint, you have to consider what you are asking for. Should a general-purpose program such as Maple be able to "beat" all programs optimised for a specific application area - say quantum mechanics, general relativity, financial option pricing, optimizattion, etc, et?.This is an unrealistic expectation on several levels
    1. Compiled code will always "beat" interpreted code
    2. A program which is seriously optimised for a specific application, will always "beat" a more "general" program
    3. Is it feasible to write "interfaces" to every specialised math program - get real, there are thousands of them
  3. The possibility of converting Maple code to "acceptable" igraph code is close to nil - free advice, odn't bother trying
  4. igraph is free. If it is the better solution for your problem, then why are you bothering with Maple at all - this make no sense

  restart;

  with(GraphTheory):
  with(SpecialGraphs):
#
# OP's original code - nothing *wrong* with it.
# just a bit verbose and difficult to tollow.
#
# (And hard-wired to find cliques with 4 vertices,
#  which might be consodered a bit "limited")
#
  NumberOfK4:= proc(g::Graph)
                    local  n,N, i, j, k,l, g1,G1;
                    uses GraphTheory;
                    n:= NumberOfVertices(g);
                    g1:= RelabelVertices(g, [$ i = 1..n]):
                    N:=0:
                    for i from 1 to n-3 do
                        for j from i+1 to n-2 do
                            for k from j+1 to n-1 do
                                for l from k+1 to n do
                                    G1:=InducedSubgraph(g1,[i,j,k,l]);
                                    if   NumberOfEdges(G1)=6
                                    then N:=N+1;
                                    fi;
                               od;
                           od;
                       od;
                   od;
                   N;
              end proc:

  g:= CompleteGraph(7):
  g1:= LineGraph(g):
  CodeTools:-Usage(NumberOfK4(g1));
#
# A smarter(?), certainly "shorter" version, which runs in
# about the same time (but also handles any size clique)
#
  fn:= (g, n)-> local j;
                add( `if`( IsClique(g, j), 1, 0),
                     j in combinat:-choose(Vertices(g), n)
                   ):
  CodeTools:-Usage( fn(g1,4) );

memory used=103.97MiB, alloc change=66.00MiB, cpu time=1.79s, real time=1.78s, gc time=78.00ms

 

105

 

memory used=127.23MiB, alloc change=70.00MiB, cpu time=2.12s, real time=2.05s, gc time=156.00ms

 

105

(1)

 

Download cliques.mw

your code so many times, and I still have no idea what you are trying to achieve.

You appear to have no concept of variable "scope" - ie what variablkes with what names are operative at what "level" in code

My best guess about what you might want ,is in the attached, which is "simpler" but of course may be absolute crap

restart; G := proc (p, q) local t, y, f, j; t := Vector[column](4); y := Vector[column](4); f := proc (a, b) options operator, arrow; a+b end proc; t[1] := p; y[1] := q; for j to 3 do y[j+1] := y[j]+(1/6)*f(t[j], y[j])+(1/6)*f(t[j], y[j]+(1/2)*f(t[j], y[j])); t[j+1] := t[1]+(1/10)*j+1/10 end do; return t, y end proc; G(0, .1)

Vector[column](%id = 36893488148113924332), Vector[column](%id = 36893488148113924452)

(1)

NULL

Download opVecs.mw

When I tried to solve this problem yesterday, using a "similar" approach, I kept bumping into the error

Error, (in Optimization:-NLPSolve) abs is not differentiable at 0

The code provided by Carl Love also generates this error??? However

  1. Increasing Digits to 16 for Carl's code makes the "Error" go away
  2. I can no longer persuade my "similar" code to produce this error - not sure what I was doing yesterday

In any case, if the above error does occur, I'm betting that it is less likely to happen for increased values of Digits.

See the attached

  restart:
  interface(displayprecision=5):
#
# Carl's code
#
  ABS:= (a,b,c)-> abs(a-b*c) <= k*c:
  (L,M,U):= seq([cat](x, __, 1..3), x= [l,m,u]):
  Cons:= ABS~(
                map( op@index, [L,M,U], [1,1,2]),
                [.67, 2.5, 1.5, 1, 3, 2, 1.5, 3.5, 2.5],
                map(op@index, [U,M,L], [2,3,3])
              )[],
         add(L) + 4*add(M) + add(U) = 6,
         (L <=~ M)[],
         (M <=~ U)[]:
  Sol:= Optimization:-NLPSolve(k, {Cons}, assume=nonnegative ):
  evalf[5](Sol);

Error, (in Optimization:-NLPSolve) abs is not differentiable at 0

 

Sol

(1)

  restart:
#
# Carls's code with Digits:=16
#
  Digits:=16:
  interface(displayprecision=5):
  ABS:= (a,b,c)-> abs(a-b*c) <= k*c:
  (L,M,U):= seq([cat](x, __, 1..3), x= [l,m,u]):
  Cons:= ABS~(
                map( op@index, [L,M,U], [1,1,2]),
                [.67, 2.5, 1.5, 1, 3, 2, 1.5, 3.5, 2.5],
                map(op@index, [U,M,L], [2,3,3])
              )[],
         add(L) + 4*add(M) + add(U) = 6,
         (L <=~ M)[],
         (M <=~ U)[]:
  Sol:= Optimization:-NLPSolve(k, {Cons}, assume=nonnegative );

[.2360679774997896, [k = .2360679774997896, l__1 = .3905068282873290, l__2 = .2992380565943533, l__3 = .1577202826350433, m__1 = .4645540877822463, m__2 = .3758321518221884, m__3 = .1680772479208869, u__1 = .5148959822871970, u__2 = .4313332881862403, u__3 = .1724516119085421]]

(2)

  restart;
  Digits;
  interface(displayprecision=5):
  with(Optimization):
  constr:= [ l__1-0.67*u__2 <=  k*u__2,
             l__1-0.67*u__2 >= -k*u__2,
             m__1-1*m__2 <=  k*m__2,
             m__1-1*m__2 >= -k*m__2,
             u__1-1.5*l__2 <=  k*l__2,
             u__1-1.5*l__2 >= -k*l__2,
             l__1-2.5*u__3 <=  k*u__3,
             l__1-2.5*u__3 >= -k*u__3,
             m__1-3*m__3 <=  k*m__3,
             m__1-3*m__3 >= -k*m__3,
             u__1-3.5*l__3 <=  k*l__3,
             u__1-3.5*l__3 >= -k*l__3,
             l__2-1.5*u__3 <=  k*u__3,
             l__2-1.5*u__3 >= -k*u__3,
             m__2-2*m__3 <=  k*m__3,
             m__2-2*m__3 >= -k*m__3,
             u__2-2.5*l__3 <=  k*l__3,
             u__2-2.5*l__3 >= -k*l__3,
             l__1/6+4*m__1/6+u__1/6+l__2/6+4*m__2/6+u__2/6+l__3/6+4*m__3/6+u__3/6=1,
             l__1 <= m__1,
             m__1 <= u__1,
             l__2 <= m__2,
             m__2 <= u__2,
             l__3 <= m__3,
             m__3 <= u__3,
             l__1 >= 0,
             l__2 >= 0,
             l__3 >= 0,
             k >= 0
           ]:
  NLPSolve(k, constr);
  

10

 

[.236067977499597043, [k = HFloat(0.23606797749959704), l__1 = HFloat(0.4074736771907326), l__2 = HFloat(0.31246476361176445), l__3 = HFloat(0.1643659979307035), m__1 = HFloat(0.4542964450909135), m__2 = HFloat(0.3675335445626036), m__3 = HFloat(0.16436599793070342), u__1 = HFloat(0.5412104597819644), u__2 = HFloat(0.44971654352795554), u__3 = HFloat(0.17998460761999122)]]

(3)

 

Download fuzzyOpt.mw

  1. The Tools - Tutors - Bascis - Back-Solver seems to be broken in Maple 2021.2
  2. It "appears" functional in earlier versions of Maple, but (to me at least) seems really painful to uses
  3. I cannot imagine anyone prefering this to "Tutor" to a few commands involving fsolve(), eval() and plot
  4. See the attached, which handles one of the examples in the "Tutor"

  restart;
#
# Example equation
#
  expr:=x+sin(y)=y;

x+sin(y) = y

(1)

#
# find 'x' for a given value of y
#
  fsolve(eval(expr, y=0.3));
#
# find y for a given value of x
#
  fsolve(eval(expr, x=0.3));

0.4479793300e-2

 

1.248515468

(2)

#
# Plot computed values of 'x' for a range of y-values
#
  plot( 'fsolve'(eval(expr, y=z)),
        z=0..1,
        labels=[y,x],
        labelfont=[times, bold, 20]
      );  
#
# Plot computed values of 'y' for a range of x-values
#
  plot( 'fsolve'(eval(expr, x=z)),
        z=0..1,
        labels=[x,y],
        labelfont=[times, bold, 20]
      );  

 

 

 

 

 

 


 

Download back.mw

Another depth-first spanning tree procedure for connected graphs - just because I like alternatives

  restart;

  doDFS:= module()
                 local disc, edges, traverse;
                 export DFS;
                 option package;
                 uses GraphTheory;

                 DFS:= proc( grph, rt)
                             disc:= DEQueue(rt);
                             edges:=DEQueue();
                             traverse( grph, rt);
                             return disc, edges;
                       end proc:

                 traverse:= proc(grph, v)
                                 local ngh:= DEQueue(Neighbors(grph, v)),
                                       ch;
                                 while not empty(ngh) do
                                       ch:= pop_front(ngh):
                                       if   not has(disc, ch)
                                       then push_back(disc, ch ):
                                            push_back(edges, {v, ch});
                                            traverse( grph, ch):
                                       fi;
                                 od;
                            end proc;
  end module:
         

  with(GraphTheory):
  with(SpecialGraphs):
  with(doDFS):
  P := PetersenGraph():
  vv, ed:=DFS(P, 1);
  DrawGraph( Graph
             ( convert(vv, list),
               convert(ed, set)
             ),
             stylesheet="legacy"
           );

module DEQueue () local num, head, tail, storage, dsp; option object; end module, module DEQueue () local num, head, tail, storage, dsp; option object; end module

 

 

 

Download DFS.mw

exactly what your problem is, so you are going to have to expalin why somthing like the attached does not fulfil you requirements

  restart;

  interface(version);
#
# define an unevaluated integral
#
  myInt := Int(1/sqrt(x), x);
#
# Evaluate the integral
#
  v:=value(myInt);
#
# Return the original unevaluated integral, and
# its value
#
  myInt;
  v;

`Standard Worksheet Interface, Maple 2021.2, Windows 7, November 23 2021 Build ID 1576349`

 

Int(1/x^(1/2), x)

 

2*x^(1/2)

 

Int(1/x^(1/2), x)

 

2*x^(1/2)

(1)

 

Download expr.mw

 nothing wrong with Kitonum's solution - I just like alternatives


 

restart

with(plots)

inf := 10^(-2)

mu := exp(-lambda*theta(y, t))

eqs := [R*(diff(w(y, t), t)) = mu*(diff(w(y, t), `$`(y, 2)))+(diff(mu, y))*(diff(w(y, t), y))-M*w(y, t)+L+Gr*theta(y, t), R*Pr*(diff(theta(y, t), t)) = diff(theta(y, t), `$`(y, 2))+Br*mu*(diff(w(y, t), y))^2+Ra*((theta(y, t)+1)^4-1)]

parameters := [Gr = 1, Pr = 7.1, R = .5, lambda = .5, Bi = 1, M = .5, Ra = .1, L = 1, Br = .1]

IBC := [w(inf, t) = 0, (D[1](w))(0, t) = 0, w(y, 0) = 0, (D[1](theta))(0, t) = 0, (D[1](theta))(inf, t)+Bi*theta(inf, t) = 0, theta(y, 0) = 0]

pds := pdsolve(eval(eqs, parameters), eval(IBC, parameters), numeric, time = t, range = 0 .. 10^(-2))

pds:-plot3d(theta(y, t), t = 0 .. 1, y = 0 .. inf, style = surface, color = red, title = typeset(theta(y, t), " versus ", y, " and ", t), titlefont = [times, bold, 20])

pds:-plot3d(w(y, t), t = 0 .. 1, y = 0 .. inf, style = surface, color = blue, title = typeset(w(y, t), " versus ", y, " and ", t), titlefont = [times, bold, 20])
NULL

 

 

``

 


 

Download pdesol.mw

you were aiming for the result shown in the attached? (BTW the gridlines in the plot don't appear in the maple worksheet, so it looks "nicer" - but for some reason, uploading toi this site "inserts" gridlines - I have no idea why!)

I would like to have this worksheet displayed "inline" here, unfortunately, I cannot "submit" such a response, because this site is broken - again, so you will have to download the attached to see that it "works"

errbars.mw

is explained by the help for DEplot() which pretty much covers the phaseportrait() command as well. This entry states (with my emphasis added) - read it very carefuly

A system of two first order differential equations produces a direction field plot, provided the system is determined to be autonomous. In addition, a single first order differential equation produces a direction field (as it can always be mapped to a system of two first order autonomous differential equations). A system is determined to be autonomous when all terms and factors, other than the differential, are free of the independent variable. For more information, see DEtools[autonomous]. For systems not meeting these criteria, no direction field is produced (only solution curves are possible in such instances). There can be only one independent variable.

You have three first order differetnial equations - so you do not "meet the criteria". Hence no direction filed is produced.

Obvioulsy you can use odeplot() to produce pretty much plots of anything versus anything - so I guess you have to ask yourself, do I really need a "direction field" ???

This appears to be one of these cases where one must have the "restart" command in a separate execution group, which (to be fair) is a requirement. Excerpt from the help

It must be executed in a separate prompt (or line) from all other commands, since all commands in a prompt are passed to the kernel at once; entering other commands in the prompt could cause unexpected results.

Consider the two files attached below: on my system, the first NEVER works, and the second ALWAYS works. (Note before testing these ensure that the test file ('fileName') has been deleted from the directory ('pathName')

  restart:
  with(plots):
  with(plottools):
  with(ColorTools):
  with(FileTools):
  interface(version);
#
# OP will have to change the following pathName
# to an appropriate directory for his/her machine.
#
# Just for test purposes on my own system, I have
# used the Windows Desktop
#
  pathName:="C:/Users/TomLeslie/Desktop":
  fileName:="Test.jpg":
  P:=display(polygonplot([[0,0],[0,1],[1,1],[1,0]], color=red,axes=none,style = polygon)):
#
# Check whether directory is writable and whether the file
# already exists - for test purposes, these two commands
# should return
#
#    true
#    false
#
# respectively
#
  IsWritable(pathName);
  Exists( JoinPath([pathName, fileName]) );
#
# If the following command "works" then all it should return
# is the number of bytes written - but it fails!
#
  exportplot( JoinPath([pathName, fileName]),P, plotoptions="height=30,width=30,quality=90");

`Standard Worksheet Interface, Maple 2021.2, Windows 7, November 23 2021 Build ID 1576349`

 

true

 

false

 

 

Error, (in plottools:-exportplot) exported file C:/Users/TomLeslie/Desktop\Test.jpg could not be created

 

 

Download testRe1.mw

  restart:

  with(plots):
  with(plottools):
  with(ColorTools):
  with(FileTools):
  interface(version);
#
# OP will have to change the following pathName
# to an appropriate directory for his/her machine.
#
# Just for test purposes on my own system, I have
# used the Windows Desktop
#
  pathName:="C:/Users/TomLeslie/Desktop":
  fileName:="Test.jpg":
  P:=display(polygonplot([[0,0],[0,1],[1,1],[1,0]], color=red,axes=none,style = polygon)):
#
# Check whether directory is writable and whether the file
# already exists - for test purposes, these two commands
# should return
#
#    true
#    false
#
# respectively
#
  IsWritable(pathName);
  Exists( JoinPath([pathName, fileName]) );
#
# If the following command "works" then all it should return
# is the number of bytes written
#
  exportplot( JoinPath([pathName, fileName]),P, plotoptions="height=30,width=30,quality=90");

`Standard Worksheet Interface, Maple 2021.2, Windows 7, November 23 2021 Build ID 1576349`

 

true

 

false

 

15988

(1)

 

Download testRe2.mw

A couple of issues with fsolve() are

  1. how many "iterations" are tried before fsolve() decides to give up
  2. little detail is available on what constitutes an "acceptable" solution. In other words what residual error in each of the equations is "acceptable". The criterion is probably(?) related to the current setting of Digits and if fsolve() cannot achieve this criteria (in a specific/unlnown number of iterations), it returns unevaluated.

One has a little more control over some of the routines in the Optimization package, where one can play with options such as iterationlimit and optimality tolerance. Sometimes this will return a solution when fsolve() does not - although the returned solution may well have higher residuals than would be accepte by fsolve().

For example the attached shows NLPSolve() obtaining a result for the case tau_l=0.8, a value for whihc fsolve() "fails". The residual error in each of the equations is between 1.0e-05 and 1.0e-09. Maybe this is acceptable?

For some reason although the attached worksheet can be displayed in-line, the resuting response cannot be "submitted", so you will have to download from the link

asOpt.mw

The most basic thing you have to realise is that 1/3 is not the same as 0.333333333 - the latter is a floating-point approximation to the former.

Maple will never make such approximations unless you ask it to.

If a calculation involves only exact (ie rational) numbers, then Maple will return an exact (ie rational) answer.

If a calculation involves a mixturre of exact (ie rational) and floating-point (ie approximate) numbers, then the answer can only be "approximate", so Maple will return the floating-point approximation

Consider the following code - which representation of ans2 do you consider to be correct? The first one is the exact answer, subsequent values are approximations to this exact answer

ans1:=(654321.987*123456.789)/2;
ans2:=(654321987*123456789)/2000;
evalf(ans2);
evalf[20](ans2);
evalf[30](ans2);
                                          10
                    ans1 := 4.039024574 10  

                           80780491487119743
                   ans2 := -----------------
                                 2000       

                                      13
                        4.039024574 10  

                                           13
                   4.0390245743559871500 10  

                                                13
              4.03902457435598715000000000000 10  

 

is shown in the attached. A similar approach should work for more complicated problems.

  restart;
  ode := diff(f(x),x) = f(x);
  initial_condition := f(0)=1;
  a := 1;

  f:=  eval
       ( f(x),
         dsolve
         ( {ode,initial_condition},
           numeric,
           output=listprocedure
         )
       );
  solution := fsolve
              ( f(x) - f(a) = 1,
                x
              );

Like Scot Gould this site will display my worksheet in-line, but then "hangs" when I try to submit it

odeProb.mw

If you want to use indexed variables then the var{something] syntax is correct. If you just want a subscript, becuase it looks pretty, then the var__something syntax is correct.

In your particular case quantities such as Psi[fd] will not cause a problem because Maple will assume that 'fd' is a name which may subsequently be assigned a meaning value for an indexed variable. However i[1t] will be misinterpreted becuase '1t' is not an acceptable 'name' in MAple. As the relevant help page points out

A name in its simplest form is a letter followed by zero or more letters, digits, underscore characters (_) and question marks (?), with lowercase and uppercase letters distinct.

So you cannot begin a name with a digit -and if i[1t] is supposedly an indexed variable but the index is not a name - then what does this index mean??

First 32 33 34 35 36 37 38 Last Page 34 of 207