Finding all paths of certain length between vertices in path graph, non-simple ones too

Question

First I will tell you what my ultimate goal is, because maybe there's another way of doing it other than what I am trying.

I need a list of all $(a_1, a_2, \ldots, a_N)$ tuples that satisfy certain conditions.

• $a_i$ can take integer values between $1$ and $s$, that is, from the set $\{1, \ldots, s\}$.
• Neighbouring entries differ by $1$, that is, $|a_i - a_{i+1}| = 1$
• A "periodic boundary condition" holds. That means that $a_N$ and $a_1$ also differ by $1$ (we effectively have a ring).
• $N$ can be assumed to be even, because otherwise it's not possible to satisfy the previous two conditions at the same time.

For $N=6, s=5$, I've put $a_1$ at the back of the list in parentheses for clarity.

$2,3,4,5,4,3(,2)$ works.
$2,3,4,5,4,5(,2)$ doesn't, because $5$ and $2$ don't differ by one.
$2,4,3,4,3,3(,2)$ obviously doesn't work, because some pairs of neighbouring entries don't differ by $1$.

OK, I've been thinking that I could take PathGraph[Range[1,s]], which is a builtin Graph with 1 being adjacent to 2, 2 being adjacent to 3 (and 1), 3 being adjacent to 4 (and 2), ..., s-1 being adjacent to s (and s-2).

Now I can easily get the number of tuples with the above properties by calculating the Nth MatrixPower of the AdjacencyMatrix and taking the Trace of that matrix, because the entries on the diagonal are Paths of length N, starting at 1, ending at 1, (plus the ones starting at 2 ending at 2, plus ...) possibly going back and forth between vertices.

However, I need to have the tuples to work with them, so I wanted to get a list of paths of length N between two vertices of that simple path graph. But I can't seem to get that with the builtin functions.

FindPath[PathGraph[Range[1, 7]], 1, 3, {2}] will give the path 1 -> 2 ->3 which has length 2.

However, FindPath[PathGraph[Range[1, 7]], 1, 3, {4}] will give an empty set. It won't give me the path 1 -> 2 -> 3 -> 4 -> 3, which visits nodes twice and goes back and forth along the same edge. Is there some way to get all of these? (the mentioned one and, say 1 -> 2 -> 3 -> 2 -> 3).

Of course I've read numerous posts on paths on these forums and a couple of entries about functions concerning graphs in the Mathematica documentation, but I couldn't seem to find what I need.

I'd expect FindPath[PathGraph[Range[1, 7]], 2, 2, {2}] (same starting and ending point) to give me 2 -> 3 -> 2 and 2 -> 1 -> 2. Which would get me the tuples $2,3(,2)$ and $2,1(,2)$, speaking the language from the beginning of the post.

Answer 1

One can use this small piece of code

f[s_, n_?EvenQ] := Flatten[Outer[Join, #[[;; , -2 ;; 2 ;; -1]], #, 1] & /@ 
 GatherBy[#, #[[{1, -1}]] &] &@ Nest[Flatten[{If[# > 1, {# - 1, ##}, ## &[]], 
     If[# < s, {# + 1, ##}, ## &[]]} & @@@ #, 1] &, Transpose@{Range@s}, n/2], 2]
f[__] = {};

f[3, 4]

{{2, 1, 2, 1}, {2, 3, 2, 1}, {1, 2, 1, 2}, {1, 2, 3, 2}, {3, 2, 1, 2}, {3, 2, 3, 2}, {2, 1, 2, 3}, {2, 3, 2, 3}}

f[3, 5] (* no paths *)

{}

f[10, 20] // Length // AbsoluteTiming

{0.526044, 983740}

It is proper number of paths and time is not so big for such parameters.

What's going on? It's not a code-golf challenge!

There is explanation of the code:

• The main idea: generate all possible paths for n/2 without periodic conditions. Then combine all possible paths with the same start and the same end.

• Add one appropriate element. I go from the right to the left

  Flatten[{If[# > 1, {# - 1, ##}, ## &[]], If[# < s, {# + 1, ##}, ## &[]]} & @@@ #, 1] &
  

• Do it n/2 times

  Nest[..., Transpose@{Range@s}, n/2]
  

• Gather by the same first and the last elements

  GatherBy[#, #[[{1, -1}]] &] &@ 
  

• Combine them (Outer is faster then the previous version with Tuples)

  Outer[Join, #[[;; , -2 ;; 2 ;; -1]], #, 1] &
  

• Flatten the result

  Flatten[..., 2]
  

Answer 2

Here's my submission in the category 'shortest code':

AllTuples[s_, n_?EvenQ] := 
  Map[
    Sequence @@ Outer[Plus, Range[1 - Min@#, s - Max@#], #]& @ Accumulate @ # &, 
    Permutations @ Flatten @ ConstantArray[{-1, 1}, n/2]
  ]

This is based on the observation that desired list of tuples without the first condition is isomorphic to all permutations of {1,1,1,-1,-1,-1} for n = 6 (and generalisations thereof for other n). Taking your first example, we namely have explicitly:

{2, 3, 4, 5, 4, 3, 2} // Differences

{1, 1, 1, -1, -1, -1}

So all we need to do is compute all permutations of this list (done with Permutations), do the inverse of Differences (that's more or less Accumulate), and then add the first condition (that all element of the tuple are members lie in a given range) by hand (the Outer bit).

As a small check, it gives the correct output for the case s = 3, n = 4:

AllTuples[3, 4]

{
  {2, 3, 2, 1}, {2, 1, 2, 1}, {3, 2, 3, 2}, {3, 2, 1, 2}, 
  {1, 2, 3, 2}, {1, 2, 1, 2}, {2, 3, 2, 3}, {2, 1, 2, 3}
}

And it gives the same (number of) results as ybeltukov's code for the case s = 10, n = 20:

AllTuples[10, 20] // Length

983740

For this case it's about 2x slower than ybeltukov's code though.

Answer 3

For smallish values of n one can use integer programming. Below is one way to code that into Reduce. I don'tpect it to be competitively fast, as written, for n>10 or so.

lacedTuples[n_] := Module[
  {a, vars, fvars, c1, c2, c3, c4, c5, constraints, adjacents, 
   wrappedvars, colvals, coldiffs, soln},
  vars = Array[a, {n, n}];
  fvars = Flatten[vars];
  c1 = Map[0 <= # <= 1 &, fvars];
  c2 = Thread[Total[vars] == 1];
  adjacents = 
   Flatten[Map[Partition[PadRight[#, n + 1, First[#]], 2, 1] &, vars],
     1];
  c3 = Map[Total[#] <= 1 &, adjacents];
  wrappedvars = PadRight[vars, {n + 1, n + 1}, vars];
  colvals = Range[n].vars;
  coldiffs = Differences[PadRight[colvals, n + 1, colvals[[1]]]];
  c4 = Map[-1 <= # <= 1 &, coldiffs];
  constraints = Flatten[Join[c1, c2, c3, c4]];
  soln = Solve[constraints, fvars, Integers];
  If[Head[soln] === Solve, {},
   colvals /. Solve[constraints, fvars, Integers]]
  ]
lacedTuples[6]

(* Out[263]= {{6, 5, 6, 5, 6, 5}, {5, 6, 5, 6, 5, 6}, {5, 6, 5, 6, 5, 
  4}, {6, 5, 6, 5, 4, 5}, {5, 6, 5, 4, 5, 6}, {5, 6, 5, 4, 5, 4}, {6, 
  5, 4, 5, 6, 5}, {6, 5, 4, 5, 4, 5}, {5, 4, 5, 6, 5, 6}, {5, 4, 5, 6,
   5, 4}, {5, 4, 5, 4, 5, 6}, {5, 4, 5, 4, 5, 4}, {4, 5, 6, 5, 6, 
  5}, {4, 5, 6, 5, 4, 5}, {4, 5, 4, 5, 6, 5}, {4, 5, 4, 5, 4, 5}, {4, 
  5, 6, 5, 4, 3}, {4, 5, 4, 5, 4, 3}, {5, 6, 5, 4, 3, 4}, {5, 4, 5, 4,
   3, 4}, {6, 5, 4, 3, 4, 5}, {4, 5, 4, 3, 4, 5}, {4, 5, 4, 3, 4, 
  3}, {5, 4, 3, 4, 5, 6}, {5, 4, 3, 4, 5, 4}, {5, 4, 3, 4, 3, 4}, {4, 
  3, 4, 5, 6, 5}, {4, 3, 4, 5, 4, 5}, {4, 3, 4, 5, 4, 3}, {4, 3, 4, 3,
   4, 5}, {4, 3, 4, 3, 4, 3}, {3, 4, 5, 6, 5, 4}, {3, 4, 5, 4, 5, 
  4}, {3, 4, 5, 4, 3, 4}, {3, 4, 3, 4, 5, 4}, {3, 4, 3, 4, 3, 4}, {3, 
  4, 5, 4, 3, 2}, {3, 4, 3, 4, 3, 2}, {4, 5, 4, 3, 2, 3}, {4, 3, 4, 3,
   2, 3}, {5, 4, 3, 2, 3, 4}, {3, 4, 3, 2, 3, 4}, {3, 4, 3, 2, 3, 
  2}, {4, 3, 2, 3, 4, 5}, {4, 3, 2, 3, 4, 3}, {4, 3, 2, 3, 2, 3}, {3, 
  2, 3, 4, 5, 4}, {3, 2, 3, 4, 3, 4}, {3, 2, 3, 4, 3, 2}, {3, 2, 3, 2,
   3, 4}, {3, 2, 3, 2, 3, 2}, {2, 3, 4, 5, 4, 3}, {2, 3, 4, 3, 4, 
  3}, {2, 3, 4, 3, 2, 3}, {2, 3, 2, 3, 4, 3}, {2, 3, 2, 3, 2, 3}, {2, 
  3, 4, 3, 2, 1}, {2, 3, 2, 3, 2, 1}, {3, 4, 3, 2, 1, 2}, {3, 2, 3, 2,
   1, 2}, {4, 3, 2, 1, 2, 3}, {2, 3, 2, 1, 2, 3}, {2, 3, 2, 1, 2, 
  1}, {3, 2, 1, 2, 3, 4}, {3, 2, 1, 2, 3, 2}, {3, 2, 1, 2, 1, 2}, {2, 
  1, 2, 3, 4, 3}, {2, 1, 2, 3, 2, 3}, {2, 1, 2, 3, 2, 1}, {2, 1, 2, 1,
   2, 3}, {2, 1, 2, 1, 2, 1}, {1, 2, 3, 4, 3, 2}, {1, 2, 3, 2, 3, 
  2}, {1, 2, 3, 2, 1, 2}, {1, 2, 1, 2, 3, 2}, {1, 2, 1, 2, 1, 2}} *)

--- edit ---

Here is a related method, inspired by the (first) approach of @belisarius. We again use integer programming but avoid disjunction constraints as well as a plethora of 0-1 variables. We use n of them just to enforce the condition that successive differences be either -1 or 1.

lacedTuples2[n_, s_] := Module[
  {a, d, vars, dvars, padvars, c1, c2, c3, c4},
  vars = Array[a, {n}];
  dvars = Array[d, {n}];
  padvars = PadRight[vars, n + 1, First[vars]];
  c1 = Map[0 <= # <= 1 &, dvars];
  c2 = Map[1 <= # <= s &, vars];
  c3 = Table[
    padvars[[j]] - padvars[[j - 1]] == -1 + 2*dvars[[j - 1]], {j, 2, 
     n + 1}];
  c4 = Table[-padvars[[j]] + padvars[[j - 1]] == 
     1 - 2*dvars[[j - 1]], {j, 2, n + 1}];
  vars /. Solve[Join[c1, c2, c3], Join[vars, dvars], Integers]
  ]

Quick test.

lacedTuples2[4, 4]

(* Out[53]= {{1, 2, 1, 2}, {1, 2, 3, 2}, {2, 1, 2, 1}, {2, 1, 2, 3}, {2, 
  3, 2, 1}, {2, 3, 2, 3}, {2, 3, 4, 3}, {3, 2, 1, 2}, {3, 2, 3, 
  2}, {3, 2, 3, 4}, {3, 4, 3, 2}, {3, 4, 3, 4}, {4, 3, 2, 3}, {4, 3, 
  4, 3}} *)

Timing[lc12 = lacedTuples2[12, 12];]

(* Out[50]= {4.108000, Null} *)

Length[lc12]

(* Out[51]= 7916 *)

--- end edit ---

--- edit 2 ---

Here is a graph theory approach. I would post it as a separate answer but it gets very slow after n=10 or so. The idea is to create a directed graph where we have a "row" for each entry from 1 to n, and "columns" 1 through s. These together index the vertices. Edges connect row i to row i+1 with the condition that they only connect contiguous column numbers. Also connect the last row to the first, with the same contiguity condition. Now get cycles of length n+1 and extract the column values.

lacedTuples3[s_, n_] := Module[
  {a, n2 = Ceiling[n/2], verts, edges, gr, cycles},
  verts = Array[a, {n, s}];
  edges = 
   Flatten[Table[{DirectedEdge[a[i, j], a[i + 1, j + 1]], 
      DirectedEdge[a[i, j + 1], a[i + 1, j]]}, {i, n - 1}, {j, 
      s - 1}]];
  edges = 
   Join[edges, 
    Flatten[Table[{DirectedEdge[a[n, j], a[1, j + 1]], 
       DirectedEdge[a[n, j + 1], a[1, j]]}, {j, s - 1}]]];
  gr = Graph[edges];
  cycles = FindCycle[gr, n + 1, All];
  cycles[[All, All, 1, 2]]
  ]

lacedTuples3[4, 4]

(* ut[83]= {{3, 4, 3, 4}, {3, 4, 3, 2}, {3, 4, 3, 4}, {3, 2, 3, 4}, {2, 
  3, 4, 3}, {2, 3, 2, 3}, {2, 3, 2, 1}, {2, 1, 2, 3}, {2, 1, 2, 
  1}, {2, 3, 4, 3}, {2, 3, 2, 3}, {2, 1, 2, 3}, {1, 2, 3, 2}, {1, 2, 1, 2}} *)

--- end edit 2 ---

--- edit 3 ---

This is mostly to convince myself that the graph cycles approach was not entirely without merit.

extendCycle[cyc_List, edges_List] := 
 Map[If[# > First[cyc] && ! MemberQ[cyc, #], Append[cyc, #], 
    Null /. Null :> Sequence[]] &, edges[[Last[cyc]]]]

allCycles[mat_, k_] := 
 Module[{n = Length[mat], m2, cyc, cyclist}, 
  m2 = Map[Last, Split[Sort[mat], First[#1] == First[#2] &], {2}];
  cyclist = 
   Flatten[Drop[MapIndexed[{#2[[1]], #1} &, m2, {2}], -k + 1], 1];
  Do[cyclist = 
    Flatten[Map[extendCycle[#, m2] &, cyclist], 1], {k - 2}];
  Map[If[MemberQ[m2[[Last[#]]], First[#]], Append[#, First[#]], 
     Null /. Null :> Sequence[]] &, cyclist]]

lacedTuples4[s_, n_] := Module[
  {a, verts, edges, vertvals, cycles},
  verts = Array[a, {n, s}];
  vertvals = Range[n^2];
  edges = 
   Table[{{a[i, j], a[i + 1, j + 1]}, {a[i, j + 1], a[i + 1, j]}}, {i,
      n - 1}, {j, s - 1}];
  edges = 
   Partition[
     Flatten[Join[edges, 
       Table[{{a[n, j], a[1, j + 1]}, {a[n, j + 1], a[1, j]}}, {j, 
         s - 1}]]], 2] /. Thread[Flatten[verts] -> Range[n^2]];
  cycles = allCycles[edges, n] /. Thread[Range[n^2] -> Flatten[verts]];
  cycles[[All, All, 2]]
  ]

The code for allCycles was cribbed from a previous MSE thread. I need to ask hereabouts as to why the new in version 14 cycles code gets slow in this example; I had thought it was fairly fast on examples I tried in past. Even so, the above becomes orders of magnitude slower than the methods of @ybeltukov and others at n=12,s=12 or thereabouts.

By the way, I'm not really beating a dead horse on this; the dead horse is beating me.

--- end edit 3 ---

Answer 4

A recursive approach. Reasonably fast, but not comparable with ybeltukov's solution:

f[n_Integer, s_Integer, head1_List] :=
 Module[{newh},
  If[Length[First@head1] >= n, Return@head1,
   newh =
    Which[
       Last@# == s,                          Append[#, s - 1],
       Last@# == 1,                          Append[#, 2],
       n - Length@# == First@# - Last@# - 1, Append[#, 1 + Last[#]],
       n - Length@# == Last@# - First@# - 1, Append[#, -1 + Last[#]],
       True, 
       Sequence @@ {Append[#, 1 + Last@#], Append[#, -1 + Last[#]]}
       ] & /@ head1];
  f[n, s, newh]]

Timing[Flatten[Table[f[20, 10, {{k}}], {k, 1, 10}], 1] // Length]

(*{20.015625, 983740}*)

Edit

An ILP approach:

findTup[m_, s_] := Module[{diffs, bounds, cyclic, x}, 
  diffs  = And @@ Table[Abs[x[n + 1] - x[n]] == 1, {n, 1, m - 1}];
  bounds = And @@ Table[1 <= x[n] <= s, {n, 1, m}];
  cyclic = (Abs[x[1] - x[m]] == 1);
  Array[x, m] /. Table[Solve[And @@ {diffs, bounds, cyclic, x[1] == k},  Array[x, m]], {k, 1, s}]
  ] 
Timing[findTup[4, 3]]

(* {0.562500, {{{1, 2, 1, 2}, {1, 2, 3, 2}}, {{2, 1, 2, 1}, {2, 1, 2, 3}, 
               {{2, 3, 2, 1}, {2, 3, 2, 3}}, {{3, 2, 1, 2}, {3, 2, 3, 2}}}} *)

Answer 5

There is a way to do it the way you want, but it will not be very speedy. The main trick is redefining plus and times, thus:

myplus[l1__] := Union @@ {l1}
mymult[l1_, l2_] := Flatten[Outer[Join[#1, Drop[#2, 1]] &, l1, l2, 1],  1]
mymult[{}, l2_] := {}
mymult[l1_, {}] := {}

Now, you define new matrix multiplication, thus:

matmult[m1_, m2_, op1_, op2_] := Outer[Inner[op1, #1, #2, op2] &, m1, Transpose[m2], 1]

Note that this is just ordinary multiplication when op1 is times, and op2 is plus.

Now, you change your adjacency matrix, which you get by: AdjacencyMatrix[PathGraph[Range[3]]]

thus:

hackadj[m_] := MapIndexed[If[#1 == 1, {#2}, {}] &, m, {2}]

Finally, the union of the diagonal elements of the $n$th power of the hacked adjacency matrix will be what you want. However, note that this will use exponential space (you will be storing ALL the paths in the graph of length up to your bound), so the above should be viewed mostly as an exercise in functional programming.

Answer 6

ok, I need to thank all of you for the effort. ybeltukov's version works and is fast. I will add another answer because i came up with another solution (for the fun of it). Needless to say it isn't as good as ybeltukov's of course (which looks like a Perl-Script to me being not very familiar with Mathematica yet, i've still got to work through it and the other examples to fully understand). ;)

My solution takes about 56 times as much time for the f[10,20] example on my machine. ;) So it's already obvious that it isn't as good.

However maybe it helps someone reading this in the future.

What i did was:

I rephrased the Problem. Instead of

a) finding paths of the PathGraph[Range[1,s]] of length $N$, which may go back and forth along edges and repeat vertices multiple times

equivalently one can do

b) For a start take a (quadratic) GridGraph[{N/2+1, N/2+1}]. Now Any shortest Path from one corner to the (diagonally) opposite corner will have length $N$ (you need to traverse $\frac{N}{2}$ horizontal plus $\frac{N}{2}$ vertical edges to get there), which is what "we" need (i guess only "I" need it though).

However, there are too many vertices in this Graph (for general s and N at least), I will explain why.

For each $k \in \{1, \ldots, s\}$ (basically for each choice of $a_1$) I will associate vertices on the diagonal connecting the two corners of having value $k$, and i will associate vertices on neighbouring diagonals with values $k+1$ (above) and $k-1$ (below) (and so on, $k+2$, $k+3$ for the other diagonals above the main diagonal).

Now I can only allow those diagonals that have a value within $\{1, \ldots, s\}$, so I need to remove the diagonals that have a value outside of this range, i.e. I remove the vertices on those diagonals and the edges connecting those vertices to the rest.

Now once I did this, if I traverse the Graph from one corner to the opposite corner (a path of $N+1$ vertices and $N$ edges) and take the values that I assigned to each vertex on the way I will get the required tuples (once I drop the very last one, which will be the same as the first one, due to the requirement of periodic boundary conditions). I have to do this $s$-times, for each $k \in \{1, \ldots, s\}$ and collect the sets of Paths into one large set.

Example: For $s=3$ and $a_1 = 2 = k$ as a starting point, I only have the main diagonal and the two adjacent diagonals. The tuples I get will look like $2,3,2,1,3,2,3,1,..(,2)$.

Putting that into practice in Mathematica, the difficult part (for me) was constructing the graph (given $k$, $s$ and $N$), and i think it's more convenient to have the vertices in a two-dimensional parametrisation.

1. I constructed the set of Vertices, diagonal by diagonal.

   vert[k_] := Flatten[Table[Table[{j, j + i}, {i, Max[1 - k, 1 - j], Min[N/2 + 1 - j, s - k]}], {j, 1, N/2 + 1}], 1]
   

i is the deviation from the diagonal, or equivalently k+i is the value and it needs to fulfill the two conditions $1 \le k+i \le s$ (allowed values) and $1 \le i + j \le N/2 + 1$ (stay within the square, i+j is the $y$-component) simultaneously.

2. I added edges where two vertices differed by 1 in one coordinate and 0 in the other (took the $p$-Norm with $p=1$).

   edg[k_] := Cases[Tuples[vert[k], 2], {x_, y_} /; (Norm[x - y, 1] == 1 && ((x - y)[[1]] > 0 || (x - y)[[2]] > 0)) -> x <-> y]
   

For each pair of vertices (element of the Cartesian product of the set of vertices with itself), I add an edge whenever the $1$-Norm of the difference of two vertices is 1, i.e. two vertices are nearest neighbours in the square lattice. The additional condition behind the && is to rule out duplicates.

3. I can now define the Graph with proper VertexLabels and VertexCoordinates, due to convenient choice of naming.

   grap[k_] := Graph[vert[k], edg[k], VertexLabels -> "Name", VertexCoordinates -> vert[k]]
   

4. Now as I said, the shortest path between two opposite corners has length $L$, so if I limit the maximum length of a path to L and tell it to give me All paths between the upper-right and lower-left (or whichever direction you're looking at it).

   FindPath[grap[k], {1, 1}, {L/2 + 1, L/2 + 1}, L, All]
   

does what I need. I just need to translate the vertex coordinates/names into the assigned value which is done by taking the difference between the two coordinates and adding k (i also drop the last vertex of each path):

    pat[k_] := (#[[2]] - #[[1]] + k) & /@ (#[[1 ;; -2]]) & /@ FindPath[grap[k], {1, 1}, {L/2 + 1, L/2 + 1}, L, All]

5. All that is left, is to Join all the paths for all the k's into one list. In fact at this point one can utilize a symmetry of this problem and calculate pat[k] only for half of the k's because you get the others by mirroring along the diagonal.

   firsthalf := Table[pat[k], {k, 1, Ceiling[r/2]}]
   secondhalf := Table[r - firsthalf[[r - k]], {k, Ceiling[r/2] + 1, r - 1}]
   Join[Flatten[firsthalf, 1], Flatten[secondhalf, 1]]
   

As I said, this is a lot slower (and probably reflects that I'm still a novice with this software). I guess the Cartesian Product is kind of slow because you generate all the unnecessary pairs only to discard those that don't meet the criteria. Then you search the resulting graph for paths, when they in fact do have a kind of simple structure, which you could "just write down", because you know what they look like.

I'm thankful for other notes, why this isn't a smart way to do things. (Although I will probably just use ybeltukov's solution.) Also I hope I got the formatting right this time, I wasn't aware in the first post, thanks for reformatting it.

The full code

vert[k_] := Flatten[Table[Table[{j, j + i}, {i, Max[1 - k, 1 - j], Min[N/2 + 1 - j, s - k]}], {j, 1, N/2 + 1}], 1]
edg[k_] := Cases[Tuples[vert[k],2], {x_, y_} /; (Norm[x - y, 1] == 1 && ((x - y)[[1]] > 0 || (x - y)[[2]] > 0)) -> x <-> y]
grap[k_] := Graph[vert[k], edg[k], VertexLabels -> "Name", VertexCoordinates -> vert[k]]
pat[k_] := (#[[2]] - #[[1]] + k) & /@ (#[[1 ;; -2]]) & /@ FindPath[grap[k], {1, 1}, {L/2 + 1, L/2 + 1}, L, All]
firsthalf := Table[pat[k], {k, 1, Ceiling[(s+1)/2]}]
secondhalf := Table[s+1 - firsthalf[[s+1 - k]], {k, Ceiling[(s+1)/2] + 1, s}]
Join[Flatten[firsthalf, 1], Flatten[secondhalf, 1]]

(small edit, as i had r in it, which is s+1)