Not knowing the constants I just picked some, and did the following with no problem:

E_horiz := sqrt(E_theta^2 + E_phi^2):

param:={W=1, h=2, k0=3, kx=3, ky=4, l=5, epsr=6, seq(JX[i]=1/(1+i),i=1..80),epsre=7}:
res:=eval(E_horiz,param):
plot(res,theta=0..2*Pi);

Not knowing the constants I just picked some, and did the following with no problem:

E_horiz := sqrt(E_theta^2 + E_phi^2):

param:={W=1, h=2, k0=3, kx=3, ky=4, l=5, epsr=6, seq(JX[i]=1/(1+i),i=1..80),epsre=7}:
res:=eval(E_horiz,param):
plot(res,theta=0..2*Pi);

There is an obvious syntactical error in the statement E_horiz := subs(theta := Pi/2, E_horiz):

The assignment operator ':=' should be replaced by '=' of course.

What do you mean by

"I am looking for E_phi(theta), E_theta(theta), E_horiz(phi)."

?

When you say

'Maple cannot solve this line: "E_phi_magnitude := sqrt((Re(E_phi))^2 + (Im(E_phi))^2):".'

what do you mean by solve? There is no equation.

Since there are lots of names involved having no numeric values there is not much simplification that can be done.

There is an obvious syntactical error in the statement E_horiz := subs(theta := Pi/2, E_horiz):

The assignment operator ':=' should be replaced by '=' of course.

What do you mean by

"I am looking for E_phi(theta), E_theta(theta), E_horiz(phi)."

?

When you say

'Maple cannot solve this line: "E_phi_magnitude := sqrt((Re(E_phi))^2 + (Im(E_phi))^2):".'

what do you mean by solve? There is no equation.

Since there are lots of names involved having no numeric values there is not much simplification that can be done.

An example would help tremendously.

An example would help tremendously.

@asfn In y appear 3 square roots: sqrt(12.25 - b^2), sqrt(2.1025-b^2), and sqrt(1-b^2). This reveals where changes are likely.

Without attempting a proof that my statement about y being purely imaginary between 1 and sqrt(2.1025) you can test the validity by evaluating y at points chosen randomly:

y:=subs(`βb`=b,y);
rnd:=RandomTools:-MersenneTwister:-GenerateFloat;
#Testing rnd
rnd();

#Now evaluating y at 10 points chosen randomly in the strip mentioned:

to 10 do evalf(eval(y,{b=1+rnd()*(sqrt(2.1025)-1),lambda=1+rnd()})) end do;

Similarly, you can test the claim that y is real for sqrt(2.1025) < b < sqrt(12.1025) by

to 10 do evalf(eval(y,{b=sqrt(2.1025)+rnd()*(sqrt(12.25)-sqrt(2.1025)),lambda=1+rnd()})) end do;

Now the disturbing question is: What did implicitplot actually do? It should have painted the whole strip 1 < b < sqrt(2.1025) red, but produced some nice looking curves!

@asfn In y appear 3 square roots: sqrt(12.25 - b^2), sqrt(2.1025-b^2), and sqrt(1-b^2). This reveals where changes are likely.

Without attempting a proof that my statement about y being purely imaginary between 1 and sqrt(2.1025) you can test the validity by evaluating y at points chosen randomly:

y:=subs(`&beta;b`=b,y);
rnd:=RandomTools:-MersenneTwister:-GenerateFloat;
#Testing rnd
rnd();

#Now evaluating y at 10 points chosen randomly in the strip mentioned:

to 10 do evalf(eval(y,{b=1+rnd()*(sqrt(2.1025)-1),lambda=1+rnd()})) end do;

Similarly, you can test the claim that y is real for sqrt(2.1025) < b < sqrt(12.1025) by

to 10 do evalf(eval(y,{b=sqrt(2.1025)+rnd()*(sqrt(12.25)-sqrt(2.1025)),lambda=1+rnd()})) end do;

Now the disturbing question is: What did implicitplot actually do? It should have painted the whole strip 1 < b < sqrt(2.1025) red, but produced some nice looking curves!

G1:=Matrix(101,93,(i,j)->evalf(sin(i/7)+cos(j/5)),datatype=float[8]):

plots:-surfdata(convert(G1,listlist));

plots:-surfdata(Array(G1));

G1:=Matrix(101,93,(i,j)->evalf(sin(i/7)+cos(j/5)),datatype=float[8]):

plots:-surfdata(convert(G1,listlist));

plots:-surfdata(Array(G1));

@zhkhan There appears a G, which apparently is just a name.

@zhkhan There appears a G, which apparently is just a name.

@zhkhan I'm not sure what you mean.

But the two larger roots start their existence at the following value for l: 

fsolve({diff((rhs-lhs)(A),sigma)=0,(rhs-lhs)(A)=0},{sigma,l=3.5});

 Output: {l = 3.531408488, sigma = -32.34143629}

I just used the default setting of Digits = 10.

Both roots can be found by RootFinding:-NextZero:

r1:='r1': r2:='r2':
for l from 3.7 by .1 to 5 do
   r1:=RootFinding:-NextZero(unapply((rhs-lhs)(A),sigma),-120);
   if r1<>FAIL then r2:=RootFinding:-NextZero(unapply((rhs-lhs)(A),sigma),r1) end if;
   print(l,r1,r2);
end do:
l:='l':

There will be numerical problems when the first of the two roots is close to the singularity of rhs(A), which happens for larger values of l.

@zhkhan I'm not sure what you mean.

But the two larger roots start their existence at the following value for l: 

fsolve({diff((rhs-lhs)(A),sigma)=0,(rhs-lhs)(A)=0},{sigma,l=3.5});

 Output: {l = 3.531408488, sigma = -32.34143629}

I just used the default setting of Digits = 10.

Both roots can be found by RootFinding:-NextZero:

r1:='r1': r2:='r2':
for l from 3.7 by .1 to 5 do
   r1:=RootFinding:-NextZero(unapply((rhs-lhs)(A),sigma),-120);
   if r1<>FAIL then r2:=RootFinding:-NextZero(unapply((rhs-lhs)(A),sigma),r1) end if;
   print(l,r1,r2);
end do:
l:='l':

There will be numerical problems when the first of the two roots is close to the singularity of rhs(A), which happens for larger values of l.

@zhkhan You may want to take a look at ExportMatrix.