Your file titles indicate that you actually already know what to do: Use cat instead of ||. The parse is to convert the strings to names.

fe11_1:= p^2;
fe11_2:= p^3;
fle11_1:= p^4;
fle11_2:= p;
for appm in ["fe", "fle"] do 
   for i in ["11_1", "11_2"] do 
      f:= parse(cat(appm, i)); 
      sp:= evalf(int(f, p = 1 .. 11/3))+evalf(int(f, p = 11/3 .. 19/3))+evalf(int(f, p = 19/3 .. 9));
      print(f, sp) 
   end do
end do;

Why do you need the Excel file? You can do the entire plot with Maple:

plot3d(
   sin(A)*exp(1-L/1000), A= 0..2*Pi, L= 0..1000, coords= cylindrical,
   grid= floor~([360/1, 1000/5]), thickness= 0
);

Note that coordinate A is measured in radians rather than degrees.

It is not allowed to assign to a local variable unless the local variable evaluates to something that can be assigned to, namely a name. Do you want to change the value of b locally? Then you need a local variable. Do you want to change the value of some variable that is passed in for b? Then you need something like this:

P:= proc(a::integer, b::evaln(integer))
   if eval(b) = 2 then b:= 3 end if
end proc:
b:= 2:
P(a,b):
b;
                               3
b:= 2.5:
P(a,b);
Error, invalid input: P expects its 2nd argument, b, to be of type evaln(integer), but received b := 2.5
b:= 4:
P(a,b);
b;
                               4

Use of this technique is considered a bad programming practice called side effects.

Use fsolve, like this:

sys:= diff(ca(t), t) = -3.600000000*10^20*exp(-15098.13790/(340-20*ca(t)))*ca(t), ca(0) = 2:
sol:= dsolve({sys}, numeric):
fsolve(T-> eval(ca(t), sol(T)) - 0.2);

3.89994719649399

A variation on the above that a lot of people seem to prefer, but which makes little difference to me, is

sol:= dsolve({sys}, numeric, output= listprocedure):
fsolve(eval(ca(t), sol)(t) = 0.2);

Acer's answers are fine, but here is one to add to his list:

R:= Matrix((1,4), (i,j)-> x-> j*x);

This is what I would use because I prefer symbol infix operators over word prefix operators.

Your ODEs refer to cot(theta). This involves division by 0 at theta=0, which is your initial point. Hence, the system is "initially singular", as the error message says.

Tom Leslie's example shows that a type declaration of M::Matrix(datatype= ...) checks only the explicitly declared datatype of M. That's probably not what you want. You probably want to check that the contents of M have a specific data type. For that, you need to use coercion (see ?coerce), like this:

P:= proc(M::~Matrix(square, datatype= rational))
   #...
end proc:

Notice the tilde in front of the word Matrix. That indicates that coercion is desired.

square is a "top-level" property of Matrices; it isn't a shape. As I used it above, square will enforce that the input Matrix be square.

B:= Array((0..n-1)$2, [ListTools:-LengthSplit(L,n)]);

Edit: Corrected error pointed out by Kitonum.

Your symbolic dsolve solution contains a complicated RootOf expression whose principal value has nonzero imaginary part for the values of t that you're plotting. You can avoid this by using a numeric dsolve with parameters, like this:

sol:= dsolve(
   {eq1b, eq2b, r(0) = r__0, theta(0) = theta__0}, 
   parameters= [r__0, theta__0], numeric
):
sol(parameters= [theta__0= Pi/4, r__0= 0.1]):
plot1:= plots:-odeplot(
   sol, r(t)*~[cos,sin](theta(t)), t= 0..10,
   axiscoordinates= polar, color= red, labels= [``$2], tickmarks= [4,3]
);

Note the explicit conversion from polar coordinates: r(t)*~[cos,sin](theta(t)).

If the matrix's entries can be expressed as a function of the indices, then you can use an initializer function like this:

U:= Matrix((5,5), (i,j)-> `if`(i=j, 1, 0));

The command solve doesn't understand domain restrictions; however, it does understand inequalities if they're put with the equations in the first argument:

solve({eq4, eq5, theta > 0}, {Fh1, theta});

The numerical solver fsolve does respect domain restrictions entered as ranges:

fsolve({eq4, eq5}, {Fh1, theta}, {theta= 0..infinity});

The following code shows all the equivalence classes of adjacency matrices.This output is much neater looking in a Maple worksheet. With a small modification, it could show a drawing of a representative of each class.

n:= 3: N:= n*(n-1):
<ListTools:-Categorize(
   (M1,M2)-> GraphTheory:-IsIsomorphic(GraphTheory:-Digraph~([M1,M2], n)[]),
   map(
      L-> Matrix((n$2), [0, seq([L[n*k+1..n*(k+1)][], 0][], k= 0..n-2)]), 
      combinat:-permute([0$N, 1$N], N)
   )
)>; 

Your model function is not differentiable at t=0, yet it appears that t=0 is one of the data points. Try removing this one data point. Of course, I can't test this without having your actual data.

If all of the parameters being estimated are coefficients, then you can and should use linear regression. The model function doesn't need to be a linear function of the data variables in order to use linear regression; it merely needs to be a linear function of the parameters being estimated.

If A is the matrix, then you can use convert(A^+, listlist). The first argument, A^+, means the transpose of A. Using this doesn't require foreknowledge of the number of columns of A.