CC BY
51
150
А. С. Шевченко
УДК 519.6
ПРОГРАММНЫЙ КОМПЛЕКС ЧИСЛЕННОГО РЕШЕНИЯ ЗАДАЧИ О ПЕРЕНОСЕ ЗАРЯДА В ТРАНЗИСТОРЕ MESFET
Шевченко Алеся Сергеевна,
кандидат физ.-мат. наук, email: [email protected]
Рубцовский институт (филиал) Алтайского государственного университета, 658225, г.Рубцовск, пр-т Ленина 200-Б
Аннотация.
Актуальность работы. В данной работе мы рассматриваем гидродинамическую МЕР модель. Эта модель представляет собой квазилинейную систему уравнений, записанных в форме законов сохранения. Эти законы сохранения получены из системы моментных соотношений для уравнения переноса Больцмана.
Цель работы: Построение и реализация эффективного вычислительного алгоритма для поиска решений задачи о переносе заряда в транзисторе MESFET.
Методы исследования: Предложенный алгоритм основан на методе установления, применении сглаживающих регуляризующих операторов и идей схем без насыщения.
Результаты: Написан программный комплекс па языке Object Pascal в среде Delphi 6, реализующий вычислительный алгоритм.
Ключевые слова: гидродинамическая модель, 2D кремниевый транзистор MESFET, уравнение Пуассона, регуляризация, метод установления, алгоритмы без насыщения, интерполяционный полином, сплайн-функция, матричная прогонка
Математическое моделирование физических процессов в полупроводниковых устройствах имеет огромное значение для технических приложений и в последнее время превратилось в быстро развивающуюся область прикладной математики. Современный уровень развития микроэлектронных технологий позволяет создавать полупроводниковые приборы столь малых размеров, что точность, обеспечиваемая упрощенными аналитическими моделями, становится неприемлемой для анализа и проектирования таких приборов. Необходимы математические модели, более полно описывающие физические процессы, происходящие в полупроводнике.
В настоящей работе рассматривается одна математическая модель, описывающая с высокой степенью достоверности процесс переноса заряда в полупроводниковых устройствах субмикронного размера. При попытке поиска численных решений этой модели мы столкнулись с рядом существенных сложностей, обусловленных нелинейностью уравнений модели и наличием в ней малых параметров. С помощью традиционной теории разностных схем эту задачу решить не удалось.
Рассматриваемая модель переноса заряда в двумерном случае в безразмерном виде имеет следующий вид [1,2]:
RT + divJ =0, ^2
J +V
— RE |=7?Q+cnJ+ cl2I,
(RE)t +divl =(J,Q)+c/?o-,
I.+V
= — REQ+c2 ,J + c22l,
A(p = (3(R-p).
Здесь и - вектор скорости электронов,
(1)
(2)
q-
поток энергии, J—Ru - электрический ток, I=Rq , R - электронная плотность, Е- энергия 2
электронов, (7 — — Е — 1, ф - электрический по-
3
тенциал, Q = V^>, р — р{х,у) - плотность легирования, Р>0 - некоторая постоянная [3]. Коэффициенты c,cu,cn,c2i,c22 являются гладкими функциями от Е [2,3].
Заметим, что для нахождения приближённых решений математической модели (1), (2) в работах [2], [4] - [5] предложены различные вычислительные алгоритмы.
Граничные условия для системы (1) и уравнения (2), соответствующие задаче для транзистора MESFET прямоугольной формы (схематическое представление рассматриваемого полупроводни-
Вестник Кузбасского государственного технического университета. 2015. № 6 151
кового прибора в безразмерных переменных изображено на рис. 1.) таковы [2]:
R =
Ng на Г0, 1
(3)
сг = б на Г0,
01
(р = \ G наГ0,
[В
(\,VR) = (\,VP) = (\,V<p) = 0na rz.J
Здесь х,у - независимые безразмерные переменные, Q = (0,1)х^0,— 1 - единичный вектор
внешней нормали к границе <9Q ,
г0 =
(х,у): у = \, 0<х<
3 о
1 2 5
-<х<-, - < х < 1 3 3 6
- часть границы дО., Г/ = \ Г0 .
Функция плотности легирования имеет вид: (\,(х,у)е&+,
Р(х,У) = -
^,(л,у)е0\0+,
(1 n fs f
-,1 X
U 3 J U J V
G =
Ge
mCl
\_ J_ 4 ’ 3
B = -*L- N -
m’Cl' 5 N*
некото-
"0 ,n ^0 рые постоянные, размерные параметры е, L, N\ С0, т приведены в табл. 1
В стационарном случае математическая модель (1), (2) сводится к системе эллиптических квазилинейных уравнений (к системе трех уравнений Пуассона):
Д<х =D(CT)(V<r, V«,Q,cr, R) =
= a, |Vaj2 +—(V/?,a2Vcr + a4Q)+ W R
+a3 ( V <x, Q) + a5 |Q|2 + beer, AR=a(R)(VcryR,Q,cr,R,p) =
= blR\Va\2+b2(Va,VR) + +biR(Va,Q) + b4{VR,Q) +
+bsR |Q|2 + —@—(R - p)R + ncRcr,
1 + <J
Acp = F w(R,p) = /](R-p), (6)
а компоненты векторов u, q находятся из следу-
(5)
1/3
1/4
source t'./s gate * 2/3 Г 1 5/6drain
(> = 1 | ; /i-i 1/4 1
1 P = * i i
1 0 1/6 1/3 1/2 1 1 2/3 5/6 1 x
Рис. 1. Схематическое представление двумерного транзистора MESFET Fig.l. Schematic representation of 2D MESFET transistor
Таблица 1. Значение физических параметров
Table 1. Value of physical properties
Величина Описание Значение
e Заряд электрона 1.6хКГ19с
L Характерный линейный размер 6 х 10-7 т
N+ Плотность легирования в области Q+ 3-1023 /т3
C0 Скорость звука 1.19170 1 05т / ^
m Эффективная масса электрона 2.9152-10-31 kg
152
А. С. Шевченко
ющих соотношении:
u = -F0(£)V(tJ,
q = G(£)|-Q + ^V« + G0(£)Vc7j.
Здесь |Vof =crj +crj ит.д.,
a\ - ~a 'F(E)F0(E) + b'G(E)G0(E), a2 =-1 + (1 + сг)|б 'G(E) - a >(£)},
аъ=а F(E)-b 'G(E)-bF(E)F0(E), a4 = -b(\ +ct)F(E), as = bF(E),
6, = -m 'F(E)F0(E) + n G(E)G0(E), b2 = (1 +ct){« 'G(E)-m >(£)},
1
(7)
63 =-
64 =
(1 + <t) 1
1 + СГ
2+m F(E)-n G(E)-nF(E)F0(E) + bF(E)F0(E), b5 = nF(E),
a = a(E) = -
2 c
21
2 c
22
5 1 + cr
m = m(is) =
-cn, b = b(E) = ^—c12>
5 l + cr
C11 c12-^
F(£) = -
1 + cr
-22 ^ ^C12
' > n - n(E) - ■
1 + cr
5
det
-Ec,
G(E) = —
Fq(E) = 1--
12
-21 “ ^ ^Cll
det — Ec]
'22 ^ jCC12
’ G0(£) = l-
11
11
Сг\~Еъ
, . > da 3 da
det — cuc22 —c2iC\2’ a ит.д.
Jcr 2 dE
Коэффициенты c,cn,c12,c21,c22 и их производные рассчитываются согласно формулам, приведённым в [2,3]. С этой целью на языках Delphi 6
[6] были написаны процедуры, осуществляющие все необходимые вычисления.
Для поиска стационарных решений задачи (4) - (6), (3) используем вычислительный алгоритм подробно описанный в [7] - [9].
Вычислительный алгоритм, предложенная в
[7] - [9], основан на применении нестационарных регуляризаций, метода прямых [10], сплайн-интерполяции [11], метода прогонки и метода установления [12]. Достоинством данного алгоритма является то, что он хорошо согласуется с идеей, так называемых, алгоритмов без насыщений, т.е. в принципе автоматически учитывает априорную гладкость решения исходной задачи.
В место уравнения Пуассона
АхуФ = фхх+ф}У= /X У% X У) е Q рассматриваются его нестационарные регуляризации:
1. параболическая регуляризация
ut=Au-f(x,y), (8)
2. регуляризация Соболева
ut-Aut=Au-f{x,y), (9)
3. гиперболическая регуляризация
utt + Kut = Au-f (х, у). (Ю)
Здесь u=u(t,x,y) - новая неизвестная
функция, К > 1 - постоянная, t > 0 - временная переменная.
При этом производную в уравнениях (8) -
(10) заменяем не разностным отношением, а используем для аппроксимации ее интерполяционный многочлен [12]
sin JC.
р(х,и)=-!-Хн)'"
N “7 cos х - cos х(. •cos(Nx)-u.(t,y), 0<х<л с узлами в нулях многочлена Чебышева: 2/-1 ;
X = X, = ■
2 N
u(t,xj,y) = uj(t,y).
N- количество узлов интерполяции В результате получаем краевую задачу ОДУ второго порядка. Для поиска приближенных решений, полученных систем ОДУ, используем ин-
терполяционного кубического сплайна класса С*2
[11].
В результате приходим к трехточечной схеме:
Г hi К [ hi Л
■и„ +
+ ХХ'
■и,
k+1
(ii)
=^-|+4Л+ЛД.* = 1.*-1-
с краевыми условиями
U q — U \,U к — LU к_\- (12)
Здесь U - вектор значений неизвестных в узлах интерполяции функции, IN - единичная матрица порядка N, hv - ширина сетки сплайн интерполяции; L - диагональная матрица порядка N (диагональные элементы равны либо 0, либо 1). Конкретные выражения для матрицы Вг и вектора Р записаны в [7] - [9].
Систему алгебраических уравнений (11) - (12) решаем методом матричной прогонки. Для того чтобы воспользоваться методом матричной прогонки, запишем систему (11) в виде
Вестник Кузбасского государственного технического университета. 2015. № 6 153
AUк_\ -CUk + BUk+l = F, к = \,К-\,(Щ
а граничные условия (12) перепишем так: ~C0U0 + B0U, = -F0, (14)
к-1 — к ~ ~Fk » (1
hi
где А = В = IN — Вг,
6
С =
hi
v 5 J
Fk -—^r{Pk-\+^Pk+Pk+\}’
Г -C = R = I Л - I F = F = ft _
1 n ’ — ^ ’ 1 о 1 к u
нулевой вектор размера N .
В соответствии с методом прогонки будем искать решение (13) - (15) в виде
ик=аш6^+рк,„ * = 0,К-1. (16)
Подставляя (16) в (13), находим рекуррентные формулы для прогоночных матричных коэффициентов ССк и векторных коэффициентов (5к,
к = \,К-\:
«*+| =(C-^at)_l5,
А+|=(с-^)-'(м+г4).
(17)
При k = 1 значения прогоночных коэффициентов Ц =Cq'B0 = IN , Д =Cq1F0 = 0 определяются из граничного условия (14). Т.о., в цикле прямой прогонки по формулам (17) находим коэффициенты GCk, рк ,к = 1, К — 1. Затем используя краевое условие (15) и формулу (16), вычисляем решение при к = К :
UК ={СК — Акак) (Ак£$к + FK).
В итоге, зная Uк и прогоночные коэффициенты, находим в цикле обратной прогонки ( к = К —1,0) по формуле (16) решение задачи на текущем временном слое.
Замечание. В [13] приведены критерии устойчивости описанного метода прогонки для схемы (13) с граничными условиями (14), (15). Они имеют вид
det(C-'S0)<l, АЛ(С~К'АК)<\,
det(C“'s) + det(C“'^) < 1
Важно заметить, что при расчетах прикладных задач мы следили за значениями приведенных выше определителей и можем утверждать, что данные критерии выполняются.
На каждом временном слое в процессе установления необходимо последовательно решать три краевые задачи для регуляризованных уравнений Пуассона (4) - (6).
Для осуществления расчётов нужно опреде-
лить начальные значения переменных (Е , R , (р ) и выбрать вид нестационарной регуляризации.
Используя идею метода установления, будем переходить с я-го временного слоя на (я+1)-й до тех пор, пока решение не установится. Схема алгоритма нахождения приближенного решения (4) - (6) на каждом временном слое следующая:
1. Сначала находим величину
F ^(R11, р), где Rn = R(nA,x,y^ значение R с предыдущего слоя. Затем находим решение (р = (рп+х = (р[{п +1)A,jc,задачи (4) и его производные на слеДУЮ1Дем
слое.
2. Рассчитываем величину
F <<г) (V сг", Vtf", Q"+l, а", R") • Далее решаем задачу для (5). Находим
о = сг"+] = сг((я + 1)Д,х,^),
3. Вычисляем
F <я)(?",/?) ■
Далее решаем задачу для (6). Находим R = R”*' =R((n + l)A,x,y),
{R'"')x,(R-')r,(R-')yy.
Выполнив 1-3 шаги, алгоритм переходит на следующий временной слой. Данный цикл работает до тех пор, пока решение не установится (т.е. пока не выполнится условие
+
/и-i
a jk ~ a jk
). После того как решение установилось, вычисляются U , q , Е , J , I Далее полученные дан-
ные сохраняем в файлах.
Написан программный комплекс на языке Object Pascal в среде Delphi 6, реализующий данный вычислительный алгоритм для поиска стационарных решений задачи о переносе заряда в 2D кремниевом транзисторе MESFET (см. Рис. 2).
В программе предусмотрено следующее:
1. Задание:
а) параметров, характеризующих вычислительный алгоритм: N - число узлов интерполяции по X, К - количество узлов сетки по у , А
- шаг разностной сетки по времени, £ - точность расчета;
б) параметров, характеризующих исходную задачу: G - напряжение затвора, Ng - плот-
ность легирования на затворе, В - напряжение смещения, 8 - безразмерная плотность легирова-
154
А. С. Шевченко
Численное решение задачи о переносе заряда в 2D кремниевом транзисторе MESFET
■ а Т^ПВЛД
Параметры исходной гадачи
Параметры вычислительного алгорнтмд
cll(E) | с 12(E) | с21(Е) { с22(Е) | ПОТЕШЩАЛ ПРОИЗВОДИВ ОТ ПОТЕНЦ_
30
N=
ЛНо.01
И~{ - число уыо* интерпелляции по х - шаг разностной сетки по т
К* j*0 rrl-к
Точность EPS= ПГ^Ы Начальные данные
У Но I
Е= [Ц I
N7K 0 I1 2 3 4 [5 *
1 0,0001 HJOOO I-0,2 2839139 -0,68569659 -1,14455530 -1,60595834 -2,0
2 0,00000000 0,22308850 -0.66981522 -1,11817969 -1,56923649 -2,0
3 0,00000000 -0,21193470 0,63638014 -1,06254622 -1,49155568 -1Д
4 0,00000000 0,19433827 0,58357367 0,97448209 -1,36816785
5 0,00000000 0,17044864 0,51181699 0,85459782 -1,19973128 >1,5-v
R= 1
количество углов сетки по v
Регулирииаии С 1. Параболическая <• '2 .Соболевская; З.Гиберболпческая
S Утвердить начальные данные РАСЧЕТ ПРЕРВАТЬ РАСЧЕТ
Сохранить данные] Считать данные |
Рис. 2. Окно расчета численного решения гидродинамической модели Fig. 2. The window calculation of the numerical solution of the hydrodynamic model
ния в области fl \ Г2+ .
в) начальных данных для искомых функций
Е, R, <р.
2. Выбор регуляризации (параболическая, Соболевская или гиперболическая).
3. Хранение физических параметров, необходимых для расчетов (см. Табл.1), а также парамет-
ров, необходимых для вычисления коэффициентов С, Cj |, С| 2, 1, С22 •
Управление расчетом осуществляется с помощью кнопок «РАСЧЕТ» «ПРЕРВАТЬ РАСЧЕТ».
Все полученные результаты отражаются на вкладках. А также есть вкладка «График зависимости энергии от t» (наибольшее значение энер-
| R | Е График зависимости энергии от t | |х 4 I ►
График зависимости энергии от t
£.\S -18-1ft Zu 19 18 17 16 15 14 13 12 11 10 9 8 7 6 5 4 3 2
- -
IO 3“ Й 12' 2* 1 ft
!
•a* IU rT 8- ft
О 4 2 0,011
— : _ - - -
JOOO 2,580000 4,760000 6,940000 9,121 1000 11,430000 13,850000 16,270000 18,700000 t
Рис. 3. Вкладка «График зависимости энергии от (» Fig. 3. The tab "Graph of energy from t»
Вестник Кузбасского государственного технического университета. 2015. № 6 155
гии на /7-ом временном слое), который позволяет нам следить за поведением энергии (см. Рис.З).
Кнопки «Сохранить данные», «Считать данные» позволяют осуществлять сохранение и загрузку результатов расчета (стационарные решения для электрического потенциала, плотности электронов, энергии электронов, электрического
тока и потока энергии по компонентам х и у). Данные сохраняются в файлах специального формата, что позволяет использовать полученные результаты в качестве начальных условий для дальнейших расчетов, а также строить стационарные решения с помощью математического пакета Mathematica [14].
СПИСОК ЛИТЕРАТУРЫ
1. Anile А.М. Non parabolic band transport in semiconductors: closure of the moment equations/ A.M. Anile and V. Romano // Cont. Mech. Thermodyn, 1999. - V.l 1. - P.307-325.
2. Romano V. 2D simulation of a silicon MESFET with a non-parabolic hydrodynamical model based on the maximum entropy principle // J. Comp. Phys, 2002. - V.l76. - P.70-92.
3. Blokhin A.M. / A.M. Blokhin, R.S. Bushmanov, A.S. Rudometova, V. Romano Linear asymptotic stability of the equilibrium state for the 2D МЕР hydrodynamical model of charge transport in semiconductors // Nonlinear Analysis, 2006. - V. 65. - P. 1018-1038.
4. Romano V. 2D Numerical Simulation of the МЕР Energy-Transport Model with a Finite Difference Scheme // J. Comp. Phys, 2007. - Vol. 221. - P. 439^68.
5. Blokhin A. M./ A. M. Blokhin, S. A. Boyarsky, В. V. Semisalov/ On an approach to the construction of difference schemes for the moment equations of charge transport in semiconductors // Le Matematiche, 2009. -Vol. LXIV. Fasc. I. - P. 77-91.
6. Белов В.В. Программирование в Delphi: процедурное, объектно-ориентированное, визуальное. Учебное пособие для вузов. / В.В. Белов, В.И. Чистяков. - М.: Горячая линия-Телеком, 2014. - 240 с.
7. Блохин А.М. Об одном варианте метода прямых для уравнения Пуассона / А.М. Блохин, А.С. Ибрагимова (А.С. Шевченко), Н.Ю. Красников // Вычислительные технологии, 2007. - Т.12, - №2. - С. 33-42.
8. Блохин А.М. Конструирование вычислительного алгоритма для системы моментных уравнений, описывающих перенос заряда в полупроводниках / А.М. Блохин, А.С. Ибрагимова (А.С. Шевченко), Б.В. Семисалов // Математическое моделирование, 2009. - Т21, - №4. - С. 15-34.
9. Шевченко А. С. Алгоритм поиска приближенных решений уравнения Пуассона // Молодой ученый, 2015.-№3.-С. 18-23.
10. Березин И.С. Методы вычислений / И.С. Березин, Н.П. Жидков. т.Н. - М.: Физматгиз, 1962. -620с.
11. Завьялов Ю.С. Методы сплайн-функций / Ю.С. Завьялов, Б.И. Квасов, В.Л. Мирошниченко. - М.: Наука, 1980.-355с.
12. Бабенко К.И. Основы численного анализа. - Москва-Ижевск: НИЦ "Регулярная и хаотическая динамика", 2002. - 848с.
13. Самарский А. А. Теория разностных схем: Учебное пособие для вузов. 3-е изд., испр. М.: Наука, -1989. -616с.
14. Давыдов Е.Г. Решение математических задач с помощью программных пакетов Scientific Work-Place, Scientific Notebook, Mathcad, Mathematica и Matlab: Учебное пособие. - M.: Книжный дом "ЛИБ-РОКОМ", 2012.-240 с.
Поступило в редакцию 20.09.2015
УДК 519.6
PROGRAM COMPLEX OF NUMERICAL SOLUTION OF THE PROBLEM OF CHARGE TRANSPORT IN A MESFET TRANSISTOR
Shevchenko Alesya S.,
C. Sc. (Phys.-math.), email: [email protected]
Rubtsovsk Institute (branch) Altai state university, 200-B Lenina prospect, Rubtsovsk, 658225, Russian Federation
Abstract.
The urgency of the discussed issue. In this paper we consider a hydrodynamical model proposed recently. This model represents a quasilinear system of equations written in the form of conservation laws. The conservation laws are obtained from a system of moment relations for the Boltzmann transport equation.
156
А. С. Шевченко
The main aim of the study: The aim of this work is the construction and the realization of effective numerical algorithm for solving of the problem of charge transport in transistor MESFET.
The methods used in the study: The proposed algorithm is based on the stabilization method, the application of regularized smoothing operators and ideas of schemes without saturation.
The results: This computational algorithm was implemented as software complex developed in Object Pascal in the Delphi 6 environment.
Keywords: hydrodynamic model, 2D silicon MESFET transistor, the Poisson equation, regularization, stabilization method, the algorithm without saturation, interpolation polynomial, spline-functions, sweep method
REFERENCES
1. Anile A.M. Non parabolic band transport in semiconductors: closure of the moment equations/ A.M. Anile and V. Romano // Cont. Mech. Thermodyn, 1999. - V.l 1. - P.307-325.
2. Romano V. 2D simulation of a silicon MESFET with a non-parabolic hydrodynamical model based on the maximum entropy principle // J. Comp. Phys, 2002. - V.l76. - P.70-92.
3. Blokhin A.M. / A.M. Blokhin, R.S. Bushmanov, A.S. Rudometova, V. Romano Linear asymptotic stability of the equilibrium state for the 2D МЕР hydrodynamical model of charge transport in semiconductors // Nonlinear Analysis, 2006. - V. 65. - P. 1018-1038.
4. Romano V. 2D Numerical Simulation of the МЕР Energy-Transport Model with a Finite Difference Scheme // J. Comp. Phys, 2007. - Vol. 221. - P. 439-468.
5. Blokhin A. M./ A. M. Blokhin, S. A. Boyarsky, В. V. Semisalov/ On an approach to the construction of difference schemes for the moment equations of charge transport in semiconductors // Le Matematiche, 2009. -Vol. LXIV. Fasc. I. - P. 77-91.
6. Belov V.V., Chistyakov V.l. Programmirovanie v Delphi: protsedumoe, ob"ektno-orientirovannoe, vi-zual'noe. Uchebnoe posobie dlya vuzov. [Programming in Delphi: procedural, object-oriented, visual. Textbook for high schools]. M.: Hot line-Telecom. 2014. 240 P.
7. Blokhin A.M., Ibragimova A.S. (Shevchenko A. S.), Krasnikov N.Yu. Ob odnom variante metoda pryamykh dlya uravneniya Puassona [On a variant of method of lines for the Poisson equation]. Vychislitel'nye tekhnologii [Computational Technologies]. 2007. Vol. 12. №2. P. 33-42.
8. Blokhin A.M., Ibragimova A.S. (Shevchenko A. S.), Semisalov B.V. Konstruirovanie vychislitel'nogo al-goritma dlya sistemy momentnykh uravneniy, opisyvayushchikh perenos zaryada v poluprovodnikakh [Design of a numerical algorithm for the system of moment equations of the charge transport in semiconductors]. Ma-tematicheskoe modelirovanie [Mathematical Modelling]. 2009. T21. №4. P. 15-34.
9.Shevchenko A. S. Algoritm poiska priblizhennykh resheniy uravneniya Puassona [The search algorithm of approximate solutions of the Poisson equation]. Molodoy uchenyy [Young scientist]. 2015. №3. P. 18-23.
10. Berezin I.S., Zhidkov N.P. Metody vychisleniy [Methods of calculations] M.: Fizmatgiz, 1962. Vol. 2. 620P.
11. Zav'yalov Yu.S, Kvasov B.I., Miroshnichenko V.L. Metody splayn-funktsiy [Methods of spline functions]. M.: Nauka, 1980. -355P.
12. Babenko K.I. Osnovy chislennogo analiza [Foundations of numerical analysis] Moscow-Izhevsk, Regular and chaotic dynamics, 2002.848P.
13.Samarskiy A. A. Teoriya raznostnykh skhem: Uchebnoe posobie dlya vuzov. 3-e izd., ispr.[The theory of difference schemes: School-book for institute of higher education. 3rd edition, revised]. M.: Nauka, 1989. 616P.
14.Davydov E.G. Reshenie matematicheskikh zadach s pomoshch'yu programmnykh paketov Scientific Workplace, Scientific Notebook, Mathcad, Mathematica i Matlab: Uchebnoe posobie [The solution of mathematical problems using software packages, Scientific Workplace, and Scientific Notebook, Mathcad, Mathematica, and Matlab]. M.: Publishing house "LIBROKOM". 2012. 240P.
Received 20 September2015
