@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.

@emendes Here's a simplifcation. This works because the operation on the 1st index commutes with the operation on the 2nd index.

(T1,T2):= (table([2,3]=~ [3,2]), table([2,3,5,6,7,9]=~ [3,2,6,5,9,7])):
T:= (T::table, x)-> `if`(assigned(T[x]), T[x], x):
subsindets(varA, specindex(A), a-> op(0,a)[T(T1, op(1,a)), T(T2, op(2,a))]);

This is also easier to understand, I think. The 2nd line says that the transformation is the identity except for indices listed in the respective table.

If your actual-use cases are simply lists of names that all begin with A, then you might as well use Acer's solution. I was assuming that the actual-use cases would have the A names mixed into other expressions. In that case, you need subsindets or evalindets.

@santamartina If mapping over a list, there will be one result for every member of the list.

Your code could be replaced by this recursive procedure:

p:= (A::Array, x::posint)-> 
    ArrayTools:-Append(`if`(x > 1, thisproc(A, x-1), A), A[-1]+1)
:
p(Array([4]), 3);