=Paper= {{Paper |id=Vol-1576/077 |storemode=property |title=Parallel implementation of stochastic cellular automata model of electron-hole recombination in a semiconductor |pdfUrl=https://ceur-ws.org/Vol-1576/077.pdf |volume=Vol-1576 |authors=Anastasiya Kireeva,Karl Sabelfeld }} ==Parallel implementation of stochastic cellular automata model of electron-hole recombination in a semiconductor== https://ceur-ws.org/Vol-1576/077.pdf
Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                               agora.guru.ru/pavt




                Параллельная реализация стохастической
               клеточно-автоматной модели рекомбинации
                                                     \ast 
                 электронов и дырок в полупроводнике

                                                  1                       1
                                   А.Е. Киреева , К.К. Сабельфельд

                                                                                                       1
   Институт вычислительной математики и математической геофизики СО РАН



               Разработаны параллельные программы, реализующие стохастические клеточно-
               автоматные (КА) модели рекомбинации электронов и дырок в неоднородном
               полупроводнике в двумерном и трехмерном случаях. С помощью разработан-
               ной КА модели рекомбинации исследовано пространственно-временное распре-
               деление частиц, обнаружено и исследовано формирование макрокластеров элек-
               тронов и дырок. Параллельная реализация программы позволила вычислить за
               приемлемое время интегральные характеристики процесса: плотности частиц и
               интенсивность фотолюминесценции, для большого числа различных начальных
               условий, а также изучить кинетику процесса рекомбинации при наличии центров
               рекомбинации и диффузии частиц в двумерном и трехмерном случаях.
               Ключевые слова:  рекомбинация электронов и дырок, полупроводник, фотолюми-
               несценция, стохастический клеточный автомат, параллельная реализация.

1. Введение

          Широкозонный полупроводник нитрид галлия (GaN ) и его твердые растворы явля-
ются одним из самых востребованных и перспективных материалов в области разработок
полевых транзисторов и современной оптоэлектроники [1]. Физические свойства                     GaN и
его соединений обеспечивают создание на их основе приборов с улучшенными значениями
мощности, напряжения и тока. Для увеличения надежности и эффективности полупровод-
никовых приборов на основе GaN , а также для расширения области их применения, необ-
ходимо исследовать процессы, протекающие в полупроводниках, и изучать характеристики
этих процессов в зависимости от параметров GaN .
          В данной работе с помощью компьютерного моделирования изучается процесс реком-
бинации электронов и дырок в неоднородном полупроводнике GaN для следующего слу-
чая: с помощью фемтосекундного лазера в полупроводнике создается избыток электро-
нов и дырок, которые аннигилируют с друг другом путем туннелирования, процесс уско-
ряется диффузией и одновременно замедляется нерадиационной рекомбинацией электро-
нов и дырок в дефектах (рекомбинационных центрах) [2]. В [3–6] разработана стохастиче-
ская модель взаимодействия электронов и дырок, основанная на стохастических интегро-
дифференциальных уравнениях типа Смолуховского с учетом пространственной неодно-
родности, диффузии, наличия центров рекомбинации и возможности флуктуаций началь-
ных концентраций. На основе [3,6] разработана клеточно-автоматная (КА) модель рекомби-
нации электронов и дырок в неоднородном полупроводнике для двумерного и трехмерного
случаев.
          Компьютерное моделирование процессов в полупроводнике требует использования зна-
чительных вычислительных мощностей, так как необходимо моделировать поведение боль-
шого количества частиц в объеме большого размера в течение длительного времени. Кроме
того, интегральные характеристики процесса рекомбинации: плотности частиц и интенсив-
ности фотолюминесценции, вычисляются с помощью осреднения по большому ансамблю на-
чальных данных. Для сокращения времени КА-моделирования и вычисления интегральных
характеристик процесса разработаны параллельные программы, реализующие КА-модель

  \ast 
          Работа выполнена при финансовой поддержке Российского научного фонда, грант № 14-11-00083.




                                                      167
Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                                                             agora.guru.ru/pavt



рекомбинации в двумерном и трехмерном пространстве. С помощью разработанных про-
грамм изучается кинетика процесса аннигиляции электронов и дырок при наличии центров
рекомбинации и диффузии частиц в двумерном и трехмерном случаях.
         Статья включает в себя четыре раздела. В разделе 2 описана математическая постанов-
ка задачи, приведена система уравнений Смолуховского и формулы, описывающие асимп-
тотическое поведение характеристик процесса рекомбинации при больших временах. Раз-
дел 3 посвящен описанию КА-модели рекомбинации, в нем объясняются основные термины
теории КА. В разделе 4 приводятся результаты распараллеливания КА-модели рекомби-
нации для трехмерного и двумерного случаев. В разделе 5 представлены результаты КА-
моделирования рекомбинации, приведены эволюции КА и значения характеристик для раз-
личных параметров моделирования.


2. Математическая модель рекомбинации электронов и дырок
      в полупроводнике

         Для описания процесса рекомбинации электронов и дырок в полупроводнике в [3, 6]
предложена стохастическая модель, основанная на системе уравнений Смолуховского (1):

                                      r
                             \partial \rho n ( ; t)
                                                                                  r r
                                                       = Dn ( )\Delta \rho n ( ; t)  -  \rho n ( ; t)
                                                                                                        \int 
                                                                                                              r B(| x| )\rho p(r+x; t)dx
                                            \partial t
                                                                                                                              \int 
                                                                  -  \beta  d (r)\rho n (r; t)\rho p (r; t)  -  \rho n (r; t) bn (| x| )\rho N (r+x; t)dx;
                                                                                                                                 n


                                       \partial \rho p (r; t)
                                                                                                                        \int 
                                                                  = Dp (r)\Delta \rho p (r; t)  -  \rho p (r; t) B(| x| )\rho n (r+x; t)dx
                                                      \partial t
                                                                                                                              \int                                 (1)

                                                                  -  \beta  (r)\rho p (r; t)\rho n (r; t)  -  \rho p (r; t) bp (| x| )\rho N (r+x; t)dx;
                                                                            d
                                                                                                                                p


             \partial \rho N (r; t)
                                                                               \int                                                  \int 
                   n
                                                     =  - \rho n (r; t) bn (| x| )\rho N (r+x; t)dx + \rho p (r; t) bp (| x| )\rho N (r+x; t)dx.
                                                                                 n                                                         p
                            \partial t
         Предполагается, что в начальный момент времени электроны, дырки и рекомбинаци-
онные центры случайно распределены в d - мерном пространстве X с плотностями \rho n (                                                                            r; t)   ,
     r                   r
\rho p ( ; t) и \rho N ( ), соответственно. Символ r обозначает пространственную координату, а t
- момент времени. Рекомбинационные центры могут быть свободны, тогда они способны
захватить электрон, либо заняты электроном, тогда они способны захватить дырку. Сум-
марная плотность рекомбинационных центров \rho N (                                                         r) = \rho N (r; t) + \rho N (r; t)
                                                                                                                      n              p          в процессе моде-
лирования остается постоянной, меняются только плотности рекомбинационных центров,
свободных для электронов \rho Nn (                               r; t)   , и рекомбинационных центров, свободных для дырок
         r
\rho Np ( ; t).
         Радиационная рекомбинация электронов и дырок осуществляется путем туннелирова-
ния со скоростью B(|                      x| ) = B0 \cdot  exp( - | x| /anp)             , где |        x|      - расстояние между электроном и дыр-
кой, anp - характерное расстояние взаимодействия электронов и дырок. При радиационной
рекомбинации происходит выделение энергии в виде фотона. Захват электронов в свободные
рекомбинационные центры происходит со скоростью bn (|                                                                 x| ) = bn0 \cdot  exp( - | x| /anN )
                                                                                                                                                       n   . Анало-
гично, скорость захвата дырок в рекомбинационные центры, содержащие электроны, равна
     x                                    x
bp (|  | ) = bp0 \cdot  exp( - |  | /apNp ). При взаимодействии электронов и дырок в рекомбинацион-
ных центрах происходит нерадиационная рекомбинация, сопровождающаяся выделением
тепловой энергии.
         Электроны и дырки могут диффундировать в области с коэффициентами диффузии
      r              r
Dn ( ) и Dp ( ). При перемещении вследствие диффузии встретившиеся электроны и дырки
                                                             3
рекомбинируют с друг другом с коэффициентом \beta  ( ) = 4\pi  \cdot  D \cdot  r0 в трехмерном про-               r
                  2             r
странстве и \beta  ( ) = 4\pi  \cdot  D в двумерном пространстве, где D( ) = Dn ( ) + Dp ( ) - это                              r               r            r
относительный коэффициент диффузии, r0 - радиус частицы, далее r0 = 1 nm.



                                                                                        168
Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                                                            agora.guru.ru/pavt



     Интенсивность фотолюминесценции в [6] вычисляется по формуле:
                                            \biggl\langle \int               \int                                                         \biggr\rangle 
                               \phi (t) =
                                                                    1
                                                                  | X| 
                                                                        dr                     x                 r
                                                                                     B(|  | )\rho n ( ; t)\rho p (            r+x; t)dx ,                                  (2)


где угловые скобки обозначают математическое ожидание по начальным распределениям
частиц, | X|  - размер области.
     В [6] с помощью корреляционного анализа выведены формулы для асимптотического
поведения плотности электронов и интенсивности фотолюминесценции при t \rightarrow  \infty  для двух
случаев:
1.    для случая чистого туннелирования, то есть радиационной рекомбинации без реком-
      бинационных центров и диффузии частиц:


                         \rho \prime n_dD (t) \sim  1/(ln(t))d/2 ;                                  \prime 
                                                                                                   IdD      (t) \sim  1/\{ t \cdot  (ln(t))d/2+1 \} ,                      (3)


2.    для случая чистой диффузии, то есть радиационной рекомбинации вследствие только
      диффузии частиц:


                                      \rho \prime n_dD (t) \sim  1/td/4 ;                                  \prime 
                                                                                                          IdD      (t) \sim  1/td/4+1 ,                                    (4)


      где d - размерность пространства.
     Для решения системы уравнений Смолуховского (1) в [6] представлен алгоритм Монте-
Карло. На основе этого метода в данной работе разработана клеточно-автоматная модель
рекомбинации электронов и дырок.


3. Клеточно-автоматная модель рекомбинации электронов и ды-
     рок в полупроводнике

     Клеточный автомат (КА) - это дискретная динамическая система, состоящая из мно-
жества клеток, плотно заполняющих d-мерное пространство [7–9]. Каждая клетка харак-
теризуется парой значений (a,                                 x)   , где a - это состояние клетки,                                        x            - координата клетки в
пространстве X . Состояние клетки a                                          \in  A, где A - алфавит состояний, который опреде-
ляется возможными состояниями моделируемой системы. Состояние клетки изменяется в
соответствии с правилами переходов \Theta  в зависимости от состояния самой клетки и от со-
стояний взаимодействующих с ней клеток, которые выбираются с помощью шаблона мо-
делирования T . Шаблон моделирования может быть фиксированным и определяться как
множество соседних клеток, расположенных вокруг центральной клетки, либо это может
быть множество случайных клеток, выбираемых в соответствии с заданным распределе-
нием из множества всех клеток пространства. Применение правил переходов \Theta (                                                                                  x)   ко всем
клеткам      x \in  X   называется итерацией. Правила переходов могут применяться к клеткам
x \in  X   в различном порядке. Этот порядок называется режимом функционирования (\mu ) КА.
Для моделирования стохастических физико-химических процессов используется асинхрон-
ный режим работы КА (\mu  = \alpha ), при котором правила переходов применяются к случайно
выбранным клеткам клеточного массива, сразу же изменяя их состояния, таким образом,
новые состояния клеток вычисляются от состояний, полученных на предыдущих и на те-
кущей итерациях.
     На основании определения КА [8] КА-модель рекомбинации электронов и дырок в по-
лупроводнике задается в следующем виде:


                                                                        \aleph  = \langle A, X d , \Theta , \alpha \rangle                                                 (5)


     В соответствии с математической моделью рекомбинации (1) и алгоритмом Монте-
Карло [6] алфавит состояний выбирается в виде A = \{ n, p, Nn , Np , \emptyset \} , где символ n обо-
значает электрон, p - дырку, Nn - рекомбинационный центр свободный для электрона, Np



                                                                                          169
 Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                                                                 agora.guru.ru/pavt



- рекомбинационный центр свободный для дырки, \emptyset  обозначает свободное место. Множе-
ство координат X
                               d в двумерном случае (d = 2) представимо в виде квадратной решетки

                 x
X 2 = \{  = (i, j), i = 1 . . . Sizei , j = 1 . . . Sizej \} , а в трехмерном (d = 3) - в виде куба
                 x
X 3 = \{  = (i, j, k) i = 1 . . . Sizei , j = 1 . . . Sizej , k = 1 . . . Sizek \} . На множестве координат
вводятся периодические граничные условия.
         Правила переходов задаются в следующем виде: \Theta  = R\{ \theta 1 , \theta 2 , \theta 3 , \theta 4 , \theta 5 \} , где символ R
обозначает вероятностный выбор одного из правил \theta l ,   l = 1, . . . , 5. Правила переходов \Theta 
моделируют процессы, происходящие в полупроводнике: \theta 1 - радиационную рекомбинацию
электрона и дырки, \theta 2 - захват электрона в свободный рекомбинационный центр, \theta 3 - захват
дырки в рекомбинационный центр, содержащий электрон, \theta 4 и \theta 5 - диффузию электронов
и дырок. Формально правила переходов записываются в виде (6):


                                                     x           x              x p \cdot  \omega  x1       1
                                                                                                                                        x
                                           \theta 1 ( ) : \{ (n, ), (p, \varphi ( )\}  \rightarrow  \{ (\emptyset , ), (\emptyset , \varphi ( )\} ,
                                           \theta 2 (x) : \{ (n, x), (Nn , \varphi (x)\}  \rightarrow  \{ (\emptyset , x), (Np , \varphi (x)\} ,
                                                                                          p \cdot  \omega 
                                                                                                        2             2



                                           \theta 3 (x) : \{ (p, x), (Np , \varphi (x)\}  \rightarrow  \{ (\emptyset , x), (Nn , \varphi (x)\} ,
                                                                                         p \cdot  \omega 
                                                                                                        3             3
                                                                                                                                                                                       (6)

                                           \theta 4 (x) : \{ (n, x), (a, \psi (x)\}  \rightarrow  \{ (\emptyset , x), (a\prime  , \psi (x)\} ,
                                                                                       p \cdot  \omega 
                                                                                                    4       4



                                           \theta 5 (x) : \{ (p, x), (b, \psi (x)\}  \rightarrow  \{ (\emptyset , x), (b\prime  , \psi (x)\} ,
                                                                                      p \cdot  \omega 
                                                                                                    5       5

                                                \left\{                                                                   \left\{ 
                                                           n, если a = \emptyset ,                                                   p, если b = \emptyset ,
                                 \prime                                                                     \prime 
                                a =                        \emptyset , если a = p,                          b =                      \emptyset , если b = n,
                                                           Np , если a = Nn ,                                                        Nn , если b = Np ,
где \varphi (  x) \psi (x)
                  и          - это координаты клеток, выбранных для взаимодействия с клеткой                                                                             x  . Клет-
ка \varphi (   x) выбирается случайным образом из множества всех клеток TSize = \{                                                                              y \in  X, y \not = x\}   .
Клетка            \psi (x)
                        выбирается случайным образом из соседних клеток клетки                                                                                  x    по шаблону
”крест”, который в двумерном случае состоит из четырех клеток                          T4 = \{ (i, j  -  1), (i +
1, j), (i, j + 1), (i  -  1, j)\} , а в трехмерном случае - из шести клеток T6 = \{ (i, j  -  1, k), (i +
1, j, k), (i, j + 1, k), (i  -  1, j, k), (i, j, k  -  1), (i, j, k + 1)\} .
      При применении \Theta  выбирается одно из правил \theta l , l = 1, . . . , 5, с вероятностью pl . Веро-
ятности применения правил переходов вычисляются по формуле (7):

                                                                                                                            5
                                                                                                                          \sum 
                                                    pl = \lambda l /\lambda , l = 1, . . . , 5, \lambda  =                               \lambda l ,
                                                                                                                          l=1
                                                                                                        min
                                                    \lambda 1 = Cn \cdot  Cp \cdot  B0 \cdot  exp\{  - rnp  /anp \} ,
                                                                                                          min
                                                    \lambda 2 = Cn \cdot  CNn \cdot  bn0 \cdot  exp\{  - rnN n
                                                                                                               /anNn \} ,                                                              (7)
                                                                                                          min
                                                    \lambda 3 = Cp \cdot  CNp \cdot  bp0 \cdot  exp\{  - rpN p
                                                                                                               /apNp \} ,
                                                    \lambda 4 = Cn \cdot  D\~nd , где D\~n2 = Dn , D\~n3 = Dn \cdot  r0 ,
                                                    \lambda 5 = Cp \cdot  D\~d , где D\~2 = Dp , D\~3 = Dp \cdot  r0 ,
                                                                             p                  p                                    p

где Cn , Cp , CNn , CNp - это количество электронов, дырок, рекомбинационных центров для
электронов и рекомбинационных центров для дырок в клеточном массиве.                                                                                                B0 - это ко-
эффициент радиационной рекомбинации электронов и дырок, bn0 - коэффициент захвата
электронов в свободные рекомбинационные центры, bp0 - коэффициент захвата дырок в ре-
комбинационные центры, содержащие электроны. Dn и Dp - это коэффициенты диффузии
электронов и дырок, r0 - радиус частицы. ruv
                                                                                             min - это минимальное расстояние между всеми

частицами типа u и v , где u \in  \{ n, p\} , v \in  \{ p, Nn , Np \} . Аналогично, auv - это характерное рас-
стояние взаимодействия частиц типа u с частицами типа v , где u \in  \{ n, p\} , v \in  \{ p, Nn , Np \} .
         После выбора одного из правил \theta l , с соответствующей вероятностью pl случайным об-
разом выбирается клетка                                    x \in  X d     и взаимодействующая с ней клетка \varphi (                                           x)   либо \psi (   x)   . В



                                                                                             170
Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                                       agora.guru.ru/pavt



выбранных клетках правила переходов \theta l ,                         l = 1, 2, 3, реализуются с вероятностью \omega l =
exp((ruv  min  -  r )/a ), где u, v соответствуют типам частиц в выбранном правиле \theta  . Ве-
                   uv  uv                                                                                                  l
роятности реализации правил \theta 4 и \theta 5 зависят от состояния соседней клетки: если клетка
   x
\psi ( ) свободна (a = \emptyset , b = \emptyset ), тогда \omega l = 1, l = 4, 5, иначе \omega l = \beta  d , где коэффициент \beta  d
вычисляется как определено в разделе 2 при r0 = 1.
       Каждой попытке применения одного из правил переходов \theta l , l = 1, . . . , 5, соответствует
локальный временной шаг \Delta \tau  =  - ln(rand1 )/\lambda , где rand1 \in  (0, 1) - случайное число. Кроме
локального временного шага, в модели используется глобальный временной шаг \Delta t = t0 \cdot q ,
                                                                                                                                 k

где t0 - это начальное время в секундах, q - коэффициент, отвечающий за длину глобального
шага tk , k - номер глобального временного шага. В данной модели глобальный временной
шаг соответствует одной итерации КА.
       Основными характеристиками, измеряемыми экспериментально, являются концентра-
ции частиц и интенсивность фотолюминесценции. В КА-модели рекомбинации \aleph  концентра-
ция частиц типа u: Cu (tk ), u \in  \{ n, p, Nn , Np \} , вычисляется как число клеток с состоянием u
в клеточном массиве в момент времени tk . Интенсивность фотолюминесценции вычисляется
по формуле (8):
                                                   I(tk ) = \sigma /(tk  -  tk - 1 ),                                            (8)

где \sigma  - число фотонов, то есть число взаимодействий электронов и дырок, произошедшее
за время (tk  -  tk - 1 ).
       В начальный момент времени электроны, дырки и рекомбинационные центры случай-
но и равномерно распределены в клеточном массиве. В связи с тем, что правила перехо-
дов и начальные данные КА-модели \aleph  являются вероятностными, значения характеристик
Cu (tk ), u \in  \{ n, p, Nn , Np \} , и I(tk ) являются случайными величинами, поэтому на основе
закона больших чисел вычисляются оценки этих характеристик как математические ожи-
дания значений характеристик, полученных в результате КА-моделирования рекомбинации
для различных начальных распределений частиц [10]. Для вычисления оценок характери-
стик с высокой точностью необходимо, чтобы объем выборки был достаточно большим
(10
    5  -  106 ), для этого требуется провести большое количество численных экспериментов.

В зависимости от размеров задачи время вычислений составляет от нескольких часов до
нескольких суток. Например, КА-моделирование рекомбинации для двумерной области раз-
мером 1000 \times  1000 клеток занимает 27.7 часов. Параллельная реализация задачи дает воз-
можность существенно сократить время вычислений.


4. Параллельная реализация клеточно-автоматной модели ре-
      комбинации электронов и дырок в полупроводнике

       Традиционный метод распараллеливания клеточных автоматов - разбиение области на
подобласти и распределение этих подобластей между процессами, для КА-модели реком-
бинации оказывается малоэффективным в связи с тем, что в данной модели парами для
взаимодействия являются не соседние клетки, а любые две клетки массива. Следователь-
но, при каждом вычислении новых состояний клеток                                   x    и \varphi (   x)   необходимо убедиться, что
их состояния не были изменены другим процессом и обеспечить передачу новых состояний
клеток при использовании их другими процессами. Это требует значительных накладных
расходов на обмен данными.
       Более эффективным для КА-модели рекомбинации является метод распределения вы-
числительных экспериментов с различными начальными данными между процессами. В ра-
боте для параллельной реализации КА-модели \aleph  используется технология OpenMP. Множе-
ство вычислительных экспериментов распределяется между доступными потоками. Каж-
дый поток генерирует свою последовательность случайных чисел, которая используется
для начального распределения частиц в клеточном массиве, а также для выбора и приме-
нения правил переходов \Theta . Все процессы независимо вычисляют значения характеристик



                                                                171
   Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                                      agora.guru.ru/pavt



на временных шагах [t0 ; tf in ], затем полученные значения суммируются главным (master)
потоком.
           Эффективность параллельной реализация КА-модели \aleph  проанализирована на приме-
ре задачи рекомбинации электронов и дырок при наличии рекомбинационных центров и
диффузии электронов: B0 = 0.04 ns
                                                     - 1 , b                         - 1 , a
                                                               n0 = bp0 = 0.02 ns              np = 4 nm, anNn = apNp =
2 nm, t0 = 0.5 ns, tf in             = 108
                                    ns, Cn (0) = Cp (0) = 10000, CNn (0) = 5000, CNp (0) = 0,
Dn = 1 nm2 \cdot  ns - 1 , Dp = 0 nm2 \cdot  ns - 1 , | X d |  = 106 клеток. Число численных эксперимен-
тов с различными начальными данными, по которым оценивались значения характеристик,
выбрано равным 1024. Вычисления проводились на кластере НКС-30Т Сибирского Супер-
компьютерного Центра СО РАН.
           Для оценки качества параллельной реализации КА-модели рекомбинации вычисляются
ускорение S(th) = T (1)/T (th) и эффективность распараллеливания Q(th) = T (1)/(T (th)\cdot  th),
где T (th) - время вычислений при использовании th потоков. Значения S(th) и Q(th), по-
лученные в результате параллельной реализации КА-модели \aleph  со значениями параметров,
указанными выше, в двумерном и трехмерном случаях представлены на Рис. 1.




                            а)                             б)
Рис. 1. Оценки качества параллельной реализации КА-модели рекомбинации: а) ускорение и б)
эффективность распараллеливания.

           Ускорение и эффективность распараллеливания программы, реализующей двумерную
КА-модель рекомбинации, выше, чем программы, реализующей трехмерную КА-модель \aleph .
На Рис. 1 видно, что ускорение S(th) в обоих случаях (d = 2 и d = 3) близко к линейному,
эффективность распараллеливания Q(th) при использовании 12 потоков выше 0.85.


5. Результаты клеточно-автоматного моделирования рекомби-
           нации электронов и дырок в полупроводнике

           Основными параметрами КА-модели \aleph , влияющими на динамику процесса рекомбина-
ции электронов и дырок, являются: константы скорости взаимодействия частиц (B0 , bn , bp ),
характерные расстояния взаимодействия (anp , anNn , apNp ), коэффициенты диффузии (Dn ,
Dp ), начальные концентрации частиц (Cn (0), Cp (0), CNn (0), CNp (0)), и размер области
(Sizei , Sizej , Sizek ).
     В ходе КА-моделирования на каждом глобальном шаге tk \in  [t0 ; tf in ] вычисляются сле-
дующие характеристики процесса:
\bullet    плотности частиц - среднее число частиц, приходящееся на одну клетку:                                \rho u (tk ) =
           Cu (tk )/| X d | , где u \in  \{ n, p, Nn , Np \} , для d = 2 :      2
                                                                             | X |  = Sizei \cdot  Sizej , и для d = 3 :
           | X 3 |  = Sizei \cdot  Sizej \cdot  Sizek .
\bullet    интенсивность фотолюминесценции I(tk ) - число фотонов, образовавшихся за время
           (tk  -  tk - 1 ) (8).
           С помощью КА-моделирования исследовалась динамика процесса рекомбинации в дву-



                                                                172
Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                          agora.guru.ru/pavt



мерном и трехмерном пространстве для следующих случаев:
1.   радиационная рекомбинация электронов и дырок вследствие чистого туннелирования,
2.   радиационная рекомбинация электронов и дырок вследствие чистой диффузии частиц,
3.   радиационная и нерадиационная рекомбинация электронов и дырок,
4.   радиационная и нерадиационная рекомбинация электронов и дырок при наличии диф-
     фузии частиц.
     Рассмотрим динамику процесса рекомбинации электронов и дырок на примере КА-
модели \aleph  со следующими значениями параметров: B0 = 0.04 ns
                                                                   - 1 , b                   - 1
                                                                          n0 = bp0 = 0.02 ns ,
                                                           2  - 1
anp = 4 nm, anNn = apNp = 2 nm, Dn = Dp = 0.5 nm \cdot  ns , t0 = 0.5 ns, Cn (0) = Cp (0) =
10000, CNn (0) = 5000, CNp (0) = 0, | X d |  = 106 клеток.
     В случае чистого туннелирования (случай 1) отсутствуют рекомбинационные центры и
диффузия частиц: CNn (0) = CNp (0) = 0, bn0 = bp0 = 0 ns
                                                                 - 1 и D                   2 \cdot  ns - 1 . Эво-
                                                                           n = Dp = 0 nm
люция \aleph , имитирующая динамику процесса рекомбинации, для приведенных выше значе-
ний параметров, представлена на Рис. 2. Для наглядности на рисунке представлены части




                      а)                        б)                            в)




                  г)                        д)                           е)
Рис. 2. Пространственное распределение электронов и дырок в случае чистого туннелирования
в двумерной области: а) равномерное распределение при t0 = 0 ns, б) формирование кластеров
электронов и дырок при t = 374 ns, в) формирование макрокластеров при t = 1.06 \cdot  105 ns;
и в трехмерной области: г) равномерное распределение при t0 = 0 ns, д) формирование кластеров
электронов и дырок при t = 374 ns, е) формирование макрокластеров при t = 1.06 \cdot  105 ns.
клеточных массивов размером        X 2 = 200 \times  200 и X 3 = 50 \times  50 \times  50 клеток. В началь-
ный момент времени электроны и дырки равномерно распределены в объеме с плотностью
\rho n (0) = \rho p (0) = 0.01 nm - 1 (Рис. 2 а, г). Чем ближе расположены электрон и дырка, тем
с большей вероятностью они взаимодействуют с друг другом. Следовательно, электроны
и дырки, расположенные на близком расстоянии, достаточно быстро аннигилируют, в ре-
зультате происходит пространственное разделение электронов и дырок (Рис. 2 б, д). Далее
происходит взаимодействие электронов и дырок, расположенных на границах кластеров, в
результате мелкие кластеры исчезают и формируются макрокластеры электронов и дырок
(Рис. 2 в, е). Причем в трехмерном случае аннигиляция электронов и дырок происходит
быстрее, чем в двумерном. Это можно объяснить тем, что в трехмерной области при оди-
наковой концентрации частиц направлений для поиска электронно-дырочных пар больше,
чем в двумерной области, следовательно, пара для взаимодействия будет найдена с большей
вероятностью в трехмерном случае.



                                                173
Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                           agora.guru.ru/pavt



    Значения плотности электронов и интенсивности фотолюминесценции, полученные при
КА-моделировании чистого туннелирования в двумерном и трехмерном случаях, представ-
лены на Рис. 3 в логарифмической шкале по обеим осям. Параллельность графиков асимп-
                                              \prime 
тотического поведения характеристик \rho 2D (t), \rho 3D (t), I2D (t),
                                                                        \prime         \prime 
                                                                         \prime  (t), вычисленных по фор-
                                                                        I3D
муле (3), и графиков характеристик \rho 2D (t), \rho 3D (t), I2D (t), I3D (t), полученных с помощью
компьютерного моделирования, свидетельствует о соответствии теоретических и модельных
результатов. Отличие значений плотности электронов \rho 3D (t) от асимптотических значений
\rho \prime 3D (t) для t > 108 ns в трехмерном случае связано с высокой скоростью аннигиляции элек-
тронов и дырок и малым количеством частиц, оставшихся при больших временах t > 10 ns.
                                                                                               8




                                                                                                 а)




                                                                                                 б)
Рис. 3. Значения характеристик, полученные при КА-моделировании чистого туннелирования: а)
плотность электронов, б) интенсивность фотолюминесценции.

   В случае чистой диффузии (случай 2) коэффициенты диффузии частиц равны Dn =
Dp = 0.5 nm2 \cdot  ns - 1 , радиационная рекомбинация вследствие туннелирования и рекомби-
национные центры отсутствуют: B0 = 0 ns
                                                  - 1 , C                                        - 1
                                                          Nn (0) = CNp (0) = 0, bn0 = bp0 = 0 ns .
На Рис. 4 показан характер эволюции \aleph , моделирующей динамику процесса рекомбинации
вследствие диффузии частиц. В этом случае рекомбинация электронов и дырок происходит
при непосредственном столкновении частиц в одной клетке, пространственного разделения
электронов и дырок не происходит из-за постоянного перемешивания частиц вследствие
диффузии. Аналогично случаю чистого туннелирования, в трехмерном пространстве ско-
рость аннигиляции частиц выше, чем в двумерном.
    Значения характеристик, полученные в результате КА-моделирования чистой диффу-
зии в двумерном и трехмерном случаях, представлены на Рис. 5. В двумерном случае гра-
фики \rho 2D (t) и I2D (t), вычисленные с помощью компьютерного моделирования, параллельны
                                                        \prime               \prime 
графикам асимптотического поведения \rho 2D (t), I2D (t), посчитанным по формуле (4) с d = 2.
Тогда как в трехмерном случае графики \rho 3D (t) и I3D (t) отличаются от асимптотического
             \prime    \prime                                        \~3D (t) = 1/t, I\~3D (t) = 1/t ,
поведения \rho 3D (t), I3D (t) (4), но хорошо согласуются с кривыми \rho 
                                                                                                       2

которые являются асимптотическими значениями характеристик для случая равномерно-
го перемешивания частиц в объеме. Полученный результат показывает, что диффузия в




                                                                  174
Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                             agora.guru.ru/pavt




          а)                      б)                                в)                                   г)
Рис. 4. Пространственное распределение электронов и дырок в случае чистой диффузии:
а) d = 2, t = 374 ns, б) d = 2, t = 1.6 \cdot  104 ns; в) d = 3, t = 374 ns, г) d = 3, t = 1.6 \cdot  104 ns.
трехмерном случае приводит к равномерному однородному распределению частиц во всем
объеме.




                                                                                              а)




                                                                                              б)
Рис. 5. Значения характеристик, полученные при КА-моделировании радиационной рекомбинации
вследствие чистой диффузии частиц в двумерном и трехмерном случаях: а) плотность электронов,
б) интенсивность фотолюминесценции.

    Характер эволюции \aleph  для случая радиационной и нерадиационной рекомбинации (слу-
чай 3), когда присутствуют рекомбинационные центры CNn (0) = 5000, CNp (0) = 0, bn0 =
bp0 = 0.02 ns - 1 , но отсутствует диффузия частиц Dn = Dp = 0 nm2 \cdot  ns - 1 , показан на
Рис. 6.
    В присутствии рекомбинационных центров так же, как и в случае чистого туннелирова-
ния, в результате аннигиляции электронов и дырок, расположенных на близком расстоянии,
происходит формирование кластеров частиц. Наличие рекомбинационных центров приво-
дит к увеличению скорости аннигиляции электронов и дырок, причем в трехмерном случае
введение центров рекомбинации более существенно влияет на процесс, чем в двумерном
случае. В отличие от Рис. 2 и Рис. 4 на Рис. 6 показаны пространственные распределения
частиц в двумерной и трехмерной области на разных временных шагах. Это связано с тем,



                                                   175
Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                         agora.guru.ru/pavt




            а)                б)                     в)                               г)
Рис. 6. Пространственное распределение электронов и дырок в случае радиационной и нерадиа-
ционной рекомбинации: а) d = 2, t = 374 ns, б) d = 2, t = 1.06 \cdot  105 ns; в) d = 3, t = 35 ns,
г) d = 3, t = 374 ns.

что в трехмерном случае к моменту времени t = 1.06 \cdot  10
                                                              5 ns все частицы аннигилировали.

    Графики значений плотности электронов и интенсивности фотолюминесценции, полу-
ченные в результате КА-моделирования рекомбинации при наличии рекомбинационных
центров для двумерного и трехмерного случаев, представлены на Рис. 7. При сравнении
Рис. 7 и Рис. 3 видно, что при добавлении рекомбинационных центров значения характе-
ристик \rho 2D (t) и I2D (t) отклонились от значений, вычисленных для случая чистого тунне-
лирования, меньше, чем значения \rho 3D (t) и I3D (t).




Рис. 7. Значения характеристик, полученные при КА-моделировании рекомбинации при наличии
рекомбинационных центров в двумерном и трехмерном случаях.

    Характер эволюции \aleph  в случае радиационной и нерадиационной рекомбинации (слу-
чай 4) при наличии рекомбинационных центров CNn (0) = 5000,     CNp (0) = 0, bn0 = bp0 =
        - 1                                 2         - 1
0.02 ns , и диффузии частиц Dn = Dp = 0.5 nm \cdot ns представлен на Рис. 8. В результате
аннигиляции плотность частиц быстро уменьшается, и, несмотря на перемешивание частиц
вследствие диффузии, на каждом временном шаге можно выделить кластеры электронов
и дырок, но эти кластеры быстро распадаются и носят случайный характер.
    На Рис. 9 представлены графики характеристик, полученные в результате КА-моделиро-
вания рекомбинации при наличии рекомбинационных центров и диффузии частиц (слу-
чай 4). Аналогично остальным случаям (1, 2, 3), плотность электронов и интенсивность
фотолюминесценции, вычисленные для трехмерной области, убывают быстрее, чем \rho 2D (t)
и I2D (t). Кроме того, при сравнении Рис. 7 и Рис. 9 видно, что диффузия оказывает большее
влияние на скорость рекомбинации в двумерном случае: значения \rho 3D (t) и I3D (t), получен-
ные при наличии рекомбинационных центров, практически не изменились при добавлении
диффузии частиц, тогда как \rho 2D (t) и I2D (t) убывают существенно быстрее по сравнению
со случаем 3.




                                               176
Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                         agora.guru.ru/pavt




            а)                      б)              в)                       г)
Рис. 8. Пространственное распределение электронов и дырок в случае радиационной и неради-
ационной рекомбинации при наличии диффузии: а) d = 2, t = 35 ns, б) d = 2, t = 374 ns;
в) d = 3, t = 35 ns, г) d = 3, t = 374 ns.




Рис. 9. Значения характеристик, полученные при КА-моделировании рекомбинации при наличии
рекомбинационных центров и диффузии частиц в двумерном и трехмерном случаях.

6. Заключение

   На основе системы уравнений Смолуховского (1) и алгоритма Монте-Карло [3, 6] раз-
работана клеточно-автоматная модель рекомбинации электронов и дырок в двумерном и
трехмерном пространстве. Для моделирования областей больших размеров и вычисления
интегральных характеристик процесса (плотности частиц и интенсивности фотолюминес-
ценции) разработаны параллельные программы, реализующие КА-модель рекомбинации в
двумерном и трехмерном случаях. Параллельные программы реализованы с помощью тех-
нологии OpenMP. Эффективность распараллеливания на кластере НКС-30Т Сибирского
Суперкомпьютерного Центра СО РАН выше 85% при использовании 12 потоков.
   С помощью разработанной параллельной программы изучена кинетика процесса реком-
бинации для четырех случаев: 1) радиационной рекомбинации электронов и дырок вслед-
ствие чистого туннелирования, 2) радиационной рекомбинации электронов и дырок вслед-
ствие чистой диффузии, 3) радиационной и нерадиационной рекомбинации при наличии
центров рекомбинации, 4) радиационной и нерадиационной рекомбинации при наличии цен-
тров рекомбинации и диффузии частиц. В результате исследования эволюции КА-модели
\aleph  обнаружено формирование кластеров электронов и дырок при отсутствии диффузии ча-
стиц. Выявлено, что во всех четырех случаях электроны и дырки аннигилируют быстрее в
трехмерной области. Кроме того, наличие центров рекомбинации оказывает большее вли-
яние на процесс рекомбинации в трехмерном случае, тогда как диффузия существеннее
влияет на скорость рекомбинации в двумерном случае.




                                               177
Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                         agora.guru.ru/pavt



Литература

1. Туркин А. Нитрид галлия как один из перспективных материалов в современной
   оптоэлектронике // Компоненты и технологии. 2011. № 5. С. 6–10.


2. Gorgis A., Flissikowski T., Brandt O., Cheze C., Geelhaar L., Riechert H., and Grahn H.T.
   Time-resolved photoluminescence spectroscopy of individual GaN nanowires // Physical
   review B. 2012. No. 86. 041302(R).


3. Sabelfeld K.K., Brandt O., Kaganer V.M. Stochastic model for the fluctuation-limited
   reaction-diffusion kinetics in inhomogeneous media based on the nonlinear Smoluchowski
   equations // J. Math. Chem. 2015. Vol. 53, Issue 2. P. 651–669.


4. Kolodko A.A. and Sabelfeld K.K. Stochastic Lagrangian model for spatially inhomogeneous
   Smoluchowski equation governing coagulating and diffusing particles // Monte Carlo
   Methods and Applications. 2001. Vol. 7, No. 3-4. P. 223–228.


5. Kolodko A., Sabelfeld K. and Wagner W. A stochastic method for solving Smoluchowski’s
   coagulation equation // Mathematics and Computers in Simulation. 1999. Vol. 49, No 1-2.
   P. 57–79.


6. Sabelfeld K.K., Levykin A.I., Kireeva A.E // Stochastic simulation of fluctuation-induced
   reaction-diffusion kinetics governed by Smoluchowski equations, Monte Carlo Methods and
   Applications. 2015. Vol. 21, No 1. P. 33–48. (DOI:10.1515/mcma-2014-0012)


7. Toffoli T., Margolus N. Cellular Automata Machines: A New Environment for Modeling //
   USA: MIT Press, 1987. 259 p.


8. Бандман О.Л. Клеточно-автоматные модели пространственной динамики //
   Системная информатика. Методы и модели современного программирования. 2006.
   № 10, С. 59–113.


9. Bandman O.L. Mapping physical phenomena onto CA-models // AUTOMATA-2008. In:
   Adamatzky A., Alonso-Sanz R., Lawniczak A., Martinez G.J., Morita K., Worsch T. (eds.)
   Theory and Applications of Cellular Automata. Luniver Press, UK. 2008. P. 381–397.


10. Ермаков С.М., Михайлов Г.А. Статистическое моделирование // М.: ФИЗМАТЛИТ,
   2-е изд., дополн., 1982. 296 C.




                                               178
Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                          agora.guru.ru/pavt




Parallel implementation of stochastic cellular automata
            model of electron-hole recombination in a
                                     semiconductor

                                                1                    1
                                A.E. Kireeva , K.K. Sabelfeld ,

                                                                                                     1
   Institute of Computational Mathematics and Mathematical Geophysics SB RAS



        Parallel programs implementing stochastic cellular automata (CA) model of electron-
        hole recombination in an inhomogeneous semiconductor for two- and three-dimensional
        cases are developed. The spatio-temporal distributions of particles are investigated by
        the CA simulation. Spatial separation of electrons and holes with clusters formation is
        found and analyzed. Parallel implementation of the CA model allows us to calculate
        integral characteristics of the recombination process (particle densities and radiative
        intensity) in acceptable time. Recombination kinetics in the vicinity of the recombination
        centers and diffusion in two- and three-dimensional space is investigated using the
        parallel program.
        Keywords:  electron-hole recombination, semiconductor, parallel implementation, stochastic
        cellular automata, radiative intensity.

References

1. Turkin А. Nitrid galliya kak odin iz perspektivnykh materialov v sovremennoy
   optoelektronike [Nitride Gaul as one of the perspective materials in modern optoelectronics]
   // Komponenty i tekhnologii [Components and technologies], 2011. No. 5. P. 6–10.


2. Gorgis A., Flissikowski T., Brandt O., Cheze C., Geelhaar L., Riechert H., and Grahn H.T.
   Time-resolved photoluminescence spectroscopy of individual GaN nanowires // Physical
   review B, 2012. No 86. 041302(R).


3. Sabelfeld K.K., Brandt O., Kaganer V.M. Stochastic model for the fluctuation-limited
   reaction-diffusion kinetics in inhomogeneous media based on the nonlinear Smoluchowski
   equations // J. Math. Chem, 2015. Vol. 53, Issue 2. P. 651–669.


4. Kolodko A.A. and Sabelfeld K.K. Stochastic Lagrangian model for spatially inhomogeneous
   Smoluchowski equation governing coagulating and diffusing particles // Monte Carlo
   Methods and Applications, 2001. Vol. 7, No. 3-4. P. 223–228.


5. Kolodko A., Sabelfeld K. and Wagner W. A stochastic method for solving Smoluchowski’s
   coagulation equation // Mathematics and Computers in Simulation, 1999. Vol. 49, No 1-2.
   P. 57–79.


6. Sabelfeld K.K., Levykin A.I., Kireeva A.E // Stochastic simulation of fluctuation-induced
   reaction-diffusion kinetics governed by Smoluchowski equations, Monte Carlo Methods and
   Applications, 2015. Vol. 21, No 1. P. 33–48. (DOI:10.1515/mcma-2014-0012)


7. Toffoli T., Margolus N. Cellular Automata Machines: A New Environment for Modeling //
   USA: MIT Press, 1987. 259 p.


8. Bandman O.L. Kletochno-avtomatnye modeli prostranstvennoy dinamiki [Cellular
   automata model of spatial dynamics] // Sistemnaya informatika. Metody i modeli
   sovremennogo programmirovaniya [System informatics. Methods and models of modern
   programming], 2006. No 10. P. 59–113.




                                                 179
Параллельные вычислительные технологии (ПаВТ’2016) || Parallel computational technologies (PCT’2016)

                                         agora.guru.ru/pavt



 9. Bandman O.L. Mapping physical phenomena onto CA-models // AUTOMATA-2008. In:
   Adamatzky A., Alonso-Sanz R., Lawniczak A., Martinez G.J., Morita K., Worsch T. (eds.)
   Theory and Applications of Cellular Automata. Luniver Press, UK. 2008. P. 381–397.


10. Ermakov S.M., Mikhaylov G.A. Statisticheskoe modelirovanie [Statistical modeling] // M.:
   FIZMATLIT, 2-e izd., dopoln. [Moscow, publishing in FIZMATLIT, the second edition
   augmented], 1982. 296 P.




                                               180