You have a space after  solve  command, which Maple in 2d math interprets as multiplication. Remove it and everything will work:

solve(arctan((2*x^2-1)/(2*x^2+1)) = 0, x);

                                         (1/2)*sqrt(2), -(1/2)*sqrt(2)

A little shorter:

e1_1:=-gamma*r*theta/(w*beta*(theta-1));
subs(1/(theta-1)=-1/(1-theta), e1_1);

or even shorter:

subs(theta-1=-``(1-theta), e1_1);

Edit.

A:=[<5,5,5>, <1,2,3>, <-5,1,2>];
LinearAlgebra:-Basis(A);

You can calculate this double integral by writing it as an iterated one:

restart;
with(Student[Calculus1]):
a := Int(Int(exp(cos(x)), y = 0 .. sin(x)), x = 0 .. (1/2)*Pi) ;
ShowSolution(a);
                

Or

sol:=dsolve(diff(y(x),x)= x/(sqrt(x^2-16))*1/(2*y(x)), implicit);
applyop(simplify, [1, 2], %);
solve(%, y(x));

or

sol:=dsolve(diff(y(x),x)= x/(sqrt(x^2-16))*1/(2*y(x)), implicit);
evalindets(%, `*`, simplify);
solve(%, y(x));

or

sol:=dsolve(diff(y(x),x)= x/(sqrt(x^2-16))*1/(2*y(x)), implicit);
solve(map(simplify, lhs(%)), y(x));
 

The procedure classifies the terms of an arbitrary multivariate polynomial with numeric coefficients according to the indicated characteristics:

ClassTerms:=proc(P)
local S, C, t, T, Terms, Tp, Tn, Tind;
S:=indets(P);
C:=[coeffs(P,S,'t')];
T:=[t];
Terms:=convert(zip(`[]`,C,T),set);
Tp:=select(p->p[1]>0 and `and`(seq(type(degree(p[2],s),even), s=S)), Terms);
Tn:=select(p->p[1]<0 and `and`(seq(type(degree(p[2],s),even), s=S)), Terms minus Tp);
Tind:=Terms minus `union`(Tp,Tn);
print(P_positive={seq(p[1]*p[2],p=Tp)});
print(P_negative={seq(p[1]*p[2],p=Tn)});
print(P_indeterminate={seq(p[1]*p[2],p=Tind)});
end proc:


Example of use:

P:=randpoly([x,y,z], dense, degree = 4);

ClassTerms(P);

ClassTerms.mw

If the function is linear, then you can delete this option in your code or do not use it at all if you make a linearity check as follows:

type(a*x+b*y+c, linear(x, y));
                              true 

type(a*x^2+b*y+c, linear(x, y));
                              false

Do

labels=[z,`&Delta;z`/`&Delta;t`]
 

See the toy example:

f:=x[1]->x[1]^2;

whattype(x[1]);

                       indexed

A workaround:

f:=x1 -> x1^2;  # Or
f:=x__1 -> x__1^2;

Two triangles  ABC  and  A'B'C'  are similar. The first triangle can be obtained from the second one by successive transformations: rotation, translation and stretching. Animation of this:
restart;
with(plots): with(plottools):
T1:=display(polygon([[-6,0],[32/3,0],[0,8]],color=cyan), textplot([[-6,0,A],[32/3,0,C],[0,8,B]], font=[times,bold,18], align=right), scaling=constrained, axes=none):
T2:=display(polygon([[14,0],[17,0],[17,-4]],color=yellow), textplot([[14,0,"A'"],[17,0,"B'"],[17,-4,"C'"]], font=[times,bold,18], align=right), scaling=constrained, axes=none):
f:=(p,phi)->rotate(p,phi,[14,0]):
An1:=animate(display,[f(T2,phi)],phi=0..arctan(4/3), background=T1):
T3:=f(T2,arctan(4/3)):
g:=(p,h)->translate(p,h,0):
An2:=animate(display,[g(T3,h)],h=0..-15.8, background=T1):
T4:=g(T3,-15.8,0):
h:=(p,k)->homothety(p,k):
An3:=animate(display,['h'(T4,k)],k=1..10/3, background=T1):
display([An1,An2,An3], insequence=true, size=[800,800]);

Example:

sol:=RootFinding:-Isolate(x^4-3*x^2+2);
eval(x, sol[1]);  # The first root

In your for loop  m  should be a numeric not a symbol.
In addition, if you want to refer to the entries of your matrix 2 by 2 in a for loop, then use a nested loop:
...
for i from 1 to 2 do
for j from 1 to 2 do
...

To find positive solutions to this system with some additional conditions, you do not need to solve it several times in a loop. Maple actually everything has already solved  for your assumption  assume(0<d, d<1). You must only substitute the desired values of the parameter _Z1 into the found formula for d. See

restart;
assume(0<d, d<1):
assume(-0.01<a, a<0):
sys:={Re((-80*Pi*I*a/((a+1)^3))*exp(4*Pi*I*d)) = -0.4, Im((-80*Pi*I*a/((a+1)^3))*exp(4*Pi*I*d)) = 0.8}:
sol:=solve(sys, {a,d}, useassumptions = true,AllSolutions=true);
about(_Z1);
seq(eval(sol[2],_Z1=z), z in [0,1]);

We found the first and second positive values  d  satisfying conditions  d>0  and  d<1

Maple correctly decides your system. See
restart;
assume(d::real, d>0):
assume(a::real, -0.01 < a, a < 0):
sys:={-800*Pi*a*cos(6.557377048*Pi*(3.470797713+d))/(a+1)^3 = -.9396060697, 800*Pi*a*sin(6.557377048*Pi*(3.470797713+d))/(a+1)^3 = -.3238482794};
solve(sys, {a,d}, useassumptions=true, AllSolutions=true);
about(_Z1);

We see that the system has an infinite number of solutions that satisfy the given conditions, and the smallest positive  d  will be for  _Z1=11

Visualization:

sol:=solve(sys, {a,d}, useassumptions=true, AllSolutions=true):
eval(rhs(sol[2]), _Z1=11);  # The smallest positive  d
plots:-implicitplot(convert(sys, list), a=-0.01..0, d=0..0.03, color=[red,blue], gridrefine=5); 

Edit.

restart;
convert(D(s)(t), diff);