Carl Love

Carl Love

28025 Reputation

25 Badges

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

MaplePrimes Activity


These are answers submitted by Carl Love

TrueIndices:= (N,F)-> [seq(seq(`if`(F(i,j), [i,j], NULL), i= 0..N), j= 0..N)];

fulfills the requirements of your exercise, although I'd still prefer to make it a proc so that I could declare i and j local.

series(sin(x)^cos(x), x= Pi) indicates that it does not converge. There is a (x - Pi)^(-1) term.

Assuming that p(x), q(x), r(x), a, b, alpha, and beta are suitably defined, then it is done like this:

restart:
p:= x-> x+2:  q:= x-> x:  r:= x-> 1:
a:= 0:  b:= 1:  alpha:= 0:  beta:= 1:
dsolve({diff(y(x),x$2)= p(x)*diff(y(x),x)+q(x)*y(x)+r(x), y(a)=alpha, y(b)=beta}, y(x));

The symbolic answer to this example is too long to copy here.

If you want a numeric, rather than symbolic, solution, then just change the last line to

dsolve({diff(y(x),x$2)= p(x)*diff(y(x),x)+q(x)*y(x)+r(x), y(a)=alpha, y(b)=beta}, y(x), numeric);

 

To do what you want, you would need to know the correspondence between the lettered column names in Excel and the numbered columns in the Maple Matrix. Once you know that, the rest is absolutely trivial. So, do you know what numbers corresponds to columns A and AB?

Let's say A=1 and AB=27. Then do

MediaRighe(M[1..15, 1..27]);
MediaColonne(M[1..15, 1..27]);

Note that your procedures duplicate the work done by the Statistics:-Mean command. Your call to MediaColonne (to get the columnwise means) could be replaced by

M:= subs("-"= undefined, M):
Statistics:-Mean(M[1..15, 1..27], ignore);

And, after that, your call to MediaRighe (to get the rowwise means) could be replaced by

Statistics:-Mean(M[1..15, 1..27]^%T, ignore);

The same is true for all the other dataset summary statistics, such as StandardDeviation.

 

To get the slope, you extract the points from the plot, then you use linear least-squares regression on the logarithms of the points. Here's an extended example:

 

restart:

Digits:= 15:

f:= (t,y)-> I*y:

Exact:= unapply(rhs(dsolve({diff(y(t),t) = f(t,y(t)), y(0)= -I})), t);

(1)

Mid:= proc(h::{positive, realcons})
local
     y0:= -I, #Initial condition
     t,
     y1:= evalf(y0+h/2*(f(0,y0)+f(h, y0+h*f(0,y0)))), #one step of Heun
     y2
;
     for t from h by h to 1 do
          y2:= y0 + 2*h*f(t,y1);
          y0:= y1;
          y1:= y2
     end do;
     #Return the error
     evalf(abs(Exact(t-h)-y0))
end proc:
          

plots:-loglogplot(Mid, 1e-5..1e-1, labels= [h, `error`]);

 

Extract points data matrix from the plot.

XY:= op([1,1], %);

(2)

To get the slope of the line, perform linear least-squares regression on the logarithms of the points data.

Statistics:-LinearFit([1,x], ln~(XY), [x]);

(3)

So the slope is the coefficient of x.

 

Download Midpoint_error.mw

The unusual form of the solution is due to gamma's special role as a constant. Although I don't know the exact reason for the mysterious appearance of I in the solution, I do know that it can be avoided by declaring gamma as local:

restart;
local gamma;
gamma = 1/sqrt(1-beta^2);
solve(%, beta
);

It is very difficult for me to read the MapleMath "font" of this website. It would best if you didn't use it and just used plaintext instead. Too bad that we don't have MathJax (a LaTeX subset used on many websites). I think that what you've entered is the arclength formula:

int(sqrt(1+diff(f(x), x)^2), x= 0..17)

Is that correct?

I presume that your function f(x) is real valued. Is that correct? Then the imaginary part is spurious, perhaps the result of evaluating an elliptic function in the symbolic antiderivative.

You'll very likely get faster, more accurate, and entirely real results if you use numeric integration:

evalf(Int(sqrt(1+diff(f(x), x)^2), x= 0..17));

Note that Int is with an uppercase I for numeric integration.

If the total number of index values (i.e., mul(k_||j+1, j= 1..n)) is small enough that the Cartesian product can easily fit in memory (say, less than 10-100 million), then a non-iterated solution is probably much faster. Here's one:

CartProdSeq := proc(L::seq(list))
local Seq,i,j;
option `Copyright (C) 2007, Joseph Riel. All rights reserved.`;
     eval([subs(Seq= seq, foldl(Seq
                              , [cat(i,1..nargs)]
                              , seq(cat(i,j)=L[j],j=nargs..1,-1)
                             ))]);
end proc:

#K:= [k1, ..., kn]:
K:= [2,3,4]:
add(f(i[]), i in CartProdSeq(seq([$0..kj], kj in K)));

While the example seems silly to me, the output seems no more jumbled than any other polynomial in Maple. Polynomials are not generally sorted by powers in Maple.

Let expr be your arbitrarily complicated expression. Then

evalindets(expr, realcons, evalf[30]);

Change float to float[8].

Your code can be sped up significantly by using evalhf. It may be sufficiently faster to just replace evalf with evalhf. But it may be made even faster than that by filling the entire Array with four for loops inside a single call to evalhf.

It's a bit tricky because the number of entries in each row is not the same. We get around this by using an internal invocation of display so that the Vector constructor < > only sees one entry in the second row. Specifically, it can be done like this:

plots:-display(< H, plots:-display(< P | Q >) >);

The < H, plots:-display(< P | Q >) > specifies a column Vector whose first row is H and whose second row is the result of another display of a Vector. The < P | Q > specifies a row Vector whose entries are P and Q.

Here's the worksheet:


restart:

H:= plots:-animate(plot3d, [x-k*y+1, x=-10..10, y=-10..10], k=-10..0, frames=4):

P:= plots:-animate(plot, [sin(x+t), x=-Pi..Pi], t=-Pi..Pi, frames=8):

Q:= plots:-animate(plot, [cos(x+t), x=-Pi..Pi], t=-Pi..Pi, frames=8):

 

plots:-display(< H, plots:-display(< P | Q >) >);

 

Download array_animation.mw

 

Sorry, it seems that MaplePrimes will not allow me to upload an array of animations. But it is in the worksheet.

Here's a hack to do it with no mouse and no GUI. This involves using a custom operator &= instead of the usual :=.

restart:
`&=`:= proc(f::symbol, expr::uneval)
local x:= eval(expr);
     print(op(1, subs(_f= f, _x= x, proc(_f:= expr=_x) end proc)));
     assign(f,x)
end proc:

a:= 10:  b:= 3:
f&= (a+b^2);


There are two caveats to using this: The expression to the left of &= must be a simple symbol rather than an indexed name (f[1]) or a function (f(1)). And the expression to the right of &= must be in parentheses (unless it only has one part, such as a function call), because it's impossible to change an operator's precedence (see ?precedence ).

Let me know if that's satisfactory because I have a few other ideas if it's not.

See ?printlevel . The default value of the environment variable printlevel is 1. Each level of loop nesting (or if statement nesting) adds one level. So, if printlevel is set to 2 or more, then you will see what's happening inside your double loop.

Local variables and parameters evaluate differently from global variables. Quoting from ?eval (seventh paragraph of Description):

 The default evaluation rules in Maple are full evaluation for global variables, and one-level evaluation for local variables and parameters.

When your eq2 (in procedure f) is evaluated to one level, that leaves its unevaluated. So, the cure is to use eval(eq2), which forces a full recursive evaluation.

restart;
f:= proc()
local
     x,y,
     eq1:= 5+3*x=0,
     eq2:= 2+7*x-3*y-5*x*y=0
;
     x:= solve(eq1,x):
     y:= solve(eval(eq2),y):
lprint('x'=x,'y'=y);
end proc:

f();
x = -5/3, y = 29/16

To address your second question, about ||: The || operator simply was not intended to concatenate arbitrary expressions. It was intended as a convenient way to create indexed symbols. What you are trying to achieve can be easily done with the formatted print statement, printf:

printf("x is %a", -5/3);

First 316 317 318 319 320 321 322 Last Page 318 of 395