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:

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 x 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);

The zip command is used to apply a two-parameter procedure to two matched arrays (or lists) of arguments.

restart:
B:= (x,y)-> fsolve(y = (1-exp(-x*b))/(1-exp(-50*b)), b):
X:= #Some array of x values
Y:= #Some array of corresponding y values
zip(B, X, Y);

Bisection is a method of last resort. Why do you want to use it? You can quickly get any number of roots of f(p) with RootFinding:-NextZero. After making Preben's unapply correction, do this:

N:= 100:  #Number of roots desired.
Root:= Array(1..N):
Root[1]:= -RootFinding:-NextZero(p-> f(-p), 0):
for n to N-1 do
     Root[n+1]:= -RootFinding:-NextZero(p-> f(-p), -Root[n])
end do:

The minus signs are because it seems that you want the roots in the negative direction.

All that you need is

series(y(x+h) - y(x-h), h=0, 3);

NumberToList:= proc(x::realcons)
uses S= StringTools;
     parse~(S:-Explode(S:-Select(S:-IsDigit, sprintf("%a", x))))
end proc:

Change the last line of your last for loop from

theta[ii](tau):=rhs(%);

to

theta[ii](tau):= value(rhs(%));

The value command makes Maple perform the integrations.

Since the desired simplification is not generally true for all bases, you have to use the symbolic option to simplify:

simplify(%, symbolic);

Using Maple 17, when I use solve followed by evalf, I get

Try doing a restart, then try the solve and evalf again.

A real answer can be obtained with fsolve or RootFinding:-NextZero in addition to the method that you mentioned.

KernelDensityPlot is a command in the Statistics package that will output the same plot in the normal inline manner:

Statistics:-KernelDensityPlot(Cdata);

AB,C:= LinearAlgebra:-GenerateMatrix(
     [F[k] $ k= 0..3],
     [u[k,n] $ k= 0..3, u[k,n-1] $ k= 0..3],
     augmented= false
);
A:= AB[.., 1..4];  B:= AB[.., 5..8];

1. In most cases, you reduce the error by increasing the value of Digits, which has the default value 10.

Digits:= 15:
fsolve(eqs2);

2. To get an exact solution to a decimal problem (if it's possible), use convert(..., rational):

solve(convert(eqs2, rational));