- Self-avoiding pathというのがある
- 1度通過した点を2度以上訪れない軌跡
- 酔歩的に行うことがある
- 今、それのモデルとして、n次元空間上の乱点を三角化したグラフ上のSelf-avoiding pathをシミュレーションすることはできる
- その軌道をk=1次元のSelf-avoiding 多様体であると考えることにする
- さて
- このようなランダムに広がるk次元多様体をその内部と辺縁とで考えることにする。k=1の場合の辺縁とは、パスの両端のノードのこと
- このノードのどちらかから、エッジを選んでk=1次元多様体を伸ばす。そのときにSelf-avoidする
- 多様体次元を上げても、内部と辺縁にするということは変わらない
- 辺縁を構成するのが、k-1次元単体の集合になるだけ
- Self-avoidingのルールはどうするか
- 内部に取り込まれたk-1次元単体が接続するk次元単体は、もう選んではいけない。これがSelf-avoiding のルール
- k=1の場合には、それは、軌跡の両端ではないノードに接続するエッジは選べないということ
- まだバグがありそうな、お試しソース
library(geometry)
library(gtools)
library(igraph)
n.pt <- 100
pc = matrix(rnorm(3*n.pt),ncol=3)
tc = delaunayn(pc)
library(rgl)
clr = rep(1,3) %o% (1:nrow(tc)+1)
rgl.viewpoint(60,fov=20)
rgl.light(270,60)
tetramesh(tc,pc,alpha=0.7,col=clr)
my.make.elist <- function(tc){
k <- ncol(tc)
r <- nrow(tc)
cmb <- combinations(k,2)
v1 <- v2 <- c()
for(i in 1:nrow(cmb)){
v1 <- c(v1,tc[,cmb[i,1]])
v2 <- c(v2,tc[,cmb[i,2]])
}
el <- cbind(v1,v2)
el <- t(apply(el,1,sort))
unique(el)
}
el <- my.make.elist(tc)
g <- graph.edgelist(el)
E <- as.matrix(get.adjacency(g))
E <- sign(E + t(E))
E[lower.tri(E)] <- 0
S.list <- list()
n <- nrow(E)
S.list[[1]] <- diag(rep(1,nrow(E)))
cnt <- 1
loop <- TRUE
while(loop){
if(length(S.list[[cnt]][,1])==0){
break
}
tmp <- S.list[[cnt]] %*% E
tmp2 <- which(tmp==cnt,arr.ind=TRUE)
if(length(tmp2[,1])>0){
cnt <- cnt
S.list[[cnt+1]] <- matrix(0,length(tmp2[,1]),n)
for(i in 1:length(tmp2[,1])){
u <- tmp2[i,1]
v <- tmp2[i,2]
pre <- which(S.list[[cnt]][u,]==1)
S.list[[cnt+1]][i,c(pre,v)] <- 1
}
}else{
loop <- FALSE
}
cnt <- cnt+1
}
tmp <- sapply(S.list,nrow)
wc <- apply(pc,2,mean)
pc.wc.L <- apply((t(pc)-wc)^2,2,sum)
ctr.id <- which(pc.wc.L==min(pc.wc.L))[1]
k <- 2
tmp.k <- S.list[[k+1]]
tmp.k_ <- S.list[[k]]
tmptmp <- tmp.k %*% t(tmp.k_)
S <- sign(tmptmp == k)
S.ori <- sign(tmptmp == k)
W <- B <- list()
W[[1]] <- sample(which(S.list[[k+1]][,ctr.id]==1),1)
B[[1]] <- which(S[W[[1]],] == 1)
S[W[[1]],] <- 0
n.step <- 500
for(i in 1:n.step){
print("----")
current.B <- B[[i]]
print(current.B)
candidates <- which(matrix(S[,current.B],ncol=length(current.B))==1,arr.ind=TRUE)
candidates.W <- candidates[,1]
candidates.B <- current.B[candidates[,2]]
print(candidates.W)
print(candidates.B)
selection <- sample(1:length(candidates.W),1)
selected.W <- candidates.W[selection]
selection2 <- which(candidates.W == selected.W)
print(selection)
print(selection2)
print(selected.W)
BofnewW <- which(S[candidates.W[selection],]==1)
print(BofnewW)
selected.B <- candidates.B[selection2]
print(selected.B)
W[[i+1]] <- c(W[[i]],selected.W)
B[[i+1]] <- unique(c(B[[i]],BofnewW))
print("pre B")
print(B[[i+1]])
B[[i+1]] <- B[[i+1]][-(which(B[[i+1]] %in% selected.B))]
print("post B")
print(B[[i+1]])
toberemoved <- which(matrix(S[,selected.B]==1,ncol=length(selected.B)),arr.ind=TRUE)
print(toberemoved)
S[toberemoved[,1],] <- 0
}
plot3d(pc)
cmb <- combinations(k+1,2)
ss <- 1:length(W)
for(i in ss){
for(j in 1:length(W[[i]])){
tmp <- W[[i]][j]
v <- which(S.list[[k+1]][tmp,]==1)
segments3d(pc[v[t(cmb)],])
}
}
ss <- length(W)
ss <- 20
for(i in ss){
plot3d(pc)
v <- c()
for(j in 1:length(W[[i]])){
tmp <- W[[i]][j]
tmp.v <- which(S.list[[k+1]][tmp,]==1)
v <- c(v,tmp.v)
}
triangles3d(pc[v,],alpha=0.8,col=1:length(v))
}