\( \newcommand{\bm}[1]{\boldsymbol{\mathbf{#1}}} \)

Model Equations

\[\begin{equation} \begin{split} &\frac{\partial N(x, t)}{\partial t} = r N(x, t) + D \frac{\partial^2 N(x, t)}{\partial x^2} \\ & D=1, ~ 0<x<L=10.0, ~ t>0 \\ \\ & Fixed ~ Boundary ~ Conditions: ~~ N(0,t)=0=N(L,t) \\ & Initial ~ Conditions: ~~if~ i=M/2=50, ~n(i,0)>0. ~~ Otherwise ~ n(i,0)=0 \end{split} \tag{1} \end{equation}\]

Discretization

\[\begin{equation} \begin{split} \frac{d N_i(t)}{d t} &= r N_i(t) + \frac{D}{(\Delta x)^2} [n_{i+1}(t) - 2 n_{i}(t) + n_{i-1}(t)] \\ \Delta x &= x_{i+1} - x_i = \frac{L}{M} \end{split} \tag{2} \end{equation}\]

Critical Value

The population goes extinct if:

\[\begin{equation} \begin{split} L &< \pi \sqrt{\frac{D}{r}}~,~~ or ~~ r &< \frac{D \pi^2}{L^2} = \frac{1 \pi^2}{10^2} = 0.098696 \end{split} \end{equation}\]

Case 1: Population Extinction

\(r=0.09, ~ N_{50}(0) = 100\)

library(ggplot2)
library(dplyr)
library(magrittr)
library(tidyr)
data_extinct <- readr::read_delim("data_0.09.txt", "\t", escape_double = FALSE, col_names = FALSE, trim_ws = TRUE) %>%
    set_names(c("time", paste0("",seq(1,100,1))))

df <- data_extinct %>% 
    gather(key=i, value=N, -time) %>%
    mutate(i = as.integer(i)) %>%
    arrange(time, i)

ggplot(df, aes(x=time, y=i, z = N))+
    stat_contour(geom = "polygon", 
                 aes(fill = ..level..),
                 bins = 6000) +
    scale_fill_distiller(palette="Spectral", 
                         trans="log10")+
    geom_contour(aes(colour = ..level..), 
                 binwidth=0.25,
                 lwd=0.3, show.legend = F)+
    scale_color_distiller(palette="Spectral", 
                         trans="log10")+
    labs(x="Time", y="Position", fill="Pop. Size",
         caption= "Diff. between 2 contours = 0.25")+
    theme_bw()