@emendes You will not be able to understand these until you understand the purpose of the selection-criterion procedure B from my previous examples!

Both of these methods produce the same list result. Both of them preserve the order in which the elements were generated.

restart:
S:= {'rand'()$9}:
B:= s-> s[1]::even:

# The Array method:
A:= Array(1..0): 
for C in Iterator:-Combination(nops(S), 6) do
    if B((s:= S[[seq(C+~1)]])) then A,= s fi
od:
L:= [seq(A)];

# The embedded-loop method: 
L:= [
    for C in Iterator:-Combination(nops(S), 6) do
        if B((s:= S[[seq(C+~1)]])) then s fi
    od
];

@fatemeh1090 Do you mean to consider those inequalities individually, or all at once?

For nondegenerate cases, there are only two possible solutions---the solutions that you already had. (Degenerate cases are when the parameters are such that the coefficient of q^2 in the minimal polynomial is 0.) All those hundreds of cases all amount to simply deciding, based on the parameters, which of those two solutions work in the original equation. The cases do not change the two solutions.

@emendes Your example doesn't work because you didn't define a procedure B. Like I said, B needs to be a procedure that returns true or false when passed a subset of S.

What's an .mm file?

@emendes The answer to your first question is an unequivocal "Yes, you can do the same thing using seq."

seq(S[[seq(C+~1)]], C= Iterator:-Combination(nops(S), 6))

The answer to your second question is more nuanced. Let's suppose that you have a "simple" boolean procedure B that returns true or false depending on whether you want to select the subset. By "simple", I mean that the procedure looks only at one subset at a time; it's not allowed to check the part of the sequence already selected. Then you can use seq to create a list like this

B1:= (s::set)-> `if`(B(s), s, NULL);
L:= [seq(B1(S[[seq(C+~1)]]), C= Iterator:-Combination(nops(S), 6))]

Now, let's admit the possibility that B is not simple. And, let's suppose that you want a list of all subsets that satisfy B, but you don't care about the order. In Maple 2017, the easiest efficient thing to do is

R:= table(): #optional, but recommended, table intialization
for C in Iterator:-Combination(nops(S), 6) do
    s:= S[[seq(C+~1)]];
    if B(s) then R[s]:= () fi; #entry is NULL; index is desired info.
od:
L:= [indices(R, 'nolist')]; #convert table indices to list

Never try to create a list, set, or sequence by adding one element at a time in a for-do loop! It is extremely inefficient. Instead, inside the loop, build a table or Array, and then convert it to a list, set, or sequence outside the loop.

Maple 2019 introduced several major syntax additions that make this easier.

Let me know if these satisfy your need.

I just Googled "GPS coordinates of London Metro stations", and the first hit had downloadable coordinates.

@acer You're right, of course. I just had trouble reading the fine print in the image, and I missed noticing the absence of the dash.

I don't see how this could be useful. The number of possible types is clearly infinite. One could check against all (non-parameterized) types for which code already exists, but what's the point of that?

Considering parameterized types opens up another dimension of complexity.

@Scythor Your error message above was caused by

G:= GT:-DrawGraph(...)

This makes G a plot, not a graph.

The vertex positions can be set explicitly with SetVertexPositions. Since you presumably already have a map of the Metro system, this shouldn't be too difficult. You could even use latitude and longitude coordinates.

@Scythor 

Perhaps it would help to remove the degree-two nodes without reconnecting their neighbors. This would change the topology, but it should open things up.

GT:= GraphTheory:  LT:= ListTools:
G:= GT:-RandomGraphs:-RandomGraph(267, 727, connected);
do
    v:= LT:-Search(2, GT:-DegreeSequence(G));
    if v=0 then break fi;
    G:= GT:-DeleteVertex(G, GT:-Vertices(G)[v]);   
od:
    
GT:-DrawGraph(G, layout= spring, dimension= 3);

The problem seems akin to the software that displays folded proteins.

Perhaps if you broke the graph into two roughly equally sized connected components using a small set of cuts, then it might be easier to visualize.

If your graph is indeed a protein, try cutting the disulfide bonds. Then using the 3D spring model should open up the visualization.

@Anthrazit Surely anything that breaks the kernel connection is a serious bug.

What happens if you put the whole filename in libname?

I said that it would be possible to construct the polynomial (G0^m)[1,1] by modular imaging and polynomial interpolation. Here I present code that does this. The problem is that this code is slower than direct computation by powering the matrix of polynomials. But, I present it in the hope that someone may be able to suggest a way to make it faster.

Domineering:= proc(
    m::posint, n::posint, x::algebraic, y::algebraic, 
    {method::identical(naive, modular):= naive}
)
uses LAM= LinearAlgebra:-Modular, CF= CurveFitting;
local G0:= <<1>>, G1:= <<0>>, Z:= <<0>>;
    to n do
        (G1,G0,Z):= (<y*G0, Z; Z, Z>, <G0+G1, x*G0; G0, Z>, <Z, Z; Z, Z>)
    od;

    # Naive  method
    if [x,y]::[name$2] implies method=naive then return (G0^m)[1,1] fi;

    #-----------------------
    # Modular Imaging Method
    #-----------------------

    local 
        #-d..d will be the interpolation points (based on degree of result in x and y).
        d:= ceil(m*iquo(n,2)/2),

        #Need primes whose product is greater than the maximum possible polynomial
        #evaluation in the [1,1] position of the mth power.
        p:= 2^25, #upperbound of float[8] primes
        bits:= ilog2(p)-1, `2^n`:= 2^n, `2^m/2`:= 2^m/2,
        primes:= [
            to ceil(ilog2((`2^n`*eval(G0[1,1], [x,y]=~ d))^`2^m/2`/`2^n`)/bits) do 
                p:= prevprime(p) 
            od
        ],
        primes1, i, j, G0ij, Pxy:= table(), #integer interpolation points
        st:= time(), ps:= 0 #efficiency info (userinfo only)
    ;
    :-`mod`:= mods; #for chrem
    userinfo(
        1, Domineering, 
        "Using d =", d, ",", (2*d+1)^2, "evaluation points, and",
        nops(primes), "primes max." 
    ); 
    #
    #This loop is what I want to make more efficient. Everything else is fine.
    #
    for i from -d to d do  for j from -d to d do
        G0ij:= eval(G0, [x,y]=~ [i,j]);
        primes1:= primes[
            ..ceil(ilog2((`2^n`*max(2,abs(G0ij[1,1])))^`2^m/2`/`2^n`)/bits)
        ];         
        userinfo(1, Domineering, NoName, proc() ps+= nops(primes1); [][] end proc());
        userinfo(2, Domineering, 'i'=i, j='j', 'nops(primes1)'=nops(primes1));
        Pxy[i,j]:= chrem(
            [seq](
                trunc(LAM:-MatrixPower(p, LAM:-Mod(p, G0ij, float[8]), m)[1,1]), 
                p= primes1 
            ),
            primes1
        )
    od od;
    userinfo(
        1, Domineering, 
        "Main loop cputime =", time()-st, "with", ps, "evaluations."
    );
    
    local Px:= table(); #interpolants wrt y
    for i from -d to d do
        Px[i]:= CF:-PolynomialInterpolation([$-d..d], [seq(Pxy[i,j], j= -d..d)], y)
    od;

    #final interpolation, wrt x
    CF:-PolynomialInterpolation([$-d..d], [seq(Px[i], i= -d..d)], x)
end proc
:         
#Example usage:
p1:= CodeTools:-Usage(Domineering(6,6,x,y));
memory used=51.45MiB, alloc change=0 bytes, cpu time=235.00ms, 
real time=237.00ms, gc time=0ns
                 [long polynomial omitted]
infolevel[Domineering]:= 1:
p2:= Domineering(6,6,x,y, method= modular);
Domineering: Using d = 9 , 361 evaluation points, and 22 primes max.
Domineering: Main loop cputime = 22.078 with 5757 evaluations.
                 [long polynomial omitted]
simplify(p1 - p2); 
                 0

@acer At first, I used the approach that you suggest. The reason that I went with op(0,a) is that it leaves only one occurence of A in the whole code, making it easier to change that to another name. Of course, the name could be made a procedure parameter also.

Yes, I realize that your way is better if we're just dealing with a list of As. I mentioned that a few Replies back in the last paragraph. 

@emendes 

The type specindex(A) means any indexed name that starts with A. The spec is short for specific. 

Let a be an indexed name, for example, A[ j1, j2 ]. Then op(0,a) is A, op(1,a) is j1, and op(2,a) is j2.

@emendes I was hoping for the desired Maple output, but it's a moot point now.