Check out this routine, for analysis of asymptotes of rational (polynomial) functions like your example.

Student[Precalculus][RationalFunctionTutor]();

To answer your other questions, where we are talking about functions y(x) with the x-axis being the horizontal axis and the y-axis being the vertical axis,

• No. The function sin(x)/x is undefined as x=0, and it does not have a vertical asymptote at that point. A function may be undefined at a point and still not have an infinite limiting value at that point. As the sin(x)/x example shows, the graph needn't be approaching a vertical slope at that point either.
• Yes, a function's graph can cross its horizontal asymptotes. The example of (x-2)/(x^2-5x) does this, by cross its horizontal asymptote (which occurs at y=0, as x tends to either -infinity or +infinity) at the point x=2. No, a function's graph may not cross its vertical asymptotes (else those would not actually be asymptotes, or it would be a graph of a a relation rather than of a function).
• A function's graph can have at most two horizontal asymptotes. Remember that horizontal asymptotes occur according to behaviour as x tends to -infinity or +infinity (which only happens once, for each of those).

acer

What sort of a fit do you want? You mentioned "least squares", but even so there is some choice.

Do you want usual linear regression (minimizing the sum of the squares of distance with one particular component such as z be taken as dependent)?

Or do you want to find a line which minimizes the sum of the Euclidean distance between the data points and that line. That is the othogonal distance regression, a four-parameter problem sometimes solved using singular value decomposition, but which can also be solved by generating an objective from all the data points'  point-to-line distance expressions.

There are some other choices, but those are the common ones. The first might be possible using Statistics[LinearFit] and some footwork. The second might be more "usual".

acer

The gif image of the desired equivalent result is messed up and corrupt, as I see it displayed here on mapleprimes. But the Properties of it say that the Alternate text for it is this,

75*Pi^5*(k*T)^4/(30*h^3)

Is that correct, or did you intend something else? For example, maybe you intended it as 7*Pi^5*(k*T)^4/(30*h^3) , with a 7 instead of 75?

acer

You could teach your students that trying to obtain a RREF for a floating-point numeric Matrix is not an good goal. There are computations appropriate for exact rational Matrices, and there are computations appropriate for floating-point Matrices, and those sets of computations are not necessarily the same for both those domains. That is a much better thing to teach students than some tricks with Maple's convert(...,rational) utility.

First, always ask the question, what is the purpose of the numerical computation in question? Are the numeric results supposed to give immediate insight, or are they intended as inputs to further calculations? Let's take one of each such purpose, which for exact rationals the RREF might serve.

The first is that the RREF might be used to gauge the rank of the Matrix. But row-reduction for the purpose of rank determinantion of a floating-point matrix is misguided, not just in Maple but in general. You might wonder whether some sort of small tolerance could be used, to identify a pivot as being "not close enough to zero to be taken as zero." But a little thought and experiment can show that when trying to use such a tolerance it is the collection of all other entries which seem to change the judgement of whether a pivot should be zero. It's not an absolute question, involving just that entry. So then we might ask, well, is there some other method, to gauge the rank, which does take into account in some way the complete set of entries where they affect each other? The answer is yes, there is a usual such method, and it is based on computing the singular values. Usually the greatest singular value is taken and the other compared to it. The rank can then be assessed as the number of singular values "close enough" to the largest. The working precision (in Maple, Digits) should also bear upon this result. You should be able to find lots of details about this on the web. It's how Maple computes the rank for a floating-point Matrix, and it's how Matlab does it too. (NB. The singular values are themselves not computed as the square-roots of the eigenvalues of the M^T.M, because that too is an algorithm appropriate for exact data but not for floating-point data.)

Perhaps an example would illustrate the difficulties inherent in trying to use the RREF for rank determination of a floating-point Matrix.

Consider Matrix M as given immediately below. Should the small entries in the second row be considered as negligible? Some people might reasonable say "no", because of their size relative to the entries in the first row.

> M := Matrix([[1e-7,2e-8],[1e-14,2e-16]]);
                              [      -6           -7 ]
                              [0.1 10       0.2 10   ]
                         M := [                      ]
                              [      -13          -15]
                              [0.1 10       0.2 10   ]

> with(LinearAlgebra):
> Rank(M);
                                       2

But now let's set M[1,1] to 1000.0 which is relatively large enough compared to the second row entries that many people might expect a rank calculation to return 1 , in a fixed Digits=10 precision environment.

> M[1,1]:=1000.0:
> Rank(M);
                                       1

> ReducedRowEchelonForm(M);
                                  [1.    0.]
                                  [        ]
                                  [0.    1.]

> Digits := 24:
> Rank(M);
                                       2

A second goal might be to solve linear systems, by further reduction above the diagonal of a Matrix that had been augmented by the identity matrix, say. The resulting inverse of the original Matrix might then be multiplied against any number of RHS vectors or matrices so as to solve the various linear systems. But this can lead to serious numerical error when attempted in the floating-point case, and should be discouraged. An appropriate method for solving repeated linear systems (same LHS Matrix) would be to compute the so-called LU factorization, where the L and the U may be used repeatedly for distinct RHS's.

Similar remarks go for computing eigenvalues, or for computing a nullspace (via singular vectors, in the floating-point case), etc, etc. Numerical (floating-point) linear algebra is often taught as a separate (advanced) course in itself, and introductory linear algebra concepts usually best presented using exact domains.

acer

There are two modes for the Standard Maple (Java) GUI. One is Document mode, and the other is Worksheet mode.

What you have described is the same as the behaviour for output in a Document. The behaviour in a Worksheet is the same as you've described, where the e:= assignment is also echoed in the output. In Document mode, only the assigned value is echoed.

A Document is a newer piece of Maple technology, and I wouldn't be surprised if the manuals you've seen were authored in the older Worksheet environment.

You should always be able to forcibly start either, Document or Worksheet, from the menubar's "File" -> "New" choices.

You should also be able to set a default for new sessions, to be either mode, by going to the menubar's choice "Tools" -> "Options" -> Interface tab, and then adjusting the "Default format for new worksheets" (dropdown). Even after changing the default for new sessions, you would still be able to open either, if you wished, using the File -> New menu choice.

If you look around this site, you can find lots of posts where people express opinions and experiences with Document vs Worksheet as well as another related choice of 1D vs 2D Math input (also with a default in the Options menus, and also forcible either way at any time via the F5 key).

acer

This usenet post might help with that question. (Casper was running through those Oxford exercises on Maple too.)

acer

The term mean is often used as another word for the average.

When used on its own as just "mean", it usually refers to the arithmetic mean.

The Maple function stats[mean] computes the arithmetic mean, which is a weighted arithmetic average. Each value is multiplied by its own weight, then those values added together, and then the sum is divided by the total weight.

Another common reference is the geometric mean, which is similar to the arithmetic mean but involves multiplication rather than addition.

You can get full definitions of these terms, in Maple, by looking at the following help-pages,

Definition/arithmeticmean

Definition/geometricmean

Definition,mean

acer

The short answer is to use add() instead of sum().

The longer answers involve explaining how sum() sees g(M,i,2) and tries to evaluate it before `i` is set with an integer value in your loop. See what happens when you enter g(M,i,2) all on its own, as input.

acer

They are not all not the same `i`. The important distinction is due to some of those `i`s occuring within operator definitions.

The operators in a[1], a[2], and a[3] are procedures with no variable `i` declared in them. Also, in those procedures the variable `i` is not assigned to. See the section "Implicit Local Variables" in the help-page ?procedure for an explanation of how the value of such names is inferred when the procedure gets run.

As a consequence of this, when you run it ( eg. a[1](t) ) you get t*4 as it find a value of 4 for i in the surrounding level. That's because, after your loop finishes, the value of the `i` at your current level is 4. If you had assigned 17 to i, afterwards, then subsequent invocation of your operators would produce 17*i.

One way to work with this behaviour of Maple is to create an operator with a dummy name in it, and then to substitute the loop variable i's current value for that dummy each time through the loop.

> for i from 1 to 3 do
>     a[i] := subs(ii=i,t -> ii*t);
> end do;
                                a[1] := t -> t
 
                               a[2] := t -> 2 t
 
                               a[3] := t -> 3 t

For your example, you could also define that template operator just once, outside of the loop. The call to eval() below is there just because operators have so-called last-name evaluation rules, see ?last_name_eval for details of this side issue.

> restart:
> myop:=t -> ii*t:
> for i from 1 to 3 do
> a[i] := subs(ii=i,eval(myop)); # or eval(myop,1)
> end do;
                                a[1] := t -> t
 
                               a[2] := t -> 2 t
 
                               a[3] := t -> 3 t

You could also do this sort of thing more directly, using an Array instead of an array.

> restart:
> myop:=t -> ii*t:
> A := Array(1..3, (i)->subs(ii=i,eval(myop,1)) );

                       A := [t -> t, t -> 2 t, t -> 3 t]

Lastly, some easier ways to do it.

> restart:
> A := Array(1..3, (i) -> t -> i*t );

                      A := [t -> 1 t, t -> 2 t, t -> 3 t]

> vector(3, (i) -> t -> i*t );

                        [t -> 1 t, t -> 2 t, t -> 3 t]
> lprint(%);
array(1 .. 3,[(1)=(t -> 1*t),(2)=(t -> 2*t),(3)=(t -> 3*t)])

What happens during the second example just above is essentially this.

> f := i -> t -> i*t:
> a := array(1..3):
> for i from 1 to 3 do
>     a[i] := f(i);
> end do;
                               a[1] := t -> 1 t
 
                               a[2] := t -> 2 t
 
                               a[3] := t -> 3 t

acer

Not only can you do this with a simple Maple commands, you can also make a procedure that can be re-used to do it for any function you want.

Here's a re-usable routine to do the job.

> T := proc(F,a,b,inc)
> local i;
> for i from a to b by inc do
> printf("%a %a \n",i,F(i));
> end do;
> end proc:

And here's an example of using it.

> f := x -> x^2:

> T( f, 1.0, 3.0, 0.5 );

1.0 1.00
1.5 2.25
2.0 4.00
2.5 6.25
3.0 9.00

Of course, if you just need to do it once you can get similar effects with calls like this.

seq( printf("%a %a \n",0.5+0.5*i,f(0.5+0.5*i)), i=1..5 );

acer

The task you are trying to do is not uncommon. Your attempt was on the right track. What should work here is the `if` operator.

Matrix(3, (i,j) -> `if`(P[i,j]=1, x[i,j], P[i,j]) );

acer

There are two ways to do that sort of thing. The first is to qualify the command by appending `assuming` on the end of it. The second is to place an assumption on the name such as `n` by calling assume(). That second method cuts down on the typing, but makes it harder to use the name `n` elsewhere without any assumptions on it.

> sin(n*Pi) assuming n::integer;
                                  0
> restart:
> assume(n::integer):
> cos(n*Pi);

                                        n~
                                    (-1)

acer

I suspect that what is happening is that the internal iteration limit of the least-squares solver is being hit before the optimality tolerance (goal) is met, or before some numerical effect takes over the computation.

I've extracted what Statistics:-NonlinearFit passes to Optimization:-LSSolve , thanks to Maple debugger.

Once the call to LSSolve has been produced, it can be modified to alter the iterationlimit and also the choice of numerical method. See ?LSSolve for more.

First I run with method=sqp and iterationlimit=10000. I also set infolevel[Optimization] to 10, to get more information displayed about what it's doing.

> fnew := proc (w, z, needfi, ruser) local i; for i to 4 do z[i] := evalf(proc ( p, r, i) p[1]*(1-p[2]*cos(Pi/(r[i,1]+1)))^(1/2)-r[i,2] end proc(w,ruser,i)) end do end proc:
>
> fjac := proc (w, z, ruser, iuser) local i; for i to 4 do z[i,1] := evalf(proc (p, r, i) (1-p[2]*cos(Pi/(r[i,1]+1)))^(1/2) end proc(w,ruser,i)); z[i,2] := eval f(proc (p, r, i) -1/2*p[1]/(1-p[2]*cos(Pi/(r[i,1]+1)))^(1/2)*cos(Pi/(r[i,1]+1)) end proc(w,ruser,i)) end do end proc:
>
> infolevel[Optimization]:=10:
>
> sol := Optimization:-LSSolve(
  [2, 4],
  fnew,
  datapoints = (Matrix(4, 2, [[8.,-1.73120199199999991],[12.,-1.9665823319999999 9],[16.,-2.12740288800000021],[20.,-2.18699629199999990]], datatype = float[8])) ,
  output = tables,
  initialpoint = (Vector(2, [-9.,.699999999999999956], datatype = float[8])),
  objectivejacobian = fjac,
  method = sqp,
  iterationlimit = 10000 );

LSSolve:   calling nonlinear LS solver
SolveGeneral:   using method=sqp
SolveGeneral:   number of problem variables   2
SolveGeneral:   number of residuals   4
SolveGeneral:   number of nonlinear inequality constraints   0
SolveGeneral:   number of nonlinear equality constraints   0
SolveGeneral:   number of general linear constraints   0
PrintSettings:   feasibility tolerance set to   .1053671213e-7
PrintSettings:   optimality tolerance set to   .3256082241e-11
PrintSettings:   iteration limit set to   10000
PrintSettings:   infinite bound set to   .10e21
SolveGeneral:   trying evalhf mode
Warning, limiting number of major iterations has been reached
E04USA:   number of major iterations taken   10000
                       sol := [settingstab, resultstab]
 

> eval(sol[2]);

table(["visible" = ["objectivevalue", "solutionpoint", "residuals",
    "evaluations", "iterations"], "iterations" = 10000,
    "objectivevalue" = 0.0492868830043316011,
                      [-0.0488744033423635441]
    "solutionpoint" = [                      ],
                      [ -1732.28244861786970 ]
    "residuals" = [-0.241298347220388010, -0.0384274262517967369,
    0.110016320471836959, 0.163611172536819893]
    ])

> eval(a*sqrt(1-b*cos(Pi/(v+1))),[a=sol[2]["solutionpoint"][1],
>                                 b=sol[2]["solutionpoint"][2]]);
                               /                             Pi   \1/2
        -0.0488744033423635441 |1 + 1732.28244861786970 cos(-----)|
                               \                            v + 1 /

I notice that the residuals and the objective function above are decreasing quite slowly as the iterationlimit is increased. You may wish to play around with this a bit more. You could set the method to modifiednewton, or adjust the toleration of iterationlimit some more. If you decide that roundoff error is preventing it from attaining the best solution then you could increase Digits. What is interesting is that, depending on how I alter those options, very different solution points can be produced which each give quite close to the same residuals.

However, with method=modifiednewton and iterationlimit=100000 I got result without seeing the warning.

With those settings, I saw userinfo messages like this,

E04FCF:   number of major iterations taken   1474
E04FCF:   number of times the objective was evaluated   15018

and a final result that included this,

    "objectivevalue" = 0.0492807204246035,


                      [-0.0192100345224023965]
    "solutionpoint" = [                      ],
                      [ -11218.7180037680664 ]

> eval(a*sqrt(1-b*cos(Pi/(v+1))),[a=sol[2]["solutionpoint"][1],
>                                 b=sol[2]["solutionpoint"][2]]);
                               /                             Pi   \1/2
        -0.0192100345224023965 |1 + 11218.7180037680664 cos(-----)|
                               \                            v + 1 /

One possible conclusion from this (not discounting Georgios' comments in a post below about the merit of the nonlinear form that you have chosen to fit to) is that Statistics:-NonlinearFit could benefit from taking a new option like iterationlimit=<value>

acer

What dharr is trying to tell you, I believe, is that the inert form uppercase-e Eval() will display like a usual notation for "y evaluated at x=0". That uppercase Eval is the so-called inert form of eval. You should be able to assign calls like Eval(y,x=0) to variables, or simply use them in your Document, and so on.

On the other hand, the lowercase-e eval() can actually do the evaluation and is not inert. The two often look good in combination, on either side of an equation for example. Try executing this in Maple

y := cos(x);
Eval(y,x=0) = eval(y,x=0);

acer

There are probably better ways. It might also be possible to get the relational symbols from one of the palettes, and subs into that beast (lprinted to see it, say).

> f := (t,s) -> convert(StringTools:-SubstituteAll(
        cat("-T <= ",s," <= T"), "T", convert(t,string)),
                    name):

> f(2,"x");

                                -2 <= x <= 2

acer