On the application of the Galerkin projection method to the nonstationary diffusion equation with a variable coefficient

Full Text

1. Введение

В [11] рассмотрен процесс нестационарной диффузии неравновесных неосновных носителей заряда (ННЗ) в однородной полупроводниковой мишени после прекращения воздействия низкоэнергетического (менее 10 кэВ) остро сфокусированного электронного пучка, электронного зонда, на полупроводниковую мишень. Рассматривался один канал рекомбинации, и потому электрофизические параметры мишени — коэффициенты уравнения диффузии — при проведении расчетов полагались постоянными. В то же время имеющиеся экспериментальные результаты (см. [5, 6, 10]) позволяют предположить наличие нескольких каналов рекомбинации в этих мишенях. В настоящей работе рассмотрена математическая модель, описывающая два механизма рекомбинации неравновесных ННЗ при выключении внешнего воздействия, при этом зависимость числа неравновесных ННЗ от времени описывается двумя экспонентами. В этом случае коэффициенты уравнения диффузии будут переменными.

Ранее метод Галеркина эффективно применялся для решения стационарного уравнения диффузии с постоянными коэффициентами [3], и была получена порядковая оценка погрешности невязки, соответствующей приближенному уравнению диффузии. Настоящая работа продолжает такие исследования. Рассмотрение проведено для нестационарного уравнения диффузии, описывающего зависимость от времени концентрации неравновесных ННЗ, генерированных электронным зондом в однородном полупроводнике, после прекращения действия внешнего низкоэнергетического источника.

2. Постановка задачи

При стационарном облучении полупроводникового материала электронным зондом в последнем генерируются неравновесные ННЗ, после чего происходит их диффузия и рекомбинация в объёме полупроводника. В математических моделях, описывающих стационарное квазиравновесное состояние, скорости генерации и рекомбинации ННЗ считаются одинаковыми [4, 8]. Но после прекращения действия электронного зонда генерации ННЗ не происходит, и концентрация ННЗ в полупроводнике  при наличии зависимости времени жизни ННЗ  от времени  может быть найдена как решение следующего нестационарного уравнения диффузии [7]:

$\frac{\partial c\left(x,y,t\right)}{\partial t}=D\text{Δ}c\left(x,y,t\right)-\frac{c\left(x,y,t\right)}{\tau \left(t\right)}$ (1)

 с граничными условиями

$\begin{array}{l}c\left(x,y\mathrm{,0}\right)=n\left(x,y\right), c\left(\pm \infty ,y,t\right)=\mathrm{0,} c\left(x,\pm \infty ,t\right)=\mathrm{0,}\\ -\infty <x<\infty , -\infty <y<\infty , 0\le t<\infty .\end{array}$ (2)

Здесь $\text{Δ}={\partial }^2/\partial x^2+{\partial }^2\partial y^2$ — двумерный оператор Лапласа, $D$ — коэффициент диффузии и $\tau \left(t\right)$ — зависимость времени жизни ННЗ от времени после выключения импульса возбуждения. Функция $n\left(x,y\right)$ удовлетворяет стационарному дифференциальному уравнению, описывающему диффузию ННЗ в состоянии квазиравновесия (до выключения электронного пучка):

$D\text{Δ}n\left(x,y\right)-\frac{n\left(x,y\right)}{\tau }=-\rho \left(x,y\right).$ (3)

В состоянии квазиравновесия $\tau =const$, а $\rho \left(x,y\right)$ — функция, которая описывает концентрацию генерированных в единицу времени неравновесных ННЗ до их диффузии в мишени. Функция $\rho \left(x,y\right)$ пропорциональна ${\rho }^{\text{*}}\left(x,y\right)$ — плотности мощности, выделяемой электронами пучка в мишени; для полупроводниковых материалов $\rho \left(x,y\right)$ может быть получена делением ${\rho }^{\text{*}}\left(x,y\right)$ на энергию образования электронно-дырочной пары или экситона. В случае узкого электронного пучка, электронного зонда, зависимость ${\rho }^{\text{*}}\left(x,y\right)$ может быть описана функцией распределения типа Гаусса, что даёт для правой части (3) соотношение, позволяющее найти $\rho \left(x,y\right)$ (см. [8, 11]):

${\rho }^{\text{*}}\left(x,y\right)=\frac{\mathrm{1,085} \left(1-\eta \right)E_0}{{\pi }^{3/2} a_1^2{z}_{ms}\left(1-\eta +\eta z_{ss}/z_{ms}\right)}\left\{\text{exp}\left\{-\left[\frac{x^2+y^2}{a_1^2}\right]\right\}+\frac{\eta a_1^2}{\left(1-\eta \right)a_2^2}\text{exp}\left\{-\left[\frac{x^2+y^2}{a_2^2}\right]\right\}\right\}.$

Здесь $E_0$ — энергия электронов зонда, рассеянная в мишени в единицу времени, $z_{ms}$ — глубина максимальных потерь энергии первичными электронами, испытавшими малоугловое рассеяние и поглощёнными в мишени; $z_{ss}$ — глубина максимальных потерь энергии обратно рассеянными электронами, вышедшими из мишени; $\eta$ — коэффициент обратного рассеяния электронов зонда. Параметры $a_1$ и $a_2$ могут быть определены из соотношений

$a_1^2=z_{ms}^2+\mathrm{0,72} d_z^2, a_2^2=\mathrm{0,25} z_{ss}^2+\mathrm{0,72} d_z^2,$

где $d_z$ — диаметр электронного зонда; для остро сфокусированного пучка электронов $d_z=0$.

При построении математической модели нестационарной диффузии с переменным коэффициентом может использоваться подход, предложенный нами в [7] для описания диффузии экситонов, генерированных низкоэнергетическим электронным зондом в монокристаллическом нитриде галлия — перспективном материале полупроводниковой микро-, оптоэлектроники и СВЧ-техники. При наличии двух независимых каналов рекомбинации профиль спада концентрации экситонов может быть описан суммой двух экспонент и их эффективное время жизни $\tau \left(t\right)$ в уравнении (1) находится по следующей формуле (см. [7]):

$\tau \left(t\right)=-t/\text{ln}\left(\alpha \text{exp}\left(-t/{\tau }_1\right)+\left(1-\alpha \right)\text{exp}\left(-t/{\tau }_2\right)\right).$

Здесь ${\tau }_1$ и ${\tau }_2$ — время жизни экситонов для первого и второго канала рекомбинации соответственно, $\alpha$ — безразмерный параметр.

Используя замену

$c\left(x,y,t\right)=\text{exp}\left(-{\int }_0^t\frac{d}{\tau \left(s\right)}\right)v\left(x,y,t\right),$

от задачи (1), (2) перейдем к уравнению

$\frac{\partial v\left(x,y,t\right)}{\partial t}=D\text{Δ}v\left(x,y,t\right)$ (4)

с начальным условием

$v\left(x,y\mathrm{,0}\right)=n\left(x,y\right)$ (5)

В [12] получено решение задачи (4), (5):

$V\left(r,t\right)=\frac{1}{2Dt}{\int }_0^{+\infty }\text{exp}\left(-\frac{r^2+{\xi }^2}{4Dt}\right)I_0\left(\frac{r\xi }{2Dt}\right)n\left(\xi \right)\xi d\xi .$

Тогда решение исходной задачи (1)–(3) примет следующий вид:

$c\left(r,t\right)=\frac{1}{2Dt}\text{exp}\left(-{\int }_0^t\frac{ds}{\tau \left(s\right)}\right){\int }_0^{+\infty }\text{exp}\left(-\frac{r^2+{\xi }^2}{4Dt}\right)I_0\left(\frac{r\xi }{2Dt}\right)n\left(\xi \right)\xi d\xi .$ (6)

Аналитическое решение задачи (1), (2) с постоянным коэффициентом $\tau$ в уравнении (1) (т.е. если имеется один канал рекомбинации) имеет вид (см. [12]):

$c\left(r,t\right)=\frac{\text{exp}\left(-t/\tau \right)}{2Dt}{\int }_0^{+\infty }\text{exp}\left(-\frac{r^2+{\xi }^2}{4Dt}\right)I_0\left(\frac{r\xi }{2Dt}\right)n\left(\xi \right)\xi d\xi .$ (7)

Здесь $r$ — полярный радиус; $I_0\left(x\right)$ — модифицированные функции Бесселя нулевого порядка первого рода.

В [3, 11] рассмотрены некоторые возможности использования проекционного метода Галеркина для моделирования двумерной диффузии ННЗ с постоянными электрофизическими параметрами в полупроводниковом материале, и получена порядковая оценка погрешности невязки, соответствующей приближенному уравнению диффузии. В настоящей работе предлагается использование этого метода для нахождения концентрации ННЗ с переменным электрофизическим параметром и в качестве примера использовано два независимых канала рекомбинации.

3. Проекционная аппроксимация исходной модели, основанная на применении метода Галеркина

Переходя к цилиндрической системе координат, от уравнения (1) перейдем к следующему уравнению:

$D\tau \left(t\right)r\frac{{\partial }^{2}c\left(r,t\right)}{\partial r^2}+D\tau \left(t\right)\frac{\partial c\left(r,t\right)}{\partial r}-\tau \left(t\right)r\frac{\partial c\left(r,t\right)}{\partial t}-rc\left(r,t\right)=0$ (8)

с граничными условиями

$c\left(r\mathrm{,0}\right)=n\left(r\right), c\left(+\infty ,t\right)=0.$ (9)

Тогда уравнение (3) примет вид

$D\tau r\frac{d^{2}n\left(r\right)}{d{r}^2}+D\tau \frac{dn\left(r\right)}{dr}-rn\left(r\right)=-\tau r\rho \left(r\right).$ (10)

 Сравнение формул (6) и (7) позволяет привести следующий алгоритм метода Галеркина.

Для нахождения приближенного решения уравнения (8)–(10) найдем сначала приближенное решение уравнения (8) с постоянным коэффициентом $\tau$. Поскольку одно из граничных условий задано на бесконечности, то для реализации проекционного метода Галеркина выберем двумерный базис из модифицированных функций Лагерра с параметрами, ускоряющими сходимость ряда (см. [9]):

${\phi }_{i,j}\left(r,t\right)={\phi }_i^{{\alpha }_1,{\gamma }_1}\left(r\right){\phi }_j^{{\alpha }_2,{\gamma }_2}\left(t\right)=\text{exp}\left(-\frac{{\gamma }_{1}r}{2}\right)L_i\left({\gamma }_{1}r;{\alpha }_1\right)\text{exp}\left(-\frac{{\gamma }_{2}t}{2}\right)L_j\left({\gamma }_{2}t;{\alpha }_2\right),$

которые определяются через многочлены Чебышева– Лагерра $L_i\left({\gamma }_{1}r;{\alpha }_1\right)$ по переменной $r$ и многочлены $L_j\left({\gamma }_{2}t;{\alpha }_2\right)$ по переменной $t$, $i,j=\mathrm{0,1,2,}\dots$. Здесь параметры ${\alpha }_1>-1$, ${\alpha }_2>-1$ и ${\gamma }_1>0$, ${\gamma }_2=2/\tau >0$ используются для оптимизации вычислительной схемы.

В методе Галеркина предполагается, что неизвестная функция $c\left(r,t\right)$ может быть достаточно точно представлена приближенным решением:

$c_{m+1}\left(r,t\right)=c_m\left(r,t\right)+c_0\left(r,t\right),$ (11)

где

$c_m\left(r,t\right)=\sum _{i=0}^{m-1}\sum _{j=0}^{m-1}{c}_{ij}{\phi }_i^{{\alpha }_1,{\gamma }_1}\left(r\right){\phi }_j^{{\alpha }_2,{\gamma }_2}\left(t\right)$

 — прямоугольная частичная сумма двойного ряда Фурье"– Лагерра порядка $m\times m$ функции $c\left(r,t\right)$,

$c_0\left(r,t\right)=c_0\left(r\right)\frac{{\phi }_m^{{\alpha }_2,{\gamma }_2}\left(t\right)}{{\phi }_m^{{\alpha }_2,{\gamma }_2}\left(0\right)}, c_0\left(r\right)=\sum _{i=0}^{m-1}{c}_i^n{\phi }_i^{{\alpha }_1,{\gamma }_1}\left(r\right)-\sum _{i=0}^{m-1}\sum _{j=0}^{m-1}{c}_{ij}{\phi }_i^{{\alpha }_1,{\gamma }_1}\left(r\right){\phi }_j^{{\alpha }_2,{\gamma }_2}\left(0\right).$

Функция $c_0\left(r,t\right)$ введена, чтобы удовлетворить граничным условиям, $c_i^n$ — коэффициенты разложения функции $n\left(r\right)$, которые находятся из решения уравнения (6), а неизвестные коэффициенты разложения $c_{{}_{ij}}$ искомой функции $c\left(r,t\right)$ согласно методу Галеркина определяются из решения следующей системы уравнений:

$\left(R, {\phi }_k^{{\alpha }_1,{\gamma }_1}\left(r\right){\phi }_l^{{\alpha }_2,{\gamma }_2}\right)=\mathrm{0,} k,l=\overline{\mathrm{0,}m-1},$ (12)

где

$R=D\tau r\frac{{\partial }^2{c}_{m+1}\left(r,t\right)}{\partial r^2}+D\tau \frac{\partial c_{m+1}\left(r,t\right)}{\partial r}-\tau r\frac{\partial c_{m+1}\left(r,t\right)}{\partial t}-r{c}_{m+1}\left(r,t\right)$

— невязка уравнения (8).

Обозначим столбцы (растянутые в столбцы матрицы) из коэффициентов разложения неизвестной функции $c\left(r,t\right)$ и функции $n\left(r\right)$ по выбранному базису через $C_{mm}$ и $C_m^n$ соответственно. Используя кронекерово произведение, введём матрицы дифференцирования по переменной $r$ и по переменной $t$:

$D_{mm}^r=D_m^{{\gamma }_1}\otimes E_m, D_{mm}^t=E_m\otimes D_m^{{\gamma }_2},$

где $E_m$ — единичная матрица, а $D_m^{\gamma }$ — матрица дифференцирования в одномерном базисе из модифицированных функций Лагерра, элементы которой находятся с помощью элементарных алгебраических операций (см. [3]):

$d_{ij}=\left\{\begin{array}{ll}-\gamma /\mathrm{2,}& i=j,\\ -\gamma ,& i<j,\\ \mathrm{0,}& i>j.\end{array}\right.$

Обозначим через $B_1$ и $B_2$ матрицы, элементы которых находятся по формулам

 $b_{ij}^1=\left\{\begin{array}{ll}{\gamma }_1\left(j-1\right)/\mathrm{4,}& i+1=j,\\ {\gamma }_1\left(1+3{\alpha }_1/4-\left(i-1\right)/4\right),& j+1=i,\\ {\gamma }_1\left(3/4+{\alpha }_1/4-i/2\right),& i=j,\\ \mathrm{0,}& i+1<j,\\ {\gamma }_1\left(1+{\alpha }_1\right),& j+1<i,\end{array}\right. b_{ij}^2=\left\{\begin{array}{ll}-i/{\gamma }_1,& i+1=j,\\ -\left(i-1+{\alpha }_1\right)/{\gamma }_1,& j+1=i,\\ \left(2i-1+{\alpha }_1\right)/{\gamma }_1,& i=j,\\ \mathrm{0,}& i+1<j,\\ \mathrm{0,}& j+1<i,\end{array}\right.$

которые устанавливаются на основании известных рекуррентных соотношений для функций Лагерра. Затем, аппроксимируя дифференциальное уравнение (10), имеем соответствующее матричное уравнение:

$\left(D\tau B_1^T+D\tau D^{{\gamma }_1}-B_2^T\right)C^n=\tau C^{\text{Φ}}.$ (13)

Здесь $C^{\text{Φ}}$ — столбец из коэффиентов разложения функции, стоящей в правой части уравнения (10). Из уравнения (13), находим столбец $C^n$.

Далее введем матрицы

$D_{mm}^{r2}=B_1^T\otimes E_m, D_{mm}^{t2}=B_2^T\otimes D_m^{{\gamma }_2}, I{i}_{mm}=B_2^T\otimes E_m.$

Перепишем систему (8) в матричном виде:

$\left(D\tau D_{mm}^{r2}-\tau D_{mm}^{t2}+D\tau D_{mm}^r+\tau {\tilde{D}}_{mm}^t-I{i}_{mm}\right)C_{mm}=\tau C_{mm}^{{\text{Φ}}_1},$ (14)

где ${\tilde{D}}_{mm}^t=B_2^T\otimes B_3$ — матрица с элементами $b_{ij}^3=-{\gamma }_2\text{Γ}\left({\alpha }_2+i\right)/\left(i-1\right)!$, а столбец $C_{mm}^{{\text{Φ}}_1}$ находится путем растяжения матрицы:

$C^{{\text{Φ}}_1}=B_2^T{C}_m^n{(-{\gamma }_2{\phi }_m^{{\gamma }_2}\left(0\right))}^T.$

Здесь ${\phi }_m^{{\gamma }_2}\left(t\right)$ — столбец из $m$ первых базисных функций.

Подстановка величин $c_{ij}$, определяемых из решения системы уравнений (14), в формулу (11), дает приближенное решение. Найденное приближенное решение уравнения (8) с постоянным коэффициентом $\tau$ фактически представляет собой разложение в ряд Фурье функции

$V\left(r,t\right)=\text{exp}\frac{t}{\tau }c\left(r,t\right)\approx V_m\left(r,t\right)=\sum _{i=0}^{m-1}\sum _{j=0}^{m-1}{c}_{ij}^v{\phi }_i^{{\alpha }_1,{\gamma }_1}\left(r\right)L_j\left(\frac{2t}{\tau };{\alpha }_2\right)+c_0\left(r\right)L_m\left(\frac{2t}{\tau }\right).$

Тогда искомое приближенное решение уравнения (8) с переменным коэффициентом в силу (6) может быть найдено по формуле

$c\left(r,t\right)\approx \text{exp}\left(-{\int }_0^t\frac{ds}{\tau \left(s\right)}\right)V_m\left(r,t\right).$

Интеграл в последней формуле может быть найден приближенно путем разложения подынтегральной функции в ряд Тейлора.

4. Условие сходимости

Используя результаты работы [3], нетрудно установить оценку невязки неоднородного уравнения (8) с постоянным коэффициентом $\tau$ и с функцией $f\left(r,t\right)\in L_2$, стоящей в правой части этого уравнения и имеющей непрерывные частные производные до порядка $2n$ по обоим пространственным направлениям. Покажем, что невязка $R$ сходится к нулю в среднем (в пределе при $m\to \infty$). Следуя [1], введем обозначения:

$D=r\frac{d^2}{d{r}^2}+\left({\alpha }_1-r+1\right)\frac{d}{dr}+t\frac{d^2}{d{t}^2}+\left({\alpha }_2-t+1\right)\frac{d}{dt};$

$L_2^n\left(D\right)$, $n=\mathrm{0,1,}\dots$, — класс таких функций $f$, что функции

$\tilde{f}\left(r,t\right)=f\left(r,t\right)\text{exp}\left(\frac{{\gamma }_{1}r}{2}\right)\text{exp}\left(\frac{{\gamma }_{2}t}{2}\right)\in L_2$

имеют обобщенные частные производные в смысле Леви ${\partial }^k\tilde{f}\left(r,t\right)/\partial r^i\partial t^j$, $i+j=k$, $k=\mathrm{0,1,}\dots$, принадлежащие пространству $L_2$, для которых $D^n\tilde{f}\in L_2$, $n=\mathrm{0,1,}\dots$, где $D^0\tilde{f}=\tilde{f}$, $D^n\tilde{f}=D\left(D^{n-1}\tilde{f}\right)$, $n=\mathrm{1,2,}\dots$, $L_2^0\left(D\right)=L_2$. Далее нам понадобится вспомогательная теорема, которой будем пользоваться ниже.

Теорема 4.1. Для любой функции $f\left(r,t\right)\in L_2^n\left(D\right)$ справедлива оценка (см. [1, 3])

$\parallel f-S_{m,m}\left(f\right)\parallel \le {[1-{(1-h)}^m]}^{-k}{m}^{-n} {\text{Ω}}_k\left(D^n\tilde{f};h\right),$ (15)

где $n=\mathrm{0,1,}\dots$, $0<h<1$. Величины ${\text{Ω}}_k\left(D^n\tilde{f};h\right)$ — обобщённые модули непрерывности функции $D^n\tilde{f}$,

$S_{m,m}\left(f\right)=\sum _{i=0}^{m-1}\sum _{j=0}^{m-1}{c}_{ij}^f{\phi }_i^{{\alpha }_1,{\gamma }_1}\left(r\right){\phi }_j^{{\alpha }_2,{\gamma }_2}\left(t\right)$

— прямоугольная частичная сумма двойного ряда Фурье функции $f\left(r,t\right)$.

Если погрешности в исходных данных и погрешности вычислений отсутстствуют, а учитываются лишь погрешности аппроксимаций, то, опираясь на результаты работы [3], в которой получены оценки для параметров ${\alpha }_1={\alpha }_2=0$, аналогично установим оценку для невязки при ${\alpha }_1>0$, ${\alpha }_2>0$.

Теорема 4.2. Пусть функция $f\left(r,t\right)\in L_2$, $n>2$, имеет непрерывные частные производные до порядка $2n$ по обоим пространственным направлениям. Тогда справедлива следующая оценка:

$\parallel R{\parallel }_{L_2}<\left(\frac{\sqrt{\text{Γ}\left({\alpha }_1+3\right)}}{{\gamma }_1}\left(3D\tau +2\tau +1\right)+2D\tau \right){[1-{(1-m^{-1})}^m]}^{-k}{m}^{-n} \text{Ω}\left(m^{-1}\right)+$

$+O\left(m^{-n+5/4}\right)\omega \left(m^{-1}\right), {\alpha }_2\ge \mathrm{0,} m\to \infty .$

Здесь $\text{Ω}\left(m^{-1}\right)$ — мажоранта обобщенных модулей непрерывности для функций из пространства $L_2$, а $\omega \left(m^{-1}\right)$ — заданная мажоранта модулей непрерывности дифференцируемых функций.

Доказательство. Для точного решения имеем

$D\tau r\frac{{\partial }^{2}c\left(r,t\right)}{\partial r^2}+D\tau \frac{\partial c\left(r,t\right)}{\partial r}-\tau r\frac{\partial c\left(r,t\right)}{\partial t}-rc\left(r,t\right)=-f\left(r,t\right).$

Тогда справедливо неравенство

$\parallel R{\parallel }_{L_2}\le D\tau \parallel r{\parallel }_{L_2,r^{{\alpha }_1}{e}^{-r}}\frac{{\partial }^2\tilde{c}\left(r,t\right)}{\partial r^2}-{\frac{{\partial }^2{\tilde{c}}_{m+1}\left(r,t\right)}{\partial r^2}}_{L_2,r^{{\alpha }_1}{e}^{-r}{t}^{{\alpha }_2}{e}^{-t}}+$

 

$+D\tau \frac{\partial \tilde{c}\left(r,t\right)}{\partial r}-{\frac{\partial {\tilde{c}}_{m+1}\left(r,t\right)}{\partial r}}_{L_2,r^{{\alpha }_1}{e}^{-r}{t}^{{\alpha }_2}{e}^{-t}}+\parallel r{\parallel }_{L_2,r^{{\alpha }_1}{e}^{-r}}\parallel \tilde{c}\left(r,t\right)-{\tilde{c}}_{m+1}\left(r,t\right){\parallel }_{L_2,r^{{\alpha }_1}{e}^{-r}{t}^{{\alpha }_2}{e}^{-t}}+$

$+\tau \parallel r{\parallel }_{L_2,r^{{\alpha }_1}{e}^{-r}}\frac{\partial \tilde{c}\left(r,t\right)}{\partial t}-{\frac{\partial {\tilde{c}}_{m+1}\left(r,t\right)}{\partial t}}_{L_2,r^{{\alpha }_1}{e}^{-r}{t}^{{\alpha }_2}{e}^{-t}}.$

Обозначим слагаемые в правой части последнего неравенства через $S_1,\dots ,S_4$ и воспользуемся оценкой (15) и возможностью дифференцирования рядов Фурье– Лагерра (см. [1-3]), т.е.

 

$S_1=D\tau \parallel r{\parallel }_{L_2,r^{{\alpha }_1}{e}^{-r}}\frac{{\partial }^2\tilde{c}\left(r,t\right)}{\partial r^2}-{\frac{{\partial }^2{\tilde{c}}_{m+1}\left(r,t\right)}{\partial r^2}}_{L_2,r^{{\alpha }_1}{e}^{-r}{t}^{{\alpha }_2}{e}^{-t}}\le$

 

$\le \frac{\sqrt{\text{Γ}\left({\alpha }_1+3\right)} D\tau }{{\gamma }_1}{[1-{(1-h)}^m]}^{-k}{m}^{-n}{\text{Ω}}_k\left(D^n\left(\widetilde{c\text{'}}{\text{'}}_{rr}\left(r,t\right)+{\gamma }_1{\widetilde{c\text{'}}}_r\left(r,t\right)+\left({\gamma }_1^2/4\right)\tilde{c}\left(r,t\right)\right);h\right)+$

$+\frac{\sqrt{\text{Γ})\left({\alpha }_1+3\right)} D\tau }{{\gamma }_1}{[1-{(1-h)}^m]}^{-k-1}{m}^{-n}{\text{Ω}}_{k+1}\left(D^n\left({\gamma }_1{\widetilde{c\text{'}}}_r\left(r,t\right)+\left({\gamma }_1^2\right)\tilde{c}\left(r,t\right)\right);h\right)+$

$+\frac{\sqrt{\text{Γ}\left({\alpha }_1+3\right)} D\tau }{{\gamma }_1}{[1-{(1-h)}^m]}^{-k-2}{m}^{-n-1}{\text{Ω}}_{k+2}\left(D^{n+1}\left(\left({\gamma }_1^2/4\right)\tilde{c}\left(r,t\right)\right);h\right),$

$S_2=D\tau \frac{\partial \tilde{c}\left(r,t\right)}{\partial r}-{\frac{\partial {\tilde{c}}_{m+1}\left(r,t\right)}{\partial r}}_{L_2,r^{{\alpha }_1}{e}^{-r}{t}^{{\alpha }_2}{e}^{-t}}\le$

$\le D\tau {[1-{(1-h)}^m]}^{-k-1}{m}^{-n}{\text{Ω}}_{k+1}\left(D^n\left({\widetilde{c\text{'}}}_r\left(r,t\right)+\left({\gamma }_1/2\right)\tilde{c}\left(r,t\right)\right);h\right)+$

$+D\tau {[1-{(1-h)}^m]}^{-k-2}{m}^{-n-1}{\text{Ω}}_{k+2}\left(D^{n+1}\left(\left({\gamma }_1/2\right)\tilde{c}\left(r,t\right)\right);h\right),$

$S_3=\parallel r{\parallel }_{L_2,r^{{\alpha }_1}{e}^{-r}}\parallel \tilde{c}\left(r,t\right)-{\tilde{c}}_{m+1}\left(r,t\right){\parallel }_{L_2,r^{{\alpha }_1}{e}^{-r}{t}^{{\alpha }_2}{e}^{-t}}\le$

$\le \frac{\sqrt{\text{Γ}\left({\alpha }_1+3\right)}}{{\gamma }_1}{[1-{(1-h)}^m]}^{-k-2}{m}^{-n-1}{\text{Ω}}_{k+2}\left(D^{n+1}\left(\tilde{c}\left(r,t\right)\right);h\right),$

$S_4=\tau \parallel r{\parallel }_{L_2,r^{{\alpha }_1}{e}^{-r}}\frac{\partial \tilde{c}\left(r,t\right)}{\partial t}-{\frac{\partial {\tilde{c}}_{m+1}\left(r,t\right)}{\partial t}}_{L_2,r^{{\alpha }_1}{e}^{-r}{t}^{{\alpha }_2}{e}^{-t}}\le$

$\le \frac{\sqrt{\text{Γ}\left({\alpha }_1+3\right)}\tau }{{\gamma }_1}{[1-{(1-h)}^m]}^{-k-1}{m}^{-n}{\text{Ω}}_{k+1}\left(D^n\left({\widetilde{c\text{'}}}_t\left(r,t\right)+\left({\gamma }_2/2\right)\tilde{c}\left(r,t\right)\right);h\right)+$

$+\frac{\sqrt{\text{Γ}\left({\alpha }_1+3\right)}\tau }{{\gamma }_1}{[1-{(1-h)}^m]}^{-k-2}{m}^{-n-1}{\text{Ω}}_{k+2}\left(D^{n+1}\left(\left({\gamma }_2/2\right)\tilde{c}\left(r,t\right)\right);h\right)+$

$+\frac{\sqrt{\text{Γ}\left({\alpha }_1+3\right)}\tau }{{\gamma }_1} \frac{\parallel {\phi }_m^{\text{'}{\alpha }_2,{\gamma }_2}\left(t\right){\parallel }_{L_2,t^{{\alpha }_2}{e}^{-t}}}{{\phi }_m^{{\alpha }_2,{\gamma }_2}\left(0\right)} \parallel \left({\tilde{c}}_0\left(r\right)\right){\parallel }_{L_2,r^{{\alpha }_1}{e}^{-r}},$                           

Оценим последнее слагаемое в $S_4$:

$\frac{\sqrt{\text{Γ}\left({\alpha }_1+3\right)}\tau }{{\gamma }_1} \frac{\parallel {\phi }_m^{\text{'}{\alpha }_2,{\gamma }_2}\left(t\right){\parallel }_{L_2,t^{{\alpha }_2}{e}^{-t}}}{{\phi }_m^{{\alpha }_2,{\gamma }_2}\left(0\right)} \parallel \left({\tilde{c}}_0\left(r\right)\right){\parallel }_{L_2,r^{{\alpha }_1}{e}^{-r}}\le$

$\le \frac{\sqrt{\text{Γ}\left({\alpha }_1+3\right)}\tau }{{\gamma }_1}\frac{\parallel {\phi }_m^{\text{'}{\alpha }_2,{\gamma }_2}\left(t\right){\parallel }_{L_2,t^{{\alpha }_2}{e}^{-t}}}{{\phi }_m^{{\alpha }_2,{\gamma }_2}\left(0\right)}\sum _{i=0}^{m-1}{c}_i^n-\sum _{i=0}^{m-1}\sum _{j=0}^{m-1}{c}_{ij}{L}_j{\left(0;{\alpha }_2\right)}_{L_2}{\widehat{\phi }}_i^{{\alpha }_1,{\gamma }_1}{\left(r\right)}_{L_2}\le$

$\le \frac{\sqrt{\text{Γ}\left({\alpha }_1+3\right)}\tau }{{\gamma }_1}\frac{\parallel {\phi }_m^{\text{'}{\alpha }_2,{\gamma }_2}\left(t\right){\parallel }_{L_2,t^{{\alpha }_2}{e}^{-t}}}{{\phi }_m^{{\alpha }_2,{\gamma }_2}\left(0\right)}\sum _{i=0}^{m-1}\left|c_i^n-\sum _{j=0}^{m-1}{c}_{ij}{L}_j\left(0;{\alpha }_2\right)\right|$

 где ${\widehat{\phi }}_i^{{\alpha }_1,{\gamma }_1}\left(r\right)$ — ортонормированные функции Лагерра. Для оценки выражения

$\sum _{i=1}^{m-1}\left|c_i^n-\sum _{j=1}^{m-1}{c}_{ij}{L}_j\left(0;{\alpha }_2\right)\right|$

воспользуемся асимптотической формулой для коэффициентов Фурье–Лагерра (см. [2]):

$c_{im}=O\left(m^{-\left({\alpha }_2+2n+2\right)/2-1/4}\right)\omega \left(m^{-1/2}\right),$ (16)

учитывая, что $c_1{m}^{{\alpha }_2}\le L_m\left(0;{\alpha }_2\right)\le c_2{m}^{{\alpha }_2}$ (см. [1]). Здесь $c_1$ и $c_2$ — некоторые фиксированные положительные постоянные, а $\omega \left(t\right)$ — заданная мажоранта модулей непрерывности. Тогда в силу (16) и оценки

$\sum _{j=m}^{\infty }{j}^{-\beta }=O\left(m^{-\beta +1}\right), \beta >\mathrm{1,}$

для отклонения сумм Фурье–Лагерра находим

$\left|c_i^n-\sum _{j=0}^{m-1}{c}_{ij}{L}_j\left(0;{\alpha }_2\right)\right|=\left|\sum _{j=m}^{\infty }{c}_{ij}{L}_j\left(0;{\alpha }_2\right)\right|=$

$=O\left(1\right)\sum _{j=m}^{\infty }{j}^{-\left(n+1\right)+{\alpha }_2/2-1/4}\omega \left(j^{-1/2}\right)=O\left(m^{-n+{\alpha }_2/2-1/4}\right)\omega \left(m^{-1/2}\right).$

Заметим, что последнее соотношение имеет место, если $n>{\alpha }_2/2-1/4$. Ясно, что если функция $c\left(r,t\right)$ бесконечное число раз непрерывно дифференцируема, то всегда будет сходиться ряд

$\sum _{j=0}^{\infty }{j}^{-\left(n+1\right)+{\alpha }_2/2-1/4}\omega \left(j^{-1/2}\right).$

Таким образом, в этом случае будем учитывать, что полученная оценка справедлива при любом $n>{\alpha }_2/2-1/4$. Итак,

$\sum _{i=0}^{m-1}\left|c_i^n-\sum _{j=0}^{m-1}{c}_{ij}{L}_j\left(0;{\alpha }_2\right)\right|=O\left(m^{-n+{\alpha }_2/2+3/4}\right)\omega \left(m^{-1/2}\right).$

Используя формулу дифференцирования для многочленов Лагерра и известную асимптотическую формулу

$\text{Γ}\left({\alpha }_2+m+1\right)=m^{{\alpha }_2+1}\left(1+O\left(m^{-1}\right)\right)\left(m-1\right)!,$

получим:

$\frac{\parallel {\phi }_m^{\text{'}{\alpha }_2,{\gamma }_2}\left(t\right){\parallel }_{L_2,t^{{\alpha }_2}{e}^{-t}}}{{\phi }_m^{{\alpha }_2,{\gamma }_2}\left(0\right)}=\frac{{\gamma }_2\parallel {\phi }_m^{{\alpha }_2,{\gamma }_2}\left(t\right){\parallel }_{L_2,t^{{\alpha }_2}{e}^{-t}}}{2{\phi }_m^{{\alpha }_2,{\gamma }_2}\left(0\right)}+{\gamma }_2\frac{мим\parallel {\phi }_{m-1}^{{\alpha }_2+\mathrm{1,}{\gamma }_2}\left(t\right){\parallel }_{L_2,t^{{\alpha }_{2+1}}{e}^{-t}}}{{\phi }_m^{{\alpha }_2,{\gamma }_2}\left(0\right)}=$

$=\frac{{\gamma }_2\sqrt{\text{Γ}\left({\alpha }_2+m+1\right)}}{2\sqrt{m!} {\phi }_m^{{\alpha }_2,{\gamma }_2}\left(0\right)}+{\gamma }_2\frac{\sqrt{\text{Γ}\left({\alpha }_2+m+1\right)}}{\sqrt{\left(m-1\right)!} {\phi }_m^{{\alpha }_2,{\gamma }_2}\left(0\right)}=O\left(1\right)m^{\left(1-{\alpha }_2\right)/2}.$

Таким образом, имеем

$\frac{\sqrt{\text{Γ}\left({\alpha }_1+3\right)}\tau }{{\gamma }_1} \frac{\parallel {\phi }_m^{\text{'}{\alpha }_2,{\gamma }_2}\left(t\right){\parallel }_{L_2,t^{{\alpha }_2}{e}^{-t}}}{{\phi }_m^{{\alpha }_2,{\gamma }_2}\left(0\right)} \parallel \left({\tilde{c}}_0\left(r\right)\right){\parallel }_{L_2,r^{{\alpha }_1}{e}^{-r}}=O\left(m^{-n+5/4}\right)\omega \left(m^{-1/2}\right).$

Последняя оценка имеет место, если $n>5/4$.

Полагая $h=m^{-1/2}$ и собирая все оценки вместе, получаем оценку для невязки:

$\parallel R{\parallel }_{L_2}<\left(\frac{\sqrt{\text{Γ}\left({\alpha }_1+3\right)}}{{\gamma }_1}\left(3D\tau +2\tau +1\right)+2D\tau \right){[1-{(1-m^{-1/2})}^m]}^{-k}{m}^{-n}\text{Ω}\left(m^{-1/2}\right)+$

$+O\left(m^{-n+5/4}\right)\omega \left(m^{-1/2}\right), {\alpha }_2\ge 0.$

 Здесь $\text{Ω}\left(m^{-1/2}\right)$ — мажоранта обобщенных модулей непрерывности (см. [1]). Теорема доказана.

5. Результаты расчетов

Расчеты проведены с помощью математического пакета Matlab (MathWorks, Inc.) для параметров, характерных для нитрида галлия. При энергии электронного пучка $E_0=8$ кэВ для $m=6$ базисных функций Лагерра число обусловленности матрицы системы (14), порожденное спектральной нормой, не превысило 60, что позволяет использовать для расчетов персональные ЭВМ. Отметим, что использование модифицированных функций Эрмита для решения поставленной задачи оказалось непригодным, поскольку в силу их рекуррентных соотношений, матрица системы (14) оказалось вырожденной.

Подынтегральная функция в экспоненциальном множителе в формуле для расчета концентрации экситонов была разложена в ряд Тейлора. Использовано 7 членов разложения, что оказалось достаточным для проведения практических расчетов.

 

Рис.1. Концентрация экситонов для случая переменного времени жизни.

 

Затраты машинного времени на расчет концентрации экситонов составили приблизительно 5 с, что говорит о вычислительной эффективности предложенного метода.

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

Метод Галеркина в комбинации с разложением в ряд Тейлора подынтегральной функции в экспоненциальном множителе, содержащей переменный коэффициент, позволяет проводить расчеты концентрации генерированных неравновесных носителей заряда с точностью достаточной для проведения практических расчетов. Также метод позволяет проводить расчеты, не используя разложения переменного коэффициента в модели диффузии в ряд Фурье, что значительно упрощает применение проекционного метода.

About the authors

E. V. Seregina

Московский государственный технический университет им. Н. Э. Баумана (национальный исследовательский университет), Калужский филиал

Author for correspondence.
Email: evfs@yandex.ru
Russian Federation

M. A. Stepovich

Калужский государственный университет им. К. Э. Циолковского

Email: m.stepovich@rambler.ru
Russian Federation

M. N. Filippov

Институт общей и неорганической химии им. Н. С. Курнакова РАН

Email: fil@igic.ras.ru
Russian Federation, Москва

References

Supplementary files