Mahrukh
Mahrukh

Reputation: 11

Energy of graphs

I am working on R for generating graph energy of rank-1 heterogeneous random graph in which edges sharing same nodes have 0.5 correlation and edges with different nodes are mutually exclusive. EN1 is the energy of graph but gives all 0 values. Where is the mistake?

library(MASS)
#library(EnvStats)
#library(energy)
library(mvtnorm)

n <- 3
Pn <- 0.1
kappa <- 1
r <- 0.5

M <- 20

EN1 <- numeric(M)


for (rep in 1:M) {
  A <- matrix(0, n, n)
  for (i in 1:(n - 1)) {
    for (j in (i + 1):n)
      for (k in (j + 1):n)
        for (l in (k + 1):n)
          A[i, j] <- rbinom(1, 1, Pn*exp(-kappa*i/n)*exp(-kappa*j/n))
    A[k, l] <- rbinom(1, 1, Pn*exp(-kappa*k/n)*exp(-kappa*l/n))
  }
  A <- A + t(A)
  A
}


for (rep in 1:M) {
  cov_matrix <- matrix(r, nrow=n, ncol=n)
  diag(cov_matrix) <- 0
  for (i in 1:(n - 1)) {
    for (j in (i + 1):n)
      for (k in (j + 1):n)
        for (l in (k + 1):n)
          B <- mvtnorm::rmvnorm(n, sigma=cov_matrix)
    B[i, j] <- rbinom(1, 1, Pn*exp(-kappa*i/n)*exp(-kappa*j/n))
    B[k, l] <- rbinom(1, 1, Pn*exp(-kappa*k/n)*exp(-kappa*l/n))
  }
  B <- B + t(B)
  B
}


C <- A + B
C <- C + t(C)
C

ei1 <- eigen(C, symmetric=TRUE, only.values=TRUE)
ei1$values
EN1[rep] <- sum(abs(ei1$values))


EN1
mean(EN1)  # report mean
sd(EN1)    # report standard deviation

Upvotes: 0

Views: 42

Answers (1)

Edward
Edward

Reputation: 19339

It looks like you want to put this code block:

C <- A + B
C <- C + t(C)

ei1 <- eigen(C, symmetric=TRUE, only.values=TRUE)
EN1[rep] <- sum(abs(ei1$values))

within the rep for loop, like so:

for (rep in 1:M) {
  cov_matrix <- matrix(r, nrow=n, ncol=n)
  diag(cov_matrix) <- 0
  for (i in 1:(n - 1)) {
    for (j in (i + 1):n)
      for (k in (j + 1):n)
        for (l in (k + 1):n)
          B <- mvtnorm::rmvnorm(n, sigma=cov_matrix)
    B[i, j] <- rbinom(1, 1, Pn*exp(-kappa*i/n)*exp(-kappa*j/n))
    B[k, l] <- rbinom(1, 1, Pn*exp(-kappa*k/n)*exp(-kappa*l/n))
  }
  B <- B + t(B)
  C <- A + B
  C <- C + t(C)
  ei1 <- eigen(C, symmetric=TRUE, only.values=TRUE)
  EN1[rep] <- sum(abs(ei1$values))
}

This gives:

> EN1
 [1]  7.076329  1.880778  7.902730  2.055223  4.731054  6.326656 10.886720  7.226352
 [9]  7.703960  8.511892  6.168880  4.670598  3.709914  5.671895  1.671177 15.059931
[17] 11.200192  4.440183  3.822154  3.661713

> mean(EN1)  # report mean
[1] 6.218917
> sd(EN1)    # report standard deviation
[1] 3.405246

Upvotes: 0

Related Questions