とにかく絵を描いていみる

• ここからの続き
• 反応拡散系から、あちこち寄り道をしているが、年の暮でもあり、絵を描いてみることにする
• こちらにあるように
• ２変量U,Vについて以下の変化を離散的に２次元平面に描いてみる
• 
• 
• 初期状態はu+v=1でそのU,V比はランダムであるもので描く

# Space
X<-100
Niter<-100


U<-V<-matrix(0,X,X)
x<-X*0.3

U[x,x]<-1
V[x,x]<-0
U[X-x,X-x]<-0
V[X-x,X-x]<-1
U<-matrix(runif(X^2),X,X)
V<--U+1

F<-0.5
k<-0.1

dif1<-0.3 # Du
dif2<-0.4 # Dv

filecounter<-1000
par(ask=TRUE)
for(i in 1:Niter){
	# Reaction
	tmpU<-tmpV<-matrix(0,X,X)
	for(j1 in 1:X){
		for(j2 in 1:X){
			tmpU[j1,j2]<--U[j1,j2]*V[j1,j2]^2+F*(1-U[j1,j2])
			tmpV[j1,j2]<-U[j1,j2]*V[j1,j2]^2-(F+k)*V[j1,j2]
		}
	}
	U<-U+tmpU
	V<-V+tmpV
	#U[x,x]<-1
	#V[x,x]<-0
	#U[X-x,X-x]<-0
	#V[X-x,X-x]<-1

	#image(U,breaks=seq(from=0,to=1,by=0.1),col=terrain.colors(10))
	# Diffusion
	tmpU<-tmpV<-matrix(0,X,X)
	for(j1 in 2:(X-1)){
		for(j2 in 2:(X-2)){
			tmpU[j1,j2]<-dif1*(U[j1,j2+1]+U[j1,j2-1]+U[j1+1,j2]+U[j1-1,j2]-4*U[j1,j2])
						tmpV[j1,j2]<-dif2*(V[j1,j2+1]+V[j1,j2-1]+V[j1+1,j2]+V[j1-1,j2]-4*V[j1,j2])
		}
	}
	U<-U+tmpU
	V<-V+tmpV
	filename<-paste("out",filecounter+i,".png")
	#png(filename)
	image(U,breaks=seq(from=0,to=1,by=0.1),col=terrain.colors(10))
	#image(U)
	#dev.off()
}