acer

32303 Reputation

29 Badges

19 years, 309 days
Ontario, Canada

Social Networks and Content at Maplesoft.com

MaplePrimes Activity


These are answers submitted by acer

If you are trying hard to work out some tricky math and maple programming technique, then the last thing you need is additional hardship due to using 2D Math input and/or Document mode. I would suggest that, until you can nail down the programming and math difficulties, you use 1D input in a Worksheet.

acer

You can do this using Maple's Re() or Im() routines, but it helps if you can also inform Maple that the other quantities H1(k), etc, are purely real.

There are a few ways to do that. One way uses the fact that Maple's evalc assumes that names are generally real.

C(k) := H1(k)/(H1(k) + I*H0(k));

RC:=evalc(Re(C(k)));

IC:=evalc(Im(C(k)));

simplify( RC+I*IC - C(k) );

A second way is to utilize assumptions on those other names via Maple's `assuming` facility.

C := H1/(H1 + I*H0);

Re(C) assuming real;

Im(C) assuming real;

acer

There were two coefficients of air resistance given: one for before deployment and one for after.

If you could alter the first, then the (terminal) velocity prior to deployment could be lower. If you could alter the latter, then the rate of change in velocity following deployment would be greater, and also the steady descent rate would be less.

Skydivers keep the former coefficient near optimal by falling spead-eagle instead of with knees tucked up in cannonball position. After deployment, the size of the parachute is the big factor in air resistance. Maybe you could show graphs of k1 and/or k2 versus impact velocity (if, say, the rip cord is pulled after terminal velocity were attained).

acer

There are a few misconceptions about the contents of the ?EffNumLA help-page apparent in your posting.

  • B:=evalf(A):
    The above does not produce a Matrix with a datatype directly suitable for external calling to the fast compiled routines for numerical linear algebra. Instead all it does it produce a Matrix with floating-point entries. See example code below.
  • DiagonalMatrix(E) won't work, as DiagonalMatrix expects a Vector or list and you have E as a scalar. I tried to guess as what you might have wanted in the code below, and you can replace in a suitable manner.
  • MatrixAdd(B+M) is likely not what you meant, if you intended to add Matrices B and M. You may instead have meant MatrixAdd(B, M).
  • C:=MatrixAdd(B, M), "Should overwrite old matrix". No, it will reassign a new Matrix to name C. The old assigned value (a Matrix, say) of C will then become garbage which takes Maple longer to manage and collect. If you really want to re-use the same memory space of the old C then use one of the float datatypes and the inplace option for MatrixAdd. See code below. (Strictly speaking it will re-use the same memory space for the data only in the hardware float[8] or complex[8] datatype case, and in the software sfloat or complex(sfloat) datatype case will re-use only the container Matrix.)
  • Numerically, it is better to do a LinearSolve than to multiply by an explicitly formed Matrix inverse. A little transposition footwork can allow this even when the eqautions look "back-to-front" from what you might expect.
  • As you can see in the comments, there are a set of progressively better ways to get DD.
  • The datatypes suitable for direct external calling are float[8], sfloat, complex[8], and complex(sfloat). Depending on Digits, the datatypes float and complex(float) will be converted to one of those accordingly.

I did not do anything about replacing the LinearAlgebra:- routines with their more efficient programming forms whose names begin with LinearAlgebra:-LA_Main . That too should be described in the ?EffNumLA help-page.

Below is a modified version of your basic code, which might help illustrate some techniques. You might want to force it through just one iteration, to check that my methods of computing `out` (or DD, as I called the final result) are equivalent.

You didn't say whether your data was symmetric. If it were, then you could get by with 'datatype'='float' throughout, and Transpose instead of HermitianTranspose.

You could likely also make LinearSolve act in-place on its second (Matrix) argument, and the Hermitian Transpose call on that argument as well. That seemed to be the best way to get the final result (`out`, or `DD` when acting in-place) , and is what I've done below (with alternatives commented out).

Digits:=trunc(evalhf(Digits)):
#Digits:=32:
N:=16:
with(LinearAlgebra):
 
A:=RandomMatrix(N,generator=-0.1..0.1,'outputoptions'=['datatype'='float']):
B:=Matrix(op(1,A),(i,j)->evalf(A[i,j]),'datatype'='float'):
 
(st,ba,bu):=time(),kernelopts(bytesalloc),kernelopts(bytesused):
for E from 0.01 by 0.0001 to 1.0 do
  C:=MatrixAdd(B,
               DiagonalMatrix(RandomVector(N,generator=-0.1..0.1,
                                           outputoptions=['datatype'='float']),
                                           outputoptions=['datatype'='float']),
               'inplace'=true):
  (l,L):=Eigenvectors(C):
  #Linv:=MatrixInverse(L):
  #DD:=L.MatrixExponential(DiagonalMatrix(l,outputoptions=['datatype'='complex'('float')])).Linv;
  # But why not this, since the MatrixExponential of a Diagonal Matrix
  # is easy.
  #DD:=L.DiagonalMatrix(map(t->exp(t),l),outputoptions=['datatype'='complex'('float')]).Linv;
  # And why not this, further eliminating a multiplication step.
  #DD := HermitianTranspose(LinearSolve(HermitianTranspose(L),
  #                                     HermitianTranspose(Matrix(N,N,(i,j)->exp(l[j])*L[i,j],
  #                                                               'datatype'='complex'('float')))));
  #out:=Eigenvalues(DD):
  DD := Matrix(N,N,(i,j)->exp(l[j])*L[i,j],'datatype'='complex'('float')):
  HermitianTranspose(DD,inplace=true): # Or just form DD with reversed indices.
  LinearSolve(HermitianTranspose(L),DD,inplace=true); # result is in DD!
  HermitianTranspose(DD,inplace=true):
end do:
time()-st,kernelopts(bytesalloc)-ba,kernelopts(bytesused)-bu;

On my machine, in Maple 11, it took 190sec and Maple claimed to have allocated about 10-14 million words, with Digits=14 and N=16, to do all 9900 loop iterations.

With Digits=32 and N=16, and running only 200 times through the loop, then Maple claimed to allocate 80 million words of memory and it was much slower. Tracking memory through the OS (top, with Linux) showed that Maple's mserver actually consumed 125MB of resident memory. You could try the above code, with Digits > evalhf(Digits) , and see if it leaks or runs out of memory while running through all the 9900 loop iterations.

acer

Using normal worked for me, in Maple 11.

> expr := (r*(r+(1/2)*s))*(sin(t)/(r+(1/\
> 2)*s)+(diff(y(x), x))/(r^2+(1/2)*r*s));
                                     /            d         \
                                     |            -- y(x)   |
                                     |sin(t)      dx        |
                 expr := r (r + s/2) |------- + ------------|
                                     |r + s/2    2          |
                                     \          r  + 1/2 r s/

> normal(expr);
                                        /d      \
                             sin(t) r + |-- y(x)|
                                        \dx     /

The simplify command also did the job. It neither of those work for you, you might want to check that you're using Maple's 2D Math's implicit multiplication (ie. no `*`) correctly.


acer

We might be able to show how it can be done, but more detail would likely be necessary for that.

For example, what is Y1? Is it the water liquid mole fraction at stage 1? And what is T1, is it the temperature at stage 1?

The key bit might be to know the full relationship between T1 and {X1,Y1,X0,Y0} with all the governing equations laid out.

As an aside, it may be unlikely that there's a chemist reading here who might have worked on this particular problem. But if the mathematical problem can be abstracted away from the discipline (here, chemistry) then the chances are quite a bit higher that the mathematical Maple programmers here might be able to show the way. So it might be better to try to remove as much of the chemistry from the description as possible, and try to couch the problem more generically. Eg. additional equations must hold -- always, perhaps, and hence at each time step in particular -- and they form an implicit relationship between some of the dependent variables or their derivatives, etc.

acer

I guess that you are already using the ':-output'=':-solutionmodule' option to NonlinearFit?

If you have residualmeansquare assigned to something, then have you tried accessing it as eq_fit:-Results("residualmeansquare") ?

Otherwise, could you post the example? It can be difficult to diagnose the problem without a complete piece of code to reproduce the behaviour.

acer

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

First 319 320 321 322 323 324 325 Last Page 321 of 336