This also works:

bound := 0:
do until SMTLIB['Satisfiable'](And(
      (154*x + 69*y) - 0 = 7**3*m1, 
      (13*x + 716*y) - 0 = 13**3*m2, 
     (23*x + 3059*y) - 0 = 23**3*m3, 
         x > 0, x < y, y = ++bound), 
                 'logic' = "QF_LIA"):
bound;
 = 
                             18837

Unfortunately, it takes several minutes to solve on my low-end PC.
That may fell a bit slow, yet the human time required to code a faster algorithm is a lot longer than several minutes and perhaps longer than one care to wait.

Edit. Of course the original problem can be reformulated as an integer linear program, but to the best of my knowledge, none of "exact" solvers are integrated into (and accordingly, are unable to take advantage of) modern computer algebra systems. (Actually, when I search for this term, I can only find exact SCIP, Exact, QSopt_ex and SoPlex, which are all irrelevant to a CAS that is capable of carrying out computations analytically.) Alas, what a pity!

If I understand right, "[x^2-x*y, x*y-y^2, x*z-y*z, -x^2+y^2, x*y+x*z]" can be obtained as follows: 

convert(expand({eval(MmaTranslator:-Mma:-Distribute([x - y]*[x, y, z], list, `*`), list = expand@``), -x^2 + y^2, x*y + x*z}), list);
 = 
      [       2    2    2   2                            ]
      [x y - y , -x  + y , x  - x y, x z - y z, x y + x z]

And "[x^3-x^2*y, x^2*y-x*y^2, x^2*z-x*y*z, x*y^2-y^3, x*y*z-y^2*z, x*z^2-y*z^2, -x^3+x*y^2, -x^2*y+y^3, -x^2*z+y^2*z, x^2*y+x^2*z, x*y^2+x*y*z, x*y*z+x*z^2, -x*y*z+z^3, x*y^2+y*z^2-z^3]" can be obtained as follows: 

convert(expand({eval(MmaTranslator:-Mma:-Distribute([x^2 - x*y, x*y - y^2, x*z - y*z, -x^2 + y^2, x*y + x*z]*[x, y, z], list, `*`), list = expand@``), -x*y*z + z^3, x*y^2 + y*z^2 - z^3}), list);
 = 
  [   2    3    2      3   2        2    3      2   3    2    
  [x y  - y , -x  y + y , x  y - x y , -x  + x y , x  - x  y, 

       2      2            3           2               2  
    x z  - y z , -x y z + z , x y z - y  z, x y z + x z , 

       2            2      2     2             2      2    
    x y  + x y z, -x  z + y  z, x  z - x y z, x  y + x  z, 

       2      2    3]
    x y  + y z  - z ]

@minhthien2016 Alternatively, 

subs('y' = 1, collect('y' * ~ [exp1, exp2], x, normal));
 = 
[ 3              2                 3            2              ]
[x  + (3 m - 1) x  - (4 m + 3) x, x  + (m + 1) x  + (4 m + 3) x]

solve('identity'(aa/f^bb = 1.234/(f^6.789), f), {aa, bb});
 = 
              {aa = 1.234000000, bb = 6.789000000}

PDETools:-Solve(aa/f^bb = 1.234/(f^6.789), {aa, bb}, 'independentof' = f);
 = 
              {aa = 1.234000000, bb = 6.789000000}

match/patmatch/typematch also works.

For some unknowable reason, Maple cannot fully simplify the L. Fortunately, a workaround exists: 
 

Download 236850-How-To-Get-A-Simpler-Expression-Of-The.mw

What about … 

nequal:=proc($)
    options hfloat;
    if procname::indexed(posint&under(k->k-2))then
        local n:=op(procname),l:='`tools/genglobal`(evaln(i))'$n-1,j,m::nonnegint:=0;
        eval(eval(foldl(uneval(seq),ifelse((n-2)*24-`+`(l)>=1,++m,NULL),seq(l[j]=ceil(((n-2)*24-`+`(l[..j-1]))/(n-(j-1)))..`if`(j>1,l[j-1],(n-2)*24-(n-1)),j=n-1..1,-1))));
        print(m)
    fi
end:

Unfortunately, it is not more efficient: 

(* a function analogous to OP's `nequal3` and `nequal4` *)
nequal5:=proc()
    option autocompile;
    local u,v,w,x,a_5;
    a_5:=0;
    for u from ceil((5-2)*24/(5-0)) to (5-2)*24-5+1 do
        for v from ceil(((5-2)*24-u)/(5-1)) to u do
            for w from ceil(((5-2)*24-u-v)/(5-2)) to v do
                for x from ceil(((5-2)*24-u-v-w)/(5-3)) to w do
                    if (5-2)*24-u-v-w-x>=1 then
                        a_5:=a_5+1;
                    end if;
                end do;
            end do;
        end do;
    end do;
    print(a_5);
end proc:

CodeTools:-Usage(nequal[5]()):
CodeTools:-Usage(nequal5()):  = 
                             10642

memory used=135.94MiB, alloc change=8.00MiB, cpu time=3.59s, real time=3.54s, gc time=187.50ms

                             10598

memory used=35.71KiB, alloc change=0 bytes, cpu time=94.00ms, real time=87.00ms, gc time=0ns

Yet the result is guaranteed to be correct.

E.g., 

decisions := "accept", "reject":
T := 2:
`?()`(`if`, [T < 1, decisions]);
 = 
                            "reject"

Note that here you cannot directly use apply.
Another way: 

eval('ifelse'(T > 1, decisions));
 = 
                            "accept"

How about

DocumentTools['Tabulate'](<"", `<,>`(F(x)) | `<|>`(C), Matrix(numelems(F), numelems(C), fill = NULL)>, width = 50, fillcolor = (T, i, jj) -> ifelse(i = 1 or jj = 1, [210, 230, 255], [255, 255, 235])):

While @acer uses the _npassed (which is equivalent to nargs), @Carl Love uses a much advanced command overload. However, if you only intend to create the desired computation sequences in seconds, the following approach will be more convenient: 

> redefine(f, [f(a::algebraic, b::algebraic, c::algebraic) = a + b*c, f(a::algebraic, b::algebraic) = a/b]);
> f(0);
                                                         f(0)

> f(1, 2);
                                                          1/2

> f(3, 4, 5);
                                                          23

> f(6, 7, 8, 9);
                                                     f(6, 7, 8, 9)

> procbody(f):
s=****************, s[1]=**************** invalid MODDEF parameter sequence
Error, (in procbody) assembled object is invalid

Which implementation is the most efficient?

In my opinion, what you need is simply the symbol nargs. (Note that this syntax differs from another CAS's one.)

You say: "But I only want the names."
But the exports command (not the Describe command) does return a "exprseq" of names. For instance,

exports(PatternMatching);
 = 
         Define, Recogniser, AlgStruct, Tools, deftypes

Describe(PackageManagement);

# Package management utilities
module PackageManagement:

    # test whether a name is exported by a package
    pmember( nom::name, pkg )

    # return a list of the exports of a package
    pexports( pname )

Did I understand right?

The reason is: the output of convert(Eval(diff(varphi(t), t), t = 0), D) becomes unevaluated (and lowercase) eval(D(varphi)(t), t = 0). Actually, in some old release of Maple, the result is correct capitalized version. Besides, convert(``(Eval(diff(varphi(t), t), t = 0)), D) gives evaluated ``(D(varphi)(0)) instead. 

Why? In my opinion, this must be a bug.

For your specific purpose ("red if s is between 0 and ell/2, and green if s is between ell/2 and ell"), the easier way is:

plots:-spacecurve([cos(theta(s)), v(s), sin(theta(s))], s = 0 .. ell, colorscheme = ["Red", "Green"]);

Another way is:

data:=<seq(<cos(theta(s))|v(s)|sin(theta(s))|s>,s=0..ell,1e-3)>:
plots:-spacecurve(data[..,..3],colorscheme=["valuesplit",data[..,4],[0..ell/2="Red",ell/2..ell="Green"]]);

@vv Strangely, I find that if I use "y" (rather than "y(x)") in the "sol", then solve works well: 

y(x) = solve(subs(y(x) = y, sol), y);
 = 
                       /                                    2 
                       |/ /                              2\\  
                       || |      /    3       2         \ ||  
              1        || \RootOf\2 _Z  + 2 _Z  - _Z + 1/ /|  
 y(x) = -------------- \\x                                 /  
                   2                                          
        (exp(c[1]))  x                                        

                                              2 
   /   /                              2     \\  
   |   |      /    3       2         \      ||  
   \exp\RootOf\2 _Z  + 2 _Z  - _Z + 1/  c[1]//  

                                    3 
   /       /    3       2         \\  
   | RootOf\2 _Z  + 2 _Z  - _Z + 1/|  
   \x                              /  

                                             3
   /   /      /    3       2         \     \\ 
   \exp\RootOf\2 _Z  + 2 _Z  - _Z + 1/ c[1]// 

                                                    \
                                                    |
                                                    |
            /    3       2         \  2            2|
    + RootOf\2 _Z  + 2 _Z  - _Z + 1/ x  (exp(c[1])) /

Why? I don't know.

This problem is due to the limited ability of subs (see showstat(isolate, 17);). Actually, "2*I*A" (as well as "2*(I*A)") is parsed as "(2*I)*A", while the subs command just does syntactical exact-match replacements.
For example: 

> dismantle(I*A);

PROD(5)
   COMPLEX(2)
      INTPOS(2): 1
   INTPOS(2): 1
   NAME(4): A
   INTPOS(2): 1


> dismantle(2*(I*A));

PROD(5)
   COMPLEX(2)
      INTPOS(2): 2
   INTPOS(2): 1
   NAME(4): A
   INTPOS(2): 1


> dismantle(i*A);

PROD(5)
   NAME(4): i
   INTPOS(2): 1
   NAME(4): A
   INTPOS(2): 1


> dismantle((2*i)*A);

SUM(3)
   PROD(5)
      NAME(4): i
      INTPOS(2): 1
      NAME(4): A
      INTPOS(2): 1
   INTPOS(2): 2

They have different structures. So,  will never work directly (I mean, without any manual steps like  and ).