2010-12-29

離散的に埋める

こちらの続き
２次元平面を正方格子にして、ある格子の隣を東西南北の４方向にするか、東西南北・北東　南東　南西　北東の８方向にするかの話がこちら
正方格子の４方向と４＋４＝８方向は、見た感じが「均一」とは言い難い(隣の定義としてはありだが、幾何的には東西南北の４つと、その中間方向の４つは、異なるということ。たとえば、ユークリッド距離にすれば、前者を１とすれば、後者はsqrt(2))
では、ユークリッド距離的に、平等にした上で、方向の数を４ではなくすることはできるだろうか
たとえば、正三角形で平面を埋めつくせば、ある正三角形の隣(辺)を共有する正三角形は３個。正６角形で埋め尽くせば、ある正六角形は６個の隣を持つ
このような３、４(正方格子)、６のそれぞれにおいて、拡散の仕方が同じようになるように拡散係数を定められるのだろうか。定められるとして、その係数の値と３，４，６という数字はどういう関係にあるのだろう
ちなみに、k角形はk-2個の三角形に分割できるので、その内角の和は $2\pi(k-2)$ 。したがって、正k角形の内角は $\frac{2\pi(k-2)}{k}$ 。正k角形をある１点で接するように並べて、平面を埋め尽くそうとすれば、 $\frac{2\pi}{\frac{2\pi(k-2)}{k}}$ が正の整数になるという条件があって、それを満足する正の整数は3,4,6のみなので、この方法で確かめられるのは、3,4,6のみ
平面を正k角形タイルで埋め尽くす方法とその拡散について考えた
シミュレーションするときに境界をどうするか、という話もあった
今、球面をタイルで埋め尽くすことを考える。サッカーボールが六角形と五角形で埋め尽くしていたように、埋め尽くし方に少し違いが出る。それは、拡散の離散的モデルである「埋め尽くし」シミュレーションにどういう影響をもたらすだろうか。そして、その影響は、現実とも関係があるだろうか。また、この場合、境界条件は「閉じている」という点で「異なる」境界条件を求められる
さらに、「無限に広がっている場合の境界条件(シミュレーションでは境界条件を設定する必要があるが、本当に無限な広がりならば、境界がないという設定になる)」で戦略を立てていた何者かがいたとする。それが、はた、と、球面のような境界条件に遭遇したとき、どんな不具合が発生するのだろうか。生物種が広まるときに、無限平面(無限立体)を仮定するべき試練にさらされてそれに対応する戦略家に育っていたとき、環境の有限条件と折り合いをつけられるのか、ということであるが。
大腸菌はシャーレの中で飽和に対する対応ができる。ヒトは地球上で飽和に対する対応戦略を生物として持っているのだろうか。
それを考えるためには、無限平面拡散と閉じた球面拡散との条件の違いから見るべきだろう。
すくなくとも粘菌は球面状の寒天培地の上でもきちんと生きていけるように想像するが。

# X,Y,Z
library(MCMCpack)
k<-4
xyz<-rdirichlet(1,rep(1,k))

Nr<-6
Nmol<-rpois(Nr,6)+1
pre<-matrix(0,Nr,k)
post<-matrix(0,Nr,k)
for(i in 1:Nr){
	prepre<-sample(1:k,Nmol,replace=TRUE)
	postpost<-sample(1:k,Nmol,replace=TRUE)
	for(j in 1:k){
		pre[i,j]<-length(which(prepre==j))
		post[i,j]<-length(which(postpost==j))
	}
}

#pre<-matrix(sample(0:3,Nr*k,replace=TRUE),Nr,k)
#post<-matrix(sample(0:3,Nr*k,replace=TRUE),Nr,k)
Rfw<-runif(Nr)
Rrev<-runif(Nr)


Niter<-1000
d<-matrix(0,Niter,k)

d[1,]<-xyz
for(i in 2:Niter){
	tmppred<-d[i-1,]
	for(j in 1:Nr){
		lnfreq<-log(tmppred)
		tmpd<-rep(0,k)

		kfw<-exp(sum(lnfreq*pre[j,]))*Rfw[j]
		krev<-exp(sum(lnfreq*post[j,]))*Rrev[j]
		tmpd<-tmpd-kfw*pre[j,]+kfw*post[j,]-krev*post[j,]+krev*pre[j,]
		tmppred<-tmppred+tmpd
	}
	d[i,]<-tmppred
}

matplot(d,type="l")

##############

X<-25
dif<-runif(k)*0.1
#dif<-rep(0,k)
A<-array(1/k,c(X,X,k))
A[((X/2)-1):(X/2),((X/2)-6):(X/2+6),1]<-0.9
A[((X/2)-1):(X/2),((X/2)-6):(X/2+6),2:k]<-0.1/(k-1)
A[((X/4)-1):(X/4),((X/4)-1):(X/4),1:(k-1)]<-0.1/(k-1)
A[((X/4)-1):(X/4),((X/4)-1):(X/4),k]<-0.9
A[((3*X/4)-1):(3*X/4),((3*X/4)-1):(3*X/4),2]<-0.9
A[((3*X/4)-1):(3*X/4),((3*X/4)-1):(3*X/4),c(1,3,4)]<-0.1/(k-1)

image(A[,,1])
Niter<-1000

for(i in 2:Niter){
	
	# Reaction
	for(j in 1:X){
		for(j2 in 1:X){
			tmppred<-A[j,j2,]

			for(j3 in 1:Nr){
				lnfreq<-log(tmppred)
				tmpd<-rep(0,k)

				kfw<-exp(sum(lnfreq*pre[j3,]))*Rfw[j3]
				krev<-exp(sum(lnfreq*post[j3,]))*Rrev[j3]
				tmpd<-tmpd-kfw*pre[j3,]+kfw*post[j3,]-krev*post[j3,]+krev*pre[j3,]
				tmppred<-tmppred+tmpd
			}
			A[j,j2,]<-tmppred

		}
	}
		# Diffusion
	for(j in 1:k){
		tmpA<-matrix(0,X,X)
		for(j2 in 2:(X-1)){
			for(j3 in 2:(X-1)){
				du<-dif[j]*(sum(A[(j2-1),j3,j],A[j2+1,j3,j],A[j2,j3-1,j],A[j2,j3+1,j])-4*A[j2,j3,j])
				tmpA[j2,j3]<-tmpA[j2,j3]+du
			}
		}
		A[,,j]<-A[,,j]+tmpA
	}
	num<-1000+i
	filename<-paste("out",num,".png")
	png(filename)
layout(matrix(c(1,2,3,4), 2, 2, byrow = TRUE))
for(j in 1:4){
	image(A[,,j],breaks=seq(from=0,to=1,by=0.1),col=terrain.colors(10))
}
dev.off()
	
A[((X/2)-1):(X/2),((X/2)-6):(X/2+6),1]<-0.9
A[((X/2)-1):(X/2),((X/2)-6):(X/2+6),2:k]<-0.1/(k-1)
A[((X/4)-1):(X/4),((X/4)-1):(X/4),1:(k-1)]<-0.1/(k-1)
A[((X/4)-1):(X/4),((X/4)-1):(X/4),k]<-0.9
A[((3*X/4)-1):(3*X/4),((3*X/4)-1):(3*X/4),2]<-0.9
A[((3*X/4)-1):(3*X/4),((3*X/4)-1):(3*X/4),c(1,3,4)]<-0.1/(k-1)

}