# Generating a matrix with specific number of 1's in every row and column

I would like to have a function f(n,k) that generates a random nxn matrix with k elements equal to 1 in each row and column and has an all-zero diagonal, with all other elements equal to 0. In other words, I would like to generate an adjacency matrix for a random regular graph. An example of such a matrix would be the 5x5 matrix with k=3:

I know there is the RegularGraph module built into Wolfram but I would like to generate the matrix to make a graph from it, not the other way around.

• What do you mean by "the RegularGraph module built into Wolfram"? Jan 10, 2022 at 12:54
• @Szabolcs This Jan 10, 2022 at 13:24
• That function is part of Combinatorica, a package that has not been supported for many years. Its sampling is far from uniform. See the Combinatorica book for details. Jan 11, 2022 at 15:19

ClearAll[f]

f[n_, k_] := Module[{r = PadRight[ConstantArray[1, k], n - 1],
ca = ConstantArray[k, n], mi},
While[Total[mi = MapIndexed[Insert[#, 0, #2[[1]]] &,
Table[RandomSample[r], n]]] != ca]; mi]


Examples:

SeedRandom[1]
Row[Grid[Through[{MatrixForm,
AdjacencyGraph[#, VertexLabels -> Automatic, ImageSize -> Small] &}@#] & /@
Table[f[##], 4]] & @@@ {{5, 2}, {5, 3}, {6, 3}},
"\[VerticalSeparator]"]


• Wow, this is very concise and it works, although I don't understand what Table[RandomSample[r], n]]] != ca]; mi] means. For instance, where does the r come from? Jan 10, 2022 at 13:34
• @Heiko, r is defined in the first argument of Module as r = PadRight[ConstantArray[1, k], n - 1] (it is an (n-1)-vector with k 1s) . In Table[RandomSample[r], n]] we take n random shufflings of r to get an (n-1)Xn matrix. We insert to this matrix a 0 diagonal to get an nXn matrix. By construction each row of this matrix has k ones. Then we loop until we get a matrix mi whose column sums (Total[mi]) is an n-vector of ks.
– kglr
Jan 10, 2022 at 13:49
• Thank you, that clears it all up Jan 10, 2022 at 15:47
• @Heiko It's concise because it's trivial rejection sampling. Try f[50,3]. The solutions I suggested handle this easily, but naive rejection sampling like this will not finish in a lifetime. Jan 11, 2022 at 14:25

My package IGraph/M has several functions for this.

But before you get started, think carefully about what you mean by "random". I am going to assume that you want to generate each such matrix (i.e. each corresponding simple labelled graph) with equal probability, i.e. you want to do uniform sampling. This is not easy to do. There are two main methods to solve this problem, the configuration model and edge-switching based MCMC sampling. You will find an easy to understand overview of the topic in section 2.1 of this paper.

If you don't require uniform sampling, think about what effects sampling bias will have on your results.

In IGraph/M, you could useIGDegreeSequenceGame[ConstantArray[k,n], Method -> "ConfigurationModelSimple"] for the configuration model. Note that the other methods available in this function do not fit your problem: "FastSimple" does not sample uniformly, "ConfigurationModel" produces multigraphs and "VigerLatapy" produces connected graphs only (also uniformity is not excellent). This method works very well for small k, but for large k it will quickly become infeasibly slow.

For MCMC sampling, first we generate one arbitrary graph with the given degrees:

g = IGRealizeDegreeSequence[ConstantArray[k, n]]


Then we do a large number of random edge switches: IGRewire[g, 10 EdgeCount[g]].

This method will be much faster for large k, but you need to be careful to perform a sufficiently large number of edge switches to get uniform sampling.

Once you have the graph, use AdjacencyMatrix to get a sparse adjacency matrix from it, then use Normal to convert it to a dense matrix.

Note that while Mathematica has a function to sample graphs with given degrees, RandomGraph[DegreeGraphDistribution[...]], the sampling is not exactly uniform, and the method that this function uses is not documented. Therefore I do not trust it.

Update: If you are looking for directed rather than undirected graphs, use the syntaxes:

With[{deg = ConstantArray[k,n]},
IGDegreeSequenceGame[deg, deg, Method -> "ConfigurationModelSimple"]]
]

IGRewire[
WithIGRealizeDegreeSequence[ConstantArray[k, n], ConstantArray[k, n]],
10*n*k
]


You can use FindInstance:

f[n_, k_, numsolutions_, seed_] := With[{
vars = Array[v, {n, n}]},
With[{
constraints = And @@ Flatten[{
Total[#] == k & /@ vars,
Total[#] == k & /@ Transpose[vars],
Total[Diagonal[vars]] == 0,
# <= 1 & /@ Flatten[vars]
}]},
vars /.
FindInstance[constraints, Flatten@vars, NonNegativeIntegers,
numsolutions, RandomSeeding -> seed]
]]

(* usage *)
(mtx = First@f[5, 3, 1, 123456]) // MatrixForm


... or more efficiently with LinearOptimization:

f[n_, k_] := With[{vars = Array[v, {n, n}]},
With[{constraints =
Join[Total[#] == k & /@ vars,
Total[#] == k & /@ Transpose[vars],
{Total[Diagonal[vars]] == 0},
0 \[VectorLessEqual] # \[VectorLessEqual] 1 & /@
Flatten[vars], # \[Element] Integers & /@ Flatten[vars]]},
vars /. LinearOptimization[0, constraints, Flatten@vars]
]
]

(* usage *)
(mtx = f[5, 3]) // MatrixForm

• OP wants a random regular graph though. Jan 10, 2022 at 12:55
• @Szabolcs as you said in your answer ... think carefully about what you mean by "random" ... it's kind of hard to make a random matrix that is so constrained and the question is a bit vague about which kind of 'random' so I'm just hoping that FindInstance/LinearOptimization is 'random enough' for OP's purposes but welcome any suggestions. Jan 10, 2022 at 13:05