격자 모델 시뮬레이션에서의 확산

격자 모델 시뮬레이션에서의 확산

시뮬레이션

SI\_diffusion.gif

이 포스트에서는 격자 공간에서 어떤 성분(Ingredient)의 확산 현상을 모방하려고 한다. 이는 그 동시에 SI 질병 확산 모델의 시뮬레이션이기도 하며, 공간이 제한되어 있다는 점에서 SIR 모델로도 볼 수 있다.

변수

파라메터

액션

(모든 칸에서, 매 턴마다)

코드 리뷰

Step 1. 초기화 및 격자공간 생성

julia> row\_size = 10
10

julia> column\_size = 10
10

julia> β = 0.5 # 확산률
0.5

julia> Random.seed!(0);

julia> stage\_lattice = rand(['S'], row\_size, column\_size)
10×10 Array{Char,2}:
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'

julia> stage\_lattice[rand(1:row\_size), rand(1:column\_size)] = 'I'; stage\_lattice
10×10 Array{Char,2}:
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'I'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'
 'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'  'S'

격자 기반 시뮬레이션 튜토리얼과 비슷하게 격자 공간을 생성하고 랜덤한 위치에 성분을 넣는다. 위의 경우에에는 $6$행 $6$열에 생성된 것을 확인할 수 있다. 이를 히트맵으로 시각화하면 다음과 같다.

06-01.png

Step 2. 계산

julia> I\_lattice = (stage\_lattice .== 'I')
10×10 BitArray{2}:
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  1  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0

julia> count\_lattice =
         vcat(I\_lattice[2:end, :], zeros(Int64, 1, column\_size)) +
         vcat(zeros(Int64, 1, column\_size), I\_lattice[1:end-1, :]) +
         hcat(I\_lattice[:, 2:end], zeros(Int64, row\_size, 1)) +
         hcat(zeros(Int64, column\_size, 1), I\_lattice[:, 1:end-1])
10×10 Array{Int64,2}:
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  1  0  0  0  0
 0  0  0  0  1  0  1  0  0  0
 0  0  0  0  0  1  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0

주변에 성분이 얼마나 있는지 세어놓은 행렬 $(n_{ij})$ 을 만든다. 예제에선 벡터 연산을 이용해 구했지만 반드시 이렇게 계산할 필요는 없다.

Step 3. 확산

julia> probability\_lattice = 1 .- (1-β).^count\_lattice
10×10 Array{Float64,2}:
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.5  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.5  0.0  0.5  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.5  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 
julia> hit\_lattice = probability\_lattice .> rand(row\_size, column\_size)
10×10 BitArray{2}:
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  1  0  0  0  0
 0  0  0  0  1  0  1  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0

각 칸 별로 확산 받을 확률을 계산해서 행렬 $\left( 1-(1-\beta)^{n_{ij}} \right)$ 를 만든다. 이렇게 계산하는 이유는 주변에 성분이 많을수록 확산률이 높아지는 것을 반영하기 위해서다. 빈 공간에서 성분을 받을 확률은 $\beta$, 받지 않을 확률은 $(1-\beta)$ 다. 빈 공간 $i$행 $j$열이 주변 $n_{ij}$개의 성분에게서 확산을 받지 않을 확률은 $(1-\beta)^{n_{ij}}$ 이고, 결국 그들 중 하나에게서라도 확산을 받았을 확률이 바로 $\left( 1-(1-\beta)^{n_{ij}} \right)$ 이다.

julia> probability\_lattice = 1 .- (1-β).^count\_lattice
10×10 Array{Float64,2}:
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.5  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.5  0.0  0.5  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.5  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0

julia> hit\_lattice = probability\_lattice .> rand(row\_size, column\_size)
10×10 BitArray{2}:
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  1  0  0  0  0
 0  0  0  0  1  0  1  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0
 0  0  0  0  0  0  0  0  0  0

이제 균등분포 $U(0,1)$에서 같은 크기의 행렬을 만들어 확산이 일어났는지 일어나지 않았는지 판정한다. 위 케이스에서는 상, 좌, 우로 성분이 퍼져나가는 것을 볼 수 있다. 실제로 히트맵으로 그려보면 다음과 같다.

06-02.png

이를 반복하면 본문과 상단과 같은 움짤을 얻을 수 있다.

댓글