tomleslie

12980 Reputation

19 Badges

12 years, 290 days

MaplePrimes Activity


These are answers submitted by tomleslie

for something like the attached, perhaps?

oddPlot.mw

As a general rule, I wouldn't do this in this manner for a couple or reasons, which are

  1. Don't generate data in a crap format, then attempt to fix it!. Generate the data in a simple-to-read format in the first place
  2. Don't use Maple as a text editor. If you have to generate data in a crap format, then it is (probably?) easier to use a decent programmer's editor to convert it to a simple-to-read format

Having got that off my chest, the attached is clunky, but it works. All data is stored in an Array. Each entry of the Array is a five-element list. The first four elements of the list are the requested info from the header for a particular dataset. The final entry in each list is a matrix containing the x-y data for that dataset.

In the attached I have used

fname:="C:/Users/TomLeslie/Desktop/2011068_2_ex.txt":

You will have to change this line so that the path corresponds to wherever you have the file stored on your machine.
readData.mw

you will have to place it "manually" using a textplot3d() command, as in the attached.

sinFun.mw

have to "assume" that t>0, as in the attached
assum.mw

the DirectSearch() package - an add-on available from the Maple Application Centre. I am however very doubtful about the actual success of any optimization package on a function which appears to have many, many local minima and maxima.. Contary to popular belief no optimisation algorithnm is capable of finding a global minimum (or maximum) for anything other than convex objective functions.

So for what iis worth you can check out the attached. NB you will not be able to execute this worksheet unless you have the DirectSearch() package installed, so treat the atached as read-only.


toughOpt.mw

@zenterix 

You state

The documentation seems like a big reference manual.

Yes the documentation is a big reference manual. However it also contains a lot of "getting started" info. In fact if you just use the menus Help->Maple Help, this will open the built-in Maple help, and the first section in the in the Table of Contents pane on the left-hand side is labelled "Getting Started" The help page displayed states (amongst other things)

Manuals
• The User Manual covers an introduction to the Maple interface and the basics of using the Maple programming language.
• The Programming Guide covers introductory, intermediate and advanced programming techniques in the Maple programming language.


Both the User Manual and Programming Guide are available as PDF downloads at the Maplesoft online documentation center.

From this starting help page, you can access the User Manual in an (executable) worksheet and work your way through Maple basics. In PDF format, this manual is about 350 pages - and by the time you have finished reading it, you will have an understanding of Maple basics. The Programming Guide is when things get serious! In pdf format about 650 pages on how to write Maple programs.

You also state

Take the documentation on the Grid package. It lists out the 18 commands in the package. Am I supposed to go through all of the 18 to find out how this package is usually used? 

No probably not! Try reading Chapter 15 Parallel Programming of the Programming Guide , which serves as an introduction to both the 'Threads' and 'Grid' models for parallel programming. This covers the basic commands in both packages in a "didactic" fashion

You further state

Have a look at this page, for Grid:-Run. This page is seemingly abandoned. The formatting of the page runs way to the right. Someone didn't care at all about the reader of this. Like, is this command so unimportant that no one cares to maintain the documentation page for it?

The equivalent page in the built-in help looks absolutely fine - why would uou use on-line help rather than built-in?

The expression you waant to plot, assigned to the upper-case Greek letter Nu (which confusingly looks like an upper-case Roman N), does not depend on your plotting variables , y or N. This expression also does not contain the variable 'x', so why substitute for it?

Thus your 2D plots will be three horizontal lines (nothing varies with 'y') and your 3D plots will be three horizontal planes (nothing varies with either 'y' or 'N' ),

Can one combine these 2D and 3D plots? - Yes, see the attached: but what is the point?

oddPlot.mw

For some reason this worksheet won't display inline here, so you will have to download to examine

logical and syntactical errors in the garbage code which you provide that I doubt anyone could either

  1. fix the garbage
  2. even figure out what you actually want

The attached may be what you are trying to produce, although I had to make so many guesses that I seriously doubt it

  restart:
  with(plots):   
  with(geometry):   
  _EnvHorizontalName := 'x':
  _EnvVerticalName := 'y':  
  EQ := proc(M, N)
             return ((y - M[2])/(x - M[1]) = (N[2] - M[2])/(N[1] - M[1])):
        end proc:

  a := 7:  b := 5:  
  c:=sqrt(a^2-b^2):
  t:= (3* Pi)/(4):    
  ellipse(el,x^2/a^2 + y^2/b^2 - 1 ):
  point(M1, a*cos(t), b*sin(t)):  
  point(F1, -c, 0):
  point(F2, c, 0):    
#  line( L1, [F1, M4]):    
  line( L1, [F1, M1]):

#  point(M1, a*cos(t),b*sin(t)):  
  EQ([a*cos(t),b*sin(t)],[-c,0]):  
  solve(`%`,y):
  m:=coeff(`%`,x):   
  line(L3,y=m*(x-c)):
  op(solve({x^2/a^2 + y^2/b^2 - 1 =0, y=m*(x-c)},{x,y},explicit)[2])[1]:
  xM2:=rhs(`%`):  
  op(solve({x^2/a^2 + y^2/b^2 - 1 =0, y=m*(x-c)},{x,y},explicit)[2])[2]:
  yM2:=rhs(`%`):  
  point(M2,xM2,yM2):
  line( tang1,x*xM2/a^2 + y*yM2/b^2 = 1):  
  intersection(C1,L1,tang1):  
  line(L4,[F1,M2]):
#  C1:=[ HorizontalCoord(C1), VerticalCoord(C1)]:
#   F1:=[-c,0]:
#  M2:= [ HorizontalCoord(M2), VerticalCoord(M2)]:
  triangle(T1, [C1,F1,M2]):
  display
  ( [ textplot
      ( [ [ -c, 0, "F1"],
          [ c,  0, "F2"] ,                         
          [ coordinates(M1)[], "M1"],                          
          [ coordinates(M2)[], "M2"],                         
          [ coordinates(C1)[], "C1"]
        ],
        align={"above",'right'}
      ),
      draw
      ( [ el(color=red),
          M1(color=black,symbol=solidcircle, symbolsize=16),
          M2(color=black,symbol=solidcircle, symbolsize=16),                      
          L1( color=black),
          L3(color=green),
          L4(color=green),                      
          tang1(color=red),
          F1(color=blue, symbol=solidcircle, symbolsize=16),                      
          F2(color=red, symbol=solidcircle, symbolsize=16),
          C1(color=black, symbol=solidcircle, symbolsize=16),
          T1(color=blue, filled=true, transparency=0.9)
        ],
        axes=none
      )
    ]
  );
 

 

 

Download ellipseStuff.mw

trivial to plot the analytical solution for the same parameters as the numerical solutions - either separately or together, see the attached.

One would have to say that the agreement is terrible - a typo somewhere maybe?

restart

ode := diff(theta(y), y, y)+Br*((-3*We*(y+1+(1/2)*x^2)*((9*(1+(1/2)*x^2)^2+9*y^2)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^3/(4*(1+(1/2)*x^2)^7)+3*y*(k+1)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^2/(1+(1/2)*x^2)^5+3*(k+1)^2*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)/(2*(1+(1/2)*x^2)^3)-(-120*lambda*(1+(1/2)*x^2)^4+(63*k^2-36*k+63)*(k-1)*(1+(1/2)*x^2)^3+288*Q*(k^2-(11/8)*k+1)*(1+(1/2)*x^2)^2+486*Q^2*(k-1)*(1+(1/2)*x^2)+324*Q^3)/(15*(1+(1/2)*x^2)^5))/(16*(1+(1/2)*x^2)^2)-3*We*(y-1-(1/2)*x^2)*((9*(1+(1/2)*x^2)^2+9*y^2)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^3/(4*(1+(1/2)*x^2)^7)+3*y*(k+1)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^2/(1+(1/2)*x^2)^5+3*(k+1)^2*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)/(2*(1+(1/2)*x^2)^3)-(-120*lambda*(1+(1/2)*x^2)^4+(63*k^2-36*k+63)*(k-1)*(1+(1/2)*x^2)^3+288*Q*(k^2-(11/8)*k+1)*(1+(1/2)*x^2)^2+486*Q^2*(k-1)*(1+(1/2)*x^2)+324*Q^3)/(15*(1+(1/2)*x^2)^5))/(16*(1+(1/2)*x^2)^2)-3*We*(y-1-(1/2)*x^2)*(y+1+(1/2)*x^2)*(9*y*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^3/(2*(1+(1/2)*x^2)^7)+(3*(k+1))*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^2/(1+(1/2)*x^2)^5)/(16*(1+(1/2)*x^2)^2)+(-3*y*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)/(1+(1/2)*x^2)^2-k-1)/(x^2+2))^2-(1/2)*We*(-3*We*(y+1+(1/2)*x^2)*((9*(1+(1/2)*x^2)^2+9*y^2)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^3/(4*(1+(1/2)*x^2)^7)+3*y*(k+1)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^2/(1+(1/2)*x^2)^5+3*(k+1)^2*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)/(2*(1+(1/2)*x^2)^3)-(-120*lambda*(1+(1/2)*x^2)^4+(63*k^2-36*k+63)*(k-1)*(1+(1/2)*x^2)^3+288*Q*(k^2-(11/8)*k+1)*(1+(1/2)*x^2)^2+486*Q^2*(k-1)*(1+(1/2)*x^2)+324*Q^3)/(15*(1+(1/2)*x^2)^5))/(16*(1+(1/2)*x^2)^2)-3*We*(y-1-(1/2)*x^2)*((9*(1+(1/2)*x^2)^2+9*y^2)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^3/(4*(1+(1/2)*x^2)^7)+3*y*(k+1)*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^2/(1+(1/2)*x^2)^5+3*(k+1)^2*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)/(2*(1+(1/2)*x^2)^3)-(-120*lambda*(1+(1/2)*x^2)^4+(63*k^2-36*k+63)*(k-1)*(1+(1/2)*x^2)^3+288*Q*(k^2-(11/8)*k+1)*(1+(1/2)*x^2)^2+486*Q^2*(k-1)*(1+(1/2)*x^2)+324*Q^3)/(15*(1+(1/2)*x^2)^5))/(16*(1+(1/2)*x^2)^2)-3*We*(y-1-(1/2)*x^2)*(y+1+(1/2)*x^2)*(9*y*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^3/(2*(1+(1/2)*x^2)^7)+(3*(k+1))*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)^2/(1+(1/2)*x^2)^5)/(16*(1+(1/2)*x^2)^2)+(-3*y*(k*(1+(1/2)*x^2)+2*Q-1-(1/2)*x^2)/(1+(1/2)*x^2)^2-k-1)/(x^2+2))^4) = 0

bc := theta(-1-(1/2)*x^2) = 0, theta(1+(1/2)*x^2) = 1

Q__vals := [.5617, .4392, .2564, .1645, 0.659e-1]

`λ__vals` := [0.96e-2, 0.10e-2, 0.93e-2, 0.75e-2, 0.31e-2]; k__vals := [.1, .3, .5, .7, .9]

We__vals := [.1, .3, .5, .7, .9]

colors := [black, red, blue, green, orange, cyan, purple, yellow, navy]; for j to numelems(k__vals) do sol := dsolve(eval([ode, bc], [k = k__vals[j], We = .1, Br = .3, x = 0, Q = Q__vals[j], lambda = `λ__vals`[j]]), numeric); pltn[j] := plots:-odeplot(sol, [y, theta(y)], y = 0 .. 1, axes = boxed, color = colors[j]) end do; plots:-display([seq(pltn[j], j = 1 .. 5)], title = "Numerical Solution")

 

Analytical_sol := -(3*y+3+3*x^2*(1/2))*((-2/3+(-(5/12)*k^2+(1/6)*k-5/12)*Br)*(1+(1/2)*x^2)^5-lambda*Br*((-(1/12)*k^2+(1/6)*k-3/4)*y+Q*(k-1))*(1+(1/2)*x^2)^4-(((1/6)*k+7/6)*y+Q)*Br*((-(1/2)*k+1/2)*y+Q)*(1+(1/2)*x^2)^3+We*Br*((1/4)*(k-1)^2*y^2-(1/3)*Q*(k-5)*y+Q^2)*y*(1+(1/2)*x^2)^2-((1-k)*y+Q)*Br*y^2*Q*(1+(1/2)*x^2)+Q^2*y^3*Br)/(4*(1+(1/2)*x^2)^6)

for j to numelems(k__vals) do sol[j] := eval(Analytical_sol, [k = k__vals[j], We = .1, Br = .3, x = 0, Q = Q__vals[j], lambda = `λ__vals`[j]]); plta[j] := plot(sol[j], y = 0 .. 1, axes = boxed, color = colors[j], linestyle = dash) end do; plots:-display([seq(plta[j], j = 1 .. numelems(k__vals))], title = "Analytical Solution"); plots:-display([seq(plots:-display([pltn[j], plta[j]]), j = 1 .. numelems(k__vals))], title = "Numerical and Analytical Solutions")

 

 

``

 

Download plots.mw

if you posted the complete worksheet for your calculation using the big green up-arrow in the Mapleprimes toolbar.

In the attached I have performed the calculation you seem to want in two completely different ways.

The first uses Maple's bult-in numapprox() and orthopoly() packages.  The second uses an "add-on" package (which you probably don't have?) called OrthogonalExpansions() which is available from the Maple Application Centre. These return the same answers: the desired series is

-0.577350269189626*T(0, x) - 0.309401076758503*T(1, x) - 0.0829037686547607*T(2, x)

The coefficients in this series disagree with those which you report - but since you do not specify how you obtained your coefficients, there isn't much I can say?

When converted to a simple polynomial in 'x', both give

-0.494446500534865 - 0.309401076758503*x - 0.165807537309521*x^2

and both return your measure of "error" as 0.02776926226

  restart;
#
# Use built-in approach with numapproz and orthopoly packages
#
# NB numapprox(0 will only give floating point coefficients
# for the series. Retyrned series is claimed to be accurate
# to 10^-(Digits) - Diigits is on default, ie 10
#
  u:=1/(x-2);
  ser0:= numapprox:-chebyshev(u, x=-1..1):
#
# Restrict to first three terms, and use orthopoly() to
# evaluate the ChebyshevT polynomials. Note that the latter
# returns the same as using the OrthogonalExpansions package
# below
#
  ser2:=add(op([1..3], ser0));
  ser3:=eval( ser2, T=orthopoly[T]);
#
# Compute OP's definition of "error"
#
  evalf(Int(abs(u - ser3), x = -1 .. 1));

1/(x-2)

 

-HFloat(0.5773502691896257)*T(0, x)-HFloat(0.309401076758503)*T(1, x)-HFloat(0.08290376865476069)*T(2, x)

 

-HFloat(0.49444650053486505)-HFloat(0.309401076758503)*x-HFloat(0.16580753730952139)*x^2

 

0.2776926226e-1

(1)

  restart;
  with(OrthogonalExpansions):
  u:=1/(x-2);
#
# Generate ChebyshevT series of 'u' above over the
# range x=-1..1
#
  ser := ChebyshevTSeries(u, x = -1 .. 1, 'n', 'Coefficients');
#
# Restrict series to th first three terms, ie, ChebyshevT(0, x)
# ChebyshevT(1, x), and ChebyshevT(2, x)
#
# Return the restricted series and the coefficients. Compare the
# latter with the values in 'ser2' above
#
  n:= 2;
  ser2:=value(ser);
  evalf(value(Coefficients));
#
# Use the orthopoly() package to convert the ChebyshevT polynomials
# Compare 'ser3' below with 'ser3' in previous execution group
#
  ser3:=evalf(eval(ser2, [ seq( ChebyshevT(i, x)=orthopoly:-T(i,x), i=0..2)]));
#
# Check OP's definition of the "error" between the series
# expansion and the orifginal function
#
  evalf(Int(abs(u - ser3), x = -1 .. 1));

u := 1/(x-2)

 

Sum(piecewise(i = 0, 1/Pi, 2/Pi)*(int(ChebyshevT(i, x)/((x-2)*(-x^2+1)^(1/2)), x = -1 .. 1))*ChebyshevT(i, x), i = 0 .. n)

 

2

 

-(1/3)*3^(1/2)*ChebyshevT(0, x)+2*(Pi-(2/3)*Pi*3^(1/2))*ChebyshevT(1, x)/Pi+2*(4*Pi-(7/3)*Pi*3^(1/2))*ChebyshevT(2, x)/Pi

 

[-.5773502693, -.3094010774, -0.8290376528e-1]

 

-.4944465040-.3094010774*x-.1658075306*x^2

 

0.2776926183e-1

(2)

 

Download cheby.mw

This illustrates how I'd probably do it.

  1. The ODE system is a simple predator-prey model.
  2. A list of solution matrices is obtained by varying one of the initial conditions
  3. The list of solution matrices is plotted as an animation by using the display() command with the option, insequence=true
  4. The list of solution matrices is exported to a (matlab) binary file - I think this is about the only option if one wants to export a list of matrices (as distinct from a single matrix). Obviously you will have to change the path name for this file to something appropriate for your machine
  5. restart;
  6. Import the list of solution matrices generated by 4 above
  7. Repeat 3 above

  restart;
  with(plots):
  odesys:= [ diff(x(t),t) = x(t) - x(t)*y(t),
             diff(y(t),t) = x(t)*y(t) - y(t)
           ] :
  tv:= Array([seq( j, j=0..10,0.1)]):
  ics:= [x(0)=k, y(0)=1]:
#
# Solve odesys, for different initial conditions
# and extract the matrix which contains the results
# to be plotted
#
# 'sol' is a list of such matrices
#  
  sol:= [ seq( dsolve( [odesys[], eval(ics, k=j)[]],
                       numeric,
                       output=tv)[2,1][..,2..3],
               j=1.5..2.5, 0.1
             )
        ]:
#
# Display the results as an animation
#
  display( [seq( listplot( sol[j]), j=1..numelems(sol))], insequence=true);
#
# Export the list of solution matrices
#
  ExportMatrix( "C:/Users/TomLeslie/Desktop/test.mat", sol, target=Matlab);

 

18095

(1)

  restart;
  with(plots):
#
# Import the list of solution matrices which was generated above
#
  ans:=[ImportMatrix( "C:/Users/TomLeslie/Desktop/test.mat", output=matrices)]:
#
# Display the results as an animation
#
  display( [seq( listplot( ans[j]), j=1..numelems(ans))], insequence=true);

 

 

Download toyAnim.mw

This is actually used as an example (with multiple plotting possibilities) on the help page for DEPlot() .

There are *a lot* of options for this command.  One possibility for your specific set of equations is shown in the attached.

If you don't want actual solution curves (for which I admit I just made up initial values), then you can just delete the lines in the attached

[ [x(0)=2,   y(0)=1],
  [x(0)=0.5, y(0)=3],
  [x(0)=1.5, y(0)=1.5]
],
linecolor=[ green, orange, yellow],

  DETools:-DEplot( [ diff(x(t),t) = x(t) - x(t)*y(t),
                     diff(y(t),t) = x(t)*y(t) - y(t)
                   ],
                   [x(t),y(t)],
                   t=0..10,
                   x=0..3,
                   y=0..3,
                   [ [x(0)=2,   y(0)=1],
                     [x(0)=0.5, y(0)=3],
                     [x(0)=1.5, y(0)=1.5]
                   ],
                   linecolor=[ green, orange, yellow],
                   arrows=mediumfill,
                   color=magnitude,
                   dirfield=[10,10],
                   axes=boxed,
                   labels=["prey/density", "predator/density"],
                   labelfont=[times, bold, 12]
                 );

 

 

Download predprey.mw

 

 

on whether you want a 2-D contour plot or a 3-D contour plot! Both come with many options. I suggest you read the help at ?contoutplot.

Examples of 2D and 3D plots for your function are shown in the attached

  restart;
  with(plots):
  f:=(x,y)-> (x - 3)*(y - 2) - 2;
  contourplot( f(x,y),
               x=-10..10,
               y=-10..10,
               contours=20,
               filledregions=true,
               coloring=[red, blue]
             );
  contourplot3d( f(x,y),
                 x=-10..10,
                 y=-10..10,
                 contours=20,
                 filledregions=true,
                 coloring=[red, blue]
               );

proc (x, y) options operator, arrow; (x-3)*(y-2)-2 end proc

 

 

 

 

Download conPlot.mw

  1. Error is not new - it occurs in Maple 2021.2. (i haven't bothered to check anything earlier
  2. "too many levels of recursion" is untrappable. see the excerpt from the Maple help below, with emphasis added

Some exceptions are untrappable. This means that a catch clause will never match such an exception, even if the catchString is missing. The finalStatSeq is also not executed if an untrappable exceptions occurs. The following exceptions are untrappable:
• Computation timed out (this can only be caught by timelimit, which raises a "time expired" exception, which can be caught).
• Computation interrupted (user pressed Ctrl+C, Break, or equivalent).
• Internal system error (which indicates a bug in Maple itself).
• ASSERT or return type or local variable type assertion failure (assertion failures are not trappable because they indicate a coding error, not an algorithmic failure).
Stack overflow (when that happens, generally stack space is insufficient for doing anything like running cleanup code). This includes the "too many levels of recursion" exception at the Maple language level.

 

you are going to have to upload an actual worksheet (using the big green up-arrow in the Mapleprimes toolbar) illustrating the problem, because the attached just works! (although for some reason it won't display inline on this site)

Download oddInt.mw

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