This seems to work:
p:=proc(u) local T;
  if type(u,`+`) then T:=[op(u)] else T:=[u] end if;
  evalindets[flat](T,`*`,op)
end proc;
u:=2*t*exp(t) + 2*t;
p(u); #Works OK
##But the following doesn't because of the product inside sin
p(sin(2*t)); #Error
## A revised version freezes functions:
q:=proc(u) local T,F;
  if type(u,`+`) then T:=[op(u)] else T:=[u] end if;
  F:=evalindets(T,function,freeze);
  thaw(evalindets[flat](F,`*`,op))
end proc;
q(u); #OK
q(sin(2*t)); OK

Assuming that your 2-D input can be interpreted to mean the same as this corrected 1-D version

ode := mu__1*E__1(x)-mu__2*(u(x)-d__1*E__1(x))/d__2+diff(gamma1*E__1(x)-gamma2*(u(x)-d__1*E__1(x))/d__2, x) = 0;

isc := E__1(0) = 0;

and assuming that u(x) is to be considered a known function, then the syntax is:

res:=dsolve({isc, ode},E__1(x));

You will notice that the result res contains an inert integral (Int, grey in print). Using value on res won't help since Maple cannot do the integral unless u is known.
Suppose u=sin, i.e. u(x) = sin(x) for all x, then
eval(res,u=sin);
value(%);
# will give you a finished result. But you can get it on a nicer form:
combine(expand(%));
resf:=collect(%,exp,factor);

tau*f(t) is not an operand in expr[2]. I think it is as simple as that.
My mistake: Apparently not quite. See bottom.
Try the similar example with I replaced by 3:
conv1 := (x) -> eval(x, tau*f(t) = f(-t)):     # Using eval
conv2 := (x) -> algsubs(tau*f(t) = f(-t),x):   # Using algsubs
expr := Vector([
     tau*f(t),
   3*tau*f(t)
]);
conv1(expr),
conv2(expr);
##
Instead you can do
conv1 := (x) -> eval(x, f(t) = f(-t)/tau):
#############
The help for eval:
Usually, the left-hand sides of the equations will be variable names being evaluated.  However, if the left-hand side x of an equation is not a simple name, a syntactic substitution may be attempted which has the same limitations as with subs.
Compare:
eval(a*b*c,a*b=4);
eval(a*b*3,a*b=4);
eval(a*b*I,a*b=4);
Similarly:
subs(a*b=4,a*b*c);
subs(a*b=4,a*b*3);
subs(a*b=4,a*b*I);

J := Int((x^2+2*x+1+(3*x+1)*sqrt(x+ln(x)))/(x*sqrt(x+ln(x))*(x+sqrt(x+ln(x)))), x);
solve(sqrt(x+ln(x))=u,x);
IntegrationTools:-Change(J,x=%);
simplify(%) assuming u>0;
value(%);
res:=subs(u=sqrt(x+ln(x)),%); #Result 2*ln(LambertW(exp(x+ln(x)))+sqrt(x+ln(x)))+2*sqrt(x+ln(x))
##Note added: #After I saw Axel's result I tried simplify and got his simpler version:
res:=simplify(subs(u=sqrt(x+ln(x)),%)) assuming x>0;
## Checking the result:
diff(res,x);
simplify(%) assuming x>0;
IntegrationTools:-GetIntegrand(J);
simplify(%%-%); #0

Note: I made the assumption u > 0. That means that x > LambertW(1) = .5671432904 to 10 digits.
solve(x+ln(x)=0,x); #LambertW(1)

Here is another version.  (Note: I briefly had a slightly different version).

restart;
with(IntegrationTools):
f:=t^2*exp(-t)/ln(t);

J:=Int(f,t=0..infinity);
Change(J,t=exp(s));
Split(%,[-eps,eps]);
J1,J2:=op(1,%),op(3,%);
J1a:=Flip(subs(t=s,Change(J1,s=-t)));
Ja:=normal(combine(J1a+J2));
J:=eval(Ja,eps=0);
## Justifying putting eps=0:
g:=GetIntegrand(J);
series(g,s=0);

Your R is just an expression in t, not a function. Compare sin(t) with sin. The first is an expression in t the other is a function.
Thus the ode should not have R(t), but just R.
restart;
R := piecewise(t < 0.1e-7, 57, t > 0.1e-7, 57-48*exp(-t+0.1e-7)/(0.1e-7));

res:=dsolve({diff(y(t), t)-R*y(t)^2 = 0, y(0) = 1}, numeric);

plots:-odeplot(res,[t,y(t)],0..10^(-7));

############
To your comment in the last question:
I know that the laplace transform is a very useful tool; shouldn't Maple also know that?

Maple obviously knows that as you are also perfectly aware since you used method=laplace in dsolve and it worked (according to you).

Try this:
eq:=11.00=11.244522435+log(x) ;
solve(eq,x);
## The answer .7830784198 should come instantaneously.

You should give us your actual example.
I made up a simple one:
restart;
abs(1-exp(I*t))<1; #The one I want to solve, t being real
u:=evalc(abs(1-exp(I*t))); # evalc assumes that variables are real
res:=solve(u<1,t); #Three ranges (there are infinitely many)
plots:-complexplot(1-exp(I*t),t=0..2*Pi); #The whole curve (circle) in the complex plane
res2:=eval([res],{Open=(x->x),RealRange=`..`}); #For use in the plots below.
plots:-complexplot(1-exp(I*t),t=res2[1],scaling=constrained);
plots:-complexplot(1-exp(I*t),t=res2[2],scaling=constrained);#Just a repeat
plots:-complexplot(1-exp(I*t),t=res2[3],scaling=constrained);#Just a repeat

You can use:
ListTools:-Search([3,5],A);

I see two obvious reasons for the failure.

1. At each loop iteration the loop variable updates. At the end of the loop the variable i thus has the value nops(datalist)+1.

2. sum(enerList[i],,i=1..2) will give you a real problem: the arguments to sum are evaluated and i now has the concrete value nops(datalist)+1. So enerList[i] will complain (as it does), but had it not done it first, sum would have complained.

The solution is to use add instead of sum. That is by the way the recommended way for finite concrete sums as yours.
So use add(enerList[i],,i=1..2) instead. add has special evaluation rules; it is not bothered by the fact that i has a value.
Of course you could just have used a different variable (like j) or used unevaluation quotes around i as in 'i', but use add!

You are creating a new matrix inside the loop a:=J[1].J[1]; (i=1,j=1)
That new matrix has the same elements as J[1], but it is a separate mutable object: Thus you can change the one without changing the other.
Try this:
restart;
A:=Matrix([[1,0],[0,1]]);
B:=Matrix([[1,0],[0,1]]);
evalb(A=B); #False
A-B; # A zero-matrix
is(A=B); #False
LinearAlgebra:-Equal(A,B); #Elemenwise check: true
## A and B reside in different locations in memory:
addressof(A);
addressof(B);
##
## So you can try
printlevel:=2:
for i to 1 do for j to 2 do a := J[i].J[j]; map2(LinearAlgebra:-Equal,a,J) end do end do;
##and
for i to 8 do for j to 8 do a := J[i].J[j]; k:=ListTools:-Search(true,map2(LinearAlgebra:-Equal,a,J)); print(i,j,k,a) end do end do:
####
## Finally about not being of type matrix inside the loop. Of course it is! What happens is illustrated here outside any loop:
whattype(A); #Answer: Matrix
print(%); #Answer: A procedure is printed, the matrix constructor procedure.
showstat(Matrix); #To see all the lines in that procedure.

There are several problems:

1. In w:=laplace(I3,t,s): you need to have inttrans[laplace] instead of just laplace.

2. In w2:=limit(sz,t=infinity): you need limit(value(sz),t=infinity):

3. If you do type(w2,range); you get the answer true.
   Trying with w2MS:=MultiSeries:-limit(value(sz),t=infinity);
    you get undefined. Thus it appears that the limit doesn't exist.
##In fact the limit doesn't seem to exist. Just try:
 plot(value(sz),t=10..10.00001);

Or look at the output of
evalf(combine(expand(value(sz))));
If you choose to ignore the wiggles (which are there) and replace w2 by -1/3 as it seems to oscillate about, then you can plot easily:

plot(Re(w1)*(-1/3),d=0..5,tickmarks=[2, 2],thickness=2,linestyle=1,color=black,axes=boxed,titlefont=[SYMBOL,14],font=[1,1,18],tickmarks=[2,3],thickness=2);

Here is an answer to the first problem.

restart;
Expr:=a*diff(x(t),t)+b*x(t)+r*diff(x(t),t,t)+a*diff(y(t),t)+b*diff(y(t),t,t)+c*y(t);
L:=convert(indets(Expr,specfunc(diff)),list);
res:=coeffs(Expr,L,k);
k;
S:=seq(res[ListTools:-Search(L[i],[k])]*L[i],i=1..4),res[ListTools:-Search(1,[k])];
Expr:='Expr';
gc();
Expr:=add(i,i=S);

##Please see my comment to Kitonum's answer for the need for the two important lines in bold face.

You should look into formatted printing. See ?printf

printf("The solution to %a = %d mod %d is x = %d\n", a^x ,b, p, 1323);

In Maple 2016:

restart;
A:=Int( arctan(x) * ln(1+1/x^2), x = 0..infinity );
value(A); #No
evalf(A);
identify(%); #Pi^2/6
IntegrationTools:-Parts(A,arctan(x));
value(%); #Pi^2/6