Carl Love

Carl Love

28055 Reputation

25 Badges

12 years, 358 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

MaplePrimes Activity


These are answers submitted by Carl Love

@mapleq2013 You asked some important questions. I'll answer them one by one.

Dt:= unapply(LinearAlgebra:-Determinant(C), (y,x)):

So the unapply here is to convert the expression in x and y into a function, is my understanding correct?

Yes.

If correct, why do we have to do this?

It isn't totally necessary. What's necessary is to compute the determinant outside the call to fsolve or else the determinant would be recomputed for every point (x,y). I could do

Dt:= LinearAlgebra:-Determinant(C):

and then in the call to fsolve replace Dt(y,x) with a call to eval. The unapply syntax is less convoluted.

solD||n:= subs(_n= n, (y,x)-> Re(fsolve(Dt(y,x), w, complex)[_n]))

I don't understand why there is a _n=n

The underscore can be used in a variable name just like any letter.

is it just to avoid using n=n?

That wouldn't work because the left n would evaluate to a number and thus wouldn't match the n (instead of _n) to the right of the ->. This latter n would still be n, because nothing on the right side of an arrow evaluates. But subs('n'= n, ...) would work.

So can I use k=n alright?

You could. However there's a convention that global variables beginning with underscore are never assigned values (except those beginning _Env, called environment variables). What I need in this situation is a global variable without an assigned value. So, by convention, these dummy variables in procedure definitions usually have names beginning with underscore.

Also you mentioned using a procedure, and the usually way I understand procedures is like:

function:= proc(parameter):: returnType
                         procedure content
   end proc

But the procedure you defined is not of the above structure, so I'm a little confused.

The three most common ways to define a procedure are with the arrow ->, with proc (as you describe), and with unapply. (There are also some much more obscure ways that are only used in low-level programming.) Procedures defined with the arrow are often colloquially called functions or operators; however, those words also have more general meanings. Only the word procedure precisely captures the meaning. After they've been created, there's no practical difference between procedures created by any of these methods.

The subs command, why do we need it?

The dummy variable _n in the procedure needs to have an integer value from 1 to 6 substituted for it.

Can I simply write 

solD||n:= (y,x)-> Re(fsolve(Dt(y,x), w, complex)[n]) 

It doesn't work because the second n, being part of a procedure definiton (because it's on the right side of an arrow), won't evaluate at the time of definition even though the n, being a for loop index, has a numeric value at the time.

And how is subs different from eval.

The difference is often subtle, and the two are often interchangeable. I recommend that you read their help pages for more details. But in this particular situation there's no subtlety: Only subs can be used on the outside of a procedure definition to replace a dummy variable inside the definition. This important and fundamental usage of subs is unfortunately not documented on its help page.

 

It isn't true that -> is a "short form" of unapply. Rather, x-> f1 is a short form of proc(x) f1 end proc. If you want a short form for unapply, you could do this:

`&->`:= (x,f)-> unapply(f,x):

Now &-> is short for unapply, and x &-> f1 will return x -> 1+x.

If you want to force the evaluation of f1 in the expression x-> f1, you could do

subs('f1'= f1, x-> f1);

This is just an example, not a recommended coding practice. Indeed, it won't work inside a larger procedure.

You need to use fsolve instead of solve. Since the equation being solved is a polynomial, fsolve (with option complex) will always give six solutions. The call to fsolve needs to be wrapped in a procedure because numeric values for x and y won't be available until they're supplied by plot3d. I made six procedures---one for each of the six eigenvalues.

At the point where you define the matrix C, get rid of the fnormal. Why reduce the precision at this point? Then include these lines:

Dt:= unapply(LinearAlgebra:-Determinant(C), (y,x)):
for n to 6 do
      solD||n:= subs(_n= n, (y,x)-> Re(fsolve(Dt(y,x), w, complex)[_n]))
end do:
plot3d(solD1, -Pi..Pi, -Pi..Pi);

Likewise, you can plot any of the six.

This code can be made faster: The same fsolve is being called six times for each pair (x,y). I can make it faster if you need that.

While waiting for your upload, here's an example of using NonlinearFit with a procedure that has an if statement:

P||(1..2):= 'randpoly(x, degree= 1)' $ 2:
F:= piecewise(x < 0, P1, P2);

X:= Vector(RandomTools:-Generate(list(float(range= -1..1, method= uniform), 99))):
Y:= unapply(F,x)~(X):
M:= proc(x, m1, m2, b1, b2)
     if x < 0 then
          m1*x+b1
     else
          m2*x+b2
     end if
end proc:
Statistics:-NonlinearFit(M, X, Y, output= parametervalues);

We see that the fit is near perfect. (But, of course, every point was exactly computed from the original function.)

It's very difficult to work with assigned names qua names in a systemmatic way because, in most contexts, the names spontaneously evaluate to their assigned values. Unevaluation quotes can be used to force individual names to be viewed as names, but it's difficult to apply this to whole batches of names. So, I wrote this module to help.

In addition to what's covered by MacDude's and Tom Leslie's Answers, this module lets you

  1. categorize the names into batches
  2. guard against using protected or already assigned names
  3. unassign batches of names
  4. use the parsed value of a string rather than the string itself as the name's value
  5. use names that are just names

This module may appear simple and short---only 21 lines without the comments---but many hours of thought and experimentation went into its creation.

 

restart:

Module SelfName : The input can be any expression e , which is converted to a name n. If e  is a string, then n won't contain the double quotes. If n is protected or already assigned, then an error is returned. Then n is given the attribute SelfNamed and any additional (optional) attribute that the user has passed in. Thus, these special names can be found and categorized. If e  is a string that can be parsed, then that parsed value is assigned to n; otherwise e  itself is assigned to n. The return value is n. The exports are procedures for listing and unassigning the variables.

SelfName:= module()
option `by Carl Love, 2015-Aug-25`;
local  
     # Just like built-in `assign`,
     # but it returns the left-hand side of the assignment.
     Assign:= proc(x,a)  assign(x,a);  'x'  end proc,
  
     (* This procedure, just like built-in `assign`, necessarily has side
        effects. Thus if this procedure is mapped over a list and errors on the
        9th entry, the first 8 entries will still have been processed. *)

     ModuleApply:= proc(e, {tag:= ()})
     local n:= nprintf(`if`(e::string, "%s", "%a"), e);
          if ormap(assigned, <n>) then
               error "name already assigned: %1", n
          end if;
          if n::protected then  error "name is protected: %1", n  end if;
          setattribute(n, ':-SelfNamed', tag);
          if e::string then try return Assign(n,parse(e)) catch: end try end if;
          Assign(n,e)
     end proc
;
export
     # Returns a list of all SelfNames, or all in a specifed batch.
     SelfNames:= proc({tag:= ()})
          select(type, [anames('alluser'), unames()], attributed(':-SelfNamed', tag))
     end proc,

     # Unassign a single name and remove its attributes.
     # Returns the name.
     UnSelfName:= proc(n::uneval)  unassign(n);  setattribute(n);  n  end proc,

     # Apply UnSelfName to all SelfNames, or to all in a specified batch.
     UnSelfNames:= proc({tag:= ()})  UnSelfName~(SelfNames(':-tag'= tag))  end proc
;
end module:

Usage examples:

FromExcel:= < "1..2", "A..B" >:

SelfName~(FromExcel, tag= Batch1);

Vector(2, {(1) = `1..2`, (2) = `A..B`})

(1)

SelfName:-SelfNames(tag= Batch1);

[`A..B`, `1..2`]

(2)

Batch2:= < 7, y, "11">:

SelfName~(Batch2, tag= "Batch2");

Vector(3, {(1) = `7`, (2) = y, (3) = `11`})

(3)

`11`;

11

(4)

Note that the returned value above is a number, not a name.

SelfName:-SelfNames();

[`A..B`, `11`, `7`, `1..2`, y]

(5)

SelfName:-SelfNames(tag= "Batch2");

[`11`, `7`, y]

(6)

SelfName:-UnSelfNames(tag= "Batch2");

[`11`, `7`, y]

(7)

`11`;

`11`

(8)

Note that the returned value above is just a name.

SelfName:-SelfNames();

[`A..B`, `1..2`]

(9)

anames(user);

`A..B`, Batch2, `1..2`, FromExcel, SelfName

(10)

 

 

Download SelfName.mw

 

The problem with your implementation is that the if statement is evaluated and executed before a value is supplied for n. Therefore, the clause if n = 1 is never executed at a time when n equals 1, so the recursion never ends. In order to delay the evaluation until a value is supplied for n, you can use a procedure. Here's a minimal procedure implementation.

B:= (m,n)-> `if`(n=1, omega(m), add(omega(j)*B(m-j,n-1), j= 1..m-1));

The following is nearly equivalent to and slightly more verbose than the above and very close to your original.

B:= (m,n)-> if n=1 then omega(m) else sum(omega(j)*B(m-j,n-1), j= 1..m-1) end if;

Here's a much more sophisticated implementation which will work for any m and n, even purely symbolic, and which makes efficient choices when either or both of the input variables aren't symbolic.

B:= proc(m, n)
option remember;
local j;
     if not n::posint then 'procname'(args)
     elif n = 1 then omega(m)
     else `if`(m::realcons,add,sum)(omega(j)*thisproc(m-j,n-1), j= 1..m-1)
     end if
end proc:

Note that a procedure can be defined with the arrow -> (see ?operators,functional) or with the keyword proc (see ?procedure).

You may have ruined your chance for a backup by opening the worksheet again; however, the following is worth a try. Go to the File menu, then Recent Documents, then Restore Backup.

If you ever again get the "Kernel connection lost" dialog box:

  1. Acknowledge the dialog box by clicking OK so that it goes away.
  2. Save the worksheet. You will not be able to do any computation in the worksheet.
  3. Go to File -> Close Document.
  4. Despite what the dialog box said, there is no need to restart Maple.
  5. Re-open the worksheet.

If your Maple session is interupted for some other reason:

  1. Start Maple.
  2. Go to File -> Recent Documents -> Restore Backup.

 

For efficiency, Maple's plot tries to use hardware floating point arithmetic (aka "double precision") when it can (via the evalhf command). But double precision can't handle numbers as big as 10^1390. On the other hand, they don't cause an error because they're just treated as infinity. So plot isn't aware that something is amiss. (It will automatically switch to software floats if it detects an error with the hardware floats.) The solution is to force plot to use software floats. There are many ways to do this. I think the easiest is to set Digits higher than 15. So,

Digits:= 16:
plot(
     .1126460951*t-0.7049523744e-1-1.090463096*10^1390*exp(-1.597924898*t) ,
     t=2000..2001
);

This produces a good plot. To maintain efficiency, remember to reset Digits to 15 or lower when you're done.  

If you want to use ArrayInterpolation, the syntax would be

f:= X-> CurveFitting:-ArrayInterpolation(x, y, X):
evalf(Int(f, xmin..xmax));

I can't tell from your post whether this is the syntax that you actually used.

However, this is wasteful and may have a higher-than-necessary error. It's wasteful because it recomputes the (local) interpolant for each value of X. It may have high error because the default interpolation is linear, so you're essentially using the trapezoid rule (except that the x-values aren't necessarily evenly spaced). You can change this with the method option to ArrayInterpolation.

A better approach may be to compute a single Spline and then integrate that.

If your x data are evenly spaced, you can (and probably should) directly apply a Newton-Cotes formula (such as Simpson's Rule, Simpson's 3/8 Rule, Boole's Rule).

I'm just expanding on Axel's Answer here.

Unless you're doing textbook problems in Calculus I or II, I'd avoid using Student:-Calculus1:-Rule. The more general command is IntegrationTools:-Parts.

What's wrong with the output as it's presented by default?

restart:
Eq1 := f(x,g(x,t)) + f(x,y):
diff(Eq1, x);

If you don't like the Ds, then do

convert(%, diff);

First, the f as you've shown it is a table, not a list. It needs to be a list or 1-D rtable (Vector or Array) to sort it. (Kitonum's Answer gets around this issue by tacitly converting the table to a list with seq.) To make it a list, define it as

f:= [[1,3,2], [2,1,3], [1,2,3], [2,3,1], [3,2,1], [3,1,2]]:

Now, you can simply do

sort(f);

This fulfills the requirements that you specified. However, if the lists have different lengths, then the default sort places priority on length over content. I don't like that. I'd like something akin to an alphabetic sort. For that you need a custom sort procedure, like this:

Sort:= proc(A::list, B::list)
local k, nA:= nops(A), nB:= nops(B);
     for k to min(nA,nB) do
          if A[k] <> B[k] then return A[k] < B[k] end if
     end do;
     nA < nB
end proc:

Then do

sort(f, Sort);

 

You have a syntax problem---the absence of multiplication operators. After correcting that, assuming that some of the constants are positive will help Maple to finish the integral.

The first problem: The coefficients C__2 and C__4 need to be followed by the multiplication operator *. So, your function becomes:

C__p:= T -> C__1 + C__2*(C__3/(T*sinh(C__3/T)))^2 + C__4*(C__5/(T*cosh(C__5/T)))^2;

The second problem can be resolved by including an assuming clause:

int(C__p(T), T= T__ref..T__sys)
     assuming T__ref > 0, T__sys > T__ref, C__3 > 0, C__5 > 0;

The output is lengthy, and you may never be able to get it simplified to the form that you suggested. But you may be able to get Maple to show that the two forms are equivalent. If that's what you want to do, let me know.

What you're doing wrong is very simple: Your Y data in the Excel file is on a percentage scale (0  to 100), but your model function is on a 0 to 1 scale. If you divide every value by 100, you'll get a near-perfect fit. You can divide it all by 100 with

Y:= Y /~ 100:

This is the fit that I got with DirectSearch:-DataFit:

[a = .522164088988118, m1 = 22.7042491482519, m2 = 9.51958388629120, s1 = 6.47780713862839, s2 = 3.34678806698068]

Here's my worksheet:

Download data_fit.mw

 (The MaplePrimes editor won't let me display this worksheet inline for some reason.)

Note that I included parameter constraints. This is not possible with Statistics:-Fit. The constraints that I used were

Constraints:= [a >= 0, a <= 1, m1 >= 0, s1 >= 0, m2 >= 0, s2 >= 0, m1 <= 60, m2 <= 60]:

Here's my final plot: The upper red line is the fitted function, the green points are the data, and the lower red line is the derivative of the fitted function, which explicitly shows the bimodality.

Try this:

Statistics:-ColumnGraph(y, datasetlabels= convert(x, list));

A histogram is something completely different: Although it may superficially appear similar to a column graph, it can only be applied to a single set of data. In Maple, a BarChart is just like a ColumnGraph with the columns being horizontal.

First 251 252 253 254 255 256 257 Last Page 253 of 395