CC BY
66
+ 2 t
q+(i-1)c
= Tq + (s - 1)maxtq+(i-1)c,
где m = - , tq+(i-i)c = tq+(i-,)c +E, q =1,c, i =1,m ,
ui, u2, ..., us-i - целые числа.
Рассмотрим далее случай, когда s = k
k>1, и Tq >
tq ах. Вычисление общего времени
ТС (P,n,s, с) с помощью функционала задачи Беллмана-Джонсона приводит к формуле ТС (P,n,s,c) =
max
1.p J-1
2 ^q+(i-i)c + 2 ^q+(i-i)c + •
km
+ 2 *q+(i-i)c
= kTq +
-11 tq .
ЗдесЬ tq+(i-1)c+km = tq+(i-!)c , q = 1С , i = 1,m .
В случае, когда s = k .P] + r, k>1, 1 ^ r <
и Tq >
Си, вычисление общего времени с
помощью функционала задачи Беллмана-Джонсона приводит к третьей формуле для вычисления
общего времени выполнения одинаково распределенных конкурирующих процессов при ограниченном числе копий ПР.
В заключение отметим, что рассмотренное обобщение математической модели с одним структурированным ПР (конвейером) на случай ограниченного числа ПР позволяет установить взаимосвязи мультиконвейерной обработки с аналогичной обработкой при одном программном конвейере, получить аналитические оценки общего времени выполнения конкурирующих процессов при ограниченном параллелизме и провести математическое исследование эффективности и оптимальности мультиконвейерной организации конкурирующих процессов, вскрыть потенциальные возможности роста ускорения вычислений, выполнить сравнительный анализ различных режимов такой обработки.
Литература
1. Воеводин В.В., Воеводин Вл.В. Параллельные вычисления. СПб: БХВ-Петербург, 2002.
2. Топорков В.В. Модели распределенных вычислений. М.: Физматлит, 2004.
3. Иванников В.П., Коваленко Н.С., Метельский В.М. О минимальном времени реализации распределенных конкурирующих процессов в синхронных режимах // Программирование. 2000. № 5. С. 44-52.
4. Павлов П.А. Оптимальность структурирования программных ресурсов при конвейерной распределенной обработке // Программные продукты и системы. 2010. N° 3. С. 76-82.
5. Kapitonova Yu.V., Kovalenko N.S., Pavlov P.A. Optimality of systems of identically distributed competing processes // Cybernetics and Systems Analysis. NY: Springer. 2006, pp. 793-799.
УДК 519.688:004.942:544.431.2
МЕТОД ОЦЕНКИ КОНСТАНТ НЕЛИНЕЙНЫХ КИНЕТИЧЕСКИХ МОДЕЛЕЙ МНОГОСТАДИЙНЫХ ХИМИЧЕСКИХ РЕАКЦИЙ
(Работа выполнялась в рамках ФЦП «Научные и научно-педагогические кадры инновационной России»
на 2009-2013 гг., гос. контракт П1550)
Е.В. Писаренко, к.т.н.
(Российский химико-технологический университет им.. Д.И. Менделеева, г. Москва,
еир1эагепко@т-аИ.ш)
Рассмотрены основные этапы проведения кинетических исследований, построения высокоточных кинетических моделей, идентификации констант моделей, планирования стартовых и прецизионных экспериментов. Разработаны комплексы алгоритмов и программ идентификации констант кинетических моделей сложных химических реакций, планирования прецизионных экспериментов и построения оптимальных планов экспериментов по критерию Б-оптимальности.
Ключевые слова: кинетика, идентификация констант, планирование эксперимента.
Рост цен на энергоносители требует разработки новых технологий производств, которые позволяют использовать как целевые продуктовые потоки, так и все многообразие побочных продук-
тов [1]. Следовательно, необходимы кинетические модели, с помощью которых можно с высокой точностью рассчитать концентрацию всех этих продуктов. Для построения подобных моделей не-
обходимо учитывать протекание в химической системе большого количества элементарных химических реакций. Значение и роль кинетических исследований [2] в общей процедуре моделирования химического процесса исключительно велики, поскольку необоснованно выбранная неполная система гипотез о механизме реакции не может привести к установлению кинетических моделей, правильно отражающих структуру химических превращений и позволяющих с высокой степенью надежности прогнозировать статические и динамические свойства изучаемой химической системы.
Из теоремы эквивалентности [3] известно, что построение высокоточных кинетических моделей, то есть моделей, чрезвычайно точно прогнозирующих изменение в ходе реакции основных показателей химического процесса, невозможно без оценки по экспериментальным данным высокоточных кинетических констант. Последнее достигается планированием кинетического эксперимента с целью получения при каждом опыте максимальной информации о численных значениях кинетических констант.
Проблема построения кинетической модели сложной химической реакции и оценки кинетических констант модели решается по следующей схеме.
• Формирование возможного механизма протекания многомаршрутной химической реакции, оценка числа маршрутов, векторов стехиометри-ческих чисел маршрутов и возможных итоговых уравнений маршрутов.
• Формирование алгоритмов построения кинетических моделей анализируемой реакционной системы, не противоречащих теории реакционно-способности химических соединений и основным положениям химической кинетики.
• Построение кинетических моделей изучаемой химической реакции.
• Планирование стартовых кинетических экспериментов.
• Качественный анализ кинетических моделей протекания реакции и оценки свойств решений уравнений кинетических моделей.
• Оценка идентифицируемости кинетических констант моделей при условии протекания химической реакции в квазистационарных условиях по промежуточным веществам и в заданной области проведения кинетического эксперимента.
• Планирование прецизионных экспериментов с целью получения оценок кинетических констант моделей с точностью, обеспечивающей заданные априори прогнозирующие возможности кинетической модели.
• Проверка соответствия построенного прецизионного плана эксперимента критерию Б-оптимальности [4].
• Построение доверительных областей для векторов кинетических констант и откликов модели, подтверждающих необходимую для последующего моделирования режимов работы реактора высокую точность расчета по кинетической модели основных показателей химического процесса.
Для возможности использования персональных компьютеров при анализе химических реагирующих систем химические реагенты и химические реакции отображаются в векторы а.,
различной размерности, компоненты которых -целые или рациональные числа. Совокупность векторов а|, 1<1<1\. образует матрицу молекулярных видов ЛКхм, где М - число структурных видов; N - число реагентов. Совокупность векторов
, 1<]<0, образует матрицу стехиометрических
коэффициентов Бдх№ где О - число химических реакций в анализируемой реакционной системе. При этом, поскольку должны выполняться основные положения химической кинетики, а именно закон действующих масс, закон сохранения массы реакционной системы и условия ее электронейтральности, имеем
Йв=Вт\У, (1)
ВА=0, (2)
где Й®х1 - вектор скоростей изменения концентраций реагентов; \¥0х1 - вектор скоростей элементарных каталитических и некаталитических химических реакций.
При условиях квазистационарного протекания химической реакции часть реагентов (боденштей-новские вещества) имеют значение скоростей Й" » 0, Ь+1<к<1\. другая часть реагентов (небо-денштейновские вещества) имеют значения величин скоростей Ф 0, 1<т<Ь. В соответствии с вышеизложенным матрица Б представима в виде двух подматриц: Бн - подматрица стехиометриче-ских коэффициентов небоденштейновских веществ; Бб - подматрица стехиометрических коэффициентов боденштейновских веществ. Таким образом, имеем Б=[Бн | Бб].
Следовательно, процедура построения кинетической модели химической реакции и оценки ее кинетических констант выполняется в следующей последовательности.
1. Матрица Б представляется в виде клеточной матрицы Б=[Бн | Бб].
2. Определяется ранг матрицы Бб, и в соответствии с правилом Хориути рассчитывается количество независимых маршрутов реакции Р, где Р=О-гапк(Вб), и количество уравнений Ь химических инвариантов для боденштейновских веществ, где Ь!=(К-Ь)-гапк(Бб).
3. Определяются уравнения химических инва-
риантов для боденштейновских веществ. Имеем
cH(t) = си(0) + В^2(в;,) ' (c61(t)-c61(0)) , (3)
где Вб=[Вб11 В2б2] и B61 s - неособенная квадратная подматрица стехиометрической матрицы В5 ранга rank(B6i)=rank(B6)=s; сб1, сб2 - подвекторы вектора концентраций боденштейновских веществ; t - независимая переменная, время, длина реактора и т.п.
4. Рассчитываются Р векторов стехиометрии -ческих чисел маршрутов v(p), 1<р<Р, решением системы линейных однородных уравнений, соответствующей матричному уравнению BТ • v = 0. Расчет независимых векторов стехиометрических чисел v(p) неоднозначен в том смысле, что их число всегда Р, но компоненты векторов могут быть различны, ибо размерность матрицы ВТ есть (N-L)xQ, причем (N-L)<Q. Следовательно, ввиду того, что коэффициенты матрицы - рациональные числа, общее количество решений счетно. Среди них выбираем только линейно независимые решения, тривиальное решение не рассматривается, так как ему соответствует отсутствие химических реакций в анализируемой системе. Каждое независимое решение определяет совокупность стехиометрических чисел стехио-метрического вектора, которому соответствует определенное выражение итоговой реакции по маршруту.
5. Рассчитываются компоненты транспонированной матрицы стехиометрических коэффициентов итоговых реакций по маршрутам ВI = ВТ • v, причем в зависимости от выбранных компонентов P векторов стехиометрических чисел (матрица v) получаются различные итоговые уравнения реакций по маршрутам.
6. Определяется ранг матрицы итоговых уравнений Bu по маршрутам. Матрица Bu представляется в виде клеточной матрицы Bu=[Bu11 Bu2], где rank(Bu)=rank(Bui). По неособенной квадратной подматрице Bu1 устанавливаются ключевые вещества, а по подматрице B^ независимые. Определяются уравнения химических инвариантов для неключевых (независимых) веществ. Последние имеют вид
*L(t) = *L(0) + Bht2(в1,) ' (cK(t)-cK(0)), (4)
где cH3(t), cK(t) - векторы концентраций независимых и ключевых веществ соответственно; снз(0), ск(0) - векторы начальных значений концентраций независимых и ключевых веществ соответственно.
Если в итоговых уравнениях маршрутов не все химические реакции независимые, число ненулевых строк подматрицы Bu1 будет меньше числа
строк матрицы Би. При этом вид уравнений (4) сохраняется.
7. По матричному соотношению рассчитываются компоненты вектора скоростей по маршрутам г(р), 1<р<Р, и компоненты вектора концентраций боденштейновских веществ как функции концентраций небоденштейновских веществ и кинетических констант:
(5) мар-
где г = (г"1,..., г(Р)) - вектор скоростей по шрутам.
Следовательно, для различных векторов стехиометрических чисел у(р) будут различные матричные выражения (5), то есть разные алгебраические уравнения для векторов скоростей по маршрутам г. Так как кинетические модели построены для условий квазистационарности по промежуточным веществам и при этом как объем эксперимента, так и область экспериментирования ограничены, не все индивидуальные кинетические константы допускают оценку, часть из них или даже все константы оцениваются в виде определенных линейных или нелинейных соотношений.
8. Проводится стартовый эксперимент, и методом наименьших квадратов или максимального правдоподобия вычисляются предварительные оценки кинетических констант.
9. Определяется количество констант (или их комбинаций), допускающих оценку. Для решения этой задачи вычисляется детерминант информационной матрицы многооткликовой модели detM(s) или его миноры в точках плана £ в выбранной области экспериментирования О, то есть
X,
11 „1
12 „1
1n „1
1<l<m,
5с( '(k,Xj)
5kf
1<i<N1, 1<f<n,
k=k„
Z"1 = {osq}, 1<s<m, 1<q<m,
(6)
(7)
(8)
(9)
где т - число откликов модели, по значениям которых оцениваются кинетические константы; N1 -число кинетических экспериментов; п - число кинетических констант в модели; Е - дисперсионно-ковариационная матрица наблюдений; - условия проведения 1-го эксперимента.
Очевидно, что detM(s)^0 только при условии, что ^т>п, ибо матрица М(е) неособенная только при условии, что число ее независимых строк и/или столбцов равно п.
k=11=1
Xif =
Вычисляется ненулевой минор Mt наибольшего порядка информационной матрицы M(e). При этом максимальное число параметров (или их комбинаций), допускающих оценку, определяется порядком подматрицы Mt.
10. Строится плотность распределения detMi(s), и проверяется гипотеза о том, что определитель выбранной подматрицы Mt статистически значимо отличается от нуля. Если это утверждение верно, по вычисленному минору устанавливаются константы модели (или их комбинации), допускающие оценку. Определим их как независимые константы.
11. Переоцениваются константы кинетической модели, уже как комплексы констант, методом наименьших квадратов. Решение дифференциальных уравнений модели осуществляется полунеявным методом Рунге-Кутта III порядка. Минимизация функционала, характеризующего соответствие модели эксперименту, осуществляется методами Маркуардта, Ньютона, Хука-Дживса, случайного поиска.
Оценки независимых констант будут определяться как результатами эксперимента, так и задаваемыми численными значениями остальных зависимых констант. В качестве зависимых констант целесообразно выбирать такие константы, объем априорной информации о которых максимален. Очевидно, что численные значения, присваиваемые зависимым константам, не могут повлиять на предсказательную силу математической модели. Мощность множества допускающих оценку констант определяется, в свою очередь, свойствами плана эксперимента. Она максимальна для оптимального плана, так как число констант со статистически значимой параметрической чувствительностью максимально. При этом для нелинейных моделей ранг матрицы M(e) и результаты анализа на параметрическую чувствительность констант зависят от численных значений констант, и, следовательно, процедура исследования максимально возможного числа оцениваемых констант должна проводиться на всех этапах последовательного прецизионного уточнения параметров.
В кинетической модели вследствие использования условий квазистационарности по боден-штейновским веществам возможно возникновение ситуаций, когда определенные совокупности констант не могут быть определены индивидуально и оцениваются только в виде некоторых функциональных соотношений (комбинаций линейных и нелинейных). Такова структура кинетической модели, и никакое планирование эксперимента не может исправить это положение. В модели существуют и константы (или комбинации констант) с низкой параметрической чувствительностью, которые при данном плане эксперимента не могут быть оценены с приемлемой точностью. Напри-
мер, константы характеризуют химическую реакцию, которая в исследуемой области практически не протекает, но с заметной скоростью идет в другой. Поэтому в качестве критерия планирования эксперимента требуется использовать и комбинированный критерий, который предусматривает максимизацию как ранга информационной матрицы плана, так и его детерминанта.
12. Точность прогноза по модели определяется дисперсионно-ковариационной матрицей d(x,k,e)
откликов cT(x,k) = (Cj(x,k),.„, cN(x, к))т вида
d(x,k,e) = FT (х,к) • (М(е))"1 • F(x,k), (10) где х - условия проведения эксперимента в области экспериментирования G; крх1 - вектор кинетических констант; е - план эксперимента; сДх, к) - значение концентрации i-ro реагента в точке х области экспериментирования; NjM(s) = =M(e):
F(x,k):
Sc Sc Sc
Skj Sk " "" Skj
9с, Sc Sc
ök2 Sk " "" Sk2
Sc
Sk"
Sc
Sc
(11)
-p 5kp 5kp y
После каждого эксперимента необходимо осуществлять проверку точности проведенных расчетов элементов информационной матрицы плана М(е) и элементов дисперсионно-ковариационной матрицы прогноза по модели d(x, k, s) по матричному уравнению
¿p„sp[2r1d(ï11,k,e)] = ii>
(12)
где хц, l<u<N - точка спектра плана; р„, l<u<N -веса плана; Е - дисперсионно-ковариационная матрица наблюдений в точке хц ;
Sp [г1 d(xu, к, б)] - след матрицы 2 1 d(xu, к,е)
Последовательное непрерывное планирование, от опыта к опыту максимизирующее detM(e), приводит к построению D-оптимального плана эксперимента, если выполняется условие
max Sp 1 d(x, k, s) J = n. (13)
Это условие по теореме эквивалентности для многооткликовых моделей эквивалентно условию M(s) = maxdetM(s), (14)
s
где s* - D-оптимальный план эксперимента.
После построения D-оптимального плана эксперимента переходят к последовательному планированию эксперимента по максимизации ранга
информационной матрицы М(е). По завершении этого этапа, при достижении увеличения ранга М(е), снова переходят к планированию прецизионного эксперимента по этапам, изложенным выше.
Если при проведении уточняющего прецизионного эксперимента сходимость к оптимальному плану эксперимента невысока и прогнозирующая способность модели не удовлетворяет априори заданным требованиям, то после выбора нового вектора констант, допускающих оценку, необходимо вернуться к этапам 1-12. При новых предположениях по результатам ранее проведенных экспериментов оцениваются отдельные константы модели или их комбинации. Определяется параметрическая чувствительность констант, и формируется новая информационная матрица М(е). Устанавливаются ее ранг и детерминант. Далее вновь проводится последовательное планирование уточняющих экспериментов с целью получения кинетической модели с заданными свойствами.
Указанная выше процедура построения кинетических моделей прошла апробацию при построении кинетических моделей многомаршрутных реакций паровой конверсии метанола и конверсии смеси метанола и диметилового эфира в олефиновые углеводороды [5] и показала себя достаточно продуктивной. Экспериментальные исследования кинетики реакции конверсии метанола и диметилового эфира в низшие олефины проводились в проточно-циркуляционном реакторе с внутренним контуром циркуляции. Реактор позволяет проводить экспериментальные исследования в диапазоне температур 350-600 оС и давлений от 0,1 МПа до 6,0 МПа. В реактор загружался мелкодисперсный высококремнистый цеолитный катализатор объемом 5-50 мл. Предложен 25-ста-дийный механизм реакции конверсии метанола и диметилового эфира в олефиновые углеводороды. Построена кинетическая модель анализируемой реакционной системы. Поскольку в моделях не все кинетические константы можно оценить по экспериментальным данным, при постановке эксперимента и обработке его результатов требовалось установить число индивидуальных констант или их комплексов, допускающих оценку.
Точность оцененных констант определяется по детерминанту дисперсионно-ковариационной матрицы оценок констант. Чем меньше детерминант, тем точнее оценки кинетических констант модели. Сущность нового метода заключается не только в определении числа констант, допускающих оценку, но и в планировании эксперимента таким образом, чтобы оценить максимальное количество кинетических констант. В противном случае кинетическая модель не будет прогнозировать результаты эксперимента с высокой точностью. По результатам эксперимента в выбранной области экспериментирования были установлены
численные значения кинетических констант и их комбинаций модели реакции получения олефи-новых углеводородов, допускающих оценку. С использованием построенной плотности распределения DetM(e) выделена критическая область значений, в которой DetM(s) статистически незначимо отличаются от нуля. Для оценки параметров кинетической модели пригодны лишь те значения DetM(e), которые не принадлежат критической области. В рассматриваемом случае из кинетической модели выделены 18 независимых комбинаций констант, которые определяет невырожденная информационная матрица плана. Следовательно, выбранные 18 комбинаций кинетических констант могут быть оценены по результатам кинетического эксперимента.
Таким образом, вследствие того, что химические реакции в реакторах протекают обычно в квазистационарном режиме, система кинетических уравнений распадается на две подсистемы -систему дифференциальных уравнений и систему алгебраических уравнений. Каждое из последних уравнений равно нулю при любом значении независимой переменной (время, длина реактора и т.п.), поэтому существуют функциональные соотношения между кинетическими константами. Причем число функциональных соотношений меньше числа кинетических констант. Следовательно, в области квазистационарного протекания реакции определить численные значения всех кинетических констант принципиально невозможно. Возникает проблема, какие кинетические константы или их комбинации можно определить по экспериментальным данным и при каких условиях проведения эксперимента константы (или их комбинации) будут оценены с наибольшей точностью для заданного объема эксперимента. На решение подобных проблем направлен изложенный выше метод оценки констант нелинейных кинетических моделей многостадийных многомаршрутных химических реакций и разработанный пакет прикладных программ «KINCON IDENTIFICATION» на Compaq Visual Fortran 6.6 и С++. Он дает возможность установить область проведения химической реакции, в которой допускается оценка максимального числа кинетических констант (или их комбинаций); для заданного объема эксперимента оценить кинетические константы с максимально возможной точностью (для D-оптимального плана эксперимента); построить кинетическую модель с высокой прогнозирующей возможностью. Предложенные методы и программы прошли проверку при решении многих практических задач.
Литература
1. Chul-Jin Lee [et al.]. Optimal gas-to-liquid selection from natural gas under uncertain price scenarios // Ind. Eng. Chem. Res., 2009. Vol. 48, pp. 794-800.
2. Wojciechowsky B.W. and Rice R.M. Experimental me-
thods in kinetic studies // Elsevier Science B.V. Amsterdam, Netherlands, 2003.
3. Федоров В.В. Теория оптимального эксперимента. М.: Наука, 1971. 312 с.
4. Berding C., Kleider W. and Gossl R. D-optimal planning & experimentation the tool for experimental design in research and
development // Computational Statistics & Data Analysis, 1996. Vol. 21. Is. 6, pp. 705-710.
5. Писаренко Е.В., Писаренко В.Н. Кинетика реакции синтеза олефинов на основе метанола и диметилового эфира // Теоретические основы химической технологии, 2008. Т. 42. № 6. С. 623-632.
УДК 004.93
АВТОМАТИЧЕСКАЯ РЕКОНСТРУКЦИЯ ТРЕХМЕРНЫХ ОБЪЕКТОВ ПО ТЕХНИЧЕСКОМУ ЧЕРТЕЖУ
(Работа выполнена при поддержке гранта МК-8838.2010.9)
А.А. Захаров, к.т.н.
(Владимирский государственный университет, aa-zaharov@ya.ru)
Рассмотрены вопросы автоматической реконструкции трехмерных объектов по техническому чертежу. Разработан алгоритм реконструкции трехмерных моделей на основе граничного представления. Алгоритм использует шаблоны для распознавания типовых конструктивных элементов модели. Приведены практические результаты работы алгоритма.
Ключевые слова: трехмерная реконструкция, технический чертеж, CAD-системы, граничное представление.
Модернизация во всех областях промышленности невозможна без широкого внедрения CALS-технологий, основанных на электронной модели изделия, одной из составляющих которой является компьютерная трехмерная модель объекта, используемая для разработки управляющих программ, инженерного анализа изделий, визуализации и т.д. Современные CAD-системы имеют широкий набор средств для создания трехмерных моделей с чистого листа: булевы операции, операции объектно-ориентированного моделирования, 2,5D-операции, операции модификации вершин, ребер и граней. Практически все они позволяют осуществлять генерацию чертежей по трехмерной модели, однако получение такой модели по чертежу вызывает у проектировщиков затруднения из-за отсутствия математического и программного обеспечения.
В архивах предприятий накоплено большое количество технических чертежей как в бумажном, так и в электронном виде. Кроме того, разработка многих объектов часто начинается не с трехмерной модели, а с чертежа. Двухмерные чертежи часто сложны для понимания, неудобны для корректировки и не могут служить основой для дальнейших разработок с использованием компьютерных систем. Создание системы, выполняющей автоматическую реконструкцию трехмерных моделей по техническому чертежу, позволит во много раз сократить время проектирования объектов различного назначения. Этим обусловлена необходимость разработки математического и программного обеспечения системы автоматической реконструкции трехмерных моделей объектов по техническому чертежу.
Все способы восстановления трехмерных моделей по чертежам можно условно разделить на две группы: CSG- и B-rep-подходы. CSG-ориенти-рованный (constructive solid geometry - конструктивная блочная геометрия) подход использует стратегию восстановления сверху вниз. Он основан на том, что каждый трехмерный объект может быть построен из определенного двухмерного примитива иерархическим способом [1]. На чертеже отыскиваются шаблоны, которые будут служить базой и использоваться для перевода в трехмерную модель. После этого сконструированные примитивы собираются в результирующую трехмерную модель при помощи булевых операций. Недостатком CSG-ориентированного подхода является то, что с его помощью трудно распознать базовые примитивы на сложных чертежах. К тому же при использовании базовых примитивов затруднительно представить поверхности сложной формы.
B-rep-ориентированный (bounding representation - граничное представление) подход использует восходящую технологию. B-rep-ориентиро-ванные алгоритмы включают обычно следующие шаги: генерирование возможных трехмерных вершин из чертежа; синтез ребер по полученным координатам вершин; конструирование граней из ребер, лежащих в одной плоскости; формирование трехмерного объекта из граней [2]. Граничное представление обеспечивает высокую точность моделирования геометрии объекта, чего нельзя достичь при использовании CSG-подхода. Однако при использовании B-rep-представления требуется больший объем памяти для хранения и обработки данных. Кроме того, создаваемая модель логиче-
