Laplacian matrix, solving a resistor mesh problem

Question

I'm trying to solve for flows on a network using a laplacian matrix. I started out by testing on this problem here: https://rosettacode.org/wiki/Resistor_mesh#Python

And it comes out with the solution perfectly, when all weights are 1, R = 1.6089whatever. I want to be able to solve for resistors that don't equal one though! then be able to get the currents on each resistor, so I tried producing a set of random values for the conductances (weights) and to check it was working see that the sum of the currents coming from the node where the current is injected is equal to the current injected as it should be. Unfortunately it was not. I've inspected the Laplacian and it looks about how I was expecting, but beyond that I'm totally lost, can anyone shed any light on whether I've even got the right laplacian here? or if I'm missing something really obvious such as some incorrect indexing?

Code is below (excuse the very amateur way I've generated the grid, any tips and tricks on that appreciated too!):

library(igraph)
library(SparseM)

# setup edgelist + grid to look at
sample_grid = function(N_x, N_y, xlim = c(-1,1), ylim = c(-1,1)) {

  
  # bounding box 
  min_x <- xlim[1]
  max_x <- xlim[2]
  min_y <- ylim[1]
  max_y <- ylim[2]

  
  x_locs <- seq(0, N_x)
  y_locs <- seq(0, N_y)
  
  N_grid <- N_x * N_y
  
  x <- rep(0, length(N_grid))
  y <- rep(0, length(N_grid))
  
  for(i in 1:N_x){
    for(j in 1:N_y) {
      x[N_y * (i-1) + j] <- x_locs[i]
      y[N_y * (i-1) + j] <- y_locs[j]
    }
  }
  
  locations <- cbind(x,y)
  
  return(locations)
}

N_x <- 10
N_y <- 10

# node locations
Vert <- sample_grid(N_x,N_y)

# edges...

# horizontally...
fromH <- c()
for(i in 1:(N_x-1)){
  fromH <- c(fromH, seq(0,(N_y*N_x-N_x), length.out = (N_y)) + i)
}
toH <- fromH + 1

# vertically...
fromV <- 1:(N_x*N_y-N_x)
toV <- fromV + N_x

# --------------------------------------------------------------------------------------------
# crux

Edges <- data.frame(from = c(fromH, fromV), to = c(toH, toV)) # , weights = rep(1,(2*N_x*N_y - (N_x+N_y))))     

# change the weights up
set.seed(1)
weight <- rlnorm(n = nrow(Edges), meanlog = 0, sdlog = 1)
Edges$weight <- weight

the_graph <- graph_from_data_frame(Edges, directed = FALSE)
lo <- layout.norm(as.matrix(Vert))
plot(the_graph, layout = lo, directed = FALSE, edge.arrow.size=0)


# solving
L <- laplacian_matrix(the_graph, weights = weight)

# boundary conditions on 68, 12: 
# draw 1 amp @ 12
# inject 1 amp @ 68

q <- matrix(rep(0, nrow(Vert)), ncol = 1)
q[68,] <- +1
q[12,] <- -1

# solve
p <- solve(L,q)


R <- p[68,] - p[12,]


# investigating why the weights aren't working! (wrong first attempt)
# neighbours <- c(78,58,67,69)  # neighbours of node 68
# neighbours_weights <- weight[neighbours]
# neighbours_potential_diffs <- p[68,] - p[neighbours,]
# neighbours_currents <- neighbours_potential_diffs * neighbours_weights

neighbours <- which(Edges$from == 68 | Edges$to == 68)
potential_diffs <- p[Edges$from,] - p[Edges$to,]
currents <- potential_diffs * Edges$weight
what <- cbind(Edges,p[Edges$from,], p[Edges$to,], currents)
what[neighbours,]



# exact solution for effective resistance between 68 and 12 with 10x10 and all 1ohm
exact <- 455859137025721/283319837425200

The problem being that the calculated currents do not sum to one as expected:

> neighbours_currents
[1] 0.05035087 0.09044874 0.03309549 0.21782845

Answer 1

Answer: the function laplacian_matrix() reorders the nodes in descending degree, which meant that all the nodes were jumbled up.

The fix was to specify the sequence of node labels so

the_graph <- graph_from_data_frame(Edges, directed = FALSE)

became

Verts <- data.frame(label = 1:(N_x*N_y))
the_graph <- graph_from_data_frame(Edges, directed = FALSE, vertices = Verts)

giving currents that sum to one:

> what[neighbours,]
    from to    weight p[Edges$from, ] p[Edges$to, ]    currents
67    67 68 0.1644813       0.2582772     0.8279230 -0.09369607
77    68 69 0.6419198       0.8279230     0.4943076  0.21415437
148   58 68 1.0175478       0.3902397     0.8279230 -0.44536367
158   68 78 0.5372635       0.8279230     0.3685843  0.24678589