There's an undocumented command inner for dot products that won't use the conjugate:

inner(vector([x,y,z]), vector([a,b,c]));

Use plottools:-sphere to plot an actual small sphere at each point.

You must be using Maple's garbage 2-D input mode. You original code probably contained the statement

ans(ind):= [op([1], f1max), roll];

If you had been using Maple's 1-D input, this statement only has one possible interpretation (in context): You are assigning an Array element. With garbage 2-D input, Maple tells you that the statement can be interpretted as either a remember table assignment or a function definition and asks you to choose which. (Note that the correct interpretation isn't even one of the choices; however, remember table assignment is effectively the same as Array assignment for the purpose of this question.) Two weeks ago, you selected "remember table assignment"; today you selected "function definition".

If you use 1-D input (aka Maple Input), you'll avoid this problem, and you won't have to answer a stupid question every time you run your code. 2-D input is garbage: No reasonable language asks the end user to provide an interpretation for the code every time a program is run.

There's a command that does the image conversion and the file writing in one step: plotsetup. You can use it like this:

currentdir(kernelopts(homedir)); #Set directory where file will go.

ExportJPGMyPlot:= proc(MyPlot, Name::string)
     plotsetup(jpeg, plotoutput= cat(Name, ".jpg"));
     Threads:-Sleep(2); #2 secs to synchronize I/O.
     print(MyPlot);
     plotsetup(default);
     Threads:-Sleep(2);
     MyPlot
end proc:

ExportJPGMyPlot(plot(x^2, x= -2..2), "XSqr");

Note that the API (call sequence, arguments, etc.) is completely different from any other Maple command. It works entirely by side effects which will redirect plots that otherwise would've been displayed in the GUI.

The list of supported image file types is documented (minimally) at ?plot,device. The list of image file types that are available on your system is returned by plotsetup(help). Unfortunately, this latter list is much longer than that documented by the former.

The Threads:-Sleep(2) isn't necessarily needed, and the 2 can possibly be reduced; it's system dependent. I think that the file I/O in this situation isn't synchronized with the Maple kernel. Maybe someone else can confirm or provide further clarification of that.

The command for filtering a set is select. The command fsolve will, by default, return all real roots (and only real roots) for a polynomial. So what you want can be done by

select(x-> g(x,a) > 0, {fsolve(f(x,a))});

To display a 2-D Array A as if it were a Matrix, use <A>.

e:= g^((2*(-sigma+k+1))/(-1+sigma))-tau^2:
ex:= op([1,2], e)/2;
thaw(factor(expand(simplify(e, {ex= freeze(ex)}))));

The problem is indeed the assignment A1:= A0, which makes the matrices identically the same.

If you want to work entirely in memory, and you have enough memory to do so, replace A1:= A0 with

A1:= copy(A0).

Think of the matrices as boxes. The statement A1:= copy(A0) creates a new box A1 and copies the contents from A0. From that point forward, A1 and A0 can be independently changed. On the other hand, the assignment A1:= A0 says that A1 and A0 are the same box, so any change to the contents of one will affect the other.

If you don't have enough primary memory to store two copies of the matrix, you can store one on disk, like this:

#Create A0.
#Do some calculations with A0.
save A0, "MyA0.m";  #Save A0 for later use.
A1:= A0:  #Uses a trivial amount of extra memory.
#Do some calculations with A1.

save A1, "MyA1.m";
A1:= 'A1': #Optional: Eliminate connection between A1 and A0.
read "MyA0.m"; #Restore A0

If you're interested in real solutions only, then using solve(P, z, parametric, real) gives a little bit more useful information than Kitonum's fsolve. In partcular, it immediately gives a floating-point approximation to the unique real value of w at which the number of real roots changes from 0 to 1 to 2, and it gives RootOf expressions that can be evalfed for those real roots.

I switched the order of your assignments:

restart:
mu:= c^(a-b)/x:
x:= (t^(1-s)+s^(1-2*s))^a:

Now compare mu with eval(mu, 1).

There are many problems:

Why are you looking for complex solutions in the fsolve? If you're expecting nonreal solutions, what do you mean by plotting them?

Are you expecting to get both roots from the fsolve? One is simply the negative of the other. I think that you intend to just use the positive root.

alpha0^2 is so small compared to the other terms in eq1 that your loop that varies alpha0 only makes the roots different after the 40th decimal place.

alpha0 is not a sequence, so alpha0[i] is nonsense. Replace alpha0[i] with 5*(i-1).

Your usage of pi is not an error, but the constant Pi (with a capital P) is already defined, so there's no need for you to set pi:= 3.14.

In the line f5:=, you have a clause Jy= J*sin (theta). You need to take away the extra space between sin and (theta).

Using LinearSolve(A,b) mod 2 won't work: It means to solve the system over the rational numbers and then take that solution mod 2.

The LinearAlgebra:-Modular subpackage is extremely fast. Use it like this:

macro(LA= LinearAlgebra, LAM= LinearAlgebra:-Modular):
vars := [indets(systems, name)[]]:
Ab:= LA:-GenerateMatrix(systems, vars, augmented):
try
     X:= LAM:-LinearSolve(2, LAM:-Mod(2, Ab, integer[]), 1, inplace= false)
catch "matrix is singular":
end try;

The variables can be automatically generated by using the symbol option to Matrix. Then the equations to solve can be automatically generated by converting a Matrix to a list.

b:= Matrix((5,1), symbol= y):
X:= Matrix((2,1), symbol= x):
eval({b,X}, solve(convert(A.b+B.X, list)));

Maple's Optimization package can solve LPs and NLPs. See ?Optimization, ?Optimization,LPSolve, and ?Optimization,NLPSolve.