The equation in question was

St[L]/exp(lambda^2)/erf(lambda)-St[S]/nu/exp(nu^2*lambda^2)/erfc(nu*lambda) 
    = lambda*Pi^(1/2)

and the solution you wanted was

lambda = 1/2*St[S]/nu/Pi^(1/2)+1/2*(2*St[L]+St[S]^2/nu^2/Pi)^(1/2)

Unfortunately, this is not a solution.  For example, with St[L] = St[S] = nu = 1, the left side of the equation with this lambda evaluates to -2.012282639 and the right side to 1.849368862.

I very much doubt that there is a closed-form solution.

Use Rename in the FileTools package.  For example:
 

> oldnames:= ["file1.jpg", "file2.jpg", "file3.jpg"]:
  newnames:= ["alpha.jpg", "beta.jpg", "gamma.jpg"];
  currentdir("c:/myfolder"):
  zip(FileTools[Rename], oldnames, newnames);

Try copying the following piece of MathML and pasting it into your MathContainer. 

<math xmlns='http://www.w3.org/1998/Math/MathML'><mtable rowalign='baseline' 
columnalign='center' groupalign='{left}' rowspacing='1.0ex'>
<mtr><mtd><mi>A</mi></mtd></mtr><mtr><mtd><mi>B</mi></mtd></mtr>
<mtr><mtd><mi>C</mi></mtd></mtr><mtr><mtd><mi>D</mi></mtd></mtr>
<mtr><mtd><mi>E</mi></mtd></mtr></mtable></math>

When Maple asks:

The clipboard contains MathML.
Do you want to convert the MathML into 2D math?

 click Yes.

Or you can do it programmatically:

> DocumentTools[SetProperty](MathContainer0,'value',
"<math xmlns='http://www.w3.org/1998/Math/MathML'><mtable rowalign='baseline' 
columnalign='center' groupalign='{left}' rowspacing='1.0ex'>
<mtr><mtd><mi>A</mi></mtd></mtr><mtr><mtd><mi>B</mi></mtd></mtr>
<mtr><mtd><mi>C</mi></mtd></mtr><mtr><mtd><mi>D</mi></mtd></mtr>
<mtr><mtd><mi>E</mi></mtd></mtr></mtable></math>");

That should be Statistics, not statistics.  The NonlinearFit function has a parameterranges option that will do this.  For example:

> NonlinearFit(((A/2)/((B-r)^2+(A/2)^2)+C1+C2*r+C3*r^2),WEE,dde,r, 
    parameterranges=[A=0..1, B=-3..3]);

frac(x/n) almost does it: the only problem is that frac(x) < 0 if x < 0.  Instead try

> `&discretemod`:= proc(x::complexcons, n::realcons)
     if type(x, realcons) then x/n - floor(x/n)
     else Re(x)/n - floor(Re(x)/n) + I*(Im(x)/n - floor(Im(x)/n))
    end if
    end proc;

> (-3.2 - 7.6 * I) &discretemod 2;

.400000000+.200000000*I

OK, So you do know which is the smallest.  Now the question is how to do the merging.  My suggestion is to write a procedure that takes two already-sorted ranges of an input list or rtable and merges them into an output rtable, starting at a given position.  For efficiency it's better to use rtables rather than lists, as assigning one element of a list requires the whole list to be copied.  Thus

Merge:= proc(a,b,c,d, In, e, Out)
# merge In[a..b] and In[c..d] into Out[e..(e+(b-a)+(d-c)+1)]

You'll use two local variables, say  i and j, to keep track of where you are in the In ranges and another, k, to keep track of where you are in Out.
So you start with i = a, j = c, k = e. 
Now in your loop, you'll first want to check if you're finished with one of the In ranges (in which case you just copy whatever is left of the other In range into the
Out range and exit).  Then compare In[i] to In[j], put the smaller one into Out, and update (i or j) and k.

There is no minimum.  Hint: try a cubic.

I hope you left out a minus sign: it should have been

 > limit(int(f(t)*exp(-s*t),t = 0..infinity),s = infinity)

Maple should not be trusted with this sort of analysis question: it often works "formally", with little regard for convergence.
The basic results having to do with interchanging limit and integral are Lebesgue's Dominated Convergence Theorem
and Monotone Convergence Theorem.  In cases such as this, involving a Laplace transform, you may also find Watson's
Lemma to be useful.
 

See the help pages

?Definition,monotoneconvergencetheorem

?Definition,dominatedconvergencetheorem

?Definition,Watsonslemma

To get one solution (or the solution that optimizes some linear objective), you could use LPSolve in the optimization package, with the option assume=binary.  For example:

> with(Optimization):
  LPSolve(b+c+e, {a+b+c=1, a+d+e=2, b+d+f=2, c+e+f=1}, assume=binary);

[0, [a = 1, b = 0, c = 0, d = 1, e = 0, f = 1]]

> LPSolve(b+c+e, {a+b+c=1, a+d+e=2, b+d+f=2, c+e+f=1}, assume=binary, maximize);

[2, [a = 0, b = 1, c = 0, d = 1, e = 1, f = 0]]

Another possibility is to use isolve, adding constraints v*(1-v) = 0 for each variable v to make them 0 or 1.

> isolve({a+b+c=1, a+d+e=2, b+d+f=2, c+e+f=1,seq(v*(1-v)=0,v=[a,b,c,d,e,f])});

{a = 0, b = 1, c = 0, d = 1, e = 1, f = 0}, {a = 1, b = 0, c = 0, d = 1, e = 0, f = 1}

pointplot in the plots package accepts a more efficient data structure, an Array of datatype float[8], but it "unpacks" that into an expression sequence of lists.  You might use quickplot in my Maple Advisor Database,
www.math.ubc.ca/~israel/advisor

The Array must start at index 0, not 1, however.
Something else in your code: sort doesn't act "in-place" on an Array, so if you want to use the sorted version of requests you must assign the result back to results.  So, after installing the database:

> genRequestDist := proc (N, t)
 local requests, i, r, pts, d; 
 r := rand(t);
 requests := Array(1 .. N,datatype=integer[8]);
 for i to N do
  requests[i] := r()
 end do;
 requests:=sort(requests);
 pts := Array(1 .. t+1, 1 .. 2,datatype=float[8]);
 for i from 0 to t do
  pts[i+1, 1] := i;
  pts[i+1, 2] := 0
 end do;
 for i to N-1 do
  d := abs(requests[i+1]-requests[i]);
  pts[d+1, 1] := d;
  pts[d+1, 2] := pts[d+1, 2]+1
 end do;
 pts
 end proc;
 A:= genRequestDist(100000,3600000);
 quickplot(A, style=point, color=red,axes=frame);

This works well in Classic: I had problems with it in the Standard interface, though.

Actually one of the solutions returned by the symbolic pdsolve works here:

> pdsolve(eq);

{i(x,t) = 0., s(x,t) = _F1(x)}, <something very complicated>

_F1 is an arbitrary function.  In your case you want _F1(x) = x, in order to satisfy the boundary conditions.  That is, the solution is i(x,t) = 0, s(x,t) = x.

I'm not sure why pdsolve(..., numeric) doesn't work here.  If it were just one PDE rather than a system, you might try the method option, but this doesn't work for systems.

Maple is quite correct: your equations 2, 4, 6, 8, 10 and 12 all contain the variable x.  This seems to come in via shape2y.  Perhaps you meant to say
 

sc := shape1xx^%T.shape1yy+shape2xx^%T.shape2yy

rather than

sc := shape1xx^%T.shape1yy+shape2xx^%T.shape2y

?

The initial condition must be y(t) for t in an interval of length 1.  I'll suppose y(t) = 1/2 for 0 <= t <= 1.  Then in each interval [n,n+1] for positive integers n, y(t) is determined by the differential equation with initial conditions at y(t) and y(t-1) given by the solution in the previous interval.  My first thought was to use

> ys[0]:= proc(t) if type(t,numeric) then 1/2 else 'procname'(t) end if end proc:;
    for n from 1 to 99 do
       ys[n]:= subs(dsolve({diff(y(t),t) = y(t) * (1 - ys[n-1](t-1)), y(n) = ys[n-1](n)}, 
                  numeric, known=[ys[n-1]], output=listprocedure),y(t))
   end do;  

But this didn't work well, because to solve the DE for ys[n] at one point dsolve must solve the DE for ys[n-1], and so on...
things slow down rapidly.  A better idea is to get an explicit expression for each ys[n].  I'll use chebpade in the numapprox package to get a good approximation by a fourth-degree polynomial.

> with(numapprox):
  ys[0]:= 1/2;
  for n from 1 to 99 do
    yd:= subs(dsolve({diff(y(t),t) = y(t) * (1 - ys[n-1](t-1)), y(n) = ys[n-1](n)}, numeric,
       output=listprocedure),y(t));
    ys[n]:= unapply(subs(T=orthopoly[T],chebpade(yd(t),t=n..n+1,4)),t);   
  end do: 

This returns an error:

Error, (in numapprox:-chebpade) for this function the maximum approximation degree is 3 when Digits = 10
 

The error occurs at n = 44, where our solution is extremely close to 1 over the whole interval.  For n from 44 I'll change the degree of the polynomial to 3.

> for n from 44 to 99 do
    yd:= subs(dsolve({diff(y(t),t) = y(t) * (1 - ys[n-1](t-1)), y(n) = ys[n-1](n)}, numeric,
       output=listprocedure),y(t));
    ys[n]:= unapply(subs(T=orthopoly[T],chebpade(yd(t),t=n..n+1,3)),t);   
  end do;

This one only goes up to n = 48 before Maple complains.  Since it's late at night, I'll quit here.

> plots[display]([seq(plot(ys[n](t), t = n .. n+1), n = 0 .. 47)]);

For example:

> evalc([solve(x^3 + x^2 - 12*x + 8)]);

[2/3*37^(1/2)*sin(1/3*arctan(6/163*669^(1/2))+1/6*Pi)-1/3, -1/3*37^(1/2)*sin(1/3*arctan(6/163*669^(1/2))+1/6*Pi)-1/3-1/3*3^(1/2)*37^(1/2)*sin(-1/3*arctan(6/163*669^(1/2))+1/3*Pi), -1/3*37^(1/2)*sin(1/3*arctan(6/163*669^(1/2))+1/6*Pi)-1/3+1/3*3^(1/2)*37^(1/2)*sin(-1/3*arctan(6/163*669^(1/2))+1/3*Pi)]

<<x+1|x|x|`...`|x>,<x|x+2|x|`...`|x>,<x|x|x+3|`...`|x>,<`...`|`...`|`...`|`...`|`...`>,
  <x|x|x|`...`|x+n>>;

              [x + 1      x        x      ...      x  ]
              [                                       ]
              [  x      x + 2      x      ...      x  ]
              [                                       ]
              [  x        x      x + 3    ...      x  ]
              [                                       ]
              [ ...      ...      ...     ...     ... ]
              [                                       ]
              [  x        x        x      ...    x + n]

Of course you can't expect Maple to understand this when you try calculating anything, it's just for display purposes.

> MM:= n -> Matrix(n,n,proc(i,j) if i=j then x+i else x end if end proc);
  MM(5);

             [x + 1      x        x        x        x  ]
             [                                         ]
             [  x      x + 2      x        x        x  ]
             [                                         ]
             [  x        x      x + 3      x        x  ]
             [                                         ]
             [  x        x        x      x + 4      x  ]
             [                                         ]
             [  x        x        x        x      x + 5]

> R:= [seq(LinearAlgebra[Determinant](MM(n)),n=1..10)];

R := [x+1, 3*x+2, 11*x+6, 50*x+24, 274*x+120, 720+1764*x, 5040+13068*x, 40320+109584*x, 362880+1026576*x, 3628800+10628640*x]

It's easy to recognize the constant term of R[n] as n!.  As for the coefficient of x:

> C:= map(coeff,R,x);

C := [1, 3, 11, 50, 274, 1764, 13068, 109584, 1026576, 10628640]

Look this up at the On-line Encyclopedia of Integer Sequences www.research.att.com/~njas/sequences/index.html: you find