CC BY
43
УДК 66.02:519.711.3
Ф. Г. Ахмадиев, Р. Ф. Гиззятов
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ И ВЫЧИСЛИТЕЛЬНЫЙ ЭКСПЕРИМЕНТ РАЗДЕЛЕНИЯ
ЗЕРНИСТЫХ МАТЕРИАЛОВ НА СИТОВЫХ КЛАССИФИКАТОРАХ
Ключевые слова: математическая модель, вероятность просеивания, зернистый материал, кинетика разделения,
оптимизация.
Построена математическая модель процесса разделения зернистых материалов на ситовых многокаскадных классификаторах с использованием теории пуассоновских процессов и проведена идентификация модели. Сформулирована и решена задача оптимального аппаратурного оформления в многокритериальной постановке.
Keywords: mathematical model, probability of screening, granular material, kinetics of separation, optimization
A mathematical model of the process of separation of granular materials on sieve multistage classifiers was constructed using the theory of Poisson processes and model identification was carried out. The problem of optimal hardware design in a multicriteria formulation is formulated and solved.
Введение
Разделение зернистых материалов по размерам является составным элементом во многих технологических процессах, в частности, в химической технологии, в производстве строительных материалов, при обогащении полезных ископаемых, в сельском хозяйстве и других областях промышленности. Для этого используются различные классификаторы, в том числе многокаскадные аппараты. Принцип работы подобного оборудования основан на различной вероятности прохождения через отверстия сит зерен различных форм и геометрических размеров. Расчет ситового классификатора сводится к определению его конструктивных и режимных параметров в зависимости от формы и размеров разделяемого материала.
Изучению разделения зернистых материалов на ситах посвящено большое количество работ, например, [1-9]. Исследованию кинетики процесса разделения посвящен целый ряд работ, например, [1-7]. При исследовании кинетики в работах [1, 3] предлагается использовать теорию марковских процессов, что позволяет определить функцию распределения частиц по размерам в аппарате. Процесс сепарации при этом представляется как диффузионный процесс с поглощающим экраном, роль которого выполняет поверхность решета и описывается уравнением Колмогорова-Фоккера-Планка (Kolmogorov - Fokker - Planck) (К.-Ф.-П.). Теория марковских процессов успешно используется и для описания других различных процессов обработки гетерогенных сред, например, [10-12]. Моделированию вероятности просеивания методом Монте-Карло посвящена работа [6]. Вероятность просеивания частицы в ячейку является основным параметром с помощью которого достигается эффект разделения, ее изучению посвящены работы [1,2]. Процессы внутрислоевого перераспределения (самосортирования) частиц изучены в работах [3,4]. Процессам тонкослойного разделения посвящена работа [1]. В работе [5] для описания кинетики процесса классификации на
вибрационных грохотах использована
стохастическая теория, при этом процесс блуждания частиц считается несимметричным и кроме действия случайных сил учитываются также действия детерминированных гравитационных сил, влияющих на вынужденный перенос частиц.
Таким образом, вероятностный характер процесса разделения и свойств частиц в смеси в целом свидетельствует о стохастическом характере процесса. Поэтому среди основных подходов при моделировании этих процессов, особое место занимают стохастические методы, в частности, теория марковских процессов [1, 11-13].
Для построения математической модели разделения на данном оборудовании можно также использовать теорию случайного блуждания частиц на плоскости [14,15]. При таком подходе моделирования процесса разделения на многокаскадных классификаторах с большим числом ярусов удается преодолеть ряд трудностей, которые возникают при использовании для этих целей уравнений (К.-Ф.-П.).
Целью работы является исследование и построение математической модели процесса разделения зернистых материалов для расчета многокаскадного классификатора на основе теории случайного блуждания частиц, проведение вычислительного эксперимента, идентификации и оптимизация процесса с целью установления оптимальных значений технологических и конструктивных параметров оборудования с использованием построенных моделей.
Математическое моделирование
Для описания положения частицы на ситах классификатора рассматривается процесс случайного блуждания частиц на плоскости, который описывается системой разностных уравнений:
р (/,У) = Р • Р (/ -1,У) + (1 - Р) • Р (/, У -1), (1)
где Р (/, У) - вероятность нахождения частицы в ] -ой ячейке / - го яруса, р - вероятность
просеивания в ячейку, i = 1, т ; у = 1, п . Здесь параметр т - количество сит, П - количество ячеек определенной формы по длине сита в направлении движения разделяемого материала. Необходимо определить распределение вероятностей Р (т, у) того, что частица окажется в ] - ой ячейке нижнего яруса. В дискретные моменты времени частица имеет вероятность р просеяться в нижний ярус, т.е.
сделать один шаг вниз и вероятность 1 - р перейти в следующую ячейку того же яруса в направлении движения частицы.
Эту систему уравнений следует решать при следующих граничных условиях, через которые учитывается физическая суть работы классификатора:
Р(0,1) = 1, Р(0, у) = 0 для
— - (2)
у = 2, п; Р(/,0) = 0, / = 1,т.
Для дальнейших расчетов по формулам (1), (2) необходимо определить вероятность р. В первом приближении вероятность просеивания можно представить в виде произведения вероятностей двух независимых событий [2]: р = рг • р^, где рг -геометрическая вероятность, зависящая от размеров и формы отверстия сита и частиц разделяемого материала и р^ - вероятность, зависящая от скорости движения частицы, которая определяется по формуле:
Я = 2 - (Ф(2) + ф(20)) , (3)
где 2 = (V - Ук) / а, 20 = Ук / а, V, а - параметры нормального закона, определяются по экспериментальным данным в процессе идентификации построенных моделей,
V - амплитуда скорости частицы относительно сита, Ф(х) - стандартная нормальная функция
распределения Ф( х) = ехр(-2 /2)с(?. Эти
вопросы рассмотрены в работе [16].
Вибрационное перемещение сыпучих материалов рассмотрено в работах, например, [1, 2,
16-21]. Следует отметить, что в зависимости от физико-механических свойств разделяемого материала, гранулометрического состава и т.д. движение дисперсных материалов по вибрирующей поверхности можно моделировать как в приближении одиночной материальной точки [1,2,
17-19], так и на основе методов механики гетерогенных сред [20-23]. Моделирование движения частиц как материальной точки можно проводить при тонкослойном движении разделяемых материалов. В сепарирующих машинах целесообразно использовать регулярный режим двухстороннего движения с мгновенными остановками без подбрасывания [2,4,17,18], что позволяет наиболее полно использовать время пребывания разделяемого материала на поверхности сита.
Таким образом, при известном значении р решение задачи (1)-(2) имеет вид:
(4)
где У-2
(т + J - 2)! = т т +1 т + у - 2
(т - 1)!(у -1)! 1 2
У -1
Распределение Р (т, у) представляет собой распределение Паскаля [15], являющееся распределением дискретной целочисленной случайной величины, равной количеству «неудач» в последовательности испытаний Бернулли с вероятностью «успеха» в одном испытании р, проводимой до т - го успеха включительно. Ее числовые характеристики могут быть выражены через производящую функцию вероятностей Р(т, у), которую можно представить в виде:
р (5) = 2 Р(т, У )5у = ртз£ ст + у• 5)у = ;
у=1 у=0 - • 5 у
- = 1-р, 15 < 1 . (5)
Тогда математическое ожидание положения частицы (порядковый номер ячейки ) на поверхности т -го сита определяется на основе (4) и (5) по формуле М[ ^ ] = Р'(1) = (р + т(1 - р)) / р, а дисперсия находится по зависимости D [Лт ] = Р ''(1) + Р '(1) - (Р '(1))2 = т(1 - р) / р2, где
знак « ' » означает дифференцирование по параметру s.
Построим функцию распределения Я/ (х) для
величин ДХ/, / = 1, т, где АХ, - расстояние, которое частица проходит при движении по поверхности I -го сита. Поскольку сита работают совместно, эти вероятности для / > 2 будут условными. С учетом нестационарности и ординарности потока событий, вероятность отсутствия события в интервале Ах в первом приближении можно оценить по
выражению Р0 (Дх) = 1 - С хДх , вероятность одного события С хДх, где С - параметр процесса, который зависит от вероятности просеивания в ячейку р . При этом математическое ожидание числа событий на участке (0; х) равно:
х)Ьх = Сх2/2 ,
(6)
Р (т, у) = ст+у - 2 рт (1 - р)у -
где /л(х) - мгновенная плотность (интенсивность) просеивания частиц в точке х . Тогда разностное уравнение для вероятности отсутствия просеивания Р0(х) имеет вид: Р0 (х + Дх) = (1 - С хДх)Р0 (х). Переходя к пределу при Дх ^ 0 , можно записать следующее
стохастическое дифференциальное уравнение: '
Р0 = -С хР0 (х). Его решение с учетом начального условия Р0 (0) = 1 (в начале сита просеивание отсутствует) имеет вид: Р0 (х) = ехр(-С х2 / 2).
Таким образом, функция распределения для первого сита имеет вид:
Р1(х) = 1 - ехр(-С х2/2).
Для остальных сит закон распределения величины ДХ(, / = 2, т , определяется из системы стохастических дифференциальных уравнений: Fy'(x) = -С xFy (х) + С xFy-1 (х), у = 1т. (7)
Прежде всего, нас интересует закон распределения величины ДХт на нижнем ярусе. С
учетом начальных условий F0 (0) = 1 и Fу (0) = 0
для всех у = 1, т, из решения системы (7) получим функцию и плотность распределения ДХт :
т-1
Fm (х) = 1 - ехр(-Сх2 /2) -X (Сх2 / 2)' / /!, (8)
/=0
/т (х) = Сх(Сх2 / 2)т-1 ехр(-Сх2 / 2) / (т -1)!.
По требованию на извлечение целевой фракции Г, зная закон распределения ДХт, можно определить координаты положения приемного бункера для целевой фракции Л и %2: Fm (х1) = (1 -г)/2 и Fm (х2) = (1 +Г*)/2 , (9) Принимая значение х2 = L, где L - длина сита, из второго уравнения (9) определяется значение параметра С, далее, из первого уравнения (9) определяется координата х1. В результате коэффициент извлечения и эффективность разделения определяются по формулам: х2
Г= |/т(х)с(х и Е = г-(1-£)• 100%, (10)
*1
где £ - чистота разделения, относительное содержание нецелевых продуктов в целевом бункере.
Расчеты по формулам (9) - (10) показывают, увеличение количества ярусов классификатора т способствует увеличению вероятности разделения (рис.1).
1 0,9 ■ Р
0,3 ■
07 ■
0,6 ■
0,5 ■ т
2 4 6 В 10 12
Рис. 1 - Зависимость вероятности разделения двух фракций от количества ярусов т при среднем квадратичном отклонении величины ДХт, равного 0,117 м для первой фракции и 0,185 м для второй
Для выбранных фракций центры распределений по длине сита числа частиц, прошедших сквозь сита, с увеличением количества ярусов расходятся (рис.2), что положительно влияет на процесс разделения.
55- М[ЛХ„],М
Ъ '
2
1.5 -
1
0.5
т
а 5 7 о 11 13
Рис. 2 - Зависимость математического ожидания распределения числа частиц, прошедших сквозь сита, от количества ярусов т при среднем квадратичном отклонении величины ДХт для мелкой (1) - 0,117 м; средней (2) - 0,185 м; крупной (3) - 0,269 м фракций
На рис.3 приводится распределение по длине сита числа частиц выбранных фракций, прошедших сквозь нижнее сито. Показано, на каком участке нижнего сита происходит просеивание частиц выбранных фракций в его ячейки.
4 Гт М,1/м
35 ■
3 ■ 1/\
25 ■
2
1,5 ■
0.5 ■ 3
0
0.4 о.е 0.8 1.2 1.4
Рис. 3 - Распределение числа частиц мелкой (1), средней (2) и крупной (3) фракций вдоль нижнего сита при т =12 и среднем квадратичном отклонении величины ДХт для мелкой (1) -0,117 м; средней (2) - 0,185 м; крупной (3) - 0,269 м фракций
Расчеты показывают, что для выбранных значений конструктивных параметров, а именно формы и размеров отверстий, длины и количества сит, зная требование на конечный продукт, можно подобрать режимные параметры работы аппарата, т.е. значения параметров А, ю, а, 3, при которых эффективность будет принимать наибольшее значение.
Идентификация построенных моделей и оптимизация
В процессе идентификации построенных моделей необходимо определить параметры скоростного режима V и а в зависимости (3). Для этого расчетные значения коэффициента извлечения сравнивают с экспериментальными значениями,
полученными при некоторых вполне определенных скоростных режимах работы аппарата.
Зная координаты приемного бункера, по экспериментальным данным для выбранных скоростных режимов определяются значения параметра С в распределении (8) и оценки числовых характеристик случайного процесса, математического ожидания и дисперсии:
х2 _ _ --
/ /„ (х; С,у )с(х = Щ , М/ = V(4т -1) / 2С/
О = 1 / 2Су. (11)
где нижний индекс / = 1,3 обозначает номер фракции ( 1 - мелкая, 2 - целевая, 3 - крупная), верхний индекс у = 1,п - номер опыта, Щ - доля
I -ой фракции в целевом бункере при скоростном режиме с порядковым номером у . Выбор скоростного режима (параметров А, ю, а, р), количества ярусов т , необходимо вести с учетом требований на извлечение и чистоту разделения:
Щ > Щ, Щ > 0 , Щ > 0 ,
(Щ + Щ) / (Щ + Щ + Щ) < &,
где параметр ¿>, - требование по чистоте разделения.
Оценки вероятности просеивания в ячейку ру определяются по выражению ру = т • 2а / М/. Эта расчетная формула получается в результате сравнения математического ожидания гамма распределения дт (х) с математическим ожиданием распределения /т (х).
На основании статистических оценок вероятности просеивания ру для выбранных скоростных режимов определяются наблюдаемые значения р - квантили скорости (безразмерное значение скорости, при котором скоростная вероятность равна р^)
иу = Ф-1(1 - ру / р8/), / = 1,3 ; у = 17п , где Ф-1 [*] - вычисление обратного значения стандартного нормального распределения с вероятностью ру. Расчетное значение р - квантили скорости можно найти по формуле: иу = Ф-1(Ф(2У) + Ф(20/) -1), / = 13 ; у = 1т, (12)
где 2у = (Ц/ - V/) / а , 20/ = V,/ / а , V/ - амплитуда относительной скорости частиц.
Неизвестные параметры зависимости (12) Ц(/, а, определяются в результате сравнения наблюдаемого и расчетного значений р - квантили скорости:
£(U;j -Uy)2 ^min , i = 1,3 .
(13)
y=1
Зная параметры скоростной вероятности V,, а по формуле р = рзр„, можно определить вероятность просеивания в ячейку. Далее по формулам (9) определяются координаты приемного
бункера для целевой фракции и по формулам (10) рассчитываются коэффициент извлечения и эффективность разделения.
Для установления оптимальных конструктивных и режимных параметров классификатора можно сформулировать задачу оптимизации в многокритериальной постановке. В качестве критериев рассматриваются производительность аппарата и эффективность разделения [16]: maxQ(A,rn,a,p,h,B) = pchBVcp,
maxE(A,a,a,p,D,L,ri, m) = i(1 -S) x 100% при условиях:
(14)
.,min ^ ^ ..max jmin ^ / / л o\ ^ a
x, < x, < x, , (A,a,a,P) <ф
y _ ,.y
min ..max
max k ,
где x;in, x;ax - наименьшее и наибольшее значение компонент вектора x = (A,m,a,p,D,L, h,|, m), фк -функциональные ограничения, связанные с выбранным скоростным режимом, pc - сыпучая плотность, B - ширина сита, L - длина сита, m -количество сит, A - амплитуда и а - частота колебаний, a - угол наклона, Р - угол вибрации, Vcp - средняя скорость перемещения зернистого
материала на сите, h - толщина слоя зернистого материала, подаваемого из загрузочного бункера, | - требование по извлечению целевой фракции.
Многокритериальная задача (14) решается поэтапно. На начальном этапе определяют оптимальные конструктивные параметры аппарата: количество ярусов m, длину и ширину сит L и B , форму и размер ячейки сита D. Длина сит L определяется с учетом коэффициента извлечения, количество ярусов m рассчитывается с учетом чистоты разделения. Далее, значения настроечных параметров, например, A,a могут быть выбраны перед эксплуатацией, а оптимальные значения управляющих параметров а,Р выбираются лицом, принимающим решение, из решения оптимизационной задачи (14), например, из некоторого множества оптимальных по Парето решений [24].
Для построения множества Парето составим линейную свертку нормированных целевых функций и найдем ее максимум для различных значений А и А :
Q (да)-Qmin E (Ра)-Emin max W = А '-^ + А ун '-^ при
Ра 1 Q - Q i ^ E - E i
^max ^min max min
(15)
А +А = 1 и о < а < 1,
где Qmax = max Q (A*,®,a*p, h*,B*) и Qmin = npinQ (A,aa*,ph* ,B') Emax = max E (A*,®,a* ,p,L*, m*,D*,|*) и E = maxE(A\®,a\p,L\m\D\|).
Значения а*, А*, Г, т*, О*, ^, В*,ц -
оптимальные значения соответствующих параметров, найденные на предыдущих этапах. В пространстве критериев строится образ
множества Парето, откуда могут быть выбраны оптимальные значения параметров , со*. Область Парето содержит, как правило, много элементов, в связи с этим возникает необходимость в дальнейшем сужении этой области. Для этого можно использовать различные процедуры. Рассмотрим процедуру выбора оптимального решения по принципу гарантированного результата:
f
max 1/Ц = max
min
Q(ß,©)-Q E(ß,©)-E
(15)
ß©e PGx (Q,E) -
Q E I
j /
множество Парето, Q, E -
характерные значения производительности и коэффициента эффективности разделения из области допустимых решений. Можно использовать и другие дополнительные процедуры [24] по усмотрению лица, принимающего решение.
Результаты и их обсуждение
На основе построенной в работе математической модели процесса разделения и поставленной задачи оптимизации были проведены численные расчеты (вычислительный эксперимент) с целью установления оптимальных параметров работы аппарата. Для проведения вычислительного эксперимента использовался комплекс программ, приведенный в работе [25]. Для разделения был использован зернистый материал на полимерной основе, частицы которого имеют цилиндрическую форму диапазоном от 0,2 мм до 1,4 мм и одинаковыми диаметрами 0,5 мм. Сыпучая плотность материала равна 1,15 г/см3, содержание целевой фракции размерами 0,5 - 0,9 мм в разделяемом материале составляет 60-75%.
На классификаторе были установлены решета длиной L = 1,6 м и шириной В = 1,0 м с размерами ячеек 2а х 2Ь = 3,5 х 3,5 мм с круглым отверстием диаметра D = 1,5 мм, количество ярусов классификатора т = 12 штук. При проведении численных расчетов были приняты: угол наклона а = -5о, амплитуда колебаний А = 5 мм, высота выходной щели загрузочного бункера h = 3 мм, требование на извлечение г« = 0,9.
В качестве критериев оптимизации выбраны: производительность оборудования (массовый расход) Q, кг/час и эффективность разделения Е , %. Параметрами оптимизации были выбраны: частота колебаний ю , сек -1 и угол вибрации 3, градус. В результате решения многокритериальной задачи оптимизации получено компромиссное решение: Е = 81,5%, Q = 1490 кг/час. Оптимальные значения параметров: ю = 54,6 сек-1 и 3 = 10,50. Средняя скорость составила V = 11,85 х 10-2 м/сек; амплитуда относительной скорости V = 0,36 м/сек. Координаты положения приемного бункера для целевой фракции равны: х1 = 0,98 м, х2 = 1,6 м, чистота разделения
(относительное содержание не целевых продуктов в целевом бункере) S = 10,4% .
Заключение
Теория случайных процессов (блуждания частиц), дополненная экспериментальными исследованиями для определения параметров модели, позволяет построить математические модели процесса разделения зернистых материалов по размерам на ситовых классификаторах, которые являются основой для оптимизации и управления процессом разделения. Для практической организации процесса классификации в зависимости от типа (формы) частиц, фракционного состава выбираются тип и размер отверстия ячейки, а длина аппарата и угол наклона сита, амплитуда колебаний в результате соответствующих расчетов. Оптимальные значения угла вибрации и частоты колебания определяются из решения соответствующей оптимизационной задачи. При найденных этих значениях параметров работы классификатора был проведен вычислительный эксперимент по изучению процесса разделения специфических полимерных зернистых материалов на многокаскадных классификаторах.
Литература
1. Ахмадиев Ф.Г., Гиззятов Р.Ф., Киямов Х.Г. Математическое моделирование процесса тонкослойного разделения зернистых материалов по размерам на ситовых классификаторах.//Теоретические основы химической технологии. 2013. Т.47. № 3. С.309-317.
2. Гортинский В.В., Демский А.Б., Борискин М.А. Процессы сепарирования на зерноперерабатывающих предприятиях. М.: Колос, 1973.
3. Непомнящий Е.А. Кинетика некоторых процессов переработки дисперсных материалов.//Теоретические основы химической технологии. 1973. Т.7. №5. С.754-763.
4. Блехман И.И., Хайнман В.Я. О теории вибрационного разделения сыпучих смесей.//Изв. АН СССР. Механика, 1965. №5. С. 22-30.
5. Огруцов В.А., Егоров А.Д. Метод расчета процесса классификации сыпучих материалов на виброгрохотах.//Ивановский химико-технологический институт, межвузовский сборник научных трудов, 1987, С. 93-96.
6. Beeckmans J., Jutan A., Can J. Monte carlo simulation of a probability screen.//Chem. Eng. - The Canadian Journal of Chemical Engineering. - 1989, V.67, Issue 2, p.329-336.
7. Мизонов В.Е., Огурцов В.А., Федосов С.В. Моделирование классификации полидисперсных материалов на виброгрохотах.//Строительные материалы, 2007, № 11.- С.26-28
8. Вайсберг Л.А. Проектирование и расчет вибрационных грохотов. - М.: Недра, 1986.
9. Акулова М.В., Алешина А.П., Огурцов А.В., Огурцов А.В. Моделирование процесса классификации сыпучих материалов на виброгрохотах с многоярусной компоновкой сит. //Вестник МГСУ. 2013. № 2. С. 80-87
10. Капранова А.Б., Лебедев А.Е., Бытев Д.О., Зайцев А.И. Стохастическое описание движения осветленной фракции суспензии. //Известия высших учебных заведений. Серия: Химия и химическая технология. 2004. Т. 47. № 6. С. 99-101
11. Ахмадиев Ф.Г., Назипов И.Т. Стохастическое моделирование кинетики процессов обработки гетерогенных систем.//Теорет. основы хим. технологии. 2013. Т.47. № 2. С.182.
12. Иванова С.А. Павский В.А. Стохастическое моделирование процесса пенообразования белковых растворов.//Теорет. основы хим. технологии. 2014. Т.48. № 6. С.701-708.
13. Тихонов В.И., Миронов М.А. Марковские процессы. М.: Сов. радио, 1977.
14. Баруча-Рид А.Т. Элементы теории марковских процессов и их приложения: пер. с англ./ под ред. А.Н. Ширяева. - М.:Наука, 1969. -512 с.
15. Феллер В. Введение в теорию вероятностей и ее приложения. В 2-х томах. Т 1: Пер. с англ.- М.: Мир, 1984.- 528 с.
16. Ахмадиев Ф. Г., Гиззятов Р.Ф., Назипов И.Т. Гидрогазодинамика и кинетика разделения дисперсных сред на ситовых классификаторах.// Инженерно-физический журнал. 2017. Т.90. №5. С.1135-1145.
17. Блехман И.И., Джанелидзе Г.Ю. Вибрационное перемещение. - М.:Наука,1964.
18. Нагаев Р.Ф. Периодические режимы вибрационного перемещения. - М.: Наука, 1964.
19. Гончаревич И.Ф., Фролов К.В. Теория вибрационной техники и технологии. М.: Наука, 1981
20. Генералов М.Б. Истечение сыпучих материалов из аппаратов. //Теоретические основы химической технологии. 1985. Т.19. № 1. С.53-58.
21. Генералов М.Б. Механика твердых дисперсных сред в процессах химической технологии. - Калуга: Издательство Бочкаревой, 2002.
22. Соу С. Гидродинамика многофазных систем: пер. с англ./под ред. М.Е. Дейча. - М.:Мир, 1971
23. Псевдоожижение. Под ред. И. Дэвидсона, Д. Харрисона. Пер. с англ. М.: Химия, 1974.
24. Ахмадиев Ф.Г. Некоторые задачи многокритериальной оптимизации технологических процессов.//Теоретические основы химической технологии. 2014. Т.48. № 5. С. 574-582.
25. Ахмадиев Ф.Г., Гиззятов Р.Ф., Назипов И.Т. Компьютерное моделирование и оптимизация кинетики процесса обработки гетерогенных сред. Свидетельство о государственной регистрации программы для ЭВМ № 2016661926, 2016.
© Ф. Г. Ахмадиев - доктор технических наук, профессор, заведующий кафедрой Прикладная математика КГАСУ, Akhmadiev@kgasu.ru; Р. Ф. Гиззятов - кандидат технических наук, доцент кафедры Прикладная математика КГАСУ, renatgrf@mail.ru.
© F. G. Akhmadiev - doctor of sciences in technology, professor, head of the department Applied Mathematics KSUAE, Akhmadiev@kgasu.ru; R. F. Gizzyatov - candidate of sciences in technology, associate professor of the department Applied Mathematics KSUAE, renatgrf@mail.ru.
