@Joe Riel Oh, I just misunderstood your last sentence to mean that the attributes trick could be done simpler, faster, and smaller in a loop. You just meant a straightforward max-finding loop. It's all clear now.

Please repost the array Data in a plaintext form so that I can cut-and-paste it.

@Joe Riel I'm having trouble understanding your last sentence---why it's simpler, faster, smaller. Could you show an example doing it in a do loop?

It just requires a slight modification of the original:

add(s^LL[1]*t^LL[2]/mul(x!, x= LL), LL= L);

@reemeaaaah I see what's wrong. When you typed into MaplePrimes, you typed B1:= ..., B2:= ..., but that should've been B[1], B[2]---indicating a subscripted variable. Change the sort to

sort([seq(B[k], k= 1..n)], (A,B)-> (A[2]-A[1]) > (B[2]-B[1]))[1];

@Kitonum Your comparison is invalid because the sums in the first plot command are being evaluated to exp(z)*w. For a valid comparison you need to use Sum (not sure it will work) or do something else to evaluate the sums numerically.

The exp(z)*w is obviously wrong. Just divide both sides by w, and then the sum depends on w while the other side does not.

I created tag "hypergeom" and applied it to your Question.

@J4James I can get portions of the above graph, but I am having trouble solving the BVP for b > 16. I get error message "Newton iteration is not converging". I am trying to adjust the numerous parameters to dsolve. Here is what I have so far (this plots a wider range than my first):

Eq:= diff(f(y),y$3) + f(y)*diff(f(y),y$2) - diff(f(y),y$1)^2 -
     c*S*(diff(f(y),y$1) + y/2*diff(f(y),y$2))=0:
bc:= D(f)(0)=1, f(0)=0, (D@@2)(f)(b)=0:

F:= proc(_b, _S)
local Sol;
     try
          Sol:= dsolve(eval({Eq,bc}, {S= _S, b= _b}), f(y), numeric, continuation= c);
     catch:
          userinfo(1, F, StringTools:-FormatMessage(lastexception[2..-1]));
          error
     end try;
     eval(f(y), Sol(_b)) - _S*_b/2
end proc:
     
infolevel[F]:= 1:

Note the inclusion of continuation parameter c in Eq.

Don't use the above for b=0 because then it's an IVP not a BVP and different coding is required.

@nm You wrote:

fouriercos(F(x),x,w);

That should be F(w), not F(x).

I believe that what Axel was saying is that fouriercos is its own inverse.

@AlStiegman 

Is your version of Maple before Maple 13? Then replace

{vars} =~ 0..1

with

map(`=`, {vars}, 0..1)

To plot just a single variable, use

plot([seq([TT(j), eval(x, Sols[j])], j= 0..10)]);

@AlStiegman 

You can do the five plots together like this:

vars:= w,v,x,y,z:
TT:= j-> 800+10*j;
for j from 0 to 10 by 1 do
     T:= TT(j);
     Sols[j]:= fsolve({fsc, fbd, fhc, fmassbal, fmb}, {vars} =~ 0..1)
end do:

plot(
     [seq([seq([TT(j), eval(V, Sols[j])], j= 0..10)], V= [vars])],
     color= [black, red, green, blue, brown], legend= [vars]
);

@kminor2578 This is not your fault, but the file didn't upload. This is a known bug in MaplePrimes.

Your post is missing the equations. Try attaching a worksheet to a Reply.

By solving for where the arguments of the three logarithms are positive, we can see that the answer returned by solve is the domain of the function, not where the function is positive. So this is a bug. After that, I solved it by the plot and fsolve method, but Thomas Richard posted it first.

@Joe Riel I believe that Acer, in this thread, applied the same asolve to the OP's function but on the domain 0..100 and showed that it didn't get all the roots.

The procedure asolve seems to be coded to include the possibility that fsolve could return more than one root per call. When does that happen for a nonpolynomial? (For a polynomial, asolve is unnecessary because fsolve returns all real roots.)