Do you want to create a single 3x3 matrix each element of which is a 3x3 matrix? Or do you want to create a single 9x9 matrix each 3x3 block of which is one of your original matrices?

The first thing is done like this:

A:= Matrix(3, 3, (i,j)-> LinearAlgebra:-RandomMatrix(3,3));

Using those as the starting 9 matrices, the second thing is done like this:

B:= <seq(`<|>`(seq(A[i,j], j= 1..3)), i= 1..3)>;

I will describe the way that works for me, which is not the only way, it's just the way that I know.

Pull down the Tools menu and select Options. Go to the second tab, Display. The first item is "Input display". Its value should be "Maple Notation". Change it if it isn't. Now go to the third tab, Interface. The fourth item is "Default format for new worksheets". The value should be "Worksheet". Change it if it isn't. At the bottom of the window, click Apply Globally.

Now open a new worksheet by typing Ctrl-N. In the "Select Math Engine" dialog, select the first item "New". Now look at your toolbar. Do you see in the picture that you inserted in your Question where in the lowest toolbar the leftmost item is a pull-down labelled "P Normal"? In your new worksheet, that should be "C Maple Input". Change it (and let me know!) if it isn't.

Now paste your code to the position that your cursor is now in. If you press Enter, the code should execute. Please let me know how it goes.

From now on, all your worksheets should open being ready to except Maple Input.

expr:= x^2 + x + 2:
subs(x= x+1, expr);

I'd recommend that you not copy output to input by cut-and-paste. In addition to the problems that you've already noticed, it tends to make your code less readable, less modifiable, and less able to be executed by someone else.

#Works like 'select' except that it returns a list of the indices.
SelectIndices:= proc(f::appliable, L::{list,set,Vector})
local k;
     [seq(`if`(f(L[k], _rest), k, [][]),  k= 1..numelems(L))]
end proc:

     
a:= Vector([2,3,4,5]):
SelectIndices(`>=`, a, 3);

Edit: Answer updated to use _rest based on the Reply by Joe Riel below.

If your data is still accessible via the %, %%, or %%% operators (i.e., if it is one the three most-recent results), then it can't be garbage collected. Observe:

Download bytesalloc.mw

@tomet The norm of a vector is nonnegative. So the plot will have all the "mountains" proceeding along the positive z-axis. Maple is making no mistake in this. What you want is either the Divergence of the vector field or perhaps the Norm times the signum of the Divergence.

plot3d(VectorCalculus:-Divergence(E), x= -10..10, y= -10..10, view= -2*10^10..2*10^10);

I got quick results (a few seconds) using method= rosenbrock, and reasonable results (20-30 seconds) with method= mebdfi. The default method= rkf45 takes a very long time; I didn't let it finish. All other methods, including all the classical methods, produce garbage results, oddly enough.

Sol:= dsolve(
     #All initial conditions set to 0.
     {eq||(1..5), seq(q[k](0)=0, k= 1..5), seq(D(q[k])(0)=0, k= 1..3)},
     numeric, maxfun= 0, method= rosenbrock
):
plots:-odeplot(Sol, [seq([T,q[k](T)], k= 1..5)], T= 0..10);

plots:-odeplot(Sol, [seq([T,q[k](T)], k= 2..5)], T= 0..0.1);

A:= Matrix([[a,b,c], [d,e,f], [g,h,i]]);

#90-degree clockwise rotation:
ArrayTools:-FlipDimension(A,1)^%T;

#Flip over horizontal axis:
ArrayTools:-FlipDimension(A,1);

All eight can be built from those two.

add(nops(select(`=`, convert(n, base, 10), 5)), n= 1..1000);

     300

Go to http://www.mapleprimes.com/users/casperyc/questions

This can be gotten around by using delay-evaluation quotes:

seq('`if`'(a[i] < b[i], c[i], d[i]), i= 1..10);

seq has special evaluation rules, and `if` has what I call very special evaluation rules. I don't think it's a bug; it turns out to be useful sometimes. Here's an example of the rule:

restart:
A:= a:  B:= 2:
F:= _a-> eval(`if`(A,B,C), a= _a):
F(true);
                               B
eval(B);
                               2

There was a discussion of this in the MaplePrimes archives, but all the Answers and Replies disappeared!!!

There is a discontinuity plainly visible on the south polar contour of the sample plot that you posted. This extends from the pole along a meridian. So, this is not the fault of my code.

Areas of "no data" are bright red. These include a large, nearly circular area at the north pole and the evenly spaced streaks roughly perpendicular to the equator. It should be possible to use some interpolation to get rid of the streaks.

Here are precisely the steps that I used to import your data:

1. Using FireFox, I clicked on the hyperlink to your data file that you posted.
2. Using FireFox, I clicked Edit -> Select All and used Control-C to copy.
3. I opened a new Notepad document on my desktop.
4. I clicked in the blank document and used Control-V to paste.
5. I positioned my cursor at the beginning of the first line of actual data.
6. I hit backspace to delete the first (blank) line of the document.
7. I saved the document to my desktop as globedata.txt.

Then I executed the following worksheet.

Download globeplot.mw

The following procedure uses foldl to build the appropriate nested add command, and then evaluates it. Thus, this procedure maintains all of the efficiency associated with add.

CartProdAdd:= proc(f::{name,appliable}, L::seq(`..`))
local A, _i, V:= _i||(2..nargs);
     eval(foldl(A, f(V), V=~ L), A= add)
end proc:

Example of its use:

CartProdAdd(f, 1..2, 3..4);

The contour plot on a sphere (or on any plottable surface) can be achieved if the data can be arranged in a grid. That is, you need to have an evenly spaced set T of values of theta and an evenly spaced set P of values of phi such that there is a value of D for every element of the cartesian product T x P.  And T x P must cover the globe. It may require a significant amount of interpolation for you to achieve this grid.

Let's suppose that the grid is stored in a matrix G:= Matrix(nops(T),nops(P)) each entry of which is a value of D. We do a linear rescaling of G to the range 0..1. Then we take a basic plot of the unit sphere in spherical coordinates with a specific coloring and replace the color matrix by G, like this:


# Linear rescaling to 0..1:
map[inplace](evalf, G):  M:= max(G):   m:= min(G):
map[inplace](x-> (x-m)/(M-m), G):

subsop(
     [1,2]= COLOR(HUE, ListTools:-Flatten(convert(G, listlist))[]),
     plot3d(
          1, theta= 0..2*Pi, phi= 0..Pi,
          grid= [nops(T),nops(P)], coords= spherical, color= theta,
          style= patchnogrid, lightmodel= none
     )
);

If you have trouble with the interpolation, then post your data.

Edit: Text and code updated with the linear rescaling idea due to a bug discovered by Markiyan Hirnyk in the Reply "Yet output" below.