Разработка математической модели пневматической системы термостабилизации

ISSN 0868-5886

НАУЧНОЕ ПРИБОРОСТРОЕНИЕ, 2006, том 16, № 1, c. 94-106

ОРИГИНАЛЬНЫЕ СТАТЬИ

УДК 536.74: 62 - 522.7

© А. В. Чернышев, А. А. Крутиков, К. Е. Демихов,

А. В. Полынков, О. В. Белова, И. А. Смирнов

РАЗРАБОТКА МАТЕМАТИЧЕСКОЙ МОДЕЛИ ПНЕВМАТИЧЕСКОЙ СИСТЕМЫ ТЕРМОСТАБИЛИЗАЦИИ

В статье представлен метод расчета пневмоэлектрического устройства для проведения полимеразной цепной реакции деления молекулы ДНК — пневматического амплификатора. Особенностями данной системы являются высокие скорости изменения температуры рабочего тела — воздуха при одновременном обеспечении точности изменения температуры биологического образца в пластиковой микропробирке. Целью исследований является моделирование тепловых процессов для создания оптимальных условий для проведения реакции деления молекулы ДНК. Представлена математическая модель расчета, основанная на методе сосредоточенных параметров. Описан экспериментальный стенд. Приведены результаты расчетно-теоретических исследований и их сравнение с результатами экспериментальных исследований. Показано, что предложенный метод позволяет оценить масштабы переходных процессов, протекающих в пневматическом амплификаторе. Сделан вывод о необходимости разработки математической модели, основанной на методах описания систем с распределенными параметрами.

ВВЕДЕНИЕ

В современной технике широко применяются системы термостабилизации, предназначенные для поддержания определенной температуры объекта с заданной точностью. В настоящее время появился новый класс технических задач, в которых термостатирующие устройства не только поддерживают заданную температуру объекта но изменяют ее по заданной программе с требуемой скоростью и точностью. Данные устройства применяются в приборостроении, электровакуумной технике, медицине и т. д. По роду используемого рабочего тела, предназначенного для подвода и отвода тепла от объекта теплового воздействия, термостабилизирующие устройства могут быть: твердотельными; жидкостными и газовыми (пневматическими). В данной статье рассматривается вопрос разработки математической модели пневматической системы на примере термостатирую-щего устройства медицинского назначения. К таким устройствам относятся так называемые ам-плификаторы ДНК [1-3], в которых температура исследуемого образца должна изменяться циклически с высокой точностью и по программе, задаваемой пользователем. В основе большинства выпускаемых сегодня амплификаторв ДНК заложена твердотельная схема [4], но в связи с тем, что к этим устройствам сегодня предъявляется очень сложное требование — повысить скорость изменения температуры исследуемого образа до 10 °С/с и выше, разработчики данного оборудования обратили внимание на пневматическую систему термостабилизации [5].

Возможны самые различные варианты реализации конструктивной схемы с точки зрения применяемого газа. Но наиболее высокие потребительские свойства будут иметь конструкции, в которых в качестве рабочего тела используется воздух из окружающей среды и без предварительной подготовки. В устройствах данного типа побудитель расхода обеспечивает перемещение рабочего тела через преобразователь энергии, в котором происходит подвод тепла, и затем оно вступает в тепловой контакт с объектом теплового воздействия для его нагрева. Для охлаждения объекта происходит замещение горячего воздуха холодным.

Разработка математической модели пневматической системы термостабилизации проводилась на примере системы для изменения и термостабилизации температуры жидкости, размещенной в стандартной, тонкостенной пластиковой микропробирке.

Расчетная термодинамическая схема данной установки приведена на рис. 1. Система состоит из набора условных проточных пневматических емкостей, разделенных между собой и отделенных от подводящих и отводящих магистралей условными дросселями. На входе в первую емкость установлен побудитель расхода (вентилятор), во второй емкости расположен преобразователь энергии — проволочный электронагреватель, а в третьей емкости размещается микропробирка с исследуемым образцом.

С позиции классической термодинамики рассматриваемое устройство представляет собой открытую термодинамическую систему [6, 7], границы которой проницаемы как для потока вещества,

Ра, Та Нагреватель Вентилятор Ра, 1а

Рис. 1. Расчетная схема пневматической системы термостабилизации.

1 — первая (входная) емкость с вентилятором, с параметрами воздуха: давлением р1, объемом Уь температурой Т1; 2 — вторая емкость (преобразователь) с электронагревателем, с параметрами р2, У2, Т2 ; 3 — третья емкость с микропробиркой и образцом в ней (жидкость) с температурой Тж, с параметрами р3, У3, Т3; Оо_3 — потоки воздуха; ра, Та — параметры атмосферы; /а1, /12, /23, /3л — дроссели

так и для теплового взаимодействия. В представленной системе протекают процессы различной физической природы: термо- и газодинамические, теплопереноса и преобразования электрической энергии в тепло. Математическое описание данных процессов может быть основано на использовании системы уравнений в частных производных (при учете распределенных параметров) или обыкновенных дифференциальных уравнений, что возможно только при описании системы с сосредоточенными параметрами. Модель с сосредоточенными параметрами применима в случае квазистатических процессов. Известно, что изменение состояния можно рассматривать как квазистатическое, если скорости протекания рабочих процессов тела малы по сравнению со скоростью распространения звука в газе при данных условиях [6]. При решении данной задачи изменения состояний можно с допустимой погрешностью считать квазистатическими, поскольку скорость потока рабочего тела и темпы изменения параметров системы как минимум на два порядка меньше указанной скорости звука. Система с распределенными параметрами по своей сущности ближе к реальным физически процессам и наиболее полно и точно описывает переходные процессы, протекающие в конкретных физических устройствах. Однако трудоемкость написания математических моделей и программ расчета данных систем, как правило, существенно превосходят трудоемкость для систем с сосредоточенными параметрами. Учитывая данный фактор, а также приемлемую точность расчетов с применением модели с сосредоточен-

ными параметрами, рабочие процессы, протекающие в пневматической системе термостабилизации, описаны системой обыкновенных дифференциальных уравнений.

ОСНОВНЫЕ ЗАВИСИМОСТИ, ИСПОЛЬЗУЕМЫЕ ДЛЯ РАСЧЕТА ПЕРЕХОДНЫХ ПРОЦЕССОВ В ПНЕВМАТИЧЕСКОЙ СИСТЕМЕ ТЕРМОСТАБИЛИЗАЦИИ

Зависимости для расчета изменения давления и температуры газа в открытой термодинамической системе при условии наличия теплообмена между рабочим телом и окружающей средой (рассчитываемая модель графически представлена рис. 2) описываются уравнениями [6, 7]:

dp _ а? _

_ У

аг

ал

X

Т , рУ' ¿у 1т,ф

Р^+У ? - -

где р — давление газа, Па; Т — температура газа, К; У — объем газа, м3; к — показатель адиабаты; Я — газовая постоянная, Дж/(кг-К); Опр —

массовый расход приходящего газа, кг/с; G}х — массовый расход уходящего газа, кг/с; Q — количество теплоты, подводимое к газу, Дж.

Зависимость для определения расхода газа в условном дросселе, соединяющем соседние емкости системы, (модель на рис. 3) для докритиче-ского режима течения [7]:

G = f 1Р1

2k 1

k -1 RT1

Pl

к+1

( Р. Л * Л

при вр < ^ < 1

Р1

Рис. 2. Рассчитываемая модель открытой термодинамической системы при условии наличия теплообмена между рабочим телом и окружающей средой

Рис. 3. Рассчитываемая модель соединения соседних емкостей системы через условный дроссель

в емкости, из которой вытекает газ, К.

Зависимость для определения расхода газа в условном дросселе, соединяющем соседние емкости системы, для сверхкритического режима течения [7]:

G = f 1Р1

к+1

к ' 2 ^к-1

,к +1J

при 0< ^<

Р1

где f — площадь поперечного сечения в условном дросселе, м2; | — коэффициент расхода; р1 — давление газа в емкости, из которой вытекает газ, Па; р2 — давление газа в емкости, в которую поступает газ, Па; Т1 — температура газа

В пневматической системе термостабилизации вентилятор обеспечивает прокачку рабочего тела через электрический нагреватель, состоящий из нихромовой нити, подключенной к источнику тока.

Рассмотрим процесс передачи теплоты от электрического нагревателя к рабочему телу. Тепловой баланс электрического нагревателя (без учета теплообмена излучением):

Q = Q + Q

эл сп ¿^ш

где Qэл — количество теплоты, выделяемое в электрическом нагревателе, Дж; Qcn — количество теплоты, расходуемое на нагрев спирали, Дж; QKoнe — количество теплоты, затрачиваемое на нагрев рабочего тела через конвективный теплообмен, Дж.

Количество теплоты, выделяемое в электрическом нагревателе:

2 и2

Qэл = и • I • х = I2 • R • х = —•х = Р • х,

R

где Р — мощность, прикладываемая к спирали, Вт; и — напряжение, прикладываемое к спирали, В; I — ток, протекаемый по спирали, А; X — время, с.

Количество теплоты, расходуемое на нагрев спирали:

Q = с • т АТ ,

сп сп сп сп '

где есп — теплоемкость материала спирали, Дж/(кг • К); тсп — масса спарили, кг; АТсп — разность между конечной и начальной температурой спирали, К.

Количество теплоты, затрачиваемое на нагрев рабочего тела через конвективный теплообмен [10]:

Q = а^Б •(Т -Т2)х

¿^конв \ сп 2 /

где а — коэффициент теплоотдачи, Вт/(м2 • К); Б — площадь поверхности теплообмена, м2; Тсп — температура спирали, К; Т2 — температура газа, К.

Дифференцирование приведенных выше зависимостей по времени дает

авэл _ о+аок _ Р dt dt dt ' dt

^ _ ссп • тсп • %; 5-(Тсп - Т2).

м си dt ' & Дифференциальное уравнение изменения тем-

пературы спирали:

аТп _ Р(Тп - Т2) dt с т

Т3 -ТСИ1 _о/(«1 • 51)_о• ;

• передачу теплоты от внешней стенки пробирки с температурой Тст1 к внутренней стенке пробирки с температурой Тст2 путем теплопроводности

Тст1 Тст2 б • ;

Яа1 + Яа2

Т - Т

13 1ж

+ -

(^•п • • 11 2 •п •!• 11 а2 • п • • 12

где Яа1 — термическое сопротивление теплоотдачи от рабочего тела к внешней стенке пробирки, Вт/К; Яя — термическое сопротивление теплопроводности от внешней стенки пробирки к внутренней стенке пробирки, Вт/К; Яа2 — термическое сопротивление теплоотдачи от внутренней

Рассмотрим процесс передачи теплоты от рабочего тела к жидкости, находящейся в пробирке. Пробирка рассматривается как тело цилиндрической формы, направление оси которого совпадает с направлением потока. Процесс переноса теплоты от рабочего газа к жидкости через разделяющую их стенку пробирки (рис. 4) включает в себя [11]:

• передачу теплоты от рабочего тела с температурой Т3 к внешней стенке пробирки с температурой Тст1 путем конвективного теплообмена

• передачу теплоты от внутренней стенки пробирки с температурой Тст2 к жидкости с температурой Тж путем конвективного теплообмена

Тст 2 - Тж _ б/

При стационарном режиме тепловой поток б во всех трех процессах одинаков и равен

Т - Т

б _ 3 ж _

Рис. 4. К расчету переноса теплоты от рабочего газа к жидкости в пробирке через разделяющую их стенку пробирки

стенки пробирки к жидкости, Вт/К; а — коэффициент теплоотдачи от рабочего тела к внешней стенке пробирки, Вт/(м2 К); а2 — коэффициент теплоотдачи от внутренней стенки пробирки к жидкости, Вт/(м2 К); 1 — коэффициент теплопроводности материала пробирки Вт/(м К); ё1 — внешний диаметр пробирки, м; d2 — внутренний диаметр пробирки, м; 11 — длина пробирки, м; 12 — высота жидкости в пробирке, м.

Температура на внешней стенке пробирки Тст1 и температура на внутренней стенке пробирки Тст 2 определяются по следующим зависимостям:

Тст1 _ Т3 - б • Яа1 , Тст2 _ Тст1 - б • .

Задача о нахождении теплового потока б и температур стенок пробирки Тст1 и Тст 2 решается методом последовательных приближений. Вначале для расчета значений а1 и а2 приходится задаваться температурами поверхностей пробирки из условия Т3 > Тст1 > Тст2 > Тж . Затем рассчитывается тепловой поток через цилиндрическую стенку пробирки по зависимости (14) и уточняются значения температур поверхностей стенок пробирки по зависимостям (15).

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

а

б

в

Рис. 5. К расчету обтекания системы цилиндрических тел (ЦТ) [12]. а — коридорный пучок ЦТ; б — шахматный пучок ЦТ; в — цилиндр, обтекаемый под углом атаки ф. Воспроизводится по [12]

шаемым задачам — обтеканию электрического нагревателя, представляющего собой пучок тел цилиндрической формы, и обтеканию пробирки, представляющей собой цилиндрическое тело, — имеется в работах [6, 7, 9-13].

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

Коридорный пучок цилиндрических тел (рис. 5, а):

Ни = 1.2Яе033 • Рг033 • РГ

'Рг

при

V < 1.25, 82А < 1.25, 10 < Яе < 150; /а ' /а

Ни = 0,26 Яе065 • Рг033 •

(л V

• I Рг/

/рг^

при

1.24 < уа < 4, 1.24 <Луа < 4, 103 < Яе < 105;

Ни = 0.02Яе084 • Рг036 • Рг

'Рг

1.3 < V < 2.5, 1.3 < 52/ < 2.3.

105 < Яе < 106, 0.7 < Рг < 500.

при

Шахматный пучок цилиндрических тел (рис. 5, б):

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

Обычно при решении большого класса задач коэффициент теплоотдачи между нагревательным элементом и потоком газа определяется через число Нуссельта:

лт а • а Ни •Л Ни = —:— и, следовательно, а =-,

Л

а

Ни = 1.8Яе-Рг

0.33пг0.33^ I Рг

Рг

при

51/ < 1.25; V; < 1.25, 10 < Яе < 200;

Ни = 0.41Яе°'6-Рга33-( Рг

'Рг.

1.3 < V < 2.6, 0.6< 52/ < 4, аа

103 < Яе < 105, 0.7 < Рг < 500,

пРи

где а — характерный определяющий размер, м; Л — коэффициент теплопроводности газа Вт/(мК).

Расчетные зависимости или численные значения числа Нуссельта широко представлены в технической литературе по тепломассообмену. Обычно это экспериментальные данные, полученные для решения различных случаев обтекания твердых тел рабочим телом. Аналог, близкий ре-

где £ = Г1'

при у < 2 и

^ = 1.12 при У > 2;

Ни = 0.021Яеа84• Рг036• I Рг

1.2 < V, < 2.5, 0.9 <"У, < 1.5, 10 < Яе < 10ь.

Обтекание цилиндра (рис. 5, в):

Ып _ 0.5Яе°'5-Рга38- (

'Рг

при

5 < Яе < 10

3

038 • I Рг

Ып _ 0.026Яе0^ Рг

103 < Яе < 2 • 105;

Ып _ 0,023Яе0'8^Рг0'37'

3405 <Яе<2406.

Ргс

при

Рг/

\ 0,25

при

е _ 1 - 0.54 -сс^2 р .

5

Ып _ с •

где:

Ып _ (0.4•Яeа5 + 0.06Яе/э )• Рг04'

( \0 Пь

при 1 < Яе < 105, 0.7 < Рг < 300,

• /

где пь, П — коэффициент вязкости при средней и текущей температуре потока соответственно.

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

Коэффициент теплоотдачи между нагревательным элементом и потоком газа определяется по следующей зависимости:

Ып ж • 1

а_-

В предварительно турбулезированном потоке, который соответствует рассматриваемому устройству (за вентилятором — побудителем расхода), теплоотдача на 50-60 % выше расчетной [12]. Если угол ф между направлением потока и осью цилиндра (угол атаки) меньше 90°, значение коэффициента а следует умножить на поправочный коэффициент

В работе [13] для расчета теплоотдачи тела цилиндрической формы потоком воздуха приводятся следующие зависимости:

где Ыпж — значение числа Нуссельта; 1 — коэффициент теплопроводности газа при определяющей температуре, Вт/(мК); й — диаметр проволоки спирали, м.

Определяющим размером является диаметр проволоки спирали. Определяющей температурой является температура Т _ (Т2 + Тсп ))2 (Т2 — средняя массовая температура газа, К; Тсп — температура спирали, К).

Число Нуссельта определяется по критериальной зависимости для поперечного обтекания пучка тел цилиндрической формы, расположенных в шахматном порядке [10]:

Ып _ 0.49Яе05 при 103 < Яе < 2 40

Ы _ 0.0245Яе06 при Яе > 103.

В работе [13] для расчета теплоотдачи проволочного нагревателя используется обобщенная формула

Ыпж _ 0.41 • Яеж0'6 • Ргж033- е,; Яеж _

_(',1/ ^ )1/6 при ^ < 2;

е, _ 1.12 при ,2 > 2,

м> •й

с _ 0.5, п _ 0 при (СгРг)_ М0-4...М0-3, с _ 0.5, п _ при (СгРг)_ 1 • 10-3...5 • 102,

с _ 0.54, п _ У4 при (СгРг)_ 5 • 102...2• 107,

Ч

с _ 0.135, п _ Уъ при (СгРг)_ 2 • 10...1 • 1013.

Обработка экспериментальных данных по обтеканию тела цилиндрической формы потоком воздуха, полученных из различных источников, позволила автору работы [12] предложить следующую приближенную расчетную формулу:

где • — средняя скорость газа в самом узком сечении пучка, вычисленная при определяющей температуре, м/с; Ргж — значение числа Прандтля для газа при определяющей температуре; Ргс — значение числа Прандтля для газа при температуре спирали; V — кинематический коэффициент вязкости газа при определяющей температуре, м2/с;

Коэффициент теплоотдачи между потоком газа и пробиркой определяется по следующей зависимости:

Ып ж • 1

где Ыпж — значение числа Нуссельта; 1 — коэффициент теплопроводности газа при определяющий температуре, Вт/(м К); й1 — наружный диаметр пробирки, м.

Определяющим размером является наружный диаметр пробирки. Определяющей температурой

является температура Т = (Т3 + Тпр1 ))2 (Т3 — средняя массовая температура газа, К; Тпр1 — температура внешней стенки пробирки, К).

Число Нуссельта определяется по критериальной зависимости для поперечного обтекания одиночной цилиндрической трубы [11]:

Ни ж =(0.43 + 0.55^ЯеЖ5-РгТ));

Яе =

где ^ — средняя скорость газа, вычисленная при определяющей температуре, м/с; Ргж — значение числа Прандтля для газа при определяющей температуре; V — коэффициент кинематической вязкости газа при определяющей температуре, м2/с;

= 0.5 — коэффициент учитывающий угол между направлением течения потока и осью трубы.

Коэффициент теплоотдачи между пробиркой и жидкостью определяется по следующей зависимости:

Ни ж Л

а2 =

а2

где Ни ж — значение числа Нуссельта; Л — коэффициент теплопроводности газа при определяющей температуре, Вт/(м К); а2 — внутренний диаметр пробирки, м.

Определяющим размером является внутренний диаметр пробирки. Определяющей температурой

является температура Т = (Тж + Тпр2)) (Тж — средняя температура жидкости, К; Тпр 2 — температура внутренней стенки пробирки, К).

Число Нуссельта определяется по критериальной зависимости для естественной конвекции [9]:

Ни ж = 1.18 • (О-Рг )8;

в\Тт 2 - Тж И

Ог =

где Ог — значение числа Грасгофа; Рг — значение числа Прандтля для жидкости при определяющей температуре; V — коэффициент кинематической вязкости газа при определяющей температуре, м2/с; g = 9.81 м/с2 — значение ускорения свободного падения; в — температурный коэффициент объемного расширения, 1/К.

Тепловой баланс пробирки (без учета теплообмена излучением):

Qконв Qпn + Qж ,

где QKонв — количество теплоты, затрачиваемое на нагрев системы "пробирка—жидкость" через конвективный теплообмен, Дж; Qпр — количество

теплоты, расходуемое на нагрев пробирки, Дж; Qж — количество теплоты, затрачиваемое на нагрев жидкости, Дж.

Количество теплоты, затрачиваемое на нагрев системы "пробирка—жидкость" через конвективный теплообмен:

Qконв = Q•t.

Количество теплоты, расходуемое на нагрев пробирки:

Q =с -т •АТ ,

пр пр пр пр'

где спр — теплоемкость материала пробирки,

Дж/(кг-К); тпр — масса пробирки, кг; АТпр —

разность между конечной и начальной температурой пробирки, К.

Количество теплоты, затрачиваемое на нагрев жидкости:

Qж = сж • тж • АТж,

где сж — теплоемкость жидкости, Дж/(кг • К); тж — масса жидкости, кг; АТж — разность между конечной и начальной температурой жидкости, К.

Дифференцируем приведенные выше зависимости по времени:

= Q;

Он^=+. ах ах ах

^ = ^ аОж =

спр тпр л „ ; сж тж

ах

ах ах

аТж

ах

Для получения упрощенной записи дифференциального уравнения изменения температуры жидкости принимается допущение о равенстве скоростей изменения температуры жидкости и температуры пробирки аТж/ах = аТпр/ах, и дифференциальное уравнение изменения температуры жидкости принимает вид

аТж

Q

ах с т + с т

пр пр ж ж

V

2

V

МАТЕМАТИЧЕСКАЯ МОДЕЛЬ РАСЧЕТА ПЕРЕХОДНЫХ ПРОЦЕССОВ,

ПРОИСХОДЯЩИХ В ПНЕВМАТИЧЕСКОЙ СИСТЕМЕ ТЕРМОСТАБИЛИЗАЦИИ

Уравнения модели

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

• рабочее тело — идеальный газ, подчиняющийся термодинамическому уравнению состояния;

• процесс истечения — адиабатический (без теплообмена);

• поток газа — установившийся;

• течение газа происходит без трения (отсутствуют диссипативные силы).

Система обыкновенных дифференциальных уравнений, описывающих переходные процессы, протекающие в системе:

^ = k- RT1G1);

^ = V^ - RT1 (G - G)); dt pV { 1 dt w J

dTc„ = P-aS(Tcn -T2) dQ

dt

dP2 k f

dt V2 {

dT2 = T2

RT1G1 - RT2G2 +

dt dQ:

= a' S-(Tcn - T2);

dt

dt dT„

P2V2

V

dp2 dt

dQ2

- RT2 ( - G2 ) + (k - 1)|;

Q

dt c m + cm

пр пр ж ж

dp1 dt

dT3

kf

V3

RT2G2 - RT3G3 -

dt dQ;

dt

= Q;

T3

f

dt p3V3

V3

dp

dQs

- RT3 (G2 - G3)--

dt n 2 з; dt

ловных дросселях:

w = 0 ; О0 = О1 = в2 = О3 = 0.

Давление рабочего тела в полостях устройства равно давлению окружающей среды:

Р = Р2 = Р3 = Ра .

Напряжение на электрическом нагревателе отсутствует:

и = 0.

Граничные условия

Тепловое взаимодействие рабочего тела с электрическим нагревателем описывается граничными условиями третьего рода:

Q = a- S - (Tcn - T2).

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

В начальный момент времени температура спирали и температуры газа в объемах равны температуре окружающей среды:

Т = Т2 = тз = Гсп = Гпр = Тж = Та.

Скорость рабочего тела равна 0 и соответственно отсутствуют расходы рабочего тела в ус-

Тепловое взаимодействие рабочего тела с пробиркой описывается граничными условиями третьего рода:

Я = а Л

Конечные условия

Конечным условием расчета является достижение заданной температуры жидкости в пробирке:

Тж = Тк.

Решение

Решение полученной системы представляет собой задачу Коши для однородных дифференциальных уравнений. Такая задача формулируется в следующем виде:

^ = {(у,X), хе [Х0,1п], у(¿0) = у0,

Шх

у (х) е К", {(у, х) е К".

В некоторых системах однородных дифференциальных уравнений (ОДУ), являющихся математическими моделями физических систем, возникает так называемое явление жесткости. Оно возникает в случае, когда в физических системах, описываемых этой системой ОДУ, протекают процессы с существенно отличающимися постоянными времени. Система ОДУ, описывающих переходные процессы, протекающие в пневматической системе термостабилизации является жесткой, т. к. постоянные времени процессов изменения давления и температуры в газовых емкостях различаются на несколько порядков.

Для решения жестких систем ОДУ использу-

ются специальные жестко-устоичивые методы (методы Гира). Введение понятия жесткоИ устойчивости метода позволило сконструировать группу методов, в которых величина шага интегрирования выбирается так, чтобы быстро затухающие и не оказывающие существенного влияния компоненты решения аппроксимировались устойчиво, тогда как для компонент с большими постоянными времени гарантировалась точность аппроксимации. Наиболее распространенным классом линейных многошаговых методов для жестких задач являются "формулы дифференцирования назад" (более общее название — методы Гира). Вычислительная формула методов Гира основывается на интерполяционной формуле Лагранжа для правой части системы ОДУ:

Табл. 1. Физические свойства сухого воздуха при давлении ра = 1.013 • 105 Па

/ (^т+1, Ут+1) , X ,=0 агУт-1 +1 ;

п,

к

п

-

г - г

а -

^т ^т-I р-1,р^I ^т-г ^т-р

I -1,2,...,к;

--X к=

к

а 1.

1 -1 1

Выражение, приведенное выше, представляет собой нелинейное уравнение относительно ут+1. Для вычисления начального приближения Ут(0+)1 по значениям ут,ут-1,...,ут-к строится интерполяционный полином к-го порядка для функции у(г) и посредством экстраполирования вычисляется значение уто+1. Следует заметить, что по этим формулам можно вести интегрирование с переменным шагом. Вместе с тем при фиксированном шаге коэффициенты аг не изменяются, что существенно упрощает метод.

РЕЗУЛЬТАТЫ РАСЧЕТНО-ТЕОРЕТИЧЕСКИХ ИССЛЕДОВАНИЙ

В пневматической системе термостабилизации в качестве рабочего тела используется сухой воздух со следующими характеристиками: Я - 287 Дж/(кгхК), к -1.41. В табл. 1 приведены физические свойства сухого воздуха, используемые в расчете, в зависимости от температуры.

Параметры окружающей среды принимаются равными нормальным условиям (температура окружающей среды Та - 300 К, давление окружающей среды ра -1.013 • 105 Па).

Термодинамическая расчетная схема пневматической системы термостабилизации (рис. 1) представлена в виде набора трех газовых емкостей, разделенных между собой, (объем первой емкости У1- 5-10-5 м3, объем второй емкости

г, °С Р, кг/м3 Я • 102, Вт/ (мхК) V • 106, м2 /с Рг

0 1.293 2.44 13.28 0.707

10 1.247 2.51 14.16 0.705

20 1.205 2.59 15.06 0.703

30 1.165 2.67 16.00 0.701

40 1.128 2.76 16.96 0.699

50 1.093 2.83 17.95 0.698

60 1.060 2.90 18.97 0.696

70 1.029 2.96 20.02 0.694

80 1.000 3.05 21.09 0.692

90 0.972 3.13 22.10 0.690

100 0.946 3.21 23.13 0.688

120 0.898 3.34 25.45 0.686

140 0.854 3.49 27.80 0.684

160 0.815 3.64 30.09 0.682

объем третьей емкости У3 -

У2 - 3.5 • 10-6 м3,

- 9-10-6 м3).

Газовые емкости отделены от подводящих и отводящих магистралей условными дросселями (площадь условного дросселя, соединяющего первую и вторую газовые емкости, /12 - 2 40-4 м2; площадь условного дросселя, соединяющего вторую и третью газовые емкости, /23 - 2 40-4 м2; площадь условного дросселя, соединяющего третью газовую емкость с окружающей средой, /За - 6.5 • 10-4 м2).

В качестве побудителя расхода используется вентилятор с максимальным расходом, равным в0 - 0.25 • 10-3 м3/с .

В качестве электрического нагревателя используется проволока из нихрома (диаметр проволоки: ё - 0.3 • 10-3 м, общая длина проволоки: Ь - 0.5 м, теплоемкость нихрома: сншром - 448 Дж/(кг*К),

плотность нихрома: Рнихром - 8800 кг/м3). Мощ-

нихром

ность, подводимая к электрическому нагревателю, составляет Р - 8 Вт.

Пробирка изготовлена из полипропилена. Коэффициент теплопроводности полипропилена Япр - 0.2 Дж/(м*К). Теплоемкость полипропилена

спр - 1927 Дж/(кг^К). Внутренний диаметр пробирки ё2 - 3-10-3 м . Внешний диаметр пробирки ё1 - 3.5 • 10-3 м . Длина пробирки I -1 • 10-2 м .

Табл. 2. Физические свойства воды на линии насыщения

X, А-102, V -106, Рг в104,

°С Вт/ (мхК) м2 / с 1/К

0 55.1 1.789 13.67 -0.63

10 57.4 1.306 9.52 +0.70

20 59.9 1.006 7.02 1.82

30 61.8 0.805 5.42 3.21

40 63.5 0.659 4.31 3.87

50 64.8 0.556 3.54 4.49

60 65.9 0.478 2.98 5.11

70 66.8 0.415 2.55 5.70

80 67.4 0.365 2.21 6.32

90 68.0 0.326 1.95 6.95

100 68.3 0.295 1.75 7.52

с

8 бо

10 20

40 50

Время, с

V 1 еп

Тг г

Тж

Г

Рис. 6. Зависимость температуры элементов расчетной схемы (рис. 1) от времени. Теп — температура спирали; Т1 — температура газа в первом объеме; Т2 — температура газа во втором объеме; Т3 — температура газа в третьем объеме; Тж — температура жидкости в пробирке

с

,ме

ч е

Р1

Р2 г~

Р3

Время, с

Рис. 7. Зависимость избыточного давления газа в объемах от времени (расчет). р1 — избыточное давление газа в первом объеме; р2 — избыточное давление газа во втором объеме; р3 — избыточное давление газа в третьем объеме

110

100

Н 50

Тс т1 \Тж

Тет2

/

Время, с

Рис. 8. Зависимость температуры элементов 3-й емкости от времени (расчет). Т3 — температура газа в третьем объеме, Тж — температура жидкости в пробирке; Тет1 — температура наружной стенки пробирки; Тет2 — температура внутренней стенки пробирки

Масса пробирки т = 3 -10-5 кг . В пробирке находится жидкий биологический препарат, физические параметры которого сходны с физическими параметрами воды. Теплоемкость воды с = 4200 Дж/(кгхК). Масса воды т = 5 • 10-5 кг .

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

На основе предложенной математической модели проведены расчетно-теоретические исследования переходных процессов, происходящих в пневматической системе термостабилизации. Результаты представлены на рис. 6, 7, 8.

Устройство сопряжения

Персонал ьный компьютер

Рис. 9. Конструктивная схема экспериментального стенда

ЭКСПЕРИМЕНТАЛЬНАЯ ПРОВЕРКА ТЕОРЕТИЧЕСКИХ ДАННЫХ, ПОЛУЧЕННЫХ С ПОМОЩЬЮ МАТЕМАТИЧЕСКОЙ МОДЕЛИ

Нагреватель установки (рис. 1, 9), обеспечивающий нагрев рабочего тела, и вентилятор, обеспечивающий прокачку рабочего тела, потребляют энергию от внешних источников питания.

Регулирование температуры внутри рабочей камеры происходит двояко:

• изменением напряжения на нагревателе, что приводит к изменению его температуры, а следовательно, и температуры рабочего тела;

• управлением расхода вентилятора путем изменения напряжения, подаваемого на вентилятор, а следовательно, и скорости вращения вентилятора.

Все изменения температуры фиксируются термодатчиками, расположенными в рабочей зоне и в пробирке, и передаются для обработки и визуализации на ЭВМ (рис. 9). Сигнал с термодатчика поступает в устройство сопряжения, представляющее собой аналого-цифровой преобразователь, где он из аналогового сигнала преобразуется в цифровой. Цифровой сигнал кодируется и по последовательному интерфейсу передается в персональный компьютер, который декодирует полученную информацию и выводит ее на монитор в виде графика изменения температуры.

Датчик 1 и датчик 3 фиксируют температуру газа в третьем объеме. Датчик 2 фиксирует температуру жидкости в пробирке. В эксперименте использовались датчики фирмы HONEYWELL Hel 700-Txx со следующими характеристиками:

• Диапазон рабочих температур: -200...+540 °С;

• сопротивление при 25 °С: 100 Ом;

• время ответа: 0.01 с;

с

s

« 40

V сп

Тз Т2 Г

Т3

Тж

Тж

-Ti г~

Время, с

Рис. 10. Зависимость температуры от времени (эксперимент).

Т1 — температура газа в первом объеме; Т2 — температура газа во втором объеме; Т3 — температура газа в третьем объеме; Тсп — температура спирали; Тж — температура жидкости в пробирке; Т3(Э) — температура газа, измеренная датчиком 1; Тж(э) — температура жидкости в пробирке, измеренная датчиком 2

• точность температуры: ±0.2 °С;

• рекомендуемый измерительный ток: 1 мА;

• чувствительность: 0.3 Ом/град.

Экспериментальные результаты представлены на рис. 10.

Фирма Roch [15] в качестве емкостей для жидкости использовала вместо стандартных полипропиленовых пробирок тонкостенные капилляры

Время, с

Рис. 11. Динамика температуры жидкости в тонкостенном капилляре в зависимости от температуры наружного воздуха. Данные фирмы Roch

диаметром 1.5 мм и максимальным объемом 20 мкл. Капилляры имеют большее отношение площади поверхности капилляра к объему капилляра по сравнению с полипропиленовыми пробирками, что позволяет достигать больших скоростей нагрева и охлаждения. Об этом свидетельствуют данные, приводимые фирмой (рис. 11).

ВЫВОДЫ

Как видно из рис. 6, температура в 1-й газовой емкости практически не изменяется и равна 20 °С, т. е. равна температуре окружающей среды. Температуры во второй и третьей газовых емкостях практически равны между собой в каждый момент времени.

Как видно из рис. 7, максимальное избыточное давление соответствует первой емкости и составляет чуть более 3 Па. Следовательно, можно пренебречь распределением давления по газовым емкостям и считать его постоянным и равным давлению окружающей среды:

p1 = p2 = p3 = ра = 1.013 • 105 Па = const.

Скорость нагрева жидкости в полипропиленовой пробирке меньше скорости нагрева воздуха, которым омывается пробирка. Этого не наблюдается при использовании тонкостенного капилляра с большим отношением площади поверхности капилляра к объему капилляра. Там температура жидкости практически не отстает от температуры газа (рис. 11). Для достижения больших скоростей нагрева и охлаждения необходимо использовать тонкостенные капилляры.

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

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

Работа выполняется при частичной поддержке госконтракта № 02.467.11.3005. с Федеральным агентством по науке и инновациям.

СПИСОК ЛИТЕРАТУРЫ

1. Fitch J.P. An engineering introduction to biotechnology // Tutorial texts in optical engineering. Includes bibliographical references and index. Published by SPIE, 2002. V. TT55. 127 p.

2. Чернышев А.В., Друца В.Л. Проблемы создания оборудования для медицинской ПЦР диагностики // Биомедицинские технологии и радиоэлектроника. 2004. № 12. С. 18.

3. Чернышев А.В. Основы теории расчета электропневмомеханического оборудования для анализа ДНК // Научное приборостроение. 2002. Т. 12, № 1. С. 53-65.

4. Patent 5038852 (USA). Apparatus and method for performing automated amplification of nucleic ACID sequences and assays using heating and cooling steps /Larry J. Johnson, San Jose, Joseph T.Widunas /Assignee: Cetus Corporation, USA. Aug. 13, 1991.

5. Patent 5455175 (USA). Rapid thermal cycling device /Carl Wittwer, David R.Hillyrd, Kirk M.Ririe /Assignee: The Dow Chemical Company, USA. Feb. 16, 1993.

6. Бэр Г.Д. Техническая термодинамика. Теоретические основы и технические приложения. М.: Мир, 1977. 518 с.

7. Техническая термодинамика: Учебник для вузов /Под ред. В.И. Крутова. М.: Высшая школа, 1981. 439 с.

8. Мамонтов М.А. Тепломеханика тела переменной массы — основа теории пневмогазоприво-дов // Пневматические приводы и системы управления. М.: Наука, 1971. С. 8-18.

9. Теплотехнический справочник. В 2-х томах / Под общ. ред. В.Н. Юренева, П.Д. Лебедева. М.: Энергия, 1975, 1976.

10. Исаченко В.П., Осипова В.А., Сукомел А.С. Теплопередача. М.: Энергия, 1981. 417 с.

11. Теплотехника / Под ред. А.П. Баскакова. М.: Энергоиздат, 1982. 261 с.

12. Тепло и массообмен. Теплотехнический эксперимент: Справочник / Под общ. редакцией В. А. Григорьева и В.М. Зорина. М.: Энерго-атомиздат, 1988. 512 с.

13. Куличенко В.Р. Справочник по теплообменным расчетам. Киев: Техника, 1990.

14. (http://www.corbettresearch.com).

15. (http://www.biochem. roche.com/lightcycler).

МГТУ им. Н.Э. Баумана

ГНЦ РФ ИМБП РАН (Смирнов И.А.) Материал поступил в редакцию 29.11.2005.

DEVELOPMENT OF A MATHEMATICAL MODEL FOR A PNEUMATIC THERMAL STABILIZATION SYSTEM

A. V. Chernyshev, A. A. Krutikov, K. E. Demikhov, A. V. Polynkov, O. V. Belova, I. A. Smirnov*

N.E. Bauman Moscow State Technical University IMBP State Science Centre RAS

The paper presents a calculational procedure used to develop a pneumoelectric device for polymerase chain reaction (PCR) splitting of DNA molecules — a pneumatic amplificator. This system is characterized by a high rate of change of the working medium (air) temperature combined with a high accuracy of change of biological sample temperature in plastic microvials. The aim of the present research is modeling thermal processes to find optimal conditions for PCR DNA splitting. A mathematical model used in calculations based on the lumped parameter method is introduced. An experimental stand is described. The results of theoretic prediction studies are given and compared with experimental results. It is shown that the proposed technique allows one to evaluate the extent of transients occurring in the pneumatic amplificator. It is inferred that a mathematical model based on distributed parameter methods of system description are needed.